A guaranteed a posteriori error estimator for certified boundary - - PowerPoint PPT Presentation

SLIDE 1

A guaranteed a posteriori error estimator for certified boundary variation algorithm

Matteo Giacomini

CMAP - Centre de Math´ ematiques Appliqu´ ees, ´ Ecole Polytechnique IPSA - Institut Polytechnique des Sciences Avanc´ ees 6th Workshop FreeFem++ Days - LJLL UPMC, December 11th 2014 Joint work with Olivier Pantz (CMAP ´ Ecole Polytechnique) and Karim Trabelsi (IPSA)

• M. Giacomini

Error estimates for shape optimization

• Dec. 11th 2014 (FreeFem++)

1 / 20

SLIDE 2

Outline

1 Shape optimization and shape identification problems ◮ A scalar model problem ◮ Differentiation with respect to the shape ◮ The boundary variation algorithm 2 A posteriori estimators of the error in the shape derivative ◮ Goal-oriented residual-type error estimator ◮ Validation of the goal-oriented estimator 3 Adaptive shape optimization procedure ◮ A first test case ◮ Improved adaptive shape optimization procedure

• M. Giacomini

Error estimates for shape optimization

• Dec. 11th 2014 (FreeFem++)

2 / 20

SLIDE 3

Outline

1

Shape optimization and shape identification problems A scalar model problem Differentiation with respect to the shape The boundary variation algorithm

2

A posteriori estimators of the error in the shape derivative Goal-oriented residual-type error estimator Validation of the goal-oriented estimator

3

Adaptive shape optimization procedure A first test case Improved adaptive shape optimization procedure

• M. Giacomini

Error estimates for shape optimization

• Dec. 11th 2014 (FreeFem++)

3 / 20

SLIDE 4

Electrical Impedance Tomography

Neumann problem (N):          − k∆uN + uN = 0 in Ω Γ uN = 0

• n Γ

k∇uN · n = 0

• n Γ

k1∇uN · n = g

• n ∂Ω

Dirichlet problem (D):          − k∆uD + uD = 0 in Ω Γ uD = 0

• n Γ

k∇uD · n = 0

• n Γ

uD = UD

• n ∂Ω

k = k0χω + k1(1 − χω)

• M. Giacomini

Error estimates for shape optimization

• Dec. 11th 2014 (FreeFem++)

4 / 20

SLIDE 5

Electrical Impedance Tomography

Neumann problem (N):          − k∆uN + uN = 0 in Ω Γ uN = 0

• n Γ

k∇uN · n = 0

• n Γ

k1∇uN · n = g

• n ∂Ω

Dirichlet problem (D):          − k∆uD + uD = 0 in Ω Γ uD = 0

• n Γ

k∇uD · n = 0

• n Γ

uD = UD

• n ∂Ω

k = k0χω + k1(1 − χω) Objective functional:

J(ω) = 1 2

• Ω

k∇(uN(ω) − uD(ω)) · ∇(uN(ω) − uD(ω))dx+1 2

• Ω

(uN(ω) − uD(ω))2dx

• M. Giacomini

Error estimates for shape optimization

• Dec. 11th 2014 (FreeFem++)

4 / 20

SLIDE 6

Electrical Impedance Tomography

a(u, v) =

• Ω
• k∇u · ∇v + uv
• dx

Neumann problem (N): uN ∈ H1(Ω) a(uN, v) = FN(v) ∀v ∈ WN = H1(Ω) FN(v) =

• ∂Ω

gv dσ Dirichlet problem (D): uD ∈ H1

UD(Ω)

a(uD, v) = FD(v) ∀v ∈ WD = H1

0(Ω)

FD(v) = 0 k = k0χω + k1(1 − χω) Objective functional:

J(ω) = 1 2a (uN(ω) − uD(ω), uN(ω) − uD(ω))

• M. Giacomini

Error estimates for shape optimization

• Dec. 11th 2014 (FreeFem++)

4 / 20

SLIDE 7

Shape optimization approach

PDE-constrained optimization problem: ω∗ = argminω J(ω) = ⇒ Optimization variable: Shape and location of the inclusion ω Shape optimization startegy:

Given the domain Ω(0), set ℓ = 0 and iterate:

• 1. Compute

the solutions u(ℓ)

N

and u(ℓ)

D ;

• 2. Compute a descent direction θ(ℓ) and an admissible

step µ(ℓ);

• 3. Move the

interface Γ(ℓ+1) = (I + µ(ℓ)θ(ℓ))Γ(ℓ);

• 4. While

stopping criterion not fulfilled , ℓ = ℓ + 1 and repeat.

• M. Giacomini

Error estimates for shape optimization

• Dec. 11th 2014 (FreeFem++)

5 / 20

SLIDE 8

Shape optimization approach

PDE-constrained optimization problem: ω∗ = argminω J(ω) = ⇒ Optimization variable: Shape and location of the inclusion ω Shape optimization startegy:

Given the domain Ω(0), set ℓ = 0 and iterate:

• 1. Compute

the solutions u(ℓ)

N

and u(ℓ)

D ;

• 2. Compute a descent direction θ(ℓ) and an admissible

step µ(ℓ);

• 3. Move the

interface Γ(ℓ+1) = (I + µ(ℓ)θ(ℓ))Γ(ℓ);

• 4. While

stopping criterion not fulfilled , ℓ = ℓ + 1 and repeat.

Classical optimization min

x∈Rn f (x)

, f : Rn → R Steepest descent direction at point x: v(x) = −∇f (x) x(ℓ+1) = x(ℓ) + µ(ℓ) v(x(ℓ)) Shape optimization min

ω∈Uad J(ω)

, J : Uad → R Uad = {Open sets in Ω ⊂ Rd} Gradient-based descent direction at ω: θ s.t. dJ(ω), θ < 0 ω(ℓ+1) = (I + µ(ℓ)θ(ℓ))ω(ℓ)

• M. Giacomini

Error estimates for shape optimization

• Dec. 11th 2014 (FreeFem++)

5 / 20

SLIDE 9

Shape derivative

Let θ ∈ W 1,∞(Ω, R2) be an admissible smooth deformation of Ω s.t. the external boundary is fixed: θ = 0 on ∂Ω = ⇒ Small perturbation of the domain: Ω(θ) = (I + tθ)Ω If the map JΩ : θ → J(Ω(θ)) is differentiable at θ = 0, we define the shape derivative: dJ(ω), θ = J ′

Ω(0), θ = lim tց0

J((I + tθ)Ω) − J(Ω) t Shape derivative of the objective functional J(ω)1:

dJ(ω), θ =1 2

• Ω

kM(θ)∇(uN(ω) + uD(ω)) · ∇(uN(ω) − uD(ω))dx −1 2

• Ω

∇ · θ (uN(ω) + uD(ω))(uN(ω) − uD(ω))dx M(θ) = ∇θ + ∇θT − (∇ · θ)I

• 1O. Pantz. Sensibilit´

e de l’´ equation de la chaleur aux sauts de conductivit´

• e. C. R. Acad. Sci.

Paris, Ser. I(341):333-337 (2005)

• M. Giacomini

Error estimates for shape optimization

• Dec. 11th 2014 (FreeFem++)

6 / 20

SLIDE 10

The boundary variation algorithm

Gradient-based strategy: θ, δθ[H1(Ω)]d + dJ(ω), δθ = 0 ∀δθ ∈ [H1(Ω)]d = ⇒ Descent direction: θ such that dJ(ω), θ < 0 There are two possible approaches to numerical shape optimization: Discretize-then-Optimize Optimize-then-Discretize

• M. Giacomini

Error estimates for shape optimization

• Dec. 11th 2014 (FreeFem++)

7 / 20

SLIDE 11

The boundary variation algorithm

Gradient-based strategy: θ, δθ[H1(Ω)]d + dJ(ω), δθ = 0 ∀δθ ∈ [H1(Ω)]d = ⇒ Descent direction: θ such that dJ(ω), θ < 0 There are two possible approaches to numerical shape optimization: Discretize-then-Optimize Optimize-then-Discretize = ⇒ Discretized gradient-based strategy: θh, δθh[X p

h ]d + dhJ, δθh ≃ 0

∀δθh ∈ [X p

h ]d

= ⇒ Certified descent direction: θh such that dhJ, θh + |E h| < 0 Numerical error in the shape derivative: E h = dJ(ωh), θh − dhJ, θh

• M. Giacomini

Error estimates for shape optimization

• Dec. 11th 2014 (FreeFem++)

7 / 20

SLIDE 12

Why coupling error estimates with shape optmization?

Initial interface

• M. Giacomini

Error estimates for shape optimization

• Dec. 11th 2014 (FreeFem++)

8 / 20

SLIDE 13

Why coupling error estimates with shape optmization?

Reconstructed interface Objective functional

• M. Giacomini

Error estimates for shape optimization

• Dec. 11th 2014 (FreeFem++)

8 / 20

SLIDE 14

The adaptive boundary variation algorithm

Given the domain Ω(0) and tol = 10−8, set the ℓ = 0 and iterate:

• 1. Construct a coarse

mesh T (ℓ)

h

⊂ Ω(ℓ);

• 2. Compute

the primal solutions uh(ℓ)

N

and uh(ℓ)

D

;

• 3. Compute a descent

direction θh(ℓ) ;

• 4. Compute the upper bound E of the numerical error |E h|;
• 5. While dhJ, θh(ℓ) + E > 0, repeat:

(a) Adapt the computational mesh T (ℓ)

h

; (b) Re -compute the primal solutions uh(ℓ)

N

and uh(ℓ)

D

; (c) Re -compute a descent direction θh(ℓ) ; (d) Re-compute the upper bound E of |E h| and dhJ, θh(ℓ) + E;

• 6. Compute an admissible

step size µh(ℓ) ;

• 7. Move the mesh T (ℓ+1)

h

= (I + µh(ℓ)θh(ℓ))T (ℓ)

h

;

• 8. While |dhJ, θh(ℓ)| + E > tol, ℓ = ℓ + 1 and repeat.
• M. Giacomini

Error estimates for shape optimization

• Dec. 11th 2014 (FreeFem++)

9 / 20

SLIDE 15

The adaptive boundary variation algorithm

Given the domain Ω(0) and tol = 10−8, set the ℓ = 0 and iterate:

• 1. Construct a coarse

mesh T (ℓ)

h

⊂ Ω(ℓ);

• 2. Compute

the primal solutions uh(ℓ)

N

and uh(ℓ)

D

;

• 3. Compute a descent

direction θh(ℓ) ;

• 4. Compute the upper bound E of the numerical error |E h|;
• 5. While dhJ, θh(ℓ) + E > 0, repeat:

(a) Adapt the computational mesh T (ℓ)

h

; (b) Re -compute the primal solutions uh(ℓ)

N

and uh(ℓ)

D

; (c) Re -compute a descent direction θh(ℓ) ; (d) Re-compute the upper bound E of |E h| and dhJ, θh(ℓ) + E;

• 6. Compute an admissible

step size µh(ℓ) ;

• 7. Move the mesh T (ℓ+1)

h

= (I + µh(ℓ)θh(ℓ))T (ℓ)

h

;

• 8. While |dhJ, θh(ℓ)| + E > tol, ℓ = ℓ + 1 and repeat.
• M. Giacomini

Error estimates for shape optimization

• Dec. 11th 2014 (FreeFem++)

9 / 20

SLIDE 16

Outline

1

Shape optimization and shape identification problems A scalar model problem Differentiation with respect to the shape The boundary variation algorithm

2

A posteriori estimators of the error in the shape derivative Goal-oriented residual-type error estimator Validation of the goal-oriented estimator

3

Adaptive shape optimization procedure A first test case Improved adaptive shape optimization procedure

• M. Giacomini

Error estimates for shape optimization

• Dec. 11th 2014 (FreeFem++)

10 / 20

SLIDE 17

Error in the Quantity of Interest

Error operator: Gθh(u, v) = 1 2

• Ω
• kM(θh)∇(u + v) · ∇(u − v) − ∇ · θh(u + v)(u − v)
• dx

= ⇒ Numerical error:

• E h

=

• Gθh(uN, uh

N) − Gθh(uD, uh D)

• M. Giacomini

Error estimates for shape optimization

• Dec. 11th 2014 (FreeFem++)

11 / 20

SLIDE 18

Error in the Quantity of Interest

Error operator: Gθh(u, v) = 1 2

• Ω
• kM(θh)∇(u + v) · ∇(u − v) − ∇ · θh(u + v)(u − v)
• dx
• Gθh(u) = 1

2

• Ω
• kM(θh)∇u · ∇u − ∇ · θhu2

dx = ⇒ Numerical error:

• E h

=

• Gθh(uN, uh

N) − Gθh(uD, uh D)

• ≤
• Gθh(uN) −

Gθh(uh

N)

• +
• Gθh(uD) −

Gθh(uh

D)

• M. Giacomini

Error estimates for shape optimization

• Dec. 11th 2014 (FreeFem++)

11 / 20

SLIDE 19

Error in the Quantity of Interest

Error operator: Gθh(u, v) = 1 2

• Ω
• kM(θh)∇(u + v) · ∇(u − v) − ∇ · θh(u + v)(u − v)
• dx
• Gθh(u) = 1

2

• Ω
• kM(θh)∇u · ∇u − ∇ · θhu2

dx = ⇒ Numerical error:

• E h

=

• Gθh(uN, uh

N) − Gθh(uD, uh D)

• ≤
• Gθh(uN) −

Gθh(uh

N)

• +
• Gθh(uD) −

Gθh(uh

D)

• E h

N

E h

D

Adjoint problems: a∗(pi, ϕ) = Fi(ϕ) ∀ϕ ∈ Wi , i = N, D Fi(ϕ) =

• Ω
• kM(θh)∇uh

i · ∇ϕ − ∇ · θhuh i ϕ

• dx
• M. Giacomini

Error estimates for shape optimization

• Dec. 11th 2014 (FreeFem++)

11 / 20

SLIDE 20

Goal-oriented residual-type error estimator

Variational formulation of the error via the adjoint problems (i = N, D): E h

i =

Gθh(ui) − Gθh(uh

i ) ≃ Fi(ei) = a∗(pi, ei) = a∗(pi − ph i , ei) = a(ei, ǫi)

Residual equations (i = N, D): a(ei, v) = Ri(v) ∀v ∈ Wi , Ri(v) = Fi(v) − a(uh

i , v)

a∗(ǫi, ϕ) = R∗

i (ϕ)

∀ϕ ∈ Wi , R∗

i (ϕ) = Fi(ϕ) − a∗(ph i , ϕ)

= ⇒ Residual error estimator in the QoI 2:

• E h

i

• ≤
• T∈Th

ηi

• Tξi
• T

ηi primal residue ξi adjoint residue Global residual-type error estimator for the shape derivative: ζ :=

• T∈Th
• ηN
• TξN
• T + ηD
• TξD
• T
• 2R. Becker, R. Rannacher. An optimal control approach to a posteriori error estimation in Finite

Element Methods. Acta Numerica, 10:1-102 (2001)

• M. Giacomini

Error estimates for shape optimization

• Dec. 11th 2014 (FreeFem++)

12 / 20

SLIDE 21

Validation of the goal-oriented mesh adaptivity strategy

• M. Giacomini

Error estimates for shape optimization

• Dec. 11th 2014 (FreeFem++)

13 / 20

SLIDE 22

Outline

1

Shape optimization and shape identification problems A scalar model problem Differentiation with respect to the shape The boundary variation algorithm

2

A posteriori estimators of the error in the shape derivative Goal-oriented residual-type error estimator Validation of the goal-oriented estimator

3

Adaptive shape optimization procedure A first test case Improved adaptive shape optimization procedure

• M. Giacomini

Error estimates for shape optimization

• Dec. 11th 2014 (FreeFem++)

14 / 20

SLIDE 23

Adaptive shape optimization procedure 3

Given the domain Ω(0) and tol = 10−8, set the ℓ = 0 and iterate:

• 1. Construct a coarse

mesh T (ℓ)

h

⊂ Ω(ℓ);

• 2. Compute

the primal solutions uh(ℓ)

N

and uh(ℓ)

D

;

• 3. Compute a descent

direction θh(ℓ) ;

• 4. Compute the adjoint solutions ph(ℓ)

N

and ph(ℓ)

D

;

• 5. Compute the error estimator ζ(ℓ);
• 6. While dhJ, θh(ℓ) + ζ(ℓ) > 0, repeat:

(a) Adapt the computational mesh T (ℓ)

h

w.r.t. ζ(ℓ); (b) Re -compute the primal solutions uh(ℓ)

N

and uh(ℓ)

D

; (c) Re -compute a descent direction θh(ℓ) ; (d) Re-compute the adjoint solutions ph(ℓ)

N

and ph(ℓ)

D

; (e) Re-compute the error estimator ζ(ℓ) and dhJ, θh(ℓ) + ζ(ℓ);

• 7. Compute an admissible

step size µh(ℓ) ;

• 8. Move the mesh T (ℓ+1)

h

= (I + µh(ℓ)θh(ℓ))T (ℓ)

h

;

• 9. While |dhJ, θh(ℓ)| + ζ(ℓ) > tol, ℓ = ℓ + 1 and repeat.
• 3M. Giacomini, O. Pantz, K. Trabelsi. An adaptive shape optimization strategy driven by fully-

computable goal-oriented error estimators. (In preparation)

• M. Giacomini

Error estimates for shape optimization

• Dec. 11th 2014 (FreeFem++)

15 / 20

SLIDE 24

A first test case

Evolution of the interface Objective functional

All the simulations were performed using FreeFem++.

• F. Hecht. New development in FreeFem++. J. Numer. Math., 20(3-4):251-265 (2012)
• M. Giacomini

Error estimates for shape optimization

• Dec. 11th 2014 (FreeFem++)

16 / 20

SLIDE 25

A first test case

Reconstructed interface Objective functional

• M. Giacomini

Error estimates for shape optimization

• Dec. 11th 2014 (FreeFem++)

16 / 20

SLIDE 26

Improved adaptive shape optimization procedure 4

Given the domain Ω(0) and tol = 10−8, set the ℓ = 0 and iterate:

• 1. Construct a coarse mesh T (ℓ)

h

⊂ Ω(ℓ) and a fine mesh S(ℓ)

h

⊂ Ω(ℓ);

• 2. Compute

the primal solutions uh(ℓ)

N

and uh(ℓ)

D

• n S(ℓ)

h ;

• 3. Compute a descent

direction θh(ℓ)

• n S(ℓ)

h ;

• 4. Compute

the adjoint solutions ph(ℓ)

N

and ph(ℓ)

D

• n S(ℓ)

h ;

• 5. Compute

the error estimator ζ(ℓ);

• 6. While dhJ, θh(ℓ) + ζ(ℓ) > 0, repeat:

(a) Adapt the computational mesh S(ℓ)

h

w.r.t. ζ(ℓ); (b) Re -compute the primal solutions uh(ℓ)

N

and uh(ℓ)

D

• n S(ℓ)

h ;

(c) Re -compute a descent direction θh(ℓ)

• n S(ℓ)

h ;

(d) Re -compute the adjoint solutions ph(ℓ)

N

and ph(ℓ)

D

• n S(ℓ)

h ;

(e) Re -compute the error estimator ζ(ℓ) and dhJ, θh(ℓ) + ζ(ℓ);

• 7. Compute an admissible

step size µh(ℓ) ;

• 8. Move the coarse

mesh T (ℓ+1)

h

= (I + µh(ℓ)θ(ℓ))T (ℓ)

h

;

• 9. While |dhJ, θh(ℓ)| + ζ(ℓ) > tol, ℓ = ℓ + 1 and repeat.
• 4G. Allaire, O. Pantz. Structural optimization with FreeFem++. Struct. Multidiscip. Optim.

32:3 (2006)

• M. Giacomini

Error estimates for shape optimization

• Dec. 11th 2014 (FreeFem++)

17 / 20

SLIDE 27

1-mesh VS 2-mesh optimization strategy

1-mesh strategy 2-mesh strategy

• M. Giacomini

Error estimates for shape optimization

• Dec. 11th 2014 (FreeFem++)

18 / 20

SLIDE 28

1-mesh VS 2-mesh optimization strategy

Reconstructed interfaces Objective functionals

• M. Giacomini

Error estimates for shape optimization

• Dec. 11th 2014 (FreeFem++)

18 / 20

SLIDE 29

The case of multiple measures

10 measures using 2-mesh strategy: Reconstructed interface Objective functional

• M. Giacomini

Error estimates for shape optimization

• Dec. 11th 2014 (FreeFem++)

19 / 20

SLIDE 30

Conclusions

Guaranteed shape optimization strategy using certified goal-oriented estimates for the error in the shape derivative. Sharp bounds are obtained using DWR estimator. A reliable stopping criterion is derived. A 2-mesh strategy improves the convergence of the algorithm. Special thanks to FreeFem++ for: Powerful management of geometries and meshes. Adaptivity strategies using user-defined metrics. Ongoing and future investigations: Application to compliance minimization in structural optimization. Goal-oriented anisotropic error estimators.

• M. Giacomini

Error estimates for shape optimization

• Dec. 11th 2014 (FreeFem++)

20 / 20