Block-structured Adaptive Mesh Refinement Methods for Conservation - PDF document

Conservation laws Finite volume methods Upwind schemes Meshes and adaptation References Examples Examples Euler equations ∂ρ ∂ t + ∂ � � ρ u n = 0 ∂ x n ∂ + ∂ � � � � ρ u k ρ u k u n + δ kn p = 0 ∂ t ∂ x n ∂ E ∂ t + ∂ � � u n ( E + p ) = 0 ∂ x n with polytrope gas equation of state E − 1 � � p = ( γ − 1) 2 ρ u n u n have structure ∂ t q ( x , t ) + ∇ · f ( q ( x , t )) = 0 Fundamentals: Used schemes and mesh adaptation 8

Conservation laws Finite volume methods Upwind schemes Meshes and adaptation References Examples Examples II Navier-Stokes equations ∂ρ ∂ ` ´ ∂ t + ρ u n = 0 ∂ x n ∂ ∂ ` ´ ` ´ ρ u k u n + δ kn p − τ kn ρ u k + = 0 ∂ t ∂ x n ∂ E ∂ ` ´ u n ( E + p ) + q n − τ nj u j ∂ t + = 0 ∂ x n with stress tensor ` ∂ u n ∂ x k + ∂ u i − 2 3 µ∂ u j ´ τ kn = µ ∂ x j δ in ∂ x n and heat conduction q n = − λ ∂ T ∂ x n have structure ∂ t q ( x , t ) + ∇ · f ( q ( x , t )) + ∇ · h ( q ( x , t ) , ∇ q ( x , t )) = 0 Type can be either hyperbolic or parabolic Fundamentals: Used schemes and mesh adaptation 9

Conservation laws Finite volume methods Upwind schemes Meshes and adaptation References Outline Conservation laws Mathematical background Examples Finite volume methods Basics of finite difference methods Splitting methods, second derivatives Upwind schemes Flux-difference splitting Flux-vector splitting High-resolution methods Meshes and adaptation Elements of adaptive algorithms Adaptivity on unstructured meshes Structured mesh refinement techniques Fundamentals: Used schemes and mesh adaptation 10

Conservation laws Finite volume methods Upwind schemes Meshes and adaptation References Basics of finite difference methods Derivation Assume ∂ t q + ∂ x f ( q ) + ∂ x h ( q ( · , ∂ x q )) = s ( q ) Time discretization t n = n ∆ t , discrete volumes I j = [ x j − 1 2 ∆ x , x j + 1 2 ∆ x [=: [ x j − 1 / 2 , x j +1 / 2 [ Using approximations Q j ( t ) ≈ 1 Z s ( Q j ( t )) ≈ 1 Z q ( x , t ) dx , s ( q ( x , t )) dx | I j | | I j | Ij Ij and numerical fluxes ` ´ ` ´ F Q j ( t ) , Q j +1 ( t ) ≈ f ( q ( x j +1 / 2 , t )) , H Q j ( t ) , Q j +1 ( t ) ≈ h ( q ( x j +1 / 2 , t ) , ∇ q ( x j +1 / 2 , t )) yields after integration (Gauss theorem) tn +1 1 Z Q j ( t n +1 ) = Q j ( t n ) − [ F ( Q j ( t ) , Q j +1 ( t )) − F ( Q j − 1 ( t ) , Q j ( t ))] dt − ∆ x tn tn +1 tn +1 1 Z Z [ H ( Q j ( t ) , Q j +1 ( t )) − H ( Q j − 1 ( t ) , Q j ( t ))] dt + s ( Q j ( t )) dt ∆ x tn tn For instance: j − ∆ t h “ ” “ ”i Q n +1 = Q n Q n j , Q n +1 Q n j − 1 , Q n F − F − j j +1 j ∆ x ∆ t h “ ” “ ”i Q n j , Q n Q n j − 1 , Q n + ∆ t s ( Q n H − H j ) dt j +1 j ∆ x Fundamentals: Used schemes and mesh adaptation 11

Conservation laws Finite volume methods Upwind schemes Meshes and adaptation References Basics of finite difference methods Some classical definitions (2 s + 1)-point difference scheme of the form Q n +1 = H (∆ t ) ( Q n j − s , . . . , Q n j + s ) j Definition (Stability) For each time τ there is a constant C S and a value n 0 ∈ N such that �H (∆ t ) ( Q n ) � ≤ C S for all n ∆ t ≤ τ , n < n 0 Definition (Consistency) If the local truncation error 1 h i L (∆ t ) ( x , t ) := q ( x , t + ∆ t ) − H (∆ t ) ( q ( · , t )) ∆ t satisfies �L (∆ t ) ( · , t ) � → 0 as ∆ t → 0 Definition (Convergence) If the global error E (∆ t ) ( x , t ) := Q ( x , t ) − q ( x , t ) satisfies �E (∆ t ) ( · , t ) � → 0 as ∆ t → 0 for all admissible initial data q 0 ( x ) Fundamentals: Used schemes and mesh adaptation 12

Conservation laws Finite volume methods Upwind schemes Meshes and adaptation References Basics of finite difference methods Some classical definitions II Definition (Order of accuracy) H ( · ) is accurate of order o if for all sufficiently smooth initial data q 0 ( x ), there is a constant C L , such that the local truncation error satisfies �L (∆ t ) ( · , t ) � ≤ C L ∆ t o for all ∆ t < ∆ t 0 , t ≤ τ Definition (Conservative form) If H ( · ) can be written in the form j − ∆ t Q n +1 = Q n F ( Q n j − s +1 , . . . , Q n j + s ) − F ( Q n j − s , . . . , Q n ` ´ j + s − 1 ) j ∆ x A conservative scheme satisfies X Q n +1 X Q n = j j j ∈ Z j ∈ Z Definition (Consistency of a conservative method) If the numerical flux satisfies F ( q , . . . , q ) = f ( q ) for all q ∈ S Fundamentals: Used schemes and mesh adaptation 13

Conservation laws Finite volume methods Upwind schemes Meshes and adaptation References Splitting methods, second derivatives Splitting methods Solve homogeneous PDE and ODE successively! H (∆ t ) : ∆ t ⇒ ˜ ∂ t q + ∇ · f ( q ) = 0 , IC: Q ( t m ) = Q S (∆ t ) : ∆ t IC: ˜ ∂ t q = s ( q ) , Q = ⇒ Q ( t m + ∆ t ) 1st-order Godunov splitting: Q ( t m + ∆ t ) = S (∆ t ) H (∆ t ) ( Q ( t m )), 2nd-order Strang splitting : Q ( t m + ∆ t ) = S ( 1 2 ∆ t ) H (∆ t ) S ( 1 2 ∆ t ) ( Q ( t m )) 1st-order dimensional splitting for H ( · ) : ∆ t X (∆ t ) ˜ Q 1 / 2 : ∂ t q + ∂ x 1 f 1 ( q ) = 0 , IC: Q ( t m ) = ⇒ 1 X (∆ t ) ∆ t IC: ˜ ˜ Q 1 / 2 : ∂ t q + ∂ x 2 f 2 ( q ) = 0 , = ⇒ Q 2 [Toro, 1999] Fundamentals: Used schemes and mesh adaptation 14

Conservation laws Finite volume methods Upwind schemes Meshes and adaptation References Splitting methods, second derivatives Conservative scheme for diffusion equation Consider ∂ t q + c ∆ q = 0, which is readily discretized as jk − c ∆ t − c ∆ t “ ” “ ” Q n +1 = Q n Q n j +1 , k − 2 Q n jk + Q n Q n j , k +1 − 2 Q n jk + Q n jk j − 1 , k j , k − 1 ∆ x 2 ∆ x 2 1 2 or conservatively jk − c ∆ t „ « − c ∆ t „ « Q n +1 = Q n H 1 2 , k − H 1 H 2 2 − H 2 j + 1 j − 1 j , k + 1 j , k − 1 jk 2 , k ∆ x 1 ∆ x 1 2 Von Neumann stability analysis: Insert single eigenmode ˆ Q ( t ) e ik 1 x 1 e ik 2 x 2 into discretization Q n +1 = ˆ Q n + ˆ Q n + ˆ “ ” “ ” ˆ Q n e ik 1 ∆ x 1 − 2 ˆ ˆ Q n e ik 2 ∆ x 2 − 2 ˆ ˆ Q n − C 1 Q n e − ik 1 ∆ x 1 Q n e − ik 2 ∆ x 2 − C 2 with C ι = − c ∆ t ι which gives after inserting e ik ι x ι = cos( k ι x ι ) + i sin( k ι x ι ) ∆ x 2 Q n +1 = ˆ Q n (1 + 2 C 1 (cos( k 1 ∆ x 1 ) − 1) + 2 C 2 (cos( k 2 ∆ x 2 ) − 1)) ˆ Stability requires | 1 + 2 C 1 (cos( k 1 ∆ x 1 ) − 1) + 2 C 2 (cos( k 2 ∆ x 2 ) − 1) | ≤ 1 i.e. | 1 − 4 C 1 − 4 C 2 | ≤ 1 from which we derive the stability condition „ ∆ t + ∆ t « ≤ 1 0 ≤ − c ∆ x 1 ∆ x 2 2 Fundamentals: Used schemes and mesh adaptation 15

Conservation laws Finite volume methods Upwind schemes Meshes and adaptation References Flux-difference splitting Linear upwind schemes t . . . . . M − 1 M X X β 1 r 1 + δ m r m β m r m + δ M r M m =2 m =1 Consider Riemann problem M M ∂ t q ( x , t )+ A ∂ ∂ X X q L = δ m r m q R = β m r m ∂ x q ( x , t ) = 0 , x ∈ R , t > 0 m =1 m =1 x Has exact solution 0 X X X X q ( x , t ) = q L + a m r m = q R − a m r m = δ m r m + β m r m λ m < x / t λ m ≥ x / t λ m ≥ x / t λ m < x / t Use Riemann problem to evaluate numerical flux F ( q L , q R ) := f ( q (0 , t )) = Aq (0 , t ) as X X X X F ( q L , q R ) = Aq L + a m λ m r m = Aq R − a m λ m r m = δ m λ m r m + β m λ m r m λ m < 0 λ m ≥ 0 λ m ≥ 0 λ m < 0 λ − λ + Use m = max( λ m , 0) , m = min( λ m , 0) Λ + := diag( λ + Λ − := diag( λ − 1 , . . . , λ + 1 , . . . , λ − to define M ) , M ) A + := R Λ + R − 1 , A − := R Λ − R − 1 and which gives F ( q L , q R ) = Aq L + A − ∆ q = Aq R − A + ∆ q = A + q L + A − q R Fundamentals: Used schemes and mesh adaptation 16

Conservation laws Finite volume methods Upwind schemes Meshes and adaptation References Flux-difference splitting Flux difference splitting Godunov-type scheme with ∆ Q n j +1 / 2 = Q n j +1 − Q n j j − ∆ t � � Q n +1 = Q n A − ∆ Q n j +1 / 2 + A + ∆ Q n j j − 1 / 2 ∆ x Use linearization ¯ q ) = ˆ f (¯ A ( q L , q R )¯ q and construct scheme for nonlinear problem as j − ∆ t � � ˆ 2 + ˆ Q n +1 = Q n A − ( Q n j , Q n j +1 )∆ Q n A + ( Q n j − 1 , Q n j )∆ Q n j + 1 j − 1 j ∆ x 2 stability condition 2 | ∆ t j ∈ Z | ˆ max λ m , j + 1 ∆ x ≤ 1 , for all m = 1 , . . . , M [LeVeque, 1992] Fundamentals: Used schemes and mesh adaptation 17

Conservation laws Finite volume methods Upwind schemes Meshes and adaptation References Flux-difference splitting Roe’s approximate Riemann solver t n +1 Choosing ˆ A ( q L , q R ) [Roe, 1981]: ˆ (i) A ( q L , q R ) has real eigenvalues A ( q L , q R ) → ∂ f ( q ) ˆ as q L , q R → q (ii) ∂ q q l q r ˆ A ( q L , q R )∆ q = f ( q R ) − f ( q L ) (iii) t n For Euler equations: √ ρ L ρ R + √ ρ R ρ L √ ρ L v L + √ ρ R v R = √ ρ L ρ R ρ = ˆ and v = ˆ for v = u n , H √ ρ L + √ ρ R √ ρ L + √ ρ R X Wave decomposition: ∆ q = q r − q l = a m ˆ r m m X ˆ X ˆ F ( q L , q R ) = f ( q L ) − λ m a m ˆ r m = f ( q R ) + λ m a m ˆ r m ˆ ˆ λ m > 0 λ m < 0 ! = 1 X | ˆ f ( q L ) + f ( q R ) − λ m | a m ˆ r m 2 m Fundamentals: Used schemes and mesh adaptation 18

Conservation laws Finite volume methods Upwind schemes Meshes and adaptation References Flux-difference splitting Harten-Lax-Van Leer (HLL) approximate Riemann solver S L t n +1 S R t n +1 t n +1 q ⋆ q L q R t n 8 f ( q L ) , 0 < s 1 , > > > s 3 f ( q L ) − s 1 f ( q R ) + s 1 s 3 ( q R − q L ) < F HLL ( q L , q R ) = s 1 ≤ 0 ≤ s 3 , , s 3 − s 1 > > > f ( q R ) , 0 > s 3 , : s 1 = min( u 1 , L − c L , u 1 , R − c R ) , s 3 = max( u 1 , L + c l , u 1 , R + c R ) 8 q L , x < s 1 t < q ⋆ , ¯ q ( x , t ) = s 1 t ≤ x ≤ s 3 t q R , x > s 3 t : [Toro, 1999], HLLC: [Toro et al., 1994] Fundamentals: Used schemes and mesh adaptation 19

Conservation laws Finite volume methods Upwind schemes Meshes and adaptation References Flux-vector splitting Flux vector splitting Splitting F ( q L , q R ) = f + ( q L ) + f − ( q R ) f ( q ) = f + ( q ) + f − ( q ) t l +1 derived under restriction ˆ f + ( q L ) f − ( q R ) f + ( q R ) λ + f − ( q L ) m ≥ 0 and ˆ λ − m ≤ 0 for all m = 1 , . . . , M for A + ( q ) = ∂ f + ( q ) A − ( q ) = ∂ f − ( q ) t l ˆ ˆ , q L q R ∂ q ∂ q plus reproduction of regular upwinding f + ( q ) f − ( q ) λ m ≥ 0 = f ( q ) , = 0 if for all m = 1 , . . . , M f + ( q ) f − ( q ) λ m ≤ 0 = 0 , = f ( q ) if for all m = 1 , . . . , M Then use F ( q L , q R ) = f + ( q L ) + f − ( q R ) Fundamentals: Used schemes and mesh adaptation 20

Conservation laws Finite volume methods Upwind schemes Meshes and adaptation References Flux-vector splitting Steger-Warming Steger-Warming: Required f ( q ) = A ( q ) q m = 1 m = 1 λ + λ − 2 ( λ m + | λ m | ) 2 ( λ m − | λ m | ) A + ( q ) := R ( q ) Λ + ( q ) R − 1 ( q ) , A − ( q ) := R ( q ) Λ − ( q ) R − 1 ( q ) Gives f ( q ) = A + ( q ) q + A − ( q ) q and the numerical flux F ( q L , q R ) = A + ( q L ) q L + A − ( q R ) q R Jacobians of the split fluxes are identical to A ± ( q ) only in linear case A ± ( q ) q ` ´ ∂ f ± ( q ) = A ± ( q ) + ∂ A ± ( q ) = ∂ q ∂ q ∂ q ∂ q Further methods: Van Leer FVS: AUSM: [Wada and Liou, 1997] Fundamentals: Used schemes and mesh adaptation 21

Conservation laws Finite volume methods Upwind schemes Meshes and adaptation References High-resolution methods High-resolution methods Objective: Higher-order accuracy in smooth solution regions but no spurious oscillations near large gradients Consistent monotone methods converge toward the entropy solution, but Theorem A monotone method is at most first order accurate. Proof: [Harten et al., 1976] Definition (TVD property)) Scheme H (∆ t ) ( Q n ; j ) TVD if TV ( Q l +1 ) ≤ TV ( Q l ) is satisfied for all discrete sequences Q n . Herein, TV ( Q l ) := P j ∈ Z | Q l j +1 − Q l j | . TVD schemes: no new extrema, local minima are non-decreasing, local maxima are non-increasing (termed monotonicity-preserving ). Monotonicity-preserving higher-order schemes are at least 5-point methods. Proofs: [Harten, 1983] TVD concept is proven [Godlewski and Raviart, 1996] for scalar schemes only but nevertheless used to construct high resolution schemes. Note: Monotonicity-preserving scheme can converge toward non-physical weak solutions. Fundamentals: Used schemes and mesh adaptation 22

Conservation laws Finite volume methods Upwind schemes Meshes and adaptation References High-resolution methods MUSCL slope limiting Monotone Upwind Schemes for Conservation Laws [van Leer, 1979] j + 1 » – L + − ˜ 2 = Q n Q (1 − ω ) Φ 2 ∆ j − 1 2 + (1 + ω ) Φ 2 ∆ j + 1 , j + 1 j − 1 j + 1 4 2 j − 1 » – R − + ˜ 2 = Q n Q (1 − ω ) Φ 2 ∆ j + 1 2 + (1 + ω ) Φ 2 ∆ j − 1 j − 1 j + 1 j − 1 4 2 with ∆ j − 1 / 2 = Q n j − Q n j − 1 , ∆ j +1 / 2 = Q n j +1 − Q n j . ∆ j + 1 ∆ j − 1 „ « „ « + − r + r − r + r − 2 2 Φ 2 := Φ , Φ 2 := Φ with := , := j − 1 j + 1 j − 1 j + 1 j − 1 j + 1 ∆ j − 1 ∆ j + 1 2 2 2 2 2 2 and slope limiters , e.g., Minmod Φ( r ) = max(0 , min( r , 1)) Using a midpoint rule for temporal integration, e.g., j − 1 ∆ t “ ” Q ⋆ j = Q n F ( Q n j +1 , Q n j ) − F ( Q n j , Q n j − 1 ) 2 ∆ x and constructing limited values from Q ⋆ to be used in FV scheme gives a TVD method if „ 1 1 » «– (1 − ω )Φ( r ) + (1 + ω ) r Φ < min(2 , 2 r ) 2 r is satisfied for r > 0. Proof: [Hirsch, 1988] Fundamentals: Used schemes and mesh adaptation 23

Conservation laws Finite volume methods Upwind schemes Meshes and adaptation References High-resolution methods Wave Propagation with flux limiting Wave Propagation Method [LeVeque, 1997] is built on the flux differencing approach A ± ∆ := ˆ A ± ( q L , q R )∆ q and the waves W m := a m ˆ r m , i.e. X ˆ X ˆ A − ∆ q = A + ∆ q = λ m W m , λ m W m ˆ ˆ λ m < 0 λ m ≥ 0 Wave Propagation 1D: j − ∆ t − ∆ t Q n +1 = Q n “ ” “ ” A − ∆ j + 1 2 + A + ∆ j − 1 ˜ 2 − ˜ F j + 1 F j − 1 ∆ x ∆ x 2 2 with M 2 = 1 „ 1 − ∆ t « 2 = 1 „ 1 − ∆ t « ˜ X | ˆ λ m | ˆ λ m 2 | ˜ W m F j + 1 2 |A| ∆ x |A| ∆ j + 1 2 | j + 1 j + 1 j + 1 2 ∆ x 2 m =1 and wave limiter W m ˜ 2 = Φ(Θ m 2 ) W m j + 1 j + 1 j + 1 2 with 8 ˆ a m / a m λ m , ≥ 0 , j − 1 j + 1 j + 1 < Θ m 2 = 2 2 2 j + 1 ˆ a m / a m λ m , < 0 j + 3 j + 1 j + 1 : 2 2 2 Fundamentals: Used schemes and mesh adaptation 24

Conservation laws Finite volume methods Upwind schemes Meshes and adaptation References High-resolution methods Wave Propagation Method in 2D Writing ˜ A ± ∆ j ± 1 / 2 := A + ∆ j ± 1 / 2 + ˜ F j ± 1 / 2 one can develop a truly two-dimensional one-step method [Langseth and LeVeque, 2000] jk − ∆ t „ 2 , k − 1 ∆ t h A − ˜ 2 + A − ˜ i Q n +1 = Q n A − ∆ j + 1 ˜ B − ∆ j +1 , k + 1 B + ∆ j +1 , k − 1 + jk ∆ x 1 2 ∆ x 2 2 2 , k − 1 ∆ t i« h A + ˜ 2 + A + ˜ A + ∆ j − 1 ˜ B − ∆ j − 1 , k + 1 B + ∆ j − 1 , k − 1 2 ∆ x 2 2 − ∆ t „ 2 − 1 ∆ t B − ˜ 2 , k +1 + B − ˜ h i B − ∆ j , k + 1 ˜ A − ∆ j + 1 A + ∆ j − 1 + 2 , k +1 ∆ x 2 2 ∆ x 1 2 − 1 ∆ t i« h B + ˜ 2 , k − 1 + B + ˜ B + ∆ j , k − 1 ˜ A − ∆ j + 1 A + ∆ j − 1 2 , k − 1 2 ∆ x 1 that is stable for  2 | ∆ t 2 | ∆ t ff j ∈ Z | ˆ k ∈ Z | ˆ max λ m , j + 1 , max λ m , k + 1 ≤ 1 , for all m = 1 , . . . , M ∆ x 1 ∆ x 2 Fundamentals: Used schemes and mesh adaptation 25

Conservation laws Finite volume methods Upwind schemes Meshes and adaptation References High-resolution methods Further high-resolution methods Some further high-resolution methods (good overview in [Laney, 1998]): ◮ FCT: 2nd order [Oran and Boris, 2001] ◮ ENO/WENO: 3rd order [Shu, 97] ◮ PPM: 3rd order [Colella and Woodward, 1984] 3rd order methods must make use of strong-stability preserving Runge-Kutta methods [Gottlieb et al., 2001] for time integration that use a multi-step update ∆ t “ ” Q υ ˜ j = α υ Q n j + β υ ˜ Q υ − 1 2 ( ˜ Q υ − 1 ) − F j − 1 2 ( ˜ Q υ − 1 ) + γ υ F j + 1 j ∆ x Q 0 := Q n , α 1 = 1, β 1 = 0; and Q n +1 := ˜ Q Υ after final stage Υ with ˜ Typical storage-efficient SSPRK(3,3): Q 2 = 3 4 Q n + 1 Q 1 + 1 Q 1 = Q n + ∆ t F ( Q n ) , ˜ ˜ ˜ 4∆ t F ( ˜ Q 1 ) , 4 Q n +1 = 1 3 Q n + 2 Q 2 + 2 ˜ 3∆ t F ( ˜ Q 2 ) 3 Fundamentals: Used schemes and mesh adaptation 26

Conservation laws Finite volume methods Upwind schemes Meshes and adaptation References Outline Conservation laws Mathematical background Examples Finite volume methods Basics of finite difference methods Splitting methods, second derivatives Upwind schemes Flux-difference splitting Flux-vector splitting High-resolution methods Meshes and adaptation Elements of adaptive algorithms Adaptivity on unstructured meshes Structured mesh refinement techniques Fundamentals: Used schemes and mesh adaptation 27

Conservation laws Finite volume methods Upwind schemes Meshes and adaptation References Elements of adaptive algorithms Elements of adaptive algorithms ◮ Base grid ◮ Solver ◮ Error indicators ◮ Grid manipulation ◮ Interpolation (restriction and prolongation) ◮ Load-balancing Fundamentals: Used schemes and mesh adaptation 28

Conservation laws Finite volume methods Upwind schemes Meshes and adaptation References Adaptivity on unstructured meshes Adaptivity on unstructured meshes ◮ Coarse cells replaced by finer ones ◮ Global time-step ◮ Cell-based data structures ◮ Neighborhoods have to stored + Geometric flexible + No hanging nodes + Easy to implement - Higher order difficult to achieve - Cell aspect ratio must be considered - Fragmented data - Cache-reuse / vectorizaton nearly impossible - Complex load-balancing - Complex synchronization Fundamentals: Used schemes and mesh adaptation 29

Conservation laws Finite volume methods Upwind schemes Meshes and adaptation References Structured mesh refinement techniques Structured mesh refinement techniques ◮ Block-based data of equal size 0 4 3 ◮ Block stored in a quad-tree 1 2 3 4 ◮ Time-step refinement 12 11 8 5 6 7 8 9 10 2 ◮ Global index coordinate system 5 6 9 10 11 12 ◮ Neighborhoods need not be stored + Numerical scheme only for single regular block necessary + Easy to implement + Simple load-balancing + Parent/Child relations according to tree +/- Cache-reuse / vectorization only in Wasted boundary space in a quad-tree data block Fundamentals: Used schemes and mesh adaptation 30

Conservation laws Finite volume methods Upwind schemes Meshes and adaptation References Structured mesh refinement techniques Block-structured adaptive mesh refinement (SAMR) ◮ Refined block overlay coarser ones ◮ Time-step refinement ◮ Block (aka patch) based data structures G 2 , 1 ◮ Global index coordinate system G 2 , 2 + Numerical scheme only for single G 1 , 3 regular block necessary G + Efficient cache-reuse / vectorization 1 , 2 possible G 1 , 1 + Simple load-balancing + Minimal synchronization overhead - Cells without mark are refined G 0 , 1 - Cluster-algorithm necessary - Hanging nodes unavoidable - Difficult to implement Fundamentals: Used schemes and mesh adaptation 31

Conservation laws Finite volume methods Upwind schemes Meshes and adaptation References References References I [Colella and Woodward, 1984] Colella, P. and Woodward, P. (1984). The piecewise parabolic method (PPM) for gas-dynamical simulations. J. Comput. Phys. , 54:174–201. [Godlewski and Raviart, 1996] Godlewski, E. and Raviart, P.-A. (1996). Numerical approximation of hyperbolic systems of conservation laws . Springer Verlag, New York. [Gottlieb et al., 2001] Gottlieb, S., Shu, C.-W., and Tadmor, E. (2001). Strong stability-preserving high-order time discretization methods. SIAM Review , 43(1):89–112. [Harten, 1983] Harten, A. (1983). High resolution schemes for hyperbolic conservation laws. J. Comput. Phys. , 49:357–393. [Harten et al., 1976] Harten, A., Hyman, J. M., and Lax, P. D. (1976). On finite-difference approximations and entropy conditions for shocks. Comm. Pure Appl. Math. , 29:297–322. [Hirsch, 1988] Hirsch, C. (1988). Numerical computation of internal and external flows . John Wiley & Sons, Chichester. Fundamentals: Used schemes and mesh adaptation 32

Conservation laws Finite volume methods Upwind schemes Meshes and adaptation References References References II [Kr¨ oner, 1997] Kr¨ oner, D. (1997). Numerical schemes for conservation laws . John Wiley & Sons and B. G. Teubner, New York, Leipzig. [Laney, 1998] Laney, C. B. (1998). Computational gasdynamics . Cambridge University Press, Cambridge. [Langseth and LeVeque, 2000] Langseth, J. and LeVeque, R. (2000). A wave propagation method for three dimensional conservation laws. J. Comput. Phys. , 165:126–166. [LeVeque, 1992] LeVeque, R. J. (1992). Numerical methods for conservation laws . Birkh¨ auser, Basel. [LeVeque, 1997] LeVeque, R. J. (1997). Wave propagation algorithms for multidimensional hyperbolic systems. J. Comput. Phys. , 131(2):327–353. [Majda, 1984] Majda, A. (1984). Compressible fluid flow and systems of conservation laws in several space variables . Applied Mathematical Sciences Vol. 53. Springer-Verlag, New York. [Oran and Boris, 2001] Oran, E. S. and Boris, J. P. (2001). Numerical simulation of reactive flow . Cambridge Univ. Press, Cambridge, 2nd edition. Fundamentals: Used schemes and mesh adaptation 33

Conservation laws Finite volume methods Upwind schemes Meshes and adaptation References References References III [Roe, 1981] Roe, P. L. (1981). Approximate Riemann solvers, parameter vectors and difference schemes. J. Comput. Phys. , 43:357–372. [Shu, 97] Shu, C.-W. (97). Essentially non-oscillatory and weigted essentially non-oscillatory schemes for hyperbolic conservation laws. Technical Report CR-97-206253, NASA. [Toro, 1999] Toro, E. F. (1999). Riemann solvers and numerical methods for fluid dynamics . Springer-Verlag, Berlin, Heidelberg, 2nd edition. [Toro et al., 1994] Toro, E. F., Spruce, M., and Speares, W. (1994). Restoration of the contact surface in the HLL-Riemann solver. Shock Waves , 4:25–34. [van Leer, 1979] van Leer, B. (1979). Towards the ultimate conservative difference scheme V. A second order sequel to Godunov’s method. J. Comput. Phys. , 32:101–136. [Wada and Liou, 1997] Wada, Y. and Liou, M.-S. (1997). An accurate and robust flux splitting scheme for shock and contact discontinuities. SIAM J. Sci. Comp. , 18(3):633–657. Fundamentals: Used schemes and mesh adaptation 34

Serial SAMR method Parallel SAMR method Examples References Lecture 2 The SAMR method for hyperbolic problems Course Block-structured Adaptive Mesh Refinement Methods for Conservation Laws Theory, Implementation and Application Ralf Deiterding Computer Science and Mathematics Division Oak Ridge National Laboratory P.O. Box 2008 MS6367, Oak Ridge, TN 37831, USA E-mail: deiterdingr@ornl.gov The SAMR method for hyperbolic problems 1

Serial SAMR method Parallel SAMR method Examples References Outline The serial Berger-Colella SAMR method Block-based data structures Numerical update Conservative flux correction Level transfer operators The basic recursive algorithm Cluster algorithm Refinement criteria Parallel SAMR method Domain decomposition A parallel SAMR algorithm Partitioning Examples Euler equations The SAMR method for hyperbolic problems 2

Serial SAMR method Parallel SAMR method Examples References Outline The serial Berger-Colella SAMR method Block-based data structures Numerical update Conservative flux correction Level transfer operators The basic recursive algorithm Cluster algorithm Refinement criteria Parallel SAMR method Domain decomposition A parallel SAMR algorithm Partitioning Examples Euler equations The SAMR method for hyperbolic problems 3

Serial SAMR method Parallel SAMR method Examples References Block-based data structures The m th refinement grid on level l µ 1 Notations: Complete grid ◮ Boundary: ∂ G l , m µ 2 with ghost ◮ Hull: cells - G σ ¯ l , m G l , m = G l , m ∪ ∂ G l , m ◮ Ghost cell region: ˜ l , m \ ¯ G σ l , m = G σ G l , m Interior grid with buffer cells - G l , m The SAMR method for hyperbolic problems 4

Serial SAMR method Parallel SAMR method Examples References Block-based data structures Refinement data ◮ Resolution: ∆ t l := ∆ t l − 1 and ∆ x n , l := ∆ x n , l − 1 r l r l ◮ Refinement factor: r l ∈ N , r l ≥ 2 for l > 0 and r 0 = 1 ◮ Integer coordinate system for internal organization [Bell et al., 1994]: l max ∆ x n , l ∼ Y r κ = κ = l +1 ◮ Computational Domain: G 0 = S M 0 m =1 G 0 , m ◮ Domain of level l : G l := S M l m =1 G l , m with G l , m ∩ G l , n = ∅ for m � = n ◮ Refinements are properly nested: G 1 l ⊂ G l − 1 ◮ Assume a FD scheme with stencil radius s . Necessary data: ◮ Vector of state: Q l := S m Q ( G s l , m ) ◮ Numerical fluxes: F n , l := S m F n (¯ G l , m ) ◮ Flux corrections: δ F n , l := S m δ F n ( ∂ G l , m ) The SAMR method for hyperbolic problems 5

Serial SAMR method Parallel SAMR method Examples References Numerical update Setting of ghost cells Synchronization with G l - ˜ l , m = ˜ S s G s l , m ∩ G l Physical boundary conditions - ˜ l , m = ˜ P s G s l , m \ G 0 Interpolation from G l − 1 - ˜ l , m = ˜ l , m \ (˜ l , m ∪ ˜ I s G s S s P s l , m ) The SAMR method for hyperbolic problems 6

Serial SAMR method Parallel SAMR method Examples References Numerical update Numerical update Time-explicit conservative finite volume scheme H (∆ t ) : Q jk ( t +∆ t ) = Q jk ( t ) − ∆ t − ∆ t “ ” “ ” F 1 2 , k − F 1 F 2 2 − F 2 j + 1 j − 1 j , k + 1 j , k − 1 2 , k ∆ x 1 ∆ x 2 2 UpdateLevel( l ) For all m = 1 To M l Do H (∆ tl ) → Q ( G l , m , t + ∆ t l ) , F n (¯ Q ( G s l , m , t ) − G l , m , t ) If level l > 0 Add F n ( ∂ G l , m , t ) to δ F n , l If level l + 1 exists Init δ F n , l +1 with F n (¯ G l , m ∩ ∂ G l +1 , t ) The SAMR method for hyperbolic problems 7

Serial SAMR method Parallel SAMR method Examples References Numerical update Conservative flux correction Example: Cell j , k 0 1 r l +1 − 1 r l +1 − 1 jk ( t ) − ∆ t l 1 Q l ˇ jk ( t + ∆ t l ) = Q l @ F 1 , l X X F 1 , l +1 2 , k − 2 , w + ι ( t + κ ∆ t l +1 ) A j + 1 v + 1 r 2 ∆ x 1 , l l +1 κ =0 ι =0 − ∆ t l “ ” F 2 , l 2 − F 2 , l j , k + 1 j , k − 1 ∆ x 2 , l 2 Correction pass: 1. δ F 1 , l +1 2 , k := − F 1 , l j − 1 j − 1 2 , k w r l +1 − 1 1 2. δ F 1 , l +1 2 , k := δ F 1 , l +1 X F 1 , l +1 2 , k + 2 , w + ι ( t + κ ∆ t l +1 ) j − 1 j − 1 v + 1 r 2 l +1 ι =0 v v +1 jk ( t + ∆ t l ) + ∆ t l 3. ˇ Q l jk ( t + ∆ t l ) := Q l ∆ x 1 , l δ F 1 , l +1 j − 1 2 , k j − 1 j The SAMR method for hyperbolic problems 8

Serial SAMR method Parallel SAMR method Examples References Conservative flux correction Conservative flux correction II ◮ Level l cells needing correction ( G r l +1 l +1 \ G l +1 ) ∩ G l ◮ Corrections δ F n , l +1 stored on level l + 1 along ∂ G l +1 (lower-dimensional data coarsened by r l +1 ) ◮ Init δ F n , l +1 with level l fluxes F n , l (¯ G l ∩ ∂ G l +1 ) ◮ Add level l + 1 fluxes F n , l +1 ( ∂ G l +1 ) to δ F n , l F n , l F n , l +1 δ F n , l +1 Cells to correct The SAMR method for hyperbolic problems 9

Serial SAMR method Parallel SAMR method Examples References Level transfer operators Level transfer operators Conservative averaging (restriction): Replace cells on level l covered by level l + 1, i.e. j − 1 j G l ∩ G l +1 , by r l +1 − 1 r l +1 − 1 k 1 ˆ X X Q l Q l +1 jk := v v + κ, w + ι ( r l +1 ) 2 κ =0 ι =0 ( x j − 1 1 , l , x k − 1 ) Bilinear interpolation (prolongation): 2 , l k − 1 Q l +1 ˇ vw := (1 − f 1 )(1 − f 2 ) Q l j − 1 , k − 1 + f 1 (1 − f 2 ) Q l j , k − 1 + w (1 − f 1 ) f 2 Q l j − 1 , k + f 1 f 2 Q l jk 1 , l +1 − x j − 1 x v x w 2 , l +1 − x k − 1 1 , l 2 , l with factors f 1 := , f 2 := derived from the spatial ∆ x 1 , l ∆ x 2 , l coordinates of the cell centers ( x j − 1 1 , l , x k − 1 2 , l ) and ( x v 1 , l +1 , x w 2 , l +1 ). For boundary conditions on ˜ I s l : linear time interpolation „ κ « Q l +1 ( t )+ κ Q l +1 ( t + κ ∆ t l +1 ) := ˜ ˇ Q l +1 ( t +∆ t l ) ˇ 1 − for κ = 0 , . . . r l +1 r l +1 r l +1 The SAMR method for hyperbolic problems 10

Serial SAMR method Parallel SAMR method Examples References The basic recursive algorithm Recursive integration order ◮ Space-time interpolation of coarse data to set I s l , l > 0 ◮ Regridding: ◮ Creation of new grids, copy existing cells on level l > 0 ◮ Spatial interpolation to initialize new cells on level l > 0 Root Level 1 r 0 = 1 Level 1 2 5 8 11 r 1 = 4 Level 2 3 4 6 7 9 10 12 13 r 2 = 2 Time Regridding of finer levels. Base level ( ) stays fixed. The SAMR method for hyperbolic problems 11

Serial SAMR method Parallel SAMR method Examples References The basic recursive algorithm The basic recursive algorithm AdvanceLevel( l ) Repeat r l times Set ghost cells of Q l ( t ) If time to regrid? Regrid( l ) UpdateLevel( l ) If level l + 1 exists? Set ghost cells of Q l ( t + ∆ t l ) AdvanceLevel( l + 1 ) Average Q l +1 ( t + ∆ t l ) onto Q l ( t + ∆ t l ) Correct Q l ( t + ∆ t l ) with δ F l +1 t := t + ∆ t l Start - Start integration on level 0 l = 0 , r 0 = 1 AdvanceLevel( l ) [Berger and Colella, 1988][Berger and Oliger, 1984] The SAMR method for hyperbolic problems 12

Serial SAMR method Parallel SAMR method Examples References The basic recursive algorithm Regridding algorithm Regrid( l ) - Regrid all levels ι > l For ι = l f Downto l Do ◮ Refinement flags: Flag N ι according to Q ι ( t ) N l := S m N ( ∂ G l , m ) If level ι + 1 exists? ◮ Activate flags below higher Flag N ι below ˘ G ι +2 levels Flag buffer zone on N ι ◮ Flag buffer cells of b > κ r cells, G ι +1 from N ι Generate ˘ κ r steps between calls of ˘ G l := G l Regrid( l ) For ι = l To l f Do ◮ Spezial cluster algorithm C ˘ G ι := G 0 \ ˘ G ι ◮ Use complement operation to 1 G ι +1 := ˘ ˘ G ι +1 \ C ˘ G ι ensure proper nesting condition Recompose( l ) The SAMR method for hyperbolic problems 13

Serial SAMR method Parallel SAMR method Examples References The basic recursive algorithm Recomposition of data Recompose( l ) - Reorganize all levels ι > l For ι = l + 1 To l f + 1 Do Interpolate Q ι − 1 ( t ) onto ˘ Q ι ( t ) Copy Q ι ( t ) onto ˘ Q ι ( t ) Set ghost cells of ˘ Q ι ( t ) Q ι ( t ) := ˘ Q ι ( t ) , G ι := ˘ G ι ◮ Creates max. 1 level above l f , but can remove multiple level if ˘ G ι empty (no coarsening!) ◮ Use spatial interpolation on entire data ˘ Q ι ( t ) ◮ Overwrite where old data exists ◮ Synchronization and physical boundary conditions The SAMR method for hyperbolic problems 14

Serial SAMR method Parallel SAMR method Examples References Cluster algorithm Clustering by signatures x x x x x x x x x x x x 6 6 x x x x x x x x x x x x 6 6 -3 x x x 3 x x x 3 3 x x x x x x 3 3 -1 x x x x 2 2 1 x x 2 x x 2 0 x x x x 2 2 0 x x 2 x x x x x x 2 Υ 6 6 2 3 2 2 2 2 2 Υ 4 4 2 3 2 2 2 2 2 ∆ -2 3 -2 1 0 0 0 Υ Flagged cells per row/column ∆ Second derivative of Υ, ∆ = Υ ν +1 − 2 Υ ν + Υ ν − 1 Technique from image detection: [Bell et al., 1994], see also [Berger and Rigoutsos, 1991], [Berger, 1986] The SAMR method for hyperbolic problems 15

Serial SAMR method Parallel SAMR method Examples References Cluster algorithm x x x x x x x x x x x x x x x x x x x x x x x x x x x x x x 3 3 x x x 3 -1 x x x 1 x x x x 2 1 x x x x 2 0 x x 2 x x x x x x x x x x Υ 4 4 2 1 1 Υ 2 1 1 ∆ -2 1 ∆ 1 Recursive generation of ˘ G l , m 1. 0 in Υ 2. Largest difference in ∆ 3. Stop if ratio between flagged and unflagged cell > η tol The SAMR method for hyperbolic problems 16

Serial SAMR method Parallel SAMR method Examples References Refinement criteria Refinement criteria Scaled gradient of scalar quantity w | w ( Q j +1 , k ) − w ( Q jk ) | > ǫ w , | w ( Q j , k +1 ) − w ( Q jk ) | > ǫ w , | w ( Q j +1 , k +1 ) − w ( Q jk ) | > ǫ w Heuristic error estimation [Berger, 1982]: Local truncation error of scheme of order o q ( x , t + ∆ t ) − H (∆ t ) ( q ( · , t )) = C ∆ t o +1 + O (∆ t o +2 ) For q smooth after 2 steps ∆ t ( q ( · , t − ∆ t )) = 2 C ∆ t o +1 + O (∆ t o +2 ) q ( x , t + ∆ t ) − H (∆ t ) 2 and after 1 step with 2∆ t q ( x , t + ∆ t ) − H (2∆ t ) ( q ( · , t − ∆ t )) = 2 o +1 C ∆ t o +1 + O (∆ t o +2 ) Gives ( q ( · , t − ∆ t )) − H (2∆ t ) ( q ( · , t − ∆ t )) = (2 o +1 − 2) C ∆ t o +1 + O (∆ t o +2 ) H (∆ t ) 2 The SAMR method for hyperbolic problems 17

Serial SAMR method Parallel SAMR method Examples References Refinement criteria Heuristic error estimation for FV methods 2. Create temporary Grid coarsened by factor 2 Initialize with fine-grid- values of preceding 3. Compare tempo- 1. Error estimation on time step rary solutions interior cells H ∆ t l Q l ( t l − ∆ t l ) H ∆ t l ( H ∆ t l Q l ( t l − ∆ t l )) H ∆ t l Q l ( t l − ∆ t l ) = 2 H 2∆ t l ¯ Q l ( t l − ∆ t l ) The SAMR method for hyperbolic problems 18

Serial SAMR method Parallel SAMR method Examples References Refinement criteria Usage of heuristic error estimation Current solution integrated tentatively 1 step with ∆ t l and coarsened � � ¯ H ∆ t l Q l ( t l − ∆ t l ) Q ( t l + ∆ t l ) := Restrict 2 Previous solution coarsened integrated 1 step with 2∆ t l Q ( t l + ∆ t l ) := H 2∆ t l Restrict Q l ( t l − ∆ t l ) � � Local error estimation of scalar quantity w | w ( ¯ Q jk ( t + ∆ t )) − w ( Q jk ( t + ∆ t )) | τ w jk := 2 o +1 − 2 In practice [Deiterding, 2003] use τ w jk max( | w ( Q jk ( t + ∆ t )) | , S w ) > η r w The SAMR method for hyperbolic problems 19

Serial SAMR method Parallel SAMR method Examples References Outline The serial Berger-Colella SAMR method Block-based data structures Numerical update Conservative flux correction Level transfer operators The basic recursive algorithm Cluster algorithm Refinement criteria Parallel SAMR method Domain decomposition A parallel SAMR algorithm Partitioning Examples Euler equations The SAMR method for hyperbolic problems 20

Serial SAMR method Parallel SAMR method Examples References Domain decomposition Parallelization strategies Decomposition of the hierarchical data ◮ Distribution of each grid ◮ Separate distribution of each level, cf. Processor 1 Processor 2 [Rendleman et al., 2000] ◮ Rigorous domain decomposition ◮ Data of all levels resides on same node ◮ Grid hierarchy defines unique ”floor-plan” ◮ Redistribution of data blocks during reorganization of hierarchical data ◮ Synchronization when setting ghost cells The SAMR method for hyperbolic problems 21

Serial SAMR method Parallel SAMR method Examples References Domain decomposition Rigorous domain decomposition formalized Parallel machine with P identical nodes. P non-overlapping portions G p 0 , p = 1 , . . . , P as P [ G p G p 0 ∩ G q 0 = ∅ for p � = q G 0 = with 0 p =1 Higher level domains G l follow decomposition of root level G p l := G l ∩ G p 0 With N l ( · ) denoting number of cells, we estimate the workload as l max " l # X Y W (Ω) = N l ( G l ∩ Ω) r κ l =0 κ =0 Equal work distribution necessitates L p := P · W ( G p 0 ) ≈ 1 for all p = 1 , . . . , P W ( G 0 ) [Deiterding, 2005] The SAMR method for hyperbolic problems 22

Serial SAMR method Parallel SAMR method Examples References Domain decomposition Ghost cell setting Local synchronization Processor 1 Processor 2 S s , p G s , p i , m ∩ G p ˜ i , m = ˜ i Parallel synchronization S s , q ˜ i , m = ˜ G s , p i , m ∩ G q i , q � = p Interpolation and physi- cal boundary conditions remain strictly local ◮ Scheme H (∆ t l ) evaluated locally ◮ Restriction and propolongation local Ghost cell values: Interpolation Parallel synchronization Local synchronization Physical boundary The SAMR method for hyperbolic problems 23

Serial SAMR method Parallel SAMR method Examples References Domain decomposition Parallel flux correction 1. Strictly local: Init δ F n , l +1 with F n (¯ G l , m ∩ ∂ G l +1 , t ) 2. Strictly local: Add F n ( ∂ G l , m , t ) to δ F n , l 3. Parallel communication: Correct Q l ( t + ∆ t l ) with δ F l +1 Node p 1 Node q v + 2 parallel exchange w k j − 1 j j − 1 j F n , l F n , l +1 δ F n , l +1 The SAMR method for hyperbolic problems 24

Serial SAMR method Parallel SAMR method Examples References A parallel SAMR algorithm The recursive algorithm in parallel AdvanceLevel( l ) Repeat r l times Set ghost cells of Q l ( t ) If time to regrid? Regrid( l ) UpdateLevel( l ) If level l + 1 exists? Set ghost cells of Q l ( t + ∆ t l ) AdvanceLevel( l + 1 ) Average Q l +1 ( t + ∆ t l ) onto Q l ( t + ∆ t l ) Correct Q l ( t + ∆ t l ) with δ F l +1 t := t + ∆ t l UpdateLevel( l ) For all m = 1 To M l Do l , m , t ) H (∆ tl ) → Q ( G l , m , t + ∆ t l ) , F n (¯ Q ( G s − G l , m , t ) If level l > 0 Add F n ( ∂ G l , m , t ) to δ F n , l If level l + 1 exists Init δ F n , l +1 with F n (¯ G l , m ∩ ∂ G l +1 , t ) The SAMR method for hyperbolic problems 25

Serial SAMR method Parallel SAMR method Examples References A parallel SAMR algorithm Regridding algorithm in parallel Regrid( l ) - Regrid all levels ι > l For ι = l f Downto l Do Flag N ι according to Q ι ( t ) If level ι + 1 exists? Flag N ι below ˘ G ι +2 Flag buffer zone on N ι G ι +1 from N ι Generate ˘ ˘ G l := G l For ι = l To l f Do C ˘ G ι := G 0 \ ˘ G ι 1 G ι +1 := ˘ ˘ G ι +1 \ C ˘ G ι Recompose( l ) The SAMR method for hyperbolic problems 26

Serial SAMR method Parallel SAMR method Examples References A parallel SAMR algorithm Recomposition algorithm in parallel Recompose( l ) - Reorganize all levels Generate G p 0 from { G 0 , ..., G l , ˘ G l +1 , ..., ˘ G l f +1 } For ι = l + 1 To l f + 1 Do If ι > i G p ˘ ι := ˘ G ι ∩ G p 0 Interpolate Q ι − 1 ( t ) onto ˘ Q ι ( t ) else G p ˘ ι := G ι ∩ G p 0 If ι > 0 Copy δ F n ,ι onto δ ˘ F n ,ι δ F n ,ι := δ ˘ F n ,ι If ι ≥ l then κ ι = 0 else κ ι = 1 For κ = 0 To κ ι Do Copy Q ι ( t ) onto ˘ Q ι ( t ) Set ghost cells of ˘ Q ι ( t ) Q ι ( t ) := ˘ Q ι ( t ) G ι := ˘ G ι The SAMR method for hyperbolic problems 27

Serial SAMR method Parallel SAMR method Examples References Partitioning Space-filling curve algorithm Necessary domain of Space-Filling Curve Calculation domain High Workload Proc. 1 Proc. 2 Medium Workload Proc. 3 Low Workload The SAMR method for hyperbolic problems 28

Serial SAMR method Parallel SAMR method Examples References Outline The serial Berger-Colella SAMR method Block-based data structures Numerical update Conservative flux correction Level transfer operators The basic recursive algorithm Cluster algorithm Refinement criteria Parallel SAMR method Domain decomposition A parallel SAMR algorithm Partitioning Examples Euler equations The SAMR method for hyperbolic problems 29

Serial SAMR method Parallel SAMR method Examples References Euler equations Benchmark run: blast wave in 2D Task [%] P =1 P =2 P =4 P =8 P =16 ◮ 2D-Wave-Propagation Update by H ( · ) 86.6 83.4 76.7 64.1 51.9 Method with Roe’s Flux correction 1.2 1.6 3.0 7.9 10.7 approximate solver Boundary setting 3.5 5.7 10.1 15.6 18.3 Recomposition 5.5 6.1 7.4 9.9 14.0 ◮ Base grid 150 × 150 Misc. 4.9 3.2 2.8 2.5 5.1 ◮ 2 levels: factor 2, 4 Time [min] 151.9 79.2 43.4 23.3 13.9 Efficiency [ % ] 100.0 95.9 87.5 81.5 68.3 After 38 time steps After 79 time steps The SAMR method for hyperbolic problems 30

Serial SAMR method Parallel SAMR method Examples References Euler equations Benchmark run 2: point-explosion in 3D l max = 4 solution ◮ Benchmark from the Chicago workshop on AMR methods, September 2003 ◮ Sedov explosion - energy deposition in sphere of radius 4 finest cells ◮ 3D-Wave-Prop. Method ◮ Base grid 32 3 l max = 5 solution ◮ Refinement factor r l = 2 ◮ Effective resolutions: 128 3 , 256 3 , 512 3 , 1024 3 ◮ Grid generation efficiency 85% ◮ Proper nesting enforced ◮ Buffer of 1 cell The SAMR method for hyperbolic problems 31

Serial SAMR method Parallel SAMR method Examples References Euler equations Benchmark run 2: visualization of refinement l = 0 l = 1 l = 2 l = 3 l = 4 l = 5 The SAMR method for hyperbolic problems 32

Serial SAMR method Parallel SAMR method Examples References Euler equations Benchmark run 2: performance results Number of grids and cells l l max = 2 l max = 3 l max = 4 l max = 5 Grids Cells Grids Cells Grids Cells Grids Cells 0 28 32768 28 32768 33 32768 34 32768 1 8 32768 14 32768 20 32768 20 32768 2 63 115408 49 116920 43 125680 50 125144 3 324 398112 420 555744 193 572768 4 1405 1487312 1498 1795048 5 5266 5871128 P 180944 580568 2234272 8429624 Breakdown of CPU time on 8 nodes SGI Altix 3000 (Linux-based shared memory system) Task [%] l max = 2 l max = 3 l max 45 l max = 5 Integration 73.7 77.2 72.9 37.8 Fixup 2.6 46 3.1 58 2.6 42 2.2 45 Boundary 10.1 79 6.3 78 5.1 56 6.9 78 Recomposition 7.4 8.0 15.1 50.4 Clustering 0.5 0.6 0.7 1.0 Output/Misc 5.7 4.0 3.6 1.7 Time [min] 0.5 5.1 73.0 2100.0 Uniform [min] 5.4 160 ∼ 5000 ∼ 180000 Factor of AMR savings 11 31 69 86 Time steps 15 27 52 115 The SAMR method for hyperbolic problems 33

Serial SAMR method Parallel SAMR method Examples References References References I [Bell et al., 1994] Bell, J., Berger, M., Saltzman, J., and Welcome, M. (1994). Three-dimensional adaptive mesh refinement for hyperbolic conservation laws. SIAM J. Sci. Comp. , 15(1):127–138. [Berger, 1982] Berger, M. (1982). Adaptive Mesh Refinement for Hyperbolic Differential Equations . PhD thesis, Stanford University. Report No. STAN-CS-82-924. [Berger, 1986] Berger, M. (1986). Data structures for adaptive grid generation. SIAM J. Sci. Stat. Comput. , 7(3):904–916. [Berger and Colella, 1988] Berger, M. and Colella, P. (1988). Local adaptive mesh refinement for shock hydrodynamics. J. Comput. Phys. , 82:64–84. [Berger and Oliger, 1984] Berger, M. and Oliger, J. (1984). Adaptive mesh refinement for hyperbolic partial differential equations. J. Comput. Phys. , 53:484–512. [Berger and Rigoutsos, 1991] Berger, M. and Rigoutsos, I. (1991). An algorithm for point clustering and grid generation. IEEE Transactions on Systems, Man, and Cybernetics , 21(5):1278–1286. The SAMR method for hyperbolic problems 34

Serial SAMR method Parallel SAMR method Examples References References References II [Deiterding, 2003] Deiterding, R. (2003). Parallel adaptive simulation of multi-dimensional detonation structures . PhD thesis, Brandenburgische Technische Universit¨ at Cottbus. [Deiterding, 2005] Deiterding, R. (2005). Construction and application of an AMR algorithm for distributed memory computers. In Plewa, T., Linde, T., and Weirs, V. G., editors, Adaptive Mesh Refinement - Theory and Applications , volume 41 of Lecture Notes in Computational Science and Engineering , pages 361–372. Springer. [Rendleman et al., 2000] Rendleman, C. A., Beckner, V. E., Lijewski, M., Crutchfield, W., and Bell, J. B. (2000). Parallelization of structured, hierarchical adaptive mesh refinement algorithms. Computing and Visualization in Science , 3:147–157. The SAMR method for hyperbolic problems 35

Complex geometry Combustion Fluid-structure interaction Turbulence References Lecture 3 Complex hyperbolic applications Course Block-structured Adaptive Mesh Refinement Methods for Conservation Laws Theory, Implementation and Application Ralf Deiterding Computer Science and Mathematics Division Oak Ridge National Laboratory P.O. Box 2008 MS6367, Oak Ridge, TN 37831, USA E-mail: deiterdingr@ornl.gov Complex hyperbolic applications 1

Complex geometry Combustion Fluid-structure interaction Turbulence References Outline Complex geometry Boundary aligned meshes Cartesian techniques Implicit geometry representation Accuracy / Verification Combustion Equations and FV schemes Shock-induced combustion examples Fluid-structure interaction Coupling to a solid mechanics solver Rigid body motion Thin elastic structures Deforming thin structures Turbulence Large-eddy simulation Complex hyperbolic applications 2

Complex geometry Combustion Fluid-structure interaction Turbulence References Outline Complex geometry Boundary aligned meshes Cartesian techniques Implicit geometry representation Accuracy / Verification Combustion Equations and FV schemes Shock-induced combustion examples Fluid-structure interaction Coupling to a solid mechanics solver Rigid body motion Thin elastic structures Deforming thin structures Turbulence Large-eddy simulation Complex hyperbolic applications 3

Complex geometry Combustion Fluid-structure interaction Turbulence References Boundary aligned meshes SAMR on boundary aligned meshes Analytic or stored geometric mapping of the coordinates (graphic from [Yamaleev and Carpenter, 2002]) ◮ Topology remains unchanged and thereby entire SAMR algorithm ◮ Patch solver and interpolation need to consider geometry transformation ◮ Handles boundary layers well Overlapping adaptive meshes [Henshaw and Schwendeman, 2003], [Meakin, 1995] ◮ Idea is to use a non-Cartesian structured grids only near boundary ◮ Very suitable for moving objects with boundary layers ◮ Interpolation between meshes is usually non-conservative Complex hyperbolic applications 4

Complex geometry Combustion Fluid-structure interaction Turbulence References Cartesian techniques Cut-cell techniques Accurate embedded boundary method “ A n +1 / 2 V n +1 Q n +1 = V n j Q n j − ∆ t j +1 / 2 F ( Q , j ) j j ” − A n +1 / 2 j − 1 / 2 F ( Q , j − 1) Methods that represent the boundary sharply: ◮ Cut-cell approach constructs appropriate finite volumes ◮ Conservative by construction. Correct boundary flux ◮ Key question: How to avoid small-cell time step restriction in explicit metyhods? ◮ Cell merging: [Quirk, 1994a] ◮ Usually explicit geometry representation used [Aftosmis, 1997], but can also be implicit, cf. [Nourgaliev et al., 2003], [Murman et al., 2003] Complex hyperbolic applications 5

Complex geometry Combustion Fluid-structure interaction Turbulence References Cartesian techniques Embedded boundary techniques Volume of fluid methods that resemble a cut-cell technique on purely Cartesian mesh ◮ Redistribution of boundary flux achieves conservation and bypasses time step restriction: [Pember et al., 1999], [Berger and Helzel, 2002] Methods that diffuse the boundary in one cell (good overview in [Mittal and Iaccarino, 2005]): ◮ Related to the immersed boundary method by Peskin, cf. [Roma et al., 1999] ◮ Boundary prescription often by internal ghost cell values, cf. [Tseng and Ferziger, 2003] ◮ Not conservative by construction but conservative correction possible ◮ Usually combined with implicit geometry representation K. J. Richards et al., On the use of the immersed boundary method for engine modeling Complex hyperbolic applications 6

ut ut ut Complex geometry Combustion Fluid-structure interaction Turbulence References Cartesian techniques Level-set method for boundary embedding ◮ Complex boundary moving with local velocity w , treat interface as moving rigid wall ◮ Implicit boundary representation via distance function ϕ , normal n = ∇ ϕ/ |∇ ϕ | ◮ Construction of values in embedded boundary cells by interpolation / extrapolation Interpolate / constant value ex- 2 w − u j trapolate values at w ˜ x = x + 2 ϕ n u j Velocity in ghost cells ρ j − 1 ρ j ρ j ρ j − 1 u ′ = (2 w · n − u · n ) n + ( u · t ) t u j − 1 u j 2 w − u j 2 w − u j − 1 p j − 1 p j p j p j − 1 = 2 (( w − u ) · n ) n + u Complex hyperbolic applications 7

Complex geometry Combustion Fluid-structure interaction Turbulence References Implicit geometry representation Closest point transform algorithm The signed distance ϕ to a surface I satisfies the eikonal equation [Sethian, 1999] ˛ |∇ ϕ | = 1 with ϕ I = 0 ˛ Solution smooth but non-diferentiable across characteristics. Distance computation trivial for non-overlapping elementary shapes but difficult to do efficiently for triangulated surface meshes: ◮ Geometric solution approach with plosest-point-transform algorithm [Mauch, 2003] b-rep Surface mesh Distance ϕ Normal to closest point Complex hyperbolic applications 8

Complex geometry Combustion Fluid-structure interaction Turbulence References Implicit geometry representation The characteristic / scan conversion algorithm Characteristic polyhedra for faces, edges, and vertices 1. Build the characteristic polyhedrons for the surface mesh 2. For each face/edge/vertex 2.1 Scan convert the polyhedron. (a) (b) 2.2 Compute distance to that primitive for the scan converted points 3. Computational complexity. (c) (d) ◮ O ( m ) to build the b-rep and the polyhedra. ◮ O ( n ) to scan convert the polyhedra and compute the Slicing and scan conversion of apolygon distance, etc. 4. Problem reduction by evaluation only within specified max. distance [Deiterding et al., 2006] Complex hyperbolic applications 9

Complex geometry Combustion Fluid-structure interaction Turbulence References Accuracy / Verification Accuracy test: stationary vortex U ( r ) Construct non-trivial radially symmetric and stationary solution by balancing hydrodynamic pressure and centripetal force per volume element, i.e. dr p ( r ) = ρ ( r ) U ( r ) 2 d r p ( r ) For ρ 0 ≡ 1 and the velocity field 8 2 r / R if 0 < r < R / 2 , < U ( r ) = α · 2(1 − r / R ) if R / 2 ≤ r ≤ R , 0 if r > R , : one gets with boundary condition p ( R ) = p 0 ≡ 2 the pressure distribution r 2 / R 2 + 1 − 2 log 2 8 if 0 < r < R / 2 , p ( r ) = p 0 + 2 ρ 0 α 2 · < r 2 / R 2 + 3 − 4 r / R + 2 log( r / R ) if R / 2 ≤ r ≤ R , 0 if r > R . : Entire solution for Euler equations reads ρ ( x , y , t ) = ρ 0 , u ( x , y , t ) = − sin φ U ( r ) , v ( x , y , t ) = cos φ U ( r ) , p ( x , y , t ) = p ( r ) ( x − x c ) 2 + ( y − y c ) 2 and φ = arctan y − y c p for all t ≥ 0 with r = x − x c Complex hyperbolic applications 10

Complex geometry Combustion Fluid-structure interaction Turbulence References Accuracy / Verification Stationary vortex: results Compute one full rotation, Roe solver, embedded slip wall boundary conditions x c = 0 . 5 , y c = 0 . 5 , R = 0 . 4 , t end = 1 , ∆ h = ∆ x = ∆ y = 1 / N , α = R π No embedded boundary Wave Propagation Godunov Splitting N Error Order Error Order 20 0.0111235 0.0182218 40 0.0037996 1.55 0.0090662 1.01 80 0.0013388 1.50 0.0046392 0.97 160 0.0005005 1.42 0.0023142 1.00 Marginal shear flow along embedded boundary, α = R π, R G = R , U W = 0 Wave Propagation Godunov Splitting N Error Order Mass loss Error Order Mass loss 20 0.0120056 0.0079236 0.0144203 0.0020241 40 0.0035074 1.78 0.0011898 0.0073070 0.98 0.0001300 80 0.0014193 1.31 0.0001588 0.0038401 0.93 -0.0001036 160 0.0005032 1.50 5.046e-05 0.0018988 1.02 -2.783e-06 Major shear flow along embedded boundary, α = R π, R G = R / 2 , U W = 0 Wave Propagation Godunov Splitting N Error Order Mass loss Error Order Mass loss 20 0.0423925 0.0423925 0.0271446 0.0271446 40 0.0358735 0.24 0.0358735 0.0242260 0.16 0.0242260 80 0.0212340 0.76 0.0212340 0.0128638 0.91 0.0128638 160 0.0121089 0.81 0.0121089 0.0070906 0.86 0.0070906 Complex hyperbolic applications 11

Complex geometry Combustion Fluid-structure interaction Turbulence References Accuracy / Verification Verification: shock reflection ◮ Reflection of a Mach 2.38 shock in nitrogen at 43 o wedge ◮ 2nd order MUSCL scheme with Roe solver, 2nd order multidimensional wave propagation method Cartesian base grid 360 × 160 cells on domain of 36 mm × 16 mm with up to 3 refinement levels with r l = 2 , 4 , 4 and ∆ x = 3 . 125 µ m , 38 h CPU GFM base grid 390 × 330 cells on domain of 26 mm × 22 mm with up to 3 refinement levels with r l = 2 , 4 , 4 and ∆ x e = 2 . 849 µ m , 200 h CPU Complex hyperbolic applications 12

Complex geometry Combustion Fluid-structure interaction Turbulence References Accuracy / Verification Shock reflection: SAMR solution for Euler equations ∆ x = 25 mm ∆ x = 12 . 5 mm ∆ x = 3 . 125 mm ∆ x e = 22 . 8 mm ∆ x e = 11 . 4 mm ∆ x e = 2 . 849 mm 2nd order MUSCL scheme with Van Leer FVS, dimen- sional splitting ∆ x = 12 . 5 mm ∆ x = 3 . 125 mm Complex hyperbolic applications 13

Complex geometry Combustion Fluid-structure interaction Turbulence References Accuracy / Verification Shock reflection: solution for Navier-Stokes equations ◮ No-slip boundary conditions enforced ◮ Conservative 2nd order centered differences to approximate stress tensor and heat flow ∆ x = 50 mm ∆ x = 25 mm ∆ x = 12 . 5 mm , SAMR ∆ x e = 45 . 6 mm ∆ x e = 22 . 8 mm ∆ x e = 11 . 4 mm , SAMR Complex hyperbolic applications 14

Complex geometry Combustion Fluid-structure interaction Turbulence References Outline Complex geometry Boundary aligned meshes Cartesian techniques Implicit geometry representation Accuracy / Verification Combustion Equations and FV schemes Shock-induced combustion examples Fluid-structure interaction Coupling to a solid mechanics solver Rigid body motion Thin elastic structures Deforming thin structures Turbulence Large-eddy simulation Complex hyperbolic applications 15

Complex geometry Combustion Fluid-structure interaction Turbulence References Equations and FV schemes Governing equations for premixed combustion Euler equations with reaction terms: ∂ t ρ i + ∇ · ( ρ u ) = ˙ ω i for i = 1 , . . . , K ∂ t ( ρ u ) + ∇ · ( ρ u ⊗ u ) + ∇ p = 0 ∂ t ( ρ E ) + ∇ · (( ρ E + p ) u ) = 0 Ideal gas law and Dalton’s law for gas-mixtures: K K K ρ i R W i T = ρ R ρ i = ρ , Y i = ρ i X X X p ( ρ 1 , . . . , ρ K , T ) = p i = W T with ρ i =1 i =1 i =1 Caloric equation: Z T K X h i ( T ) = h 0 h ( Y 1 , . . . , Y K , T ) = Y i h i ( T ) with i + c pi ( s ) ds 0 i =1 Computation of T = T ( ρ 1 , . . . , ρ K , e ) from implicit equation K K ρ i X X ρ i h i ( T ) − R T W i − ρ e = 0 i =1 i =1 for thermally perfect gases with γ i ( T ) = c pi ( T ) / c vi ( T ) Complex hyperbolic applications 16

Complex geometry Combustion Fluid-structure interaction Turbulence References Equations and FV schemes Chemistry Arrhenius-Kinetics: “ ρ n “ ρ n M K K ” ν f ” ν r » jn − k r jn – X ( ν r ji − ν f k f Y Y ω i = ˙ ji ) i = 1 , . . . , K j j W n W n n =1 n =1 j =1 ◮ Parsing of mechanisms with Chemkin-II ◮ Evalutation of ˙ ω i with automatically generated optimized Fortran-77 functions in the line of Chemkin-II Integration of reaction rates: ODE integration in S ( · ) for Euler equations with chemical reaction ◮ Standard implicit or semi-implicit ODE-solver subcycles within each cell ◮ ρ , e , u remain unchanged! ∂ t ρ i = W i ˙ ω i ( ρ 1 , . . . , ρ K , T ) i = 1 , . . . , K Use Newton or bisection method to compute T iteratively. Complex hyperbolic applications 17

Complex geometry Combustion Fluid-structure interaction Turbulence References Equations and FV schemes Non-equilibrium mechanism for hydrogen-oxgen combustion Eact A β [ cal mol − 1 ] [ cm , mol , s ] 1 . 86 × 1014 1. H + O2 − → O + OH 0 . 00 16790 . 1 . 48 × 1013 2. O + OH − → H + O2 0 . 00 680 . 1 . 82 × 1010 3. H2 + O − → H + OH 1 . 00 8900 . 8 . 32 × 1009 4. H + OH − → H2 + O 1 . 00 6950 . 3 . 39 × 1013 5. H2O + O − → OH + OH 0 . 00 18350 . 3 . 16 × 1012 6. OH + OH − → H2O + O 0 . 00 1100 . 9 . 55 × 1013 7. H2O + H − → H2 + OH 0 . 00 20300 . 2 . 19 × 1013 8. H2 + OH − → H2O + H 0 . 00 5150 . 1 . 00 × 1013 9. H2O2 + OH − → H2O + HO2 0 . 00 1800 . 2 . 82 × 1013 10. H2O + HO2 − → H2O2 + OH 0 . 00 32790 . . . . . . . . . . . . . . . . . . . 7 . 94 × 1019 30. OH + M − → O + H + M − 1 . 00 103720 . 5 . 13 × 1015 31. O2 + M − → O + O + M 0 . 00 115000 . 4 . 68 × 1015 32. O + O + M − → O2 + M − 0 . 28 0 . 2 . 19 × 1014 33. H2 + M − → H + H + M 0 . 00 96000 . 3 . 02 × 1015 34. H + H + M − → H2 + M 0 . 00 0 . Third body efficiencies: f ( O2 ) = 0 . 40, f ( H2O ) = 6 . 50 C. K. Westbrook. Chemical kinetics of hydrocarbon oxidation in gaseous detonations. J. Combustion and Flame , 46:191–210, 1982. Complex hyperbolic applications 18

Complex geometry Combustion Fluid-structure interaction Turbulence References Equations and FV schemes Riemann solver for combustion u n , ˆ H , ˆ Y i , ˆ (S1) Calculate standard Roe-averages ˆ ρ , ˆ T . Z Tr 1 (S2) Compute ˆ γ := ˆ c p / ˆ c v with ˆ c { p / v } i = c { p , v } i ( τ ) d τ . T r − T l Tl “ u 2 ” (S3) Calculate ˆ ˆ 2 − ˆ γ R i ˆ e i or ˆ φ i := (ˆ γ − 1) h i + ˆ T with standard Roe-averages ˆ h i . ” 1 / 2 . “P K u 2 + (ˆ i =1 ˆ Y i ˆ γ − 1)ˆ (S4) Calculate ˆ c := φ i − (ˆ γ − 1)ˆ H (S5) Use ∆ Q = Q r − Q l and ∆ p to compute the wave strengths a m . K + d X (S6) Calculate W 1 = a 1 ˆ r 1 , W 2 = a ι ˆ r ι , W 3 = a K + d +1 ˆ r K + d +1 . ι =2 (S7) Evaluate s 1 = ˆ u 1 − ˆ c , s 2 = ˆ u 1 , s 3 = ˆ u 1 + ˆ c . (S8) Evaluate ρ ⋆ l / r , u ⋆ 1 , l / r , e ⋆ l / r , c ⋆ 1 , l / r from Q ⋆ l = Q l + W 1 and Q ⋆ r = Q r − W 3 . (S9) If ρ ⋆ l / r ≤ 0 or e ⋆ l / r ≤ 0 use F HLL ( Q l , Q r ) and go to (S12). (S10) Entropy correction: Evaluate | ˜ s ι | . F Roe ( Q l , Q r ) = 1 f ( Q l ) + f ( Q r ) − P 3 ` ´ ι =1 | ˜ s ι | W ι 2 (S11) Positivity correction: Replace F i by  Y l i , F ρ ≥ 0 , F ⋆ i = F ρ · Y r i , F ρ < 0 . (S12) Evaluate maximal signal speed by S = max( | s 1 | , | s 3 | ). Complex hyperbolic applications 19

Complex geometry Combustion Fluid-structure interaction Turbulence References Equations and FV schemes Riemann solver for combustion: carbuncle fix 2D modification of entropy correction Entropy corrections [Harten, 1983] [Sanders et al., 1998]: [Harten and Hyman, 1983] ( | s ι | if | s ι | ≥ 2 η i , j + 1 i + 1 , j + 1 2 2 1. | ˜ s ι | = | s 2 ι | 4 η + η otherwise i + 1 η = 1 ˘ ¯ 2 , j 2 max ι | s ι ( Q r ) − s ι ( Q l ) | 2. Replace | s ι | by | ˜ s ι | only if i , j − 1 i + 1 , j − 1 s ι ( Q l ) < 0 < s ι ( Q r ) 2 2 ˘ ¯ η i +1 / 2 , j = max ˜ η i +1 / 2 , j , η i , j − 1 / 2 , η i , j +1 / 2 , η i +1 , j − 1 / 2 , η i +1 , j +1 / 2 Roe + EC 1. Exact Riemann solver Carbuncle phenomenon ◮ [Quirk, 1994b] ◮ Test from Roe + EC 2. SW FVS, VL FVS, HLL, Roe + EC 2.+2D [Deiterding, 2003] Complex hyperbolic applications 20

Complex geometry Combustion Fluid-structure interaction Turbulence References Shock-induced combustion examples Detonations - motivation for SAMR ◮ Extremly high spatial resolution in reaction zone necessary. ◮ Minimal spatial resolution: 7 − 8 Pts / l ig − → ∆ x ≈ 0 . 2 − 0 . 175 mm ◮ Uniform grids for typical geometries: > 10 7 Pts in 2D, > 10 9 Pts in 3D − → Self-adaptive finite volume method (AMR) 3000 3000 0.014 0.014 0.012 0.012 2500 2500 Temperature [K] Temperature [K] Mass fraction [-] Mass fraction [-] 0.01 0.01 0.008 0.008 2000 T 2000 T Y H Y H 0.006 0.006 Y O Y O Y H 2 Y H 2 0.004 0.004 1500 1500 0.002 0.002 1000 0 1000 0 6.28 6.32 6.36 6.4 5.9 5.94 5.98 6.02 x 1 [cm] x 1 [cm] Approximation of H 2 : O 2 detonation at ∼ 1 . 5 Pts / l ig (left) and ∼ 24 Pts / l ig (right) Complex hyperbolic applications 21

Complex geometry Combustion Fluid-structure interaction Turbulence References Shock-induced combustion examples Detonation ignition in a shock tube ◮ Shock-induced detonation ignition of H 2 : O 2 : Ar mixture at molar ratios 2:1:7 in closed 1d shock tube ◮ Insufficient resolution leads to inaccurate results ◮ Reflected shock is captured correctly by FV scheme, detonation is resolution dependent ◮ Fine mesh necessary in the induction zone at the head of the detonation 1.0 1.0 Pressure ∆ x 1 =6.25 µ m 0.9 0.9 Level ∆ x 1 =100 µ m 0.8 ∆ x 1 =200 µ m 0.8 0.7 0.7 Pressure [MPa] Pressure [MPa] 0.6 0.6 0.5 0.5 0.4 0.4 0.3 0.3 0.2 0.2 0.1 0.1 4.5 5 5.5 6 6.5 7 4 4.5 5 5.5 6 6.5 7 7.5 8 x 1 [cm] x 1 [cm] Left: Comparison of pressure distribution t = 170 µ s after shock reflection. Right: Domains of refinement levels Complex hyperbolic applications 22

Complex geometry Combustion Fluid-structure interaction Turbulence References Shock-induced combustion examples Detonation ignition in 1d - adaptive vs. uniform Uniformly refined vs. dynamic adaptive simulations (Intel Xeon 3 . 4 GHz CPU) Uniform Adaptive ∆ x 1 [ µ m ] Cells t m [ µ s ] Time [s] l max r l t m [ µ s ] Time [s] 400 300 166.1 31 200 600 172.6 90 2 2 172.6 99 100 1200 175.5 277 3 2,2 175.8 167 50 2400 176.9 858 4 2,2,2 177.3 287 25 4800 177.8 2713 4 2,2,4 177.9 393 12.5 9600 178.3 9472 5 2,2,2,4 178.3 696 6.25 19200 178.6 35712 5 2,2,4,4 178.6 1370 ∼ 12 Pts / l ig Refinement criteria: 1.0 Uniform 0.014 0.9 Adaptive S Y i · 10 − 4 η r Y i · 10 − 3 Y i 0.8 0.012 10.0 2.0 O 2 0.7 Pressure [MPa] 0.01 0.6 H 2 O 7.8 8.0 Y O [-] 0.008 0.5 H 0.16 5.0 0.4 0.006 O 1.0 5.0 0.3 0.004 1.8 5.0 OH 0.2 0.002 0.1 1.3 2.0 H 2 0 ǫ ρ = 0 . 07 kg m − 3 , ǫ p = 50 kPa 5.95 5.97 5.99 6.01 6.03 6.05 x 1 [cm] Complex hyperbolic applications 23

Complex geometry Combustion Fluid-structure interaction Turbulence References Shock-induced combustion examples Shock-induced combustion around a sphere ◮ Spherical projectile of radius 1 . 5 mm travels with constant velocity v I = 2170 . 6 m / s through H 2 : O 2 : Ar mixture (molar ratios 2:1:7) at 6 . 67 kPa and T = 298 K ◮ Cylindrical symmetric simulation on AMR base mesh of 70 × 40 cells ◮ Comparison of 3-level computation with refinement factors 2,2 ( ∼ 5 Pts / l ig ) and a 4-level computation with refinement factors 2,2,4 ( ∼ 19 Pts / l ig ) at t = 350 µ s ◮ Higher resolved computation captures combustion zone visibly better and at slightly different position (see below) Iso-contours of p (black) and Y H 2 (white) on refinement domains for 3-level (left) and 4-level computation (right) Complex hyperbolic applications 24

Complex geometry Combustion Fluid-structure interaction Turbulence References Shock-induced combustion examples Combustion around a sphere - adaptation Refinement indicators on l = 2 at t = 350 µ s . Refinement criteria: Blue: ǫ ρ , light blue: ǫ p , green shades: η r Yi , S Y i · 10 − 4 η r Y i · 10 − 4 Y i red: embedded boundary O 2 10.0 4.0 H 2 O 5.8 3.0 H 0.2 10.0 1.4 10.0 O 2.3 10.0 OH H 2 1.3 4.0 ǫ ρ = 0 . 02 kg m − 3 , ǫ p = 16 kPa Scaling of different code portions 10 Parallel performance 100 Total time Ideal time [s] 1 time [s] 10 Fluid dynamics Chemical kinetics Boundary setting Embedded boundary Recomposition 0.1 1 1 2 4 8 16 32 64 128 1 2 4 8 16 32 64 128 CPU CPU Complex hyperbolic applications 25

Complex geometry Combustion Fluid-structure interaction Turbulence References Shock-induced combustion examples Detonation diffraction ◮ CJ detonation for H 2 : O 2 : Ar /2 : 1 : 7 at T 0 = 298 K and p 0 = 10 kPa . Cell width λ c = 1 . 6 cm ◮ Adaption criteria (similar as before): 1. Scaled gradients of ρ and p E. Schultz. Detonation diffraction through an abrupt area expan- sion . PhD thesis, California Institute of Technology, Pasadena, 2. Error estimation in Y i by California, April 2000. Richardson extrapolation ◮ 25 Pts / l ig . 5 refinement levels (2,2,2,4). ◮ Adaptive computations use up to ∼ 2 . 2 M instead of ∼ 150 M cells (uniform grid) ◮ ∼ 3850 h CPU ( ∼ 80 h real time) on 48 nodes Athlon 1.4GHz Complex hyperbolic applications 26