Efficient Simulation of Convection Diffusion Equations Mario - - PowerPoint PPT Presentation

efficient simulation of convection diffusion equations
SMART_READER_LITE
LIVE PREVIEW

Efficient Simulation of Convection Diffusion Equations Mario - - PowerPoint PPT Presentation

Westflische Institut fr Numerische Wilhelms-Universitt und Angewandte Mnster Mathematik Efficient Simulation of Convection Diffusion Equations Mario Ohlberger, Universit at M unster Computational Methods with Applications,


slide-1
SLIDE 1

Westfälische Wilhelms-Universität Münster Institut für Numerische und Angewandte Mathematik mario.ohlberger@uni-muesnter.de www.uni-muenster.de/math/u/ohlberger

Efficient Simulation of Convection Diffusion Equations

Mario Ohlberger, Universit¨ at M¨ unster

Computational Methods with Applications, Harrachov 2007

slide-2
SLIDE 2

Westfälische Wilhelms-Universität Münster Institut für Numerische und Angewandte Mathematik mario.ohlberger@uni-muesnter.de www.uni-muenster.de/math/u/ohlberger

Outline

  • Introduction
  • General concept for obtaining error control
  • Higher order DG for conservation laws
  • DUNE – adaptive and parallel programming
  • Application: Simulation of PEM fuel cells
slide-3
SLIDE 3

Westfälische Wilhelms-Universität Münster Institut für Numerische und Angewandte Mathematik mario.ohlberger@uni-muesnter.de www.uni-muenster.de/math/u/ohlberger

Basic model problem: convection-diffusion equation ∂tu + ∇ · F(u) − ∆D(u) = 0.

accumulation convection diffusion

Special interest: Convection dominated flow (D′ < < |F′|). Goal: A posteriori error control and adaptivity!

slide-4
SLIDE 4

Westfälische Wilhelms-Universität Münster Institut für Numerische und Angewandte Mathematik mario.ohlberger@uni-muesnter.de www.uni-muenster.de/math/u/ohlberger

Goal: A posteriori error estimates and adaptivity

Situation: u exact solution, uh approximate solution.

slide-5
SLIDE 5

Westfälische Wilhelms-Universität Münster Institut für Numerische und Angewandte Mathematik mario.ohlberger@uni-muesnter.de www.uni-muenster.de/math/u/ohlberger

Goal: A posteriori error estimates and adaptivity

Situation: u exact solution, uh approximate solution. First step: A posterirori error estimate.

||u − uh||1 ≤ η(uh).

slide-6
SLIDE 6

Westfälische Wilhelms-Universität Münster Institut für Numerische und Angewandte Mathematik mario.ohlberger@uni-muesnter.de www.uni-muenster.de/math/u/ohlberger

Goal: A posteriori error estimates and adaptivity

Situation: u exact solution, uh approximate solution. First step: A posterirori error estimate.

||u − uh||1 ≤ η(uh).

Second step: Definition of local error indicators.

η(uh) =

  • j

ηj(uh).

slide-7
SLIDE 7

Westfälische Wilhelms-Universität Münster Institut für Numerische und Angewandte Mathematik mario.ohlberger@uni-muesnter.de www.uni-muenster.de/math/u/ohlberger

Goal: A posteriori error estimates and adaptivity

Third step: Equidistribution strategy. Choose local mesh size such that all ηj(uh) are approximately of the same size, and η(uh) ≤ TOL!

slide-8
SLIDE 8

Westfälische Wilhelms-Universität Münster Institut für Numerische und Angewandte Mathematik mario.ohlberger@uni-muesnter.de www.uni-muenster.de/math/u/ohlberger

Goal: A posteriori error estimates and adaptivity

Third step: Equidistribution strategy. Choose local mesh size such that all ηj(uh) are approximately of the same size, and η(uh) ≤ TOL! This is done by the estimate–mark–adapt algorithm: ”estimate” ηj

slide-9
SLIDE 9

Westfälische Wilhelms-Universität Münster Institut für Numerische und Angewandte Mathematik mario.ohlberger@uni-muesnter.de www.uni-muenster.de/math/u/ohlberger

Goal: A posteriori error estimates and adaptivity

Third step: Equidistribution strategy Choose local mesh size such that all ηj(uh) are approximately of the same size, and η(uh) ≤ TOL! This is done by the estimate–mark–adapt algorithm: ”mark”

slide-10
SLIDE 10

Westfälische Wilhelms-Universität Münster Institut für Numerische und Angewandte Mathematik mario.ohlberger@uni-muesnter.de www.uni-muenster.de/math/u/ohlberger

Goal: A posteriori error estimates and adaptivity

Third step: Equidistribution strategy Choose local mesh size such that all ηj(uh) are approximately of the same size, and η(uh) ≤ TOL! This is done by the estimate–mark–adapt algorithm: ”adapt”

slide-11
SLIDE 11

Westfälische Wilhelms-Universität Münster Institut für Numerische und Angewandte Mathematik mario.ohlberger@uni-muesnter.de www.uni-muenster.de/math/u/ohlberger

General concept for obtaining error control ... (The hyperbolic case (D ≡ 0)!)

slide-12
SLIDE 12

Westfälische Wilhelms-Universität Münster Institut für Numerische und Angewandte Mathematik mario.ohlberger@uni-muesnter.de www.uni-muenster.de/math/u/ohlberger

Analytical framework: Entropy weak solution

u is called an entropy weak solution of the conservation law, if u satisfies for all entropy pairs (S, FS), and for all φ ∈ C1

0(I

Rd × I R+, I R+):

  • I

Rd

  • I

R+

(S(u)∂tφ + FS(u) · ∇φ) dt dx +

  • I

Rd

S(u0)φ(x, 0) dx ≥ 0.

Recall that (S, FS) is called an entropy - entropy flux pair, iff S is convex and F ′

S = S′f ′

.

slide-13
SLIDE 13

Westfälische Wilhelms-Universität Münster Institut für Numerische und Angewandte Mathematik mario.ohlberger@uni-muesnter.de www.uni-muenster.de/math/u/ohlberger

General concept for obtaining error control Entropy Residual RS(v) for any given approximation v: RS(v), φ :=

  • I

R2×R+S(v)∂tφ + FS(v) · ∇φ +

  • I

R2S(u0)φ(·, 0).

Fundamental error estimate:

[Eymard, Gallou¨ et, Ghilani, Herbin ’98], [Chainais-Hillairet ’99] Let S(u) := |u − κ| be the Kruzkov entropy. Suppose that for v there exist measures µv ∈ M(I Rd × I R+) and νv ∈ M(I Rd) independent of κ such that RS(v), φ ≥ −(|∂tφ| + |∇φ|, µv + |φ(·, 0)|, νv).

slide-14
SLIDE 14

Westfälische Wilhelms-Universität Münster Institut für Numerische und Angewandte Mathematik mario.ohlberger@uni-muesnter.de www.uni-muenster.de/math/u/ohlberger

General concept for obtaining error control Entropy Residual RS(v) for any given approximation v: RS(v), φ :=

  • I

R2×R+S(v)∂tφ + FS(v) · ∇φ +

  • I

R2S(u0)φ(·, 0).

Fundamental error estimate:

[Eymard, Gallou¨ et, Ghilani, Herbin ’98], [Chainais-Hillairet ’99] Let S(u) := |u − κ| be the Kruzkov entropy. Suppose that for v there exist measures µv ∈ M(I Rd × I R+) and νv ∈ M(I Rd) independent of κ such that RS(v), φ ≥ −(|∂tφ| + |∇φ|, µv + |φ(·, 0)|, νv). Then the following error estimate holds:

||u − v||L1(K) ≤ T(νv(BR+δ(x0)) + C1µv(Dδ) + C2

  • µv(Dδ)).
slide-15
SLIDE 15

Westfälische Wilhelms-Universität Münster Institut für Numerische und Angewandte Mathematik mario.ohlberger@uni-muesnter.de www.uni-muenster.de/math/u/ohlberger

General concept for obtaining error control ||u − v||L1(K) ≤ T(νv(BR+δ(x0)) + C1µv(Dδ) + C2

  • µv(Dδ)).

Cone of dependence Dδ:

for details see also [Kr¨

  • ner, Ohlberger ’00]

x0 R T K slope ω y x t

slide-16
SLIDE 16

Westfälische Wilhelms-Universität Münster Institut für Numerische und Angewandte Mathematik mario.ohlberger@uni-muesnter.de www.uni-muenster.de/math/u/ohlberger

A posteriori results in this context 1) Hyperbolic conservation laws:

1995 Cockburn, Gau: 2000 Kr¨

  • ner, Ohlberger:

2000 Gosse, Makridakis: 2003 K¨ uther, Ohlberger: 2006 Ohlberger, Vovelle: 2007 Dedner, Makridakis, Ohlberger: Finite volume schemes; Relaxation schemes; Central staggered schemes; Boundary value problems; Discontinuous Galerkin;

slide-17
SLIDE 17

Westfälische Wilhelms-Universität Münster Institut für Numerische und Angewandte Mathematik mario.ohlberger@uni-muesnter.de www.uni-muenster.de/math/u/ohlberger

A posteriori results in this context 1) Hyperbolic conservation laws:

1995 Cockburn, Gau: 2000 Kr¨

  • ner, Ohlberger:

2000 Gosse, Makridakis: 2003 K¨ uther, Ohlberger: 2006 Ohlberger, Vovelle: 2007 Dedner, Makridakis, Ohlberger: Finite volume schemes; Relaxation schemes; Central staggered schemes; Boundary value problems; Discontinuous Galerkin;

2) Degenerate parabolic problems:

2001 Ohlberger: 2002 Ohlberger, Rohde: 2002 Herbin, Ohlberger: 2004 Ohlberger: 2004 Chen, Ji: Vertex centered FV; Weakly coupled systems; Cell centered FV Higher order FV; Finite element schemes;

slide-18
SLIDE 18

Westfälische Wilhelms-Universität Münster Institut für Numerische und Angewandte Mathematik mario.ohlberger@uni-muesnter.de www.uni-muenster.de/math/u/ohlberger

A posteriori results in this context 1) Hyperbolic conservation laws:

1995 Cockburn, Gau: 2000 Kr¨

  • ner, Ohlberger:

2000 Gosse, Makridakis: 2003 K¨ uther, Ohlberger: 2006 Ohlberger, Vovelle: 2007 Dedner, Makridakis, Ohlberger: Finite volume schemes; Relaxation schemes; Central staggered scheme; Boundary value problems; Discontinuous Galerkin;

2) Degenerate parabolic problems:

2001 Ohlberger: 2002 Ohlberger, Rohde: 2002 Herbin, Ohlberger: 2004 Ohlberger: 2004 Chen, Ji: Vertex centered FV; Weakly coupled systems; Cell centered FV Higher order FV; Finite element schemes;

slide-19
SLIDE 19

Westfälische Wilhelms-Universität Münster Institut für Numerische und Angewandte Mathematik mario.ohlberger@uni-muesnter.de www.uni-muenster.de/math/u/ohlberger

Error control for Discontinuous Galerkin approximations

  • f nonlinear conservation laws

[Dedner, Makridakis, Ohlberger ’07]

slide-20
SLIDE 20

Westfälische Wilhelms-Universität Münster Institut für Numerische und Angewandte Mathematik mario.ohlberger@uni-muesnter.de www.uni-muenster.de/math/u/ohlberger

DG for nonlinear conservation laws ∂tu + ∇ · f(u) = 0 in Rd × R+ , u(·, 0) = u0 in Rd.

Sought: u ∈ BV (Rd × R+) Given: Nonlinear flux: f ∈ C1(R) Initial data: u0 ∈ BV (Rd)

slide-21
SLIDE 21

Westfälische Wilhelms-Universität Münster Institut für Numerische und Angewandte Mathematik mario.ohlberger@uni-muesnter.de www.uni-muenster.de/math/u/ohlberger

Notation (fixed mesh for all times)

Decomposition of Rd: T with control volumes Tj ∈ T , j ∈ J and faces Sjl = ¯ Tj ∩ ¯ Tl. DG approximation space: V p

h := {vh ∈ BV (Rd)| vh|Tj ∈ Pp ∀j ∈ J}.

slide-22
SLIDE 22

Westfälische Wilhelms-Universität Münster Institut für Numerische und Angewandte Mathematik mario.ohlberger@uni-muesnter.de www.uni-muenster.de/math/u/ohlberger

Notation (fixed mesh for all times)

Decomposition of Rd: T with control volumes Tj ∈ T , j ∈ J and faces Sjl = ¯ Tj ∩ ¯ Tl. DG approximation space: V p

h := {vh ∈ BV (Rd)| vh|Tj ∈ Pp ∀j ∈ J}.

Semi-discrete DG approximation without stabilization

d dt(uj(t), vj)Tj − (f(uj(t)), ∇vj)Tj +

  • l∈N(j)

(fjl(uj(t), ul(t)), vj)Sjl = 0, for all vj ∈ Pp, Tj ∈ T .

slide-23
SLIDE 23

Westfälische Wilhelms-Universität Münster Institut für Numerische und Angewandte Mathematik mario.ohlberger@uni-muesnter.de www.uni-muenster.de/math/u/ohlberger

Notation (adaptive meshes)

Decomposition of [0, T]: I = {t0, ..., tN}, In := (tn, tn+1], ∆tn := |In|. Decomposition of Rd in In : T n with control volumes Tj ∈ T n, j ∈ Jn. DG space in In: V p

h,n := {vh ∈ BV (Rd)| vh|Tj ∈ Pp ∀j ∈ Jn}.

slide-24
SLIDE 24

Westfälische Wilhelms-Universität Münster Institut für Numerische und Angewandte Mathematik mario.ohlberger@uni-muesnter.de www.uni-muenster.de/math/u/ohlberger

Notation (adaptive meshes)

Decomposition of [0, T]: I = {t0, ..., tN}, In := (tn, tn+1], ∆tn := |In|. Decomposition of Rd in In : T n with control volumes Tj ∈ T n, j ∈ Jn. DG space in In: V p

h,n := {vh ∈ BV (Rd)| vh|Tj ∈ Pp ∀j ∈ Jn}.

Limiting projection operators on In

For t ∈ In let Λn,t

h

: V p

h,n → V p h,n

and Λn,tn

h

: V p

h,n−1 → V p h,n

be two projection operators that are mass conservative. For uh ∈ V p

h,n define the projected approximation

uh through:

  • uh(t) = Λn,t

h (uh(t))

for t ∈ In, n = 0, . . . , N − 1 .

slide-25
SLIDE 25

Westfälische Wilhelms-Universität Münster Institut für Numerische und Angewandte Mathematik mario.ohlberger@uni-muesnter.de www.uni-muenster.de/math/u/ohlberger

Semi-discrete DG approximation on adaptive meshes

Set u−1

h := ΠV p

h,0(u0).

For n = 0, ..., N − 1, un

h ∈ C1(In; V p h,n) is defined through

un

h(tn) := Λn,tn h

(un−1

h

(tn)) , d dt(un

j (t), vj)Tj − (f(˜

un

j (t)), ∇vj)Tj +

  • l∈N(j)

(fjl(˜ un

j (t), ˜

un

l (t)), vj)Sjl = 0

for all vj ∈ Pp, j ∈ Jn, t ∈ In . The global approximation uh ∈ L∞(0, T; V p

h,n) is defined through uh|In := un h.

slide-26
SLIDE 26

Westfälische Wilhelms-Universität Münster Institut für Numerische und Angewandte Mathematik mario.ohlberger@uni-muesnter.de www.uni-muenster.de/math/u/ohlberger

A posteriori error estimate [Dedner, Makridakis, Ohlberger ’05]

||(u − uh)(T)||L1(BR(x0)) ≤ ||( uh − uh)(T)||L1(BR(x0)) + η0 +

  • K1η1 +
  • K2η2

with

η0 =

  • j∈J0
  • Tj

|u0 − u0

j(0)|,

η1 =

  • n
  • j∈Jn

tn+1

  • tn

hj Rn

T,j +1

2

tn+1

  • tn
  • l∈N(j)

hjl Rn

S,jl+

hj Rn

Λ,j

  • ,

η2 =

  • n
  • j∈Jn

tn+1

  • tn

|| un

j −

un

j ||∞Rn T,j +1

2

tn+1

  • tn
  • l∈N(j)

max

k∈{j,l} ||

un

k −

un

k||∞Rn S,jl+||

un−1(tn)− un−1(tn)||∞Rn

Λ,j

  • ,

Rn

T,j =

  • Tj
  • ∂t

uj + ∇ · f( uj)

  • ,

Rn

S,jl =

  • Sjl

Qjl| uj − ul|, Rn

Λ,j =

  • Tj

| un(tn) − un−1(tn)| element residual jump residual projection residual

slide-27
SLIDE 27

Westfälische Wilhelms-Universität Münster Institut für Numerische und Angewandte Mathematik mario.ohlberger@uni-muesnter.de www.uni-muenster.de/math/u/ohlberger

Choice of the projection operators Condition from the a posteriori error estimate:

Construct projection operator, such that

|| un

j (·, t) −

un

j (·, t)||L∞(Tj) ≤ λn j (t).

with

λn

j (t) ≤

hj near discontinuities, large in smooth regions.

slide-28
SLIDE 28

Westfälische Wilhelms-Universität Münster Institut für Numerische und Angewandte Mathematik mario.ohlberger@uni-muesnter.de www.uni-muenster.de/math/u/ohlberger

P–adaptive method:

slide-29
SLIDE 29

Westfälische Wilhelms-Universität Münster Institut für Numerische und Angewandte Mathematik mario.ohlberger@uni-muesnter.de www.uni-muenster.de/math/u/ohlberger

P–adaptive method:

slide-30
SLIDE 30

Westfälische Wilhelms-Universität Münster Institut für Numerische und Angewandte Mathematik mario.ohlberger@uni-muesnter.de www.uni-muenster.de/math/u/ohlberger

P–adaptive method:

slide-31
SLIDE 31

Westfälische Wilhelms-Universität Münster Institut für Numerische und Angewandte Mathematik mario.ohlberger@uni-muesnter.de www.uni-muenster.de/math/u/ohlberger

Numerical experiment .....

slide-32
SLIDE 32

Westfälische Wilhelms-Universität Münster Institut für Numerische und Angewandte Mathematik mario.ohlberger@uni-muesnter.de www.uni-muenster.de/math/u/ohlberger

Buckley–Leverett problem ut + ∂xf(u) = 0,

  • n (−1, 1) × (0, 0.4),

u(·, 0) = u0, on (−1, 1),

with f(u) = u2 u2 + 1

2(1 − u)2

and initial data u0(x) :=    1, for x < −0.6, 0, for − 0.6 ≤ x < 0.2, 1, for 0.2 ≤ x.

slide-33
SLIDE 33

Westfälische Wilhelms-Universität Münster Institut für Numerische und Angewandte Mathematik mario.ohlberger@uni-muesnter.de www.uni-muenster.de/math/u/ohlberger

TOL = 0.25 TOL = 0.125 p=1

0.2 0.4 0.6 0.8 1

  • 1
  • 0.5

0.5 1 approximate solution 1 2 p

  • 12
  • 11
  • 8
  • 6
  • 4
  • 2

log2(h) 0.2 0.4 0.6 0.8 1

  • 1
  • 0.5

0.5 1 approximate solution 1 2 p

  • 12
  • 11
  • 8
  • 6
  • 4
  • 2

log2(h)

p=2

0.2 0.4 0.6 0.8 1

  • 1
  • 0.5

0.5 1 approximate solution 1 2 p

  • 12
  • 11
  • 8
  • 6
  • 4
  • 2

log2(h) 0.2 0.4 0.6 0.8 1

  • 1
  • 0.5

0.5 1 approximate solution 1 2 p

  • 12
  • 11
  • 8
  • 6
  • 4
  • 2

log2(h)

slide-34
SLIDE 34

Westfälische Wilhelms-Universität Münster Institut für Numerische und Angewandte Mathematik mario.ohlberger@uni-muesnter.de www.uni-muenster.de/math/u/ohlberger

Results in 2D: Linear advection [Dedner, Ohlberger ’06] ∂tc + ∇ · (bc) = 0, c(·, 0) = c0.

t = 0 t = T

slide-35
SLIDE 35

Westfälische Wilhelms-Universität Münster Institut für Numerische und Angewandte Mathematik mario.ohlberger@uni-muesnter.de www.uni-muenster.de/math/u/ohlberger

Results in 2D: Linear advection

1e-06 1e-05 1e-04 0.001 0.01 0.1 1 10 100 1000 10000 100000 error cpu-time p=0 p=1 p=2

slide-36
SLIDE 36

Westfälische Wilhelms-Universität Münster Institut für Numerische und Angewandte Mathematik mario.ohlberger@uni-muesnter.de www.uni-muenster.de/math/u/ohlberger

Results in 2D: Buckley–Leverett with advection ∂tu + ∂xf(u) + a∂yu = 0, u(·, 0) = u0,

with f(u) =

u2 u2+0.5(1−u)2,

a = 1.3.

Initial data:

0.2 0.4 0.6 0.8 1

  • 0.45
  • 0.25
  • 0.05

0.15 0.35 0.55 0.75 0.95 1.15 v(x0,ζ) ζ

c0 c0(x = 0, ·)

slide-37
SLIDE 37

Westfälische Wilhelms-Universität Münster Institut für Numerische und Angewandte Mathematik mario.ohlberger@uni-muesnter.de www.uni-muenster.de/math/u/ohlberger

Results in 2D: Buckley–Leverett with advection

solution uh polynomial at t = 0.35 degree grid level of refinement

slide-38
SLIDE 38

Westfälische Wilhelms-Universität Münster Institut für Numerische und Angewandte Mathematik mario.ohlberger@uni-muesnter.de www.uni-muenster.de/math/u/ohlberger

Results in 2D: Buckley–Leverett with advection

||(u − uh)(T)||L1(Ω)

0.01 0.1 10 100 1000 10000 100000 p=1, global refine p=2, global refine p=1, local refine p=2, local refine

cpu-time

slide-39
SLIDE 39

Westfälische Wilhelms-Universität Münster Institut für Numerische und Angewandte Mathematik mario.ohlberger@uni-muesnter.de www.uni-muenster.de/math/u/ohlberger

DUNE–Distributed and Unified Numerics Environment

  • Interface concept

= ⇒ efficient reuse of existing software

slide-40
SLIDE 40

Westfälische Wilhelms-Universität Münster Institut für Numerische und Angewandte Mathematik mario.ohlberger@uni-muesnter.de www.uni-muenster.de/math/u/ohlberger

DUNE–Distributed and Unified Numerics Environment

  • Design principles

– Modularity through use of interfaces – Dimension and data structure independend programming – Efficiency through static polymorphism

  • Developer

– P. Bastian, M. Blatt, C. Engwer (Stuttgart), O. Sander (Berlin) – A. Dedner, R. Kl¨

  • fkorn (Freiburg), M. Ohlberger (M¨

unster)

  • Availability

– Homepage http://www.dune-project.org/ – Programming language C++ – Portability via ISO standard conformmity – Open source software ’GNU with linking exceptions’

slide-41
SLIDE 41

Westfälische Wilhelms-Universität Münster Institut für Numerische und Angewandte Mathematik mario.ohlberger@uni-muesnter.de www.uni-muenster.de/math/u/ohlberger

DUNE – Parallel efficiency test with loaad balancing

  • Benchmark: 3D Forward Facing Step

[Burri, Dedner, Kl¨

  • fkorn, Ohlberger ’05]
  • Speedup and efficiency comparison
  • riginal code

K CPU time S4→K E4→K 4 0.0089062 8 0.0046045 1.93424 0.96712 16 0.0023943 3.71978 0.92995 32 0.0012710 7.00712 0.87589 DUNE K CPU time S4→K E4→K 4 0.0101473 8 0.0052195 1.94411 0.97205 16 0.0026859 3.77795 0.94449 32 0.0013971 7.26325 0.90791

slide-42
SLIDE 42

Westfälische Wilhelms-Universität Münster Institut für Numerische und Angewandte Mathematik mario.ohlberger@uni-muesnter.de www.uni-muenster.de/math/u/ohlberger

Application: Simulation of PEM - fuel cells

Joint BMBF–project with:

  • Robert Kl¨
  • fkorn, Dietmar Kr¨
  • ner, AAM, Freiburg

urgen Schumacher, Fraunhofer ISE, Freiburg

  • Willy J¨

ager, Heidelberg

  • Ben Schweizer, Basel
  • Proton Motor Fuel Cell GmbH, Starnberg
  • Freudenberg FCCT OHG, Weinheim
slide-43
SLIDE 43

Westfälische Wilhelms-Universität Münster Institut für Numerische und Angewandte Mathematik mario.ohlberger@uni-muesnter.de www.uni-muenster.de/math/u/ohlberger

PEM fuel cells for small electronic devices

7 cm Prototype Fuel Cell System Powering a Camcorder 9 W m

  • ax. power

Sam

e size and sam e energy than largest rechargeable battery pack

Increase of power and energy

density by m

iniaturization of

functional units

slide-44
SLIDE 44

Westfälische Wilhelms-Universität Münster Institut für Numerische und Angewandte Mathematik mario.ohlberger@uni-muesnter.de www.uni-muenster.de/math/u/ohlberger

PEM fuel cell components

Functional units of PEM fuel cells

Cell frame Diffusion layer Hydrogen

electrode

Ion conducting

membrane

Air electrode H2 2H

+ + 2e−

1/ 2 O2 +2H

+ + 2e− 2H2O

Hydrogen Diffusion layer Air Water

− +

An

  • d

e

Cath

  • d

e

Mem

brane

Flow Field / Bipolar Plate Diffusion layer

slide-45
SLIDE 45

Westfälische Wilhelms-Universität Münster Institut für Numerische und Angewandte Mathematik mario.ohlberger@uni-muesnter.de www.uni-muenster.de/math/u/ohlberger

Sketch of a fuel cell

cathodic catalyst layer: O2 + 4e− + 4H+ → 2H2O anodic catalyst layer: 2H2 → 4H+ + 4e−

membrane cathodic anodic diffusion layer diffusion layer

+ _ +

cathodic anodic gas channel gas channel anodic catalyst layer cathodic catalyst layer free gas flow of O

2 N2 H2O

free gas flow of H

2 H2O

✁✁✁✁ ✁✁✁✁ ✁✁✁✁ ✁✁✁✁ ✁✁✁✁ ✁✁✁✁ ✁✁✁✁ ✁✁✁✁ ✁✁✁✁ ✁✁✁✁ ✁✁✁✁ ✁✁✁✁ ✁✁✁✁ ✁✁✁✁ ✁✁✁✁ ✁✁✁✁ ✁✁✁✁ ✂✁✂✁✂✁✂✁✂ ✂✁✂✁✂✁✂✁✂ ✂✁✂✁✂✁✂✁✂ ✂✁✂✁✂✁✂✁✂ ✂✁✂✁✂✁✂✁✂ ✂✁✂✁✂✁✂✁✂ ✂✁✂✁✂✁✂✁✂ ✂✁✂✁✂✁✂✁✂ ✂✁✂✁✂✁✂✁✂ ✂✁✂✁✂✁✂✁✂ ✂✁✂✁✂✁✂✁✂ ✂✁✂✁✂✁✂✁✂ ✂✁✂✁✂✁✂✁✂ ✂✁✂✁✂✁✂✁✂ ✂✁✂✁✂✁✂✁✂ ✂✁✂✁✂✁✂✁✂ ✂✁✂✁✂✁✂✁✂ ✄☎✄☎✄☎✄ ✄☎✄☎✄☎✄ ✄☎✄☎✄☎✄ ✄☎✄☎✄☎✄ ✄☎✄☎✄☎✄ ✄☎✄☎✄☎✄ ✄☎✄☎✄☎✄ ✄☎✄☎✄☎✄ ✄☎✄☎✄☎✄ ✄☎✄☎✄☎✄ ✄☎✄☎✄☎✄ ✄☎✄☎✄☎✄ ✄☎✄☎✄☎✄ ✄☎✄☎✄☎✄ ✄☎✄☎✄☎✄ ✄☎✄☎✄☎✄ ✄☎✄☎✄☎✄ ✆☎✆☎✆ ✆☎✆☎✆ ✆☎✆☎✆ ✆☎✆☎✆ ✆☎✆☎✆ ✆☎✆☎✆ ✆☎✆☎✆ ✆☎✆☎✆ ✆☎✆☎✆ ✆☎✆☎✆ ✆☎✆☎✆ ✆☎✆☎✆ ✆☎✆☎✆ ✆☎✆☎✆ ✆☎✆☎✆ ✆☎✆☎✆ ✆☎✆☎✆ ✝☎✝☎✝☎✝ ✝☎✝☎✝☎✝ ✝☎✝☎✝☎✝ ✝☎✝☎✝☎✝ ✝☎✝☎✝☎✝ ✝☎✝☎✝☎✝ ✝☎✝☎✝☎✝ ✝☎✝☎✝☎✝ ✝☎✝☎✝☎✝ ✝☎✝☎✝☎✝ ✝☎✝☎✝☎✝ ✝☎✝☎✝☎✝ ✝☎✝☎✝☎✝ ✝☎✝☎✝☎✝ ✝☎✝☎✝☎✝ ✝☎✝☎✝☎✝ ✝☎✝☎✝☎✝ ✞☎✞☎✞ ✞☎✞☎✞ ✞☎✞☎✞ ✞☎✞☎✞ ✞☎✞☎✞ ✞☎✞☎✞ ✞☎✞☎✞ ✞☎✞☎✞ ✞☎✞☎✞ ✞☎✞☎✞ ✞☎✞☎✞ ✞☎✞☎✞ ✞☎✞☎✞ ✞☎✞☎✞ ✞☎✞☎✞ ✞☎✞☎✞ ✞☎✞☎✞ ✟✁✟✁✟✁✟✁✟✁✟✁✟ ✟✁✟✁✟✁✟✁✟✁✟✁✟ ✟✁✟✁✟✁✟✁✟✁✟✁✟ ✟✁✟✁✟✁✟✁✟✁✟✁✟ ✟✁✟✁✟✁✟✁✟✁✟✁✟ ✟✁✟✁✟✁✟✁✟✁✟✁✟ ✟✁✟✁✟✁✟✁✟✁✟✁✟ ✟✁✟✁✟✁✟✁✟✁✟✁✟ ✟✁✟✁✟✁✟✁✟✁✟✁✟ ✟✁✟✁✟✁✟✁✟✁✟✁✟ ✟✁✟✁✟✁✟✁✟✁✟✁✟ ✟✁✟✁✟✁✟✁✟✁✟✁✟ ✟✁✟✁✟✁✟✁✟✁✟✁✟ ✟✁✟✁✟✁✟✁✟✁✟✁✟ ✟✁✟✁✟✁✟✁✟✁✟✁✟ ✟✁✟✁✟✁✟✁✟✁✟✁✟ ✟✁✟✁✟✁✟✁✟✁✟✁✟ ✠✁✠✁✠✁✠✁✠✁✠✁✠ ✠✁✠✁✠✁✠✁✠✁✠✁✠ ✠✁✠✁✠✁✠✁✠✁✠✁✠ ✠✁✠✁✠✁✠✁✠✁✠✁✠ ✠✁✠✁✠✁✠✁✠✁✠✁✠ ✠✁✠✁✠✁✠✁✠✁✠✁✠ ✠✁✠✁✠✁✠✁✠✁✠✁✠ ✠✁✠✁✠✁✠✁✠✁✠✁✠ ✠✁✠✁✠✁✠✁✠✁✠✁✠ ✠✁✠✁✠✁✠✁✠✁✠✁✠ ✠✁✠✁✠✁✠✁✠✁✠✁✠ ✠✁✠✁✠✁✠✁✠✁✠✁✠ ✠✁✠✁✠✁✠✁✠✁✠✁✠ ✠✁✠✁✠✁✠✁✠✁✠✁✠ ✠✁✠✁✠✁✠✁✠✁✠✁✠ ✠✁✠✁✠✁✠✁✠✁✠✁✠ ✠✁✠✁✠✁✠✁✠✁✠✁✠ ✡✁✡✁✡✁✡✁✡✁✡✁✡ ✡✁✡✁✡✁✡✁✡✁✡✁✡ ✡✁✡✁✡✁✡✁✡✁✡✁✡ ✡✁✡✁✡✁✡✁✡✁✡✁✡ ✡✁✡✁✡✁✡✁✡✁✡✁✡ ✡✁✡✁✡✁✡✁✡✁✡✁✡ ✡✁✡✁✡✁✡✁✡✁✡✁✡ ✡✁✡✁✡✁✡✁✡✁✡✁✡ ✡✁✡✁✡✁✡✁✡✁✡✁✡ ✡✁✡✁✡✁✡✁✡✁✡✁✡ ✡✁✡✁✡✁✡✁✡✁✡✁✡ ✡✁✡✁✡✁✡✁✡✁✡✁✡ ✡✁✡✁✡✁✡✁✡✁✡✁✡ ✡✁✡✁✡✁✡✁✡✁✡✁✡ ✡✁✡✁✡✁✡✁✡✁✡✁✡ ✡✁✡✁✡✁✡✁✡✁✡✁✡ ✡✁✡✁✡✁✡✁✡✁✡✁✡ ☛✁☛✁☛✁☛✁☛✁☛✁☛ ☛✁☛✁☛✁☛✁☛✁☛✁☛ ☛✁☛✁☛✁☛✁☛✁☛✁☛ ☛✁☛✁☛✁☛✁☛✁☛✁☛ ☛✁☛✁☛✁☛✁☛✁☛✁☛ ☛✁☛✁☛✁☛✁☛✁☛✁☛ ☛✁☛✁☛✁☛✁☛✁☛✁☛ ☛✁☛✁☛✁☛✁☛✁☛✁☛ ☛✁☛✁☛✁☛✁☛✁☛✁☛ ☛✁☛✁☛✁☛✁☛✁☛✁☛ ☛✁☛✁☛✁☛✁☛✁☛✁☛ ☛✁☛✁☛✁☛✁☛✁☛✁☛ ☛✁☛✁☛✁☛✁☛✁☛✁☛ ☛✁☛✁☛✁☛✁☛✁☛✁☛ ☛✁☛✁☛✁☛✁☛✁☛✁☛ ☛✁☛✁☛✁☛✁☛✁☛✁☛ ☛✁☛✁☛✁☛✁☛✁☛✁☛ ☞✁☞✁☞ ☞✁☞✁☞ ☞✁☞✁☞ ☞✁☞✁☞ ☞✁☞✁☞ ☞✁☞✁☞ ☞✁☞✁☞ ☞✁☞✁☞ ☞✁☞✁☞ ☞✁☞✁☞ ☞✁☞✁☞ ☞✁☞✁☞ ☞✁☞✁☞ ☞✁☞✁☞ ☞✁☞✁☞ ☞✁☞✁☞ ☞✁☞✁☞ ✌✁✌✁✌ ✌✁✌✁✌ ✌✁✌✁✌ ✌✁✌✁✌ ✌✁✌✁✌ ✌✁✌✁✌ ✌✁✌✁✌ ✌✁✌✁✌ ✌✁✌✁✌ ✌✁✌✁✌ ✌✁✌✁✌ ✌✁✌✁✌ ✌✁✌✁✌ ✌✁✌✁✌ ✌✁✌✁✌ ✌✁✌✁✌ ✌✁✌✁✌ ✍✁✍✁✍ ✍✁✍✁✍ ✍✁✍✁✍ ✍✁✍✁✍ ✍✁✍✁✍ ✍✁✍✁✍ ✍✁✍✁✍ ✍✁✍✁✍ ✍✁✍✁✍ ✍✁✍✁✍ ✍✁✍✁✍ ✍✁✍✁✍ ✍✁✍✁✍ ✍✁✍✁✍ ✍✁✍✁✍ ✍✁✍✁✍ ✍✁✍✁✍ ✎✁✎✁✎ ✎✁✎✁✎ ✎✁✎✁✎ ✎✁✎✁✎ ✎✁✎✁✎ ✎✁✎✁✎ ✎✁✎✁✎ ✎✁✎✁✎ ✎✁✎✁✎ ✎✁✎✁✎ ✎✁✎✁✎ ✎✁✎✁✎ ✎✁✎✁✎ ✎✁✎✁✎ ✎✁✎✁✎ ✎✁✎✁✎ ✎✁✎✁✎

H2O

2

H2 H2O H2O O H+ _

slide-46
SLIDE 46

Westfälische Wilhelms-Universität Münster Institut für Numerische und Angewandte Mathematik mario.ohlberger@uni-muesnter.de www.uni-muenster.de/math/u/ohlberger

Modeling concept:

Porous layers:

  • Two phase - multi component flow in porous media with

phase transition

  • Potential flow of electrons and protons

Gas channels:

  • Free multi component gas flow

Coupling through interface conditions For details see:

  • [K¨

uhn, Ohlberger, Schumacher, Ziegler, Kl¨

  • fkorn ’03]
  • [Steinkamp, Schumacher, Goldsmith, Ohlberger, Ziegler ’07]
slide-47
SLIDE 47

Westfälische Wilhelms-Universität Münster Institut für Numerische und Angewandte Mathematik mario.ohlberger@uni-muesnter.de www.uni-muenster.de/math/u/ohlberger

Reduced model problem in three space dimensions

[Kl¨

  • fkorn, Kr¨
  • ner, Ohlberger ’07]

Two phase flow in global pressure formulation (solve for p and u) −∇ · (Kλ(sw)∇p) = 0, u = −Kλ(sw)∇p, Balance equation for liquid water saturation (solve for sw) ∂t(nsw) + ∇ · (fw(sw)(u + λg(sw)K∇pc(sw))) = rphase. Transport of species in the gas phase (solve for c = (cH2O, cO2)T) ∂t(nsgc) + ∇ · (vgc) − ∇ · (nsgDg∇c) = qg.

slide-48
SLIDE 48

Westfälische Wilhelms-Universität Münster Institut für Numerische und Angewandte Mathematik mario.ohlberger@uni-muesnter.de www.uni-muenster.de/math/u/ohlberger

Initial values and boundary conditions

Initial Values sw(., 0) = 0.1, cO2(., 0) = 0.8, cH2O(., 0) = 0.2 Boundary Values

Γ2 Γ3 Γ1 Γ1 Γ4 Γ1

m µ 200 m µ 600

GDL

(x−axis) (y−axis)

Γ1 no flow Γ2 pw = 100500 Pa, sw = 0, cO2 = 0.8, cH2O = 0.2 Γ3 pw = 100000 Pa, sw = 0, cO2 = 0.8, cH2O = 0.2 Γ4 pw no flow, sw = 1, cO2, cH2O no flow

slide-49
SLIDE 49

Westfälische Wilhelms-Universität Münster Institut für Numerische und Angewandte Mathematik mario.ohlberger@uni-muesnter.de www.uni-muenster.de/math/u/ohlberger

Numerical results Pressure, adaptation level and partitioning

pw grid level partitions

slide-50
SLIDE 50

Westfälische Wilhelms-Universität Münster Institut für Numerische und Angewandte Mathematik mario.ohlberger@uni-muesnter.de www.uni-muenster.de/math/u/ohlberger

Velocity components

ux uy uz

slide-51
SLIDE 51

Westfälische Wilhelms-Universität Münster Institut für Numerische und Angewandte Mathematik mario.ohlberger@uni-muesnter.de www.uni-muenster.de/math/u/ohlberger

Saturation, and mass concentrations

sw cH2O cO2

slide-52
SLIDE 52

Westfälische Wilhelms-Universität Münster Institut für Numerische und Angewandte Mathematik mario.ohlberger@uni-muesnter.de www.uni-muenster.de/math/u/ohlberger

Complexity of the DG discretization

  • Complexity per time step:

– Number of elements: ∼ 150.000 – Degrees of freedom per element: 43 – Degrees of freedom per time step: ∼ 6.450.000 – Computational time per time step: ∼ 22 seconds – Linear solver: preconditioned BiCG-stab from the dune-istl library [Blatt, Bastian ’06]

  • Overall complexity:

– Number of time steps: ∼ 15.000 – Computational time: ∼ 3 days on 16 processors (XC4000 linux cluster) – Time evolution with first order ODE solver

slide-53
SLIDE 53

Westfälische Wilhelms-Universität Münster Institut für Numerische und Angewandte Mathematik mario.ohlberger@uni-muesnter.de www.uni-muenster.de/math/u/ohlberger

Thank you for your attention! www.uni-muenster.de/math/u/ohlberger

Software: DUNE www.dune-project.org ALUGrid www.mathematik.uni-freiburg.de/IAM/Research/alugrid GRAPE www.iam.uni-bonn.de/grape