Development of robust and efficient solution strategies for coupled problems

Det er mange modeller i moderne vitenskap hvor sammenkoblingen mellom forskjellige fysiske prosesser er svært viktig. Disse finner man for eksempel i forbindelse med lagring av karbondioksid i undervannsreservoarer, flyt i kroppsvev, kreftsvulstvekst og geotermisk energiutvinning. Denne avhandlingen har to fokusområder som er knyttet til sammenkoblede modeller. Det første er å utvikle pålitelige og effektive tilnærmingsmetoder, og det andre er utviklingen av en ny modell som tar for seg flyt i et porøst medium som består av to forskjellige materialer. For tilnærmingsmetodene har det vært et spesielt fokus på splittemetoder. Dette er metoder hvor hver av de sammenkoblede modellene håndteres separat, og så itererer man mellom dem. Dette gjøres i hovedsak fordi man kan utnytte tilgjengelig teori og programvare for å løse hver undermodell svært effektivt. Ulempen er at man kan ende opp med løsningsalgoritmer for den sammenkoblede modellen som er trege, eller ikke kommer frem til noen løsning i det hele tatt. I denne avhandlingen har tre forskjellige metoder for å forbedre splittemetoder blitt utviklet for tre forskjellige sammenkoblede modeller. Den første modellen beskriver flyt gjennom deformerbart porøst medium og er kjent som Biot ligningene. For å anvende en splittemetode på denne modellen har et stabiliseringsledd blitt tilført. Dette sikrer at metoden konvergerer (kommer frem til en løsning), men dersom man ikke skalerer stabiliseringsleddet riktig kan det ta veldig lang tid. Derfor har et intervall hvor den optimale skaleringen av stabiliseringsleddet befinner seg blitt identifisert, og utfra dette presenteres det en måte å praktisk velge den riktige skaleringen på. Den andre modellen er en fasefeltmodell for sprekkpropagering. Denne modellen løses vanligvis med en splittemetode som er veldig treg, men konvergent. For å forbedre dette har en ny akselerasjonsmetode har blitt utviklet. Denne anvendes som et postprosesseringssteg til den klassiske splittemetoden, og utnytter både overrelaksering og Anderson akselerasjon. Disse to forskjellige akselerasjonsmetodene har kompatible styrker i at overrelaksering akselererer når man er langt fra løsningen (som er tilfellet når sprekken propagerer), og Anderson akselerasjon fungerer bra når man er nærme løsningen. For å veksle mellom de to metodene har et kriterium basert på residualfeilen blitt brukt. Resultatet er en pålitelig akselerasjonsmetode som alltid akselererer og ofte er svært effektiv. Det siste modellen kalles Cahn-Larché ligningene og er også en fasefeltmodell, men denne beskriver elastisitet i et medium bestående av to elastiske materialer som kan bevege seg basert på overflatespenningen mellom dem. Dette problemet er spesielt utfordrende å løse da det verken er lineært eller konvekst. For å håndtere dette har en ny måte å behandle tidsavhengigheten til det underliggende koblede problemet på blitt utviklet. Dette leder til et diskret system som er ekvivalent med et konvekst minimeringsproblem, som derfor er velegnet til å løses med de fleste numeriske optimeringsmetoder, også splittemetoder. Den nye modellen som har blitt utviklet er en utvidelse av Cahn-Larché ligningene og har fått navnet Cahn-Hilliard-Biot. Dette er fordi ligningene utgjør en fasefelt modell som beskriver flyt i et deformerbart porøst medium med to poroelastiske materialer. Disse kan forflytte seg basert på overflatespenning, elastisk spenning, og poretrykk, og det er tenkt at modellen kan anvendes i forbindelse med kreftsvulstmodellering.There are many applications where the study of coupled physical processes is of great importance. These range from the life sciences with flow in deformable human tissue to structural engineering with fracture propagation in elastic solids. In this doctoral dissertation, there is a twofold focus on coupled problems. Firstly, robust and efficient solution strategies, with a focus on iterative decoupling methods, have been applied to several coupled systems of equations. Secondly, a new thermodynamically consistent coupled system of equations is proposed. Solution strategies are developed for three different coupled problems; the quasi-static linearized Biot equations that couples flow through porous materials and elastic deformation of the solid medium, variational phase-field models for brittle fracture that couple a phase-field equation for fracture evolution with linearized elasticity, and the Cahn-Larché equations that model elastic effects in a two-phase elastic material and couples an extended Cahn-Hilliard phase-field equation and linearized elasticity. Finally, the new system of equations that is proposed models flow through a two-phase deformable porous material where the solid phase evolution is governed by interfacial forces as well as effects from both the fluid and elastic properties of the material. In the work that concerns the quasi-static linearized Biot equations, the focus is on the fixed-stress splitting scheme, which is a popular method for sequentially solving the flow and elasticity subsystems of the full model. Using such a method is beneficial as it allows for the use of readily available solvers for the subproblems; however, a stabilizing term is required for the scheme to converge. It is well known that the convergence properties of the method strongly depend on how this term is chosen, and here, the optimal choice of it is addressed both theoretically and practically. An interval where the optimal stabilization parameter lies is provided, depending on the material parameters. In addition, two different ways of optimizing the parameter are proposed. The first is a brute-force method that relies on the mesh independence of the scheme's optimal stabilization parameter, and the second is valid for low-permeable media and utilizes an equivalence between the fixed-stress splitting scheme and the modified Richardson iteration. Regarding the variational phase-field model for brittle fracture propagation, the focus is on improving the convergence properties of the most commonly used solution strategy with an acceleration method. This solution strategy relies on a staggered scheme that alternates between solving the elasticity and phase-field subproblems in an iterative way. This is known to be a robust method compared to the monolithic Newton method. However, the staggered scheme often requires many iterations to converge to satisfactory precision. The contribution of this work is to accelerate the solver through a new acceleration method that combines Anderson acceleration and over-relaxation, dynamically switching back and forth between them depending on a criterion that takes the residual evolution into account. The acceleration scheme takes advantage of the strengths of both Anderson acceleration and over-relaxation, and the fact that they are complementary when applied to this problem, resulting in a significant speed-up of the convergence. Moreover, the method is applied as a post-processing technique to the increments of the solver, and can thus be implemented with minor modifications to readily available software. The final contribution toward solution strategies for coupled problems focuses on the Cahn-Larché equations. This is a model for linearized elasticity in a medium with two elastic phases that evolve with respect to interfacial forces and elastic effects. The system couples linearized elasticity and an extended Cahn-Hilliard phase-field equation. There are several challenging features with regards to solution strategies for this system including nonlinear coupling terms, and the fourth-order term that comes from the Cahn-Hilliard subsystem. Moreover, the system is nonlinear and non-convex with respect to both the phase-field and the displacement. In this work, a new semi-implicit time discretization that extends the standard convex-concave splitting method applied to the double-well potential from the Cahn-Hilliard subsystem is proposed. The extension includes special treatment for the elastic energy, and it is shown that the resulting discrete system is equivalent to a convex minimization problem. Furthermore, an alternating minimization solver is proposed for the fully discrete system, together with a convergence proof that includes convergence rates. Through numerical experiments, it becomes evident that the newly proposed discretization method leads to a system that is far better conditioned for linearization methods than standard time discretizations. Finally, a new model for flow through a two-phase deformable porous material is proposed. The two poroelastic phases have distinct material properties, and their interface evolves according to a generalized Ginzburg–Landau energy functional. As a result, a model that extends the Cahn-Larché equations to poroelasticity is proposed, and essential coupling terms for several applications are highlighted. These include solid tumor growth, biogrout, and wood growth. Moreover, the coupled set of equations is shown to be a generalized gradient flow. This implies that the system is thermodynamically consistent and makes a toolbox of analysis and solvers available for further study of the model.Doktorgradsavhandlin

On the optimization of iterative schemes for solving non-linear and/or coupled PDEs

In this thesis we study the optimization of iterative schemes as both linearization methods, and as splitting methods for solving non-linear and coupled partial differential equations (PDEs). We consider two equations that are describing processes in porous media; Richards’ equation, a possibly degenerate, non-linear and elliptic/parabolic equation that models flow of water in saturated/unsaturated porous media, and Biot’s equations, a coupled system of equations that models flow in deformable porous media. For Richards’ equation we compare the numerical properties of several linearization schemes, including the Newton-Raphson method, the modified Picard method and the L-scheme. Additionally, we prove convergence of the linearly and globally convergent L-scheme and discuss theoretically and practically how to choose its stabilization parameter optimally in the sense that convergence is obtained in the least amount of iterations. The second aim of the thesis is to effectively solve the quasi-static, linear Biot model. We consider the fixed-stress splitting scheme, which is a popular method for iteratively solving Biot’s equations. It is well-known that the convergence of the method is strongly dependent on the applied stabilization parameter. We propose a new approach to optimize this parameter, and show theoretically that it does not only depend on the mechanical properties and the coupling coefficient, but also on the fluid’s flow properties. The type of analysis presented in this thesis is not restricted to a particular spatial discretization, but we require it to be inf-sup stable. The convergence proof also applies to low-compressible or incompressible fluids, and low-permeable porous media. We perform illustrative numerical examples, including a well-known benchmark problem, Mandel’s problem. The results largely agree with the theoretical findings. Furthermore, we show numerically that for conditionally inf-sup stable discretizations, the performance of the fixed-stress splitting scheme behaves in a manner which contradicts the theory provided for inf-sup stable discretizations.Masteroppgave i anvendt og beregningsorientert matematikkMAMN-MABMAB39

Well-posedness analysis of the Cahn-Hilliard-Biot model

We investigate the well-posedness of the recently proposed Cahn-Hilliard-Biot model. The model is a three-way coupled PDE of elliptic-parabolic nature, with several nonlinearities and the fourth order term known to the Cahn-Hilliard system. We show existence of weak solutions to the variational form of the equations and uniqueness under certain conditions of the material parameters and secondary consolidation, adding regularizing effects. Existence is shown by discretizing in space and applying ODE-theory (the Peano-Cauchy theorem) to prove existence of the discrete system, followed by compactness arguments to retain solutions of the continuous system. In addition, the continuous dependence of solutions on the data is established, in particular implying uniqueness. Both results build strongly on the inherent gradient flow structure of the model

The Fixed-Stress splitting scheme for Biot's equations as a modified Richardson iteration: Implications for optimal convergence

The fixed-stress splitting scheme is a popular method for iteratively solving the Biot equations. The method successively solves the flow and mechanics subproblems while adding a stabilizing term to the flow equation, which includes a parameter that can be chosen freely. However, the convergence properties of the scheme depend significantly on this parameter and choosing it carelessly might lead to a very slow, or even diverging, method. In this paper, we present a way to exploit the matrix structure arising from discretizing the equations in the regime of impermeable porous media in order to obtain a priori knowledge of the optimal choice of this tuning/stabilization parameter.acceptedVersio

A Cahn-Hilliard-Biot system and its generalized gradient flow structure

In this work, we propose a new model for flow through deformable porous media, where the solid material has two phases with distinct material properties. The two phases of the porous material evolve according to a generalized Ginzburg–Landau energy functional, with additional impact from both elastic and fluid effects, and the coupling between flow and deformation is governed by Biot’s theory. This results in a three-way coupled system which can be seen as an extension of the Cahn–Larché equations with the inclusion of a fluid flowing through the medium. The model covers essential coupling terms for several relevant applications, including solid tumor growth, biogrout, and wood growth simulation. Moreover, we show that this coupled set of equations follow a generalized gradient flow framework. This opens a toolbox of analysis and solvers which can be used for further study of the model. Additionally, we provide a numerical example showing the impact of the flow on the solid phase evolution in comparison to the Cahn–Larché system.publishedVersio

DarSIA: An Open-Source Python Toolbox for Two-Scale Image Processing of Dynamics in Porous Media

Understanding porous media flow is inherently a multi-scale challenge, where at the core lies the aggregation of pore-level processes to a continuum, or Darcy-scale, description. This challenge is directly mirrored in image processing, where pore-scale grains and interfaces may be clearly visible in the image, yet continuous Darcy-scale parameters may be what are desirable to quantify. Classical image processing is poorly adapted to this setting, as most techniques do not explicitly utilize the fact that the image contains explicit physical processes. Here, we extend classical image processing concepts to what we define as “physical images” of porous materials and processes within them. This is realized through the development of a new open-source image analysis toolbox specifically adapted to time-series of images of porous materials.publishedVersio

Splitting schemes for coupled differential equations: Block Schur-based approaches & Partial Jacobi approximation

Coupled multi-physics problems are encountered in countless applications and pose significant numerical challenges. In a broad sense, one can categorise the numerical solution strategies for coupled problems into two classes: monolithic approaches and sequential (also known as split, decoupled, partitioned or segregated) approaches. The monolithic approaches treat the entire problem as one, whereas the sequential approaches are iterative decoupling techniques where the different sub-problems are treated separately. Although the monolithic approaches often offer the most robust solution strategies, they tend to require ad-hoc preconditioners and numerical implementations. Sequential methods, on the other hand, offer the possibility to add and remove equations from the model flexibly and rely on existing black-box solvers for each specific equation. Furthermore, when problems are non-linear, inner iterations need to be performed even in monolithic solvers, making the sequential approaches an even more viable alternative. The cost of running inner iterations to recover the multi-physics coupling could, however, easily become prohibitive. Moreover, the sequential approaches might not converge at all. In this work, we present a general formulation of splitting schemes for continuous operators with arbitrary implicit/explicit splitting, like in standard iterative methods for linear systems. By introducing a generic relaxation operator, we find the conditions for the convergence of the iterative schemes. We show how the relaxation operator can be thought of as a preconditioner and constructed based on an approximate Schur complement. We propose a Schur-based Partial Jacobi relaxation operator to stabilise the coupling and show its effectiveness. Although we mainly focus on scalar-scalar linear problems, most results are easily extended to non-linear and higher-dimensional problems. The schemes presented are not explicitly dependent on any particular discretisation methodologies. Numerical tests (1D and 2D) for two PDE systems, namely the Dual-Porosity model and a Quad-Laplacian operator, are carried out to investigate the practical implications of the theoretical results

An accelerated staggered scheme for variational phase-field models of brittle fracture

There is currently an increasing interest in developing efficient solvers for variational phase-field models of brittle fracture. The governing equations for this problem originate from a constrained minimization of a non-convex energy functional, and the most commonly used solver is a staggered solution scheme. This is known to be robust compared to the monolithic Newton method, however, the staggered scheme often requires many iterations to converge when cracks are evolving. The focus of our work is to accelerate the solver through a scheme that sequentially applies Anderson acceleration and over-relaxation, switching back and forth depending on the residual evolution, and thereby ensuring a decreasing tendency. The resulting scheme takes advantage of the complementary strengths of Anderson acceleration and over-relaxation to make a robust and accelerating method for this problem. The new method is applied as a post-processing technique to the increments of the solver. Hence, the implementation merely requires minor modifications to already available software. Moreover, the cost of the acceleration scheme is negligible. The robustness and efficiency of the method are demonstrated through numerical examples.publishedVersio