Mixed and Mixed Hybrid Finite Element Methods: Theory, Implementation and Applications

(1)

F.A. Radu

Max-Planck Institute for Mathematics in the Sciences Leipzig, Germany

mailto:[email protected]

(2)

Objective

Offer the skills to be able to set, analyse and implement a mixed or mixed hybrid finite element method

for second order partial differential equations.

(3)

Lecture organization

• Two hours weekly.

• Desirable: Some of you to take on from me some tasks in the frame of the lecture, study and present them in the form of a paper.

• Exercises: I will always let you some exercises, it would be very good if you will solve them and come and discuss the solution (if necessary) with me.

(4)

Knowledge requirements

• Numerical mathematics,

→ A. Quarteroni et al., Numerical Mathematics, 2000.

• Numerical methods for partial differential equations,

→ P. Knabner and L. Angermann, Numerical methods for elliptic and parabolic partial differential equations, 2003.

• Functional Analysis,

→ B. Rynne and M. Youngson, Linear functional analysis, 2000.

(5)

What are Mixed Finite Element Methods ?

• They are FE methods founded on a variational principle

expressing an equilibrium (saddle point) condition and not a minimization principle (conforming FE).

• MFEM approximate both a scalar variable (e.g. pressure) and a vector variable (its gradient, the flux) simultaneously;

from here comes also the name mixed.

• They are nonconforming methods in the sense that the

primal variable is not necessary continuous, as in the case

(6)

History

It was first introduced by engineers in the 1960’s to solve problems in solid continua:

• Fraeijs de Veubeke, 1965

• Hellan, 1967

• Hermann, 1967

(7)

Why MFEM

Local conservation.

→ Very important when the equation to be discretized

corresponds to a conservation law (usually mass). Example:

∂_tc + ∇ · q = f mass conservation, (1) q = −D∇c Fick’s law of diffusion. (2)

→ MFEM: The equation (1) holds not only global, but also locally on each simplex.

→ The same property we have by the finite volume (FV), not on

(8)

Why MFEM

An intrinsic and accurate approximation of the flux.

→ Transport equation (convection-diffusion equation):

∂_t(Θc) − ∇ · (D∇c) + ∇ · (qc) = R

The flux q is much important than the pressure, a good approximation of it is of great interest.

→ MFEM: the normal component of the discrete flux is

continuous over edges, because the flux is in H(div; Ω). By FE, the flux is only in L²(Ω)^d, so no continuity.

(9)

Why MFEM

→ Error estimates (Poisson equation, enough regularity for the domain, homogeneuos boundary conditions, regular

triangulation, smooth data;):

FE : ku − u_hk₀ + k∇u − ∇u_hk₀ ≤ Ch MFEM(RT₀) : ku − u_hk₀ + kq − q_hk_div ≤ Ch MFEM(BDM1) : ku − u_hk₀ + kq − q_hk_div ≤ Ch

kq − q_hk₀ ≤ Ch²

(10)

Richy (left pressure, right flux)

(11)

Feflow (left pressure, right flux)

Simulations performed by Ch. Kohlepp, Univ. Erlangen-Nuernberg.

(12)

Why MFEM

Applicable for equations with jumps in the coefficients and irregular geometries.

• heterogeneous soil or materials.

• anisotrop soil or materials.

• sophisticated domains.

(13)

Reactive transport in porous medium

(14)

To describe

• Water flow, including the unsaturated zone near the subsurface.

• Advective and dispersive transport of multiple contaminants.

• Non-equilibrium and equilibrium sorption.

• Biodegradation.

(15)

An appropiate model for the water flow in porous media is Richards’ equation (here in the pressure formulation):

∂_tΘ(ψ) − ∇ · K(ψ)∇(ψ + z) = 0 in J × Ω

Water content: θ(ψ) ∈ [0,1]

Pressure head: ψ

Unsaturated hydraulic conductivity: K(ψ) Height against the gravitational direction: z Time interval: J = (0,T)

Domain: Ω in IR^d (d = 1,2 or 3)

(16)

The equation results from

• mass conservation

∂_tΘ(ψ) + ∇ · q = 0

• Darcy’s law

q = −K(ψ)∇(ψ + z)

(17)

The equation results from

• mass conservation

∂_tΘ(ψ) + ∇ · q = 0

• Darcy’s law

q = −K(ψ)∇(ψ + z)

Nonlinearities Θ, K: strictly monotone increasing for ψ ≤ 0, constant for ψ ≥ 0 (saturated region)

=⇒ elliptic - parabolic equation.

(18)

• the soil-water retention Θ(ψ),

• the unsaturated hydraulic conductivity K(Θ),

Gardner Θ_exp(ψ) = Θ_r + (Θ_s − Θ_r)e^αψ Kexp(ψ) = Kse^αψ

Haverkamp ΘHav(ψ) = Θr + ^(Θ_1+(αψ)^s^−Θ^rn⁾

K_Hav(ψ) = _1+(βψ)^K^s _p

van Genuchten- Θ_vG(ψ) = Θ_r + (Θ_s − Θ_r)Φ(ψ) Mualem Φ(ψ) = _(1+(αψ)¹ n)^m , m = 1 − _n¹

K_vG(ψ) = K_sp

Φ(ψ)(1 − (1 − Φ(ψ)^m¹ )^m)²

(19)

General model with multicomponent organic transport and biodegradation

N mobile species, M immobile species

∂t(Θci) + ρb∂tsi − ∇ · (Di∇c_i − qci) = −R_i ,

∂tsi = ki(φ(ci) − si) or si = φ(ci), i ∈ 1, ..., N

∂tci + k_dici =

“

1 − γi c_i c_imax

”

Ri, i ∈ N + 1, ..., N + M.

c_i concentration of the species, s_i concentration of the absorbed species , D_i diffusion coefficient, ρ_b bulk density, R_i degradation rate, φ sorption isotherm, k_di death rate, cimax a maximal realistic concentration, γi ∈ {0,1}.

(20)

Boundary Conditions

c_i = g_Di on J × Γ_Di,

−D_i∇c_i · n = g_{N i} on J × Γ_{N i},

−D_i∇c_i · n + c_iq · n

| {z }

q_i·n

= g_{F i} on J × Γ_{F i},

Remark. Γ_Di,Γ_{N i},Γ_{F i} are species depending.

(21)

Benzene Biodegradation

???

...

Ω₁

Ω₂ Γ₁

(0,0) (2,0)

(0,3) (2,3)

•Water Flow : F3 days rain, 4 days dry

Fvan Genuchten-Mualem Model

•Biodegradation : F2 mobile species, 1 biomass

Fno sorption

FMonod Model

(22)

∂_t(Θc_D) − ∇ · (D_D∇c_D − qc_D) = −R ,

∂_t(Θc_A) − ∇ · (D_A∇c_A − qc_A) = −α_A/DR ,

∂_tc_X + k_dc_X = Y

1 − γ_X _c ^c^X

Xmax

R

with donator/contaminant cD, acceptor cA, biomass cX.

Reactive term:

R = µ_max c_X

c_D

K_M_D + c_D







c_A

K_M_A + c_A + c²_A K_I_A





 .

(23)

Benzene concentration at T = 30, 60, 90, 120, 150, 160 days

(24)

Oxygen concentration at T = 30, 60, 90, 120, 150, 160 days

(25)

Biomass concentration at T = 30, 60, 90, 120, 150, 160 days

(26)

•Water Flow : F stationary flow

F variable permeability

•Biodegradation : F 2 mobile species, 1 biomass

F without sorption

F Monod Model

(27)

Xylene degradation: variable permeability

Domain Pressure Profile Flux

(28)

Concentration [mg/l] profiles at T = 1 [year]

(without additional delivery of contaminant)

Xylene Oxygen Biomass

(29)

Concentration [mg/l] profiles at T = 3 [years]

(30)

Concentration [mg/l] profiles at T = 5 [years]

(31)

Modelling drug release from collagen matrices

• Department of Pharmacy, Pharmaceutical Technology and Biopharmacy, University of Munich, Germany

(32)

• Matrix systems of insoluble collagen are a promising and advantageous drug delivery system for prolonged protein release over several days.

Applications

• Collagen implants have been evaluated for tumor treatment, bone, and nerve regeneration as well as therapy of infections.

Goal

• Optimizing the controlled release from degradable collagen matrices.

(33)

To describe

PHASE I:

• Swelling (very short, 20 - 30 min),

→ free boundary problem, front-tracking method:

(34)

PHASE II:

• diffusion of the enzyme in the matrix,

• adsorption of the enzyme from the fluid to the collagen fibers,

• enzymatic degradation of the polymer,

• enzyme activity (death),

• drug release.

(35)

The general behaviour of an enzymatically catalyzed degradation process can be summarized by the equations:

E + S →

^k¹

ES ES →

^k²

P + E

k₁ = the constant rate of formation of the enzyme-substrate complex.

k₂ = the catalysis rate.

(36)

∂_tC_E − ∇ · (D_E(C_K)∇C_E) + k_aktC_E = − k₁(C_E)^α

1 + _{M axSorp}^k¹ (C_E)^αC_K + k₂C_ES^γ ,

∂_tC_ES = k₁(C_E)^α

1 + _{M axSorp}^k¹ (C_E)^αC_K − k₂C_ES^γ ,

∂_tC_K = −k₁(C_E)^αC_K ,

∂_tC_P − ∇ · (D_P∇C_P) = k₂C_ES^γ .

(37)

The model for the enzymatic degradation is validated by comparison of the experimentally obtained data with the numerically simulated data for the collagen.

(38)

The model for the enzymatic degradation is validated by comparison of the experimentally obtained data with the numerically simulated data for the collagen.

Next step: Drug Release

The release of the active agent is governed by a diffusion

equation with a source term due to liberation of the immobilized active agent by matrix degradation:

∂_tC_A − ∇ · (D_A(C_K)∇C_A) = −∂_t(C_A_i) ,

where C_A, C_A denote the concentrations of free and immobilized drug.

(39)

We assume

C_A_i = f(C_A, C_K)

• the simplest approach: C_A_i = σC_K ( Tzafriri 2000)

• C_A_i = σC_K²

• C_A_i = σ√ C_K

•• σ can be obtained experimentally (if one considers release from a collagen matrices, without enzyme (i.e. no degradation), the

concentration of the collagen remaining in the matrix gives us C_Ai⁰ , and therefore also σ.)

•• the form of the function f is determinated by fitting a set of release data and then validated for an other one.

(40)

enzyme

transport sorption polymer degradation

coupled solver (Newton)

?

drug release

solver

• the algorithm was implemented in ug.

P. Bastian et al., UG-a flexible toolbox for solving partial differential equation, Comput. Visualization in Science 1, pp. 27-40, 1997.

(41)

Enzymatic degradation of collagen (left) and drug release from an insoluble collagen matrix: comparison of numerical and experimental results (points).

(42)

Concentration [µmol/cm³] profiles of collagen and drug at T = 30,60,120 [min].

(43)

References

[1] R. A. ADAMS, Sobolev Spaces, Academic Press, New York, 1975.

[2] T. ARBOGAST, M. F. WHEELER AND N. Y. ZHANG, A nonlinear mixed finite element method for a degenerate parabolic equation arising in flow in porous media, SIAM J. Numer. Anal. 33, pp. 1669–1687, 1996.

[3] J. BARANGER, J-F. MAITRE AND F. OUDIN, Connection between finite volume and mixed finite element methods, RAIRO Model. Math. Anal.

Numer., Vol. 30 No. 4, pp. 445–465, 1996.

[4] P. BASTIAN, K. BIRKEN, K. JOHANSSEN, S. LANG, N. NEUSS, H.

RENTZ-REICHERT AND C. WIENERS, UG–a flexible toolbox for solving

(44)

Springer Verlag, New York, 1991.

[6] Z. CHEN, Finite Element Methods and Their Applications, Springer Verlag, 2005.

[7] P. KNABNER AND L. ANGERMANN, Numerical methods for elliptic and parabolic partial differential equations, Springer Verlag, 2003.

[8] S. MICHELETTI, R. SACCO, F. SALERI, On Some Mixed Finite Element Methods with Numerical Integration, SIAM J. Sci. Comput. 23, No.1, 245–270, 2001.

[9] A. QUARTERONI AND A. VALLI, Numerical approximations of partial differential equations, Springer-Verlag, 1994.

[10] A. QUARTERONI, R. SACCO AND F. SALERI, Numerical mathematics, Springer-Verlag, New York, 2000.

[11] F. RADU, I. S. POP AND P. KNABNER, Order of convergence estimates

(45)

for an Euler implicit, mixed finite element discretization of Richards’

equation, SIAM J. Numer. Anal. 42, No. 4, pp. 1452–1478, 2004.

[12] B. RYNNE AND M. YOUNGSON, Linear functional analysis, Springer-Verlag, 2000.

[13] C. WOODWARD AND C. DAWSON, Analysis of expanded mixed finite element methods for a nonlinear parabolic equation modeling flow into variably saturated porous media, SIAM J. Numer. Anal. 37, pp. 701–724, 2000.

(46)

Software

• Richy1D.

• UG.

• Matlab or octave.

(47)

Homepage

http://personal-homepages.

mis.mpg.de/fradu/mfem course.html mailto:[email protected]

(48)

Lecture developing

(49)

Lecture developing

0 Preliminary topics.

(50)

Lecture developing

1 Theory of MFEM exemplified on the Poisson equation.

(51)

Lecture developing

1.1 Conforming variational formulation and equivalence with a minimization problem for the Poisson equation.

(52)

Lecture developing

1.2 Mixed variational formulation and equivalence with the conforming method.

(53)

Lecture developing

1.2 Mixed variational formulation and equivalence with the conforming method.

1.3 Abstract formulation of the continuous mixed problem.

Equivalence with a saddle point problem.

(54)

1.4 Existence and uniqueness for the continuous mixed problem.

(55)

1.5 The discrete mixed variational problem. Existence and Uniqueness. The inf-sup condition.

(56)

1.6 Error estimates.

(57)

1.7 Criteria for checking the inf-sup.

(58)

1.8 Extensions of the theory.

(59)

1.8.1 More complicated equations (parabolic problems).

(60)

1.8.1 More complicated equations (parabolic problems).

1.8.2 Error estimates through duality techniques.

(61)

2 The discrete problem.

(62)

2.1 Mixed FE spaces.

(63)

2.2 Interpolation (projection) operators.

(64)

2.2.1 Construction.

(65)

2.2.1 Construction.

2.2.2 Technical lemmas.

(66)

2.2.1 Construction.

2.2 MHFEM.

(67)

2.2.1 Construction.

2.2 MHFEM.

2.3 Implementation.

(68)

2.2.1 Construction.

2.2 MHFEM.

2.3 Implementation.

2.3.1 MFEM and multigrid.

(69)

2.2.1 Construction.

2.2 MHFEM.

2.3 Implementation.

3 Applications of MFEM and MHFEM.

(70)

2.2.1 Construction.

2.2 MHFEM.

2.3 Implementation.

3.1 Reactive flow in porous media.

(71)

2.2.1 Construction.

2.2 MHFEM.

2.3 Implementation.

3.1 Reactive flow in porous media.

(72)

(73)

4 Connection between MFEM and other numerical schemes.

(74)

4.1 Cell centered FV.

(75)

4.2 Multi point flux approximation (MPFA) method.

(76)

5 MFEM and adaptivity.

(77)

5 MFEM and adaptivity.

Hopefully we will enjoy it !!!