Advances in the mathematical theory of the finite element immersed - - PowerPoint PPT Presentation

Advances in the mathematical theory of the finite element immersed boundary method

Daniele Boffi

Dipartimento di Matematica “F. Casorati”, Universit` a di Pavia http://www-dimat.unipv.it/boffi

May 12, 2014

Outline

1

Immersed boundary method

2

Mass conservation

3

IBM with Lagrange multiplier Main collaborators: Lucia Gastaldi, Nicola Cavallini, Francesca Gardini

Outline

1

Immersed boundary method The model FE approximation CFL condition Numerical results

2

Mass conservation

3

IBM with Lagrange multiplier

Immersed boundary method Mass conservation IBM with Lagrange multiplier

IBM – Immersed boundary method

Introduced by Peskin for the simulation of the blood flow in the heart. Peskin ’72–’77 McQueen–Peskin ’83– Peskin ’02 Successfully applied to many biological problems, where a fluid interacts with a flexible structure.

Daniele Boffi IBM and Finite Elements May 12, 2014 page 3

Immersed boundary method Mass conservation IBM with Lagrange multiplier

IBM – Immersed boundary method

Introduced by Peskin for the simulation of the blood flow in the heart. Peskin ’72–’77 McQueen–Peskin ’83– Peskin ’02 Successfully applied to many biological problems, where a fluid interacts with a flexible structure. The main feature is that the structure is considered as a part of the fluid by introducing suitable additional forces and masses. The Navier–Stokes equations are solved in the whole domain (fluid + solid) by finite differences and the interaction with the structure is obtained by means of singular force and mass terms defined by a Dirac delta function localized in the solid domain.

Daniele Boffi IBM and Finite Elements May 12, 2014 page 3

Immersed boundary method Mass conservation IBM with Lagrange multiplier

Finite elements for IBM

At the beginning, we used finite elements mainly because we thought this would simplify the mathematical analysis. Indeed, it turned out that this was a good choice also from the practical point of view. B.–Gastaldi ’03 B.–Gastaldi–Heltai ’04–’07 B.–Gastaldi–Heltai–Peskin ’08 No need to approximating the Dirac delta functions, since the variational formulation takes care of it in a natural way Better interface approximation (less diffusion, sharp pressure jump) The fluid equations can be approximated with standard mixed schemes (Q2 − P1, Hood–Taylor, P1isoP2 − Pc

1, . . . )

Daniele Boffi IBM and Finite Elements May 12, 2014 page 4

Immersed boundary method Mass conservation IBM with Lagrange multiplier

Immersed elastic bodies

Fluid Ω Elastic body Bt

Immersed body of codimension 0 the fluid domain and the immersed body have the same dimension

Elastic boundary Fluid Ω Bt

Immersed body of codimension 1 the immersed body can be either a curve in 2D or a surface in 3D

Daniele Boffi IBM and Finite Elements May 12, 2014 page 5

Immersed boundary method Mass conservation IBM with Lagrange multiplier

Notation

ω Ω X(t) B Bt

Ω fluid + solid Bt deformable structure domain Ω ⊂ Rd, d = 2, 3 Bt ⊂ Rm, m = d, d − 1 x Euler. var. in Ω s Lagrangian var. in B B reference domain u(x, t) fluid velocity X(·, t) : B → Bt position of the solid p(x, t) fluid pressure F = ∂X ∂s deformation grad. (det F > 0) u(x, t) = ∂X ∂t (s, t) where x = X(s, t)

Daniele Boffi IBM and Finite Elements May 12, 2014 page 6

Immersed boundary method Mass conservation IBM with Lagrange multiplier

From conservation of momenta, in absence of external forces, it holds ρ˙ u = ρ ∂u ∂t + u · ∇ u

• = ∇ ·σ

in Ω In our case the Cauchy stress tensor has the following form σ =

• σf

in Ω \ Bt σf + σs in Bt

Daniele Boffi IBM and Finite Elements May 12, 2014 page 7

Immersed boundary method Mass conservation IBM with Lagrange multiplier

From conservation of momenta, in absence of external forces, it holds ρ˙ u = ρ ∂u ∂t + u · ∇ u

• = ∇ ·σ

in Ω In our case the Cauchy stress tensor has the following form σ =

• σf

in Ω \ Bt σf + σs in Bt Incompressible fluid: σ = σf = −pI + µ(∇ u + (∇ u)T) Visco-elastic material: σ = σf + σs with σs elastic part of the stress

Daniele Boffi IBM and Finite Elements May 12, 2014 page 7

Immersed boundary method Mass conservation IBM with Lagrange multiplier

From conservation of momenta, in absence of external forces, it holds ρ˙ u = ρ ∂u ∂t + u · ∇ u

• = ∇ ·σ

in Ω In our case the Cauchy stress tensor has the following form σ =

• σf

in Ω \ Bt σf + σs in Bt Incompressible fluid: σ = σf = −pI + µ(∇ u + (∇ u)T) Visco-elastic material: σ = σf + σs with σs elastic part of the stress Moreover, if the structural material has a density ρs different from the fluid density ρf , we have ρ =

• ρf

in Ω \ Bt ρs in Bt

Daniele Boffi IBM and Finite Elements May 12, 2014 page 7

Immersed boundary method Mass conservation IBM with Lagrange multiplier

We get the following variational formulation Navier–Stokes ρf d dt (u(t), v) + a(u(t), v) + b(u(t), u(t), v) − (div v, p(t)) = d(t), v + F(t), v ∀v ∈ H1

0(Ω)d

(div u(t), q) = 0 ∀q ∈ L2

0(Ω)

Excess Lagrangian mass density d(t), v = −δρ

• B

∂2X ∂t2 v(X(s, t)) ds Load F(t), v = −

• B

˜ P(F(s, t)) : ∇s v(X(s, t)) ds ∀v ∈ H1

0(Ω)d

Body movement ∂X ∂t (s, t) = u(X(s, t), t) ∀s ∈ B Initial conditions u(x, 0) = u0(x) ∀x ∈ Ω, X(s, 0) = X0(s) ∀s ∈ B

Daniele Boffi IBM and Finite Elements May 12, 2014 page 8

Immersed boundary method Mass conservation IBM with Lagrange multiplier

The excess Lagrangian mass density is defined as follows δρ =

• ρs − ρf

codimension 0 ts(ρs − ρf ) codimension 1 The following definition ˜ P =

• P

codimension 0 tsP codimension 1 makes use of the Piola–Kirchoff tensor and takes into account the change

• f variable

P(s, t) = |F(s, t)|σs(X(s, t), t)F−T(s, t) in order to have

• ∂Pt

σsn da =

• ∂P

PN dA ∀Pt

Daniele Boffi IBM and Finite Elements May 12, 2014 page 9

Immersed boundary method Mass conservation IBM with Lagrange multiplier

Stability

B.–Cavallini–Gastaldi ’11 Recalling that ∂X ∂t (s, t) = u(X(s, t), t) ∀s ∈ B it holds ρf 2 d dt ||u(t)||2

0 + µ|| ∇ u(t)||2 0 + d

dt E(X(t)) + 1 2(ρs − ρf ) d dt

• ∂X

∂t

• 2

B

= 0 where E is the total elastic potential energy E (X(t)) =

• B

W (F(s, t)) ds

Daniele Boffi IBM and Finite Elements May 12, 2014 page 10

Immersed boundary method Mass conservation IBM with Lagrange multiplier

Finite element approximation

Uniform background grid Th for the domain Ω (meshsize hx) Inf-sup stable finite element pair

Vh ⊂ H1

0(Ω)d

Qh ⊂ L2

0(Ω)

Grid Sh for B (meshsize hs) Piecewise linear finite element space for X Sh = {Y ∈ C 0(B; Ω) : Y ∈ P1} Notation Tk, k = 1, . . . , Me elements of Sh sj, j = 1, . . . , M vertices of Sh Eh set of the edges e of Sh

Daniele Boffi IBM and Finite Elements May 12, 2014 page 11

Immersed boundary method Mass conservation IBM with Lagrange multiplier

Discrete source term

Source term: F(t), v = −

• B P(Fh(s, t)) : ∇s v(Xh(s, t)) ds

∀v ∈ Vh Xh p.w. linear ⇒ Fh, Ph p.w. constant By integration by parts Fh(t), vh = −

Me

• k=1
• Tk

Ph : ∇s v(X(s, t)) ds = −

Me

• k=1
• ∂Tk

PhNv(X(s, t)) dA that is Fh(t), vh = −

• e∈Eh
• e

[ [Ph] ] · v(X(s, t)) dA [ [P] ] = P+N+ + P−N− jump of P across e for internal edges [ [P] ] = PN jump when e ⊂ ∂B

Daniele Boffi IBM and Finite Elements May 12, 2014 page 12

Immersed boundary method Mass conservation IBM with Lagrange multiplier

The semidiscrete problem becomes: find (uh, ph) : ]0, T[ → Vh × Qh and Xh : [0, T] → Sh such that                        ρf d dt (uh(t), v) + a(uh(t), v) + b(uh(t), uh(t), v) −(div v, ph(t)) = −

• B

(ρs − ρf )∂2Xh ∂t2 v(Xh(s, t))ds −

• e∈Eh
• e

[ [Ph] ] · v(Xh(s, t))dA ∀v ∈ Vh (div uh(t), q) = 0 ∀q ∈ Qh dXhi dt (t) = uh(Xhi(t), t) ∀i = 1, . . . , M uh(0) = u0h in Ω Xhi(0) = X0(si) ∀i = 1, . . . , M

Daniele Boffi IBM and Finite Elements May 12, 2014 page 13

Immersed boundary method Mass conservation IBM with Lagrange multiplier

Fully discrete problem

Backward Euler – BE

Find (un+1

h

, pn+1

h

) ∈ Vh × Qh e Xn+1

h

∈ Sh such that Fn+1

h

, vh = −

• e∈Eh
• e

[ [Ph] ]n+1 · v(Xn+1

h

(s))dA ∀v ∈ Vh NS                          ρf

• un+1

h

− un

h

∆t , v

• + a(un+1

h

, v) + b(un+1

h

, un+1

h

, v) −(div v, pn+1

h

) = −

• B

(ρs − ρf )Xn+1

h

− 2Xn

h + Xn−1 h

∆t2 · v(Xn+1

h

(s))ds + < Fn+1

h

, v >h ∀v ∈ Vh (div un+1

h

, q) = 0 ∀q ∈ Qh; Xn+1

hi

− Xn

hi

∆t = un+1

h

(Xn+1

hi

) ∀i = 1, . . . , M.

Daniele Boffi IBM and Finite Elements May 12, 2014 page 14

Immersed boundary method Mass conservation IBM with Lagrange multiplier

Fully discrete problem

Modified backward Euler – MBE

Step 1. Fn

h, vh = −

• e∈Eh
• e

[ [Ph] ]n · v(Xn

h(s, t)) dA

∀v ∈ Vh Step 2. find (un+1

h

, pn+1

h

) ∈ Vh × Qh such that NS                          ρf

• un+1

h

− un

h

∆t , v

• + a(un+1

h

, v) + b(un+1

h

, un+1

h

, v) −(div v, pn+1

h

) = −

• B

(ρs − ρf )Xn+1

h

− 2Xn

h + Xn−1 h

∆t2 · v(Xn

h(s))ds

+Fn

h, vh

∀v ∈ Vh (div un+1

h

, q) = 0 ∀q ∈ Qh; Step 3. Xn+1

hi

− Xn

hi

∆t = un+1

h

(Xn

hi)

∀i = 1, . . . , M.

Daniele Boffi IBM and Finite Elements May 12, 2014 page 15

Immersed boundary method Mass conservation IBM with Lagrange multiplier

Using Step 3 in Step 2 we get: Step 1. Fn

h, vh = −

• e∈Eh
• e

[ [Ph] ]n · v(Xn

h(s, t)) dA

∀v ∈ Vh Step 2. find (un+1

h

, pn+1

h

) ∈ Vh × Qh such that NS                          ρf

• un+1

h

− un

h

∆t , v

• + a(un+1

h

, v) + b(un+1

h

, un+1

h

, v) −(div v, pn+1

h

) = −

• B

(ρs − ρf )un+1

h

(Xn

h(s)) − un h(Xn−1 h

(s)) ∆t · v(Xn

h(s))ds

+Fn

h, vh

∀v ∈ Vh (div un+1

h

, q) = 0 ∀q ∈ Qh; Step 3. Xn+1

hi

− Xn

hi

∆t = un+1

h

(Xn

hi)

∀i = 1, . . . , M.

Daniele Boffi IBM and Finite Elements May 12, 2014 page 16

Immersed boundary method Mass conservation IBM with Lagrange multiplier

Discrete Energy Estimate

B.–Cavallini–Gastaldi ’11 Artificial Viscosity Theorem Let un

h, pn h and Xn h be a solution to the FE-IBM, then

ρf 2∆t

• un+1

h

2

0 − un h2

• + (µ + µa) ∇ un+1

h

2 + 1 ∆t

• E
• Xn+1

h

• − E [Xn

h]

• +

1 2∆t (ρs − ρf )

• un+1

h

(Xn

h)2 0,B − un h(Xn−1 h

2

0,B

• ≤ 0

CFL Conditions: µ + µa ≥ 0, ρs ≥ ρf (might be relaxed)

Daniele Boffi IBM and Finite Elements May 12, 2014 page 17

Immersed boundary method Mass conservation IBM with Lagrange multiplier

CFL condition

BE is unconditionally stable, while MBE requires the term µa to be not too large µa = −κmaxC h(m−2)

s

∆t h(d−1)

x

Ln Ln := max

Tk∈Sh

• max

sj,si∈V (Tk) |Xn hj − Xn hi|

• space dim.

solid dim. CFL condition 2 1 Ln∆t ≤ Chxhs 2 2 Ln∆t ≤ Chx 3 2 Ln∆t ≤ Ch2

x

3 3 Ln∆t ≤ Ch2

x/hs

Daniele Boffi IBM and Finite Elements May 12, 2014 page 18

Immersed boundary method Mass conservation IBM with Lagrange multiplier

Some numerical results

Original 2D code in Fortran 77, ported to DEAL.II (c++) (www.dealii.org) by L. Heltai (Q2 − P1) 2D Codimension 1 Codimension 0 3D Codimension 1

Daniele Boffi IBM and Finite Elements May 12, 2014 page 19

Immersed boundary method Mass conservation IBM with Lagrange multiplier

More numerical results

Fortran 90 code written by N. Cavallini (P1isoP2 − Pc

1)

Densities: ρs = 21 and ρf = 1 κ = 1 κ = 0.1 κ = 0.1 κ = 1 κ = 0.1 Heart valve (Auricchio–B.–Cavallini–Gastaldi–Lefieux)

Daniele Boffi IBM and Finite Elements May 12, 2014 page 20

Outline

1

Immersed boundary method

2

Mass conservation Inf-sup condition Main theorem Numerical results Mass conservation and IBM

3

IBM with Lagrange multiplier

Immersed boundary method Mass conservation IBM with Lagrange multiplier

Mass conservation of the IBM

B.–Cavallini–Gardini–Gastaldi ’12 Well-known and studied problem The discrete divergence free condition is imposed in a weak sense

• Ω

div uhqh dx = 0 ∀qh ∈ Qh which is not exact unless div(Vh) ⊂ Qh

Daniele Boffi IBM and Finite Elements May 12, 2014 page 22

Immersed boundary method Mass conservation IBM with Lagrange multiplier

Mass conservation of the IBM

B.–Cavallini–Gardini–Gastaldi ’12 Well-known and studied problem The discrete divergence free condition is imposed in a weak sense

• Ω

div uhqh dx = 0 ∀qh ∈ Qh which is not exact unless div(Vh) ⊂ Qh Basic remark Discontinuous pressure schemes enjoy local mass conservation properties (average of divergence is zero element by element)

Daniele Boffi IBM and Finite Elements May 12, 2014 page 22

Immersed boundary method Mass conservation IBM with Lagrange multiplier

Our elements

Hood–Taylor P1isoP2 − Pc

1

Enhanced Hood–Taylor Enhanced P1isoP2 − Pc

1

We actually considered generalized Hood–Taylor in two and three dimensions Pk+1 − Pc

k (k ≥ 1)

Local mass conservation is guaranteed by extra degree of freedom: add piecewise constant pressures

Daniele Boffi IBM and Finite Elements May 12, 2014 page 23

Immersed boundary method Mass conservation IBM with Lagrange multiplier

Analysis of our elements

Known facts Hood–Taylor Introduced in 1973 Hood–Taylor ’73 First analysis Bercovier–Pironneau ’79, Verf¨ urth ’84 Full analysis with some restrictions on boundary elements Scott–Vogelius ’85, Brezzi–Falk ’91 General analysis for the Pk+1 − Pc

k element with no restrictions (mesh

contains at least 3 elements)

• B. ’94

Daniele Boffi IBM and Finite Elements May 12, 2014 page 24

Immersed boundary method Mass conservation IBM with Lagrange multiplier

Analysis of our elements

Known facts Hood–Taylor Introduced in 1973 Hood–Taylor ’73 First analysis Bercovier–Pironneau ’79, Verf¨ urth ’84 Full analysis with some restrictions on boundary elements Scott–Vogelius ’85, Brezzi–Falk ’91 General analysis for the Pk+1 − Pc

k element with no restrictions (mesh

contains at least 3 elements)

• B. ’94

P1isoP2 − Pc

1

Same analysis as for the Hood-Taylor element can be carried on Bercovier–Pironneau ’79, Brezzi–Fortin ’91 Error estimates are suboptimal (unbalanced spaces); ease of implementation makes it appealing, in particular in 3D

Daniele Boffi IBM and Finite Elements May 12, 2014 page 24

Immersed boundary method Mass conservation IBM with Lagrange multiplier

Analysis of our elements (cont’ed)

Pressure enhancement Numerical evidence for lowest order Hood-Taylor (triangles and squares) Gresho–Lee–Chan–Leone ’80 Griffiths ’82 Tidd–Thatcher–Kaye ’88 Proof of inf-sup for lowest order Hood-Taylor (triangles and squares) Thatcher ’90, Pierre ’94, Quin–Zhang ’05

Daniele Boffi IBM and Finite Elements May 12, 2014 page 25

Immersed boundary method Mass conservation IBM with Lagrange multiplier

Analysis of our elements (cont’ed)

Theorem (B.–Cavallini–Gardini–Gastaldi ’12) The generalized enhanced Hood-Taylor scheme Pk+1 − (Pc

k + P0)

in two (k ≥ 1) and three (k ≥ 2) dimensions and the enhanced P1isoP2 − (Pc

1 + P0)

in two dimensions satisfy the inf-sup condition

Daniele Boffi IBM and Finite Elements May 12, 2014 page 26

Immersed boundary method Mass conservation IBM with Lagrange multiplier

Analysis of our elements (cont’ed)

Theorem (B.–Cavallini–Gardini–Gastaldi ’12) The generalized enhanced Hood-Taylor scheme Pk+1 − (Pc

k + P0)

in two (k ≥ 1) and three (k ≥ 2) dimensions and the enhanced P1isoP2 − (Pc

1 + P0)

in two dimensions satisfy the inf-sup condition Minimal restriction on the mesh: each element has at least one internal vertex.

Daniele Boffi IBM and Finite Elements May 12, 2014 page 26

Immersed boundary method Mass conservation IBM with Lagrange multiplier

Mesh restrictions

2D: let us understand the restrictions Standard schemes: the mesh needs at least three elements Enhanced schemes: each element needs at least an internal vertex Uniform mesh Symmetric mesh

Daniele Boffi IBM and Finite Elements May 12, 2014 page 27

Immersed boundary method Mass conservation IBM with Lagrange multiplier

Numerical results

Ω =]0, 1[×]0, 1[ f chosen such that exact solution is u(x, y) = rot ϕ(x, y) ϕ(x, y) = x2(x − 1)2y2(1 − y)2 p(x, y) = x Solution computed with the four different schemes on uniform and symmetric meshes, successively refined

Daniele Boffi IBM and Finite Elements May 12, 2014 page 28

Immersed boundary method Mass conservation IBM with Lagrange multiplier

Hood-Taylor

0.5 1 0.5 1 0.5 1 x y pressure analytical numerical

Uniform mesh

0.5 1 0.5 1 0.5 1 x y pressure analytical numerical

Symmetric mesh

Daniele Boffi IBM and Finite Elements May 12, 2014 page 29

Immersed boundary method Mass conservation IBM with Lagrange multiplier

Enhanced Hood-Taylor

0.5 1 0.5 1 0.5 1 x y pressure analytical numerical

Uniform mesh

0.5 1 0.5 1 0.5 1 x y pressure analytical numerical

Symmetric mesh

Daniele Boffi IBM and Finite Elements May 12, 2014 page 30

Immersed boundary method Mass conservation IBM with Lagrange multiplier

P1isoP2 − Pc

1

0.5 1 0.5 1 0.5 1 x y pressure analytical numerical

Uniform mesh

0.5 1 0.5 1 0.5 1 x y pressure analytical numerical

Symmetric mesh

Daniele Boffi IBM and Finite Elements May 12, 2014 page 31

Immersed boundary method Mass conservation IBM with Lagrange multiplier

Enhanced P1isoP2 − (Pc

1 + P0)

0.5 1 0.5 1 0.5 1 x y pressure analytical numerical

Uniform mesh

0.5 1 0.5 1 0.5 1 x y pressure analytical numerical

Symmetric mesh

Daniele Boffi IBM and Finite Elements May 12, 2014 page 32

Immersed boundary method Mass conservation IBM with Lagrange multiplier

Divergence convergence in L2

Hood–Taylor h HT rate HT+ rate 1/4 6.373e-03

• 6.385e-03
• 1/8

1.797e-03 1.8 1.809e-03 1.8 1/16 4.710e-04 1.9 4.724e-04 1.9 1/32 1.197e-04 2.0 1.198e-04 2.0 1/64 3.006e-05 2.0 3.007e-05 2.0 Bercovier–Pironneau h BP rate BP+ rate 1/4 1.038e-02

• 1.106e-02
• 1/8

5.435e-03 0.9 6.054e-03 0.9 1/16 2.737e-03 1.0 3.118e-03 1.0 1/32 1.368e-03 1.0 1.571e-03 1.0 1/64 6.832e-04 1.0 7.871e-04 1.0

Daniele Boffi IBM and Finite Elements May 12, 2014 page 33

Immersed boundary method Mass conservation IBM with Lagrange multiplier

A more significant example

Consider a jump in the pressure p =

• y(1 − y) exp (x − 1/2)2 + 1/2,

x ≥ 1/2 y(1 − y) exp (x − 1/2)2 − 1/2, x < 1/2

Daniele Boffi IBM and Finite Elements May 12, 2014 page 34

Immersed boundary method Mass conservation IBM with Lagrange multiplier

A more significant example

Consider a jump in the pressure p =

• y(1 − y) exp (x − 1/2)2 + 1/2,

x ≥ 1/2 y(1 − y) exp (x − 1/2)2 − 1/2, x < 1/2 Theorerically, we cannot have more than O(h1/2) when approximating the pressure with continuous finite elements, while discontinuous elements can achieve optimal convergence N.B.: the mesh fits with the jump

Daniele Boffi IBM and Finite Elements May 12, 2014 page 34

Immersed boundary method Mass conservation IBM with Lagrange multiplier

The pressure error

HT

0.2 0.4 0.6 0.8 1 0.5 1 −0.5 0.5 y x p−ph

HT+

0.2 0.4 0.6 0.8 1 0.5 1 −0.5 0.5 y x p−ph

BP

0.5 1 0.5 p−ph

BP+

1 0.5 p−ph

Daniele Boffi IBM and Finite Elements May 12, 2014 page 35

Immersed boundary method Mass conservation IBM with Lagrange multiplier

Divergence convergence in L2

Hood–Taylor h HT rate HT+ rate 1/4 1.623e-01

• 6.486e-03
• 1/8

1.251e-01 0.4 1.832e-03 1.9 1/16 9.178e-02 0.4 4.752e-04 1.9 1/32 6.604e-02 0.5 1.201e-04 2.0 1/64 4.710e-02 0.5 3.011e-05 2.0 Bercovier–Pironneau h BP rate BP+ rate 1/4 1.421e-01

• 1.107e-02
• 1/8

1.137e-01 0.3 6.055e-03 0.9 1/16 8.449e-02 0.4 3.118e-03 1.0 1/32 6.114e-02 0.5 1.571e-03 1.0 1/64 4.371e-02 0.5 7.871e-04 1.0

Daniele Boffi IBM and Finite Elements May 12, 2014 page 36

Immersed boundary method Mass conservation IBM with Lagrange multiplier

Iterative solver

Number of iterations needed to reach convergence when using conjugate gradient ` a la Glowinski Element type Iterations hp = 1/8 hp = 1/16 hp = 1/32 P2 − Pc

1

130 169 172 P2 − (Pc

1 + P0)

25 29 29 P1isoP2 − Pc

1

19 24 24 P1isoP2 − (Pc

1 + P0)

30 35 35

Daniele Boffi IBM and Finite Elements May 12, 2014 page 37

Immersed boundary method Mass conservation IBM with Lagrange multiplier

Mass conservation and IBM

Inflated balloon test case X0(s) = R cos(s/R) + 0.5 R sin(s/R) + 0.5

• ,

s ∈ [0, 2πR] F(t), v = −κ 2πR ∂X(s, t) ∂s ∂v(X(s, t)) ∂s ]

Daniele Boffi IBM and Finite Elements May 12, 2014 page 38

Immersed boundary method Mass conservation IBM with Lagrange multiplier

Mass conservation and IBM

Inflated balloon test case X0(s) = R cos(s/R) + 0.5 R sin(s/R) + 0.5

• ,

s ∈ [0, 2πR] F(t), v = −κ 2πR ∂X(s, t) ∂s ∂v(X(s, t)) ∂s ] p(x, t) = κ(1/R − πR), |x| ≤ R −κπR, |x| > R ∀t ∈]0, T[ T = 10−1 ρf = ρs = 1 µ = 1 κ = 1 hx = 1/32 hs = 2πR/1024

Daniele Boffi IBM and Finite Elements May 12, 2014 page 38

Immersed boundary method Mass conservation IBM with Lagrange multiplier

Area loss w.r.t. time

0.0001 0.04 0.08 0.5 1 1.5 2 2.5 x 10

−3

time 1−A/A0 P1isoP2/P1 P2/P1 P1isoP2/(P1+P0) P2/(P1+P0)

Area loss

0.0001 0.04 0.08 0.05 0.1 0.15 0.2 0.25 0.3 0.35 time ||∇ ⋅ u||L

2

Divergence norm

Daniele Boffi IBM and Finite Elements May 12, 2014 page 39

Outline

1

Immersed boundary method

2

Mass conservation

3

IBM with Lagrange multiplier An interface problem IBM and DLM/FD Mass conservation and stability Analysis of IBM-DLM/FD

Immersed boundary method Mass conservation IBM with Lagrange multiplier

Fictitious domain

There are strong connections between IBM and fictitious domain (FD) with distributed Lagrangian multiplier (DLM). There is a rich literature concerning DLM for FD. Glowinski–Pan–Periaux ’94-’98 Baaijens ’01 van Loon–Anderson–de Hart–Baaijens ’04 Yu ’05 van Loon–Anderson–van de Vosse–Sherwin ’07 Yu–Shao ’07 Vos–van Loon–Sherwin ’08 Pati–Ladipo–Paniagua–Glowinski ’11

Daniele Boffi IBM and Finite Elements May 12, 2014 page 41

Immersed boundary method Mass conservation IBM with Lagrange multiplier

An interface problem

Auricchio–B.–Gastaldi–Lefieux–Reali ’13 B.–Gastaldi–Ruggeri ’13 Let us consider a standard interface problem − div(β1 ∇ u1) = f1 in Ω1 − div(β2 ∇ u2) = f2 in Ω2 u1 = u2

• n Γ

β1 ∇ u1 · n1 + β2 ∇ u2 · n2 = 0

• n Γ

u1 = 0

• n ∂Ω1 \ Γ

u2 = 0

• n ∂Ω2 \ Γ

Daniele Boffi IBM and Finite Elements May 12, 2014 page 42

Immersed boundary method Mass conservation IBM with Lagrange multiplier

A mixed formulation

Find u ∈ H1

0(Ω), u2 ∈ H1(Ω2), and λ ∈ Λ = [H1(Ω2)]∗ such that

• Ω

β ∇ u · ∇ v dx + λ, v|Ω2 =

• Ω

fv dx ∀v ∈ H1

0(Ω)

• Ω2

(β2 − β) ∇ u2 · ∇ v2 dx − λ, v2 =

• Ω2

(f2 − f )v2 dx ∀v2 ∈ H1(Ω2) µ, u|Ω2 − u2 = 0 ∀µ ∈ Λ This is equivalent to the interface problem with u|Ω1 = u1 if β|Ω1 = β1 and f |Ω1 = f1 Notation: Ω = Ω1 ∪ Ω2

Daniele Boffi IBM and Finite Elements May 12, 2014 page 43

Immersed boundary method Mass conservation IBM with Lagrange multiplier

An alternative mixed formulation

Find u ∈ H1

0(Ω), u2 ∈ H1(Ω2), and ψ ∈ H1(Ω2) such that

• Ω

β ∇ u · ∇ v dx + ((ψ, v|Ω2)) =

• Ω

fv dx ∀v ∈ H1

0(Ω)

• Ω2

(β2 − β) ∇ u2 · ∇ v2 dx − ((ψ, v2)) =

• Ω2

(f2 − f )v2 dx ∀v2 ∈ H1(Ω2) ((ϕ, u|Ω2 − u2)) = 0 ∀ϕ ∈ H1(Ω2) where ((·, ·)) denotes the scalar product in H1(Ω2)

Daniele Boffi IBM and Finite Elements May 12, 2014 page 44

Immersed boundary method Mass conservation IBM with Lagrange multiplier

An alternative mixed formulation

Find u ∈ H1

0(Ω), u2 ∈ H1(Ω2), and ψ ∈ H1(Ω2) such that

• Ω

β ∇ u · ∇ v dx + ((ψ, v|Ω2)) =

• Ω

fv dx ∀v ∈ H1

0(Ω)

• Ω2

(β2 − β) ∇ u2 · ∇ v2 dx − ((ψ, v2)) =

• Ω2

(f2 − f )v2 dx ∀v2 ∈ H1(Ω2) ((ϕ, u|Ω2 − u2)) = 0 ∀ϕ ∈ H1(Ω2) where ((·, ·)) denotes the scalar product in H1(Ω2) Remark The two mixed formulations are equivalent but give rise to different discrete schemes

Daniele Boffi IBM and Finite Elements May 12, 2014 page 44

Immersed boundary method Mass conservation IBM with Lagrange multiplier

Approximation of the mixed formulations

Two meshes: Th for Ω and T2,h for Ω2 Three finite element spaces: Vh continuous p/w linears based on Th V2,h continuous p/w linears based on T2,h Λh = V2,h

Daniele Boffi IBM and Finite Elements May 12, 2014 page 45

Immersed boundary method Mass conservation IBM with Lagrange multiplier

Approximation of the mixed formulations

Two meshes: Th for Ω and T2,h for Ω2 Three finite element spaces: Vh continuous p/w linears based on Th V2,h continuous p/w linears based on T2,h Λh = V2,h Of course, several other choices are possible Remark First mixed formulation makes use of Vh, V2,h, and Λh (duality is represented by scalar product in L2(Ω2)), while second mixed formulation makes use of Vh, V2,h, and V2,h

Daniele Boffi IBM and Finite Elements May 12, 2014 page 45

Immersed boundary method Mass conservation IBM with Lagrange multiplier

Stability of the approximation

We need to show the ellipticity in the kernel and the inf-sup condition ELKER

• Ω

β|∇v|2 dx +

• Ω2

(β2 − β)| ∇ v2|2 dx ≥ κ1

• v2

H1(Ω) + v22 H1(Ω2)

• ∀(v, v2) ∈ Kh

INFSUP1 sup

(v,v2)∈Vh×V2,h

(µ, v|Ω2 − v2)

• v2

H1(Ω) + v22 H1(Ω2)

1/2 ≥ κ2µΛ ∀µ ∈ Λh INFSUP2 sup

(v,v2)∈Vh×V2,h

(ϕ, v|Ω2 − v2)H1(Ω2) (v2

H1(Ω) + v22 H1(Ω2))1/2 ≥ κ2ϕH1(Ω2)

∀ϕ ∈ V2,h

Daniele Boffi IBM and Finite Elements May 12, 2014 page 46

Immersed boundary method Mass conservation IBM with Lagrange multiplier

Theorem If β2 − β|Ω2 ≥ η0 > 0 then ELKER holds true for both formulations, uniformly in h and h2 Theorem If the mesh sequence T2,h is quasi-uniform, then INFSUP1 holds true, uniformly in h and h2 Theorem INFSUP2 holds true, uniformly in h and h2 without any additional assumptions on the mesh sequence Remark For the second mixed formulation, ELKER holds true without assumptions

• n β if h2/hd/2 is small enough and Th is quasi-uniform

Daniele Boffi IBM and Finite Elements May 12, 2014 page 47

Immersed boundary method Mass conservation IBM with Lagrange multiplier

Back to the IBM: fully variational approach

Body motion ∂X ∂t (s, t) = u(X(s, t), t) ∀s ∈ B We introduce a variational formulation:

• µ, u(X(·, t), t) − ∂X(t)

∂t

• B

= 0 ∀µ ∈ Λ where Λ is the dual space of H1(B)d.

Daniele Boffi IBM and Finite Elements May 12, 2014 page 48

Immersed boundary method Mass conservation IBM with Lagrange multiplier

Back to the IBM: fully variational approach

Body motion ∂X ∂t (s, t) = u(X(s, t), t) ∀s ∈ B We introduce a variational formulation:

• µ, u(X(·, t), t) − ∂X(t)

∂t

• B

= 0 ∀µ ∈ Λ where Λ is the dual space of H1(B)d. This equation plays the role of a constraint on the difference between the velocities in the fluid and in the solid, hence we introduce the Lagrange multiplier associated to it.

Daniele Boffi IBM and Finite Elements May 12, 2014 page 48

Immersed boundary method Mass conservation IBM with Lagrange multiplier

The IBM as a DLM/FD method

B.–Gastaldi, in progress For t ∈ [0, T], find u(t) ∈ H1

0(Ω)d, p(t) ∈ L2 0(Ω), X(t) ∈ W 1,∞(B)d, and

λ(t) ∈ Λ such that ρ d dt (u(t), v) + a(u(t), v) + b(u(t), u(t), v) − (div v, p(t)) + λ, v(X(·, t))B = 0 ∀v ∈ H1

0(Ω)d

(div u(t), q) = 0 ∀q ∈ L2

0(Ω)

δρ

• B

∂2X ∂t2 Yds +

• B

P(F(s, t)) : ∇s Y ds − λ, YB = 0 ∀Y ∈ H1(B)d

• µ, u(X(·, t), t) − ∂X(t)

∂t

• B

= 0 ∀µ ∈ Λ u(0) = u0 in Ω, X(0) = X0 in B

Daniele Boffi IBM and Finite Elements May 12, 2014 page 49

Immersed boundary method Mass conservation IBM with Lagrange multiplier

Time advancing scheme

Implicit scheme

Given u0 ∈ H1

0(Ω)d and X0 ∈ W 1,∞(B)d, for n = 1, . . . , N find

un, pn ∈ H1

0(Ω)d × L2 0(Ω), Xn ∈ W 1,∞(B)d, and λn ∈ Λ, such that

ρf un+1 − un ∆t , v

• + a(un+1, v) − (div v, pn+1)

+ λn+1, v(Xn+1B = 0 ∀v ∈ H1

0(Ω)d

(div un+1, q) = 0 ∀q ∈ L2

0(Ω)

δρ Xn+1 − 2Xn + Xn−1 ∆t2 , Y

• B

+ (P(Fn+1), ∇s Y)B − λn+1, YB = 0 ∀Y ∈ H1(B)d

• µ, un+1(Xn+1) − Xn+1 − Xn

∆t

• B

= 0 ∀µ ∈ Λ.

Daniele Boffi IBM and Finite Elements May 12, 2014 page 50

Immersed boundary method Mass conservation IBM with Lagrange multiplier

Time advancing scheme

Semi-implicit scheme

Given u0 ∈ H1

0(Ω)d and X0 ∈ W 1,∞(B)d, for n = 1, . . . , N find

un, pn ∈ H1

0(Ω)d × L2 0(Ω), Xn ∈ W 1,∞(B)d, and λn ∈ Λ, such that

ρf un+1 − un ∆t , v

• + a(un+1, v) − (div v, pn+1)

+ λn+1, v(Xn(·))B = 0 ∀v ∈ H1

0(Ω)d

(div un+1, q) = 0 ∀q ∈ L2

0(Ω)

δρ Xn+1 − 2Xn + Xn−1 ∆t2 , Y

• B

+ (P(Fn+1), ∇s Y)B − λn+1, YB = 0 ∀Y ∈ H1(B)

• µ, un+1(Xn(·)) − Xn+1 − Xn

∆t

• B

= 0 ∀µ ∈ Λ.

Daniele Boffi IBM and Finite Elements May 12, 2014 page 51

Immersed boundary method Mass conservation IBM with Lagrange multiplier

Energy estimate for the time discretization

Proposition Unconditional stability Assume that W (energy density) is convex and δρ > 0. The following estimate holds true for all n = 1, . . . , N ρf 2∆t

• un+12

0 − un2

• + µ ∇ un+12

+ δρ 2∆t

• Xn+1 − Xn

∆t

• 2

0,B

−

• Xn − Xn−1

∆t

• 2

0,B

• + 1

∆t (E(Xn+1) − E(Xn)) ≤ 0 Recall the potential energy E (X(t)) =

• B

W (F(s, t)) ds

Daniele Boffi IBM and Finite Elements May 12, 2014 page 52

Immersed boundary method Mass conservation IBM with Lagrange multiplier

IBM versus DLM-IBM

Time advancing scheme for IBM formulation Compute F FSI using Xn Solve at time tn+1 ANS(un) B⊤ B un+1 pn+1

• =

F FSI

• Update pointwise the structure by

Xn+1 − Xn ∆t (s) = un+1(Xn(s))

Daniele Boffi IBM and Finite Elements May 12, 2014 page 53

Immersed boundary method Mass conservation IBM with Lagrange multiplier

IBM versus DLM-IBM

Time advancing scheme for IBM formulation Compute F FSI using Xn Solve at time tn+1 ANS(un) B⊤ B un+1 pn+1

• =

F FSI

• Update pointwise the structure by

Xn+1 − Xn ∆t (s) = un+1(Xn(s)) Time advancing scheme for DLM-IBM formulation      ANS(un) B⊤ L⊤

f (Xn)

B As −L⊤

s

Lf (Xn) −Ls           u p X λ      =      f(un) g(Xn) d(Xn)     

Daniele Boffi IBM and Finite Elements May 12, 2014 page 53

Immersed boundary method Mass conservation IBM with Lagrange multiplier

Ellipse immersed in a static fluid

˜ P = κF Fluid initially at rest: u0h = 0

X0(s) = 0.2 cos(2πs) + 0.45 0.1 sin(2πs) + 0.45

• s ∈ [0, 1],

hx = 1/32, hs = 1/32, ∆t = 10−2, µ = 1, κ = 5

0.5 0.6 0.7 0.8 0.9 0.5 0.55 0.6 0.65 0.7 0.75 0.8 0.85 0.9 x y 0.5 0.6 0.7 0.8 0.9 0.5 0.55 0.6 0.65 0.7 0.75 0.8 0.85 0.9 x y

Standard IBM with PW update of the immersed boundary IBM with DLM

Daniele Boffi IBM and Finite Elements May 12, 2014 page 54

Immersed boundary method Mass conservation IBM with Lagrange multiplier

Stability analysis

ρ = 1, µ = 1, κ = 5, hx = 1/64, hs = 1/128

0.5 1 1.5 2 10

−1

10 10

1

time DLM PW

∆t = 0.01

0.5 1 1.5 2 10

−1

10 10

1

time DLM PW

∆t = 0.05

0.5 1 1.5 2 10

−1

10 10

1

time DLM PW

∆t = 0.1 Daniele Boffi IBM and Finite Elements May 12, 2014 page 55

Immersed boundary method Mass conservation IBM with Lagrange multiplier

Stability analysis

ρ = 1, µ = 1, κ = 5, hx = 1/64, hs = 1/128

0.5 1 1.5 2 10

−1

10 10

1

time DLM PW

∆t = 0.01

0.5 1 1.5 2 10

−1

10 10

1

time DLM PW

∆t = 0.05

0.5 1 1.5 2 10

−1

10 10

1

time DLM PW

∆t = 0.1

ρ = 1, µ = 1, κ = 5, hx = 1/64, ∆t = 0.01

0.5 1 1.5 2 10

−1

10 10

1

time DLM PW

hs = 1/128

0.5 1 1.5 2 10

−1

10 10

1

time DLM PW

hs = 1/256

0.5 1 1.5 2 10

−1

10 10

1

time DLM PW

hs = 1/512 Daniele Boffi IBM and Finite Elements May 12, 2014 page 55

Mass conservation

hs ∆t

1 64 1 48 1 40 1 32 1 24 1 16 1 8

ibm codimension one

1 · 10−2 0.943 0.942 0.942 0.941 0.941 0.941 0.940 2 · 10−2 cfl 0.942 0.941 0.941 0.940 0.940 0.940 3 · 10−2 cfl cfl cfl 0.941 0.940 0.940 0.939 5 · 10−2 cfl cfl cfl 0.795 0.940 0.939 0.938 1 · 10−1 cfl cfl cfl cfl cfl 0.574 0.936

dlm codimension one

1 · 10−2 inf-sup inf-sup inf-sup 1.023 1.022 1.022 1.023 2 · 10−2 inf-sup inf-sup inf-sup 1.022 1.022 1.022 1.022 3 · 10−2 inf-sup inf-sup 1.023 1.022 1.022 1.022 1.022 5 · 10−2 inf-sup inf-sup inf-sup 1.021 1.021 1.021 1.022 1 · 10−1 inf-sup inf-sup 1.021 1.020 1.020 1.020 1.020

Immersed boundary method Mass conservation IBM with Lagrange multiplier

Examples with thick body

Daniele Boffi IBM and Finite Elements May 12, 2014 page 57

Immersed boundary method Mass conservation IBM with Lagrange multiplier

Examples with thick body (cont’ed)

Daniele Boffi IBM and Finite Elements May 12, 2014 page 58

Immersed boundary method Mass conservation IBM with Lagrange multiplier

Saddle point problem

Let X ∈ W 1,∞(B) be invertible with Lipschitz inverse. Given f ∈ (L2(Ω))d, g ∈ (L2(B))m and d ∈ (L2(B))m, find u ∈ H1

0(Ω),

p ∈ L2

0(Ω), X ∈ H1(B) and λ ∈ (H1(B))∗ such that

af (u, v) − (div v, p) + λ, v(X(·))B = (f, v) ∀v ∈ H1

0(Ω)

(div u, q) = 0 ∀q ∈ L2

0(Ω)

as(X, Y) − λ, YB = (g, Y)B ∀Y ∈ H1(B) µ, u(X(·)) − XB = µ, dB ∀µ ∈ (H1(B))∗. where af (u, v) = α(u, v) + a(u, v) ∀u, v ∈ H1

0(Ω)

as(X, Y) = β(X, Y) + γ(P(∇s X), ∇s Y)B ∀X, Y ∈ H1(B)

Daniele Boffi IBM and Finite Elements May 12, 2014 page 59

Immersed boundary method Mass conservation IBM with Lagrange multiplier

Conclusions

1 The finite element Immersed Boundary Method provides interesting

results for the approximation of fluid-structure interaction problems. Rigorous proof of a CFL condition shows that modified BE scheme can be successfully used in this framework

2 We performed a rigorous analysis of locally mass preserving Stokes

element in a general setting

3 Adding p/w constant pressures significantly enhances the

performance of the scheme

4 We are now investigating the use of a Lagrange multiplier for

imposing the body motion in a weak sense. The new formulation looks promising both for its mass conservation properties and for its stability

Daniele Boffi IBM and Finite Elements May 12, 2014 page 60