Parallelization of stencil-based methods: step-69 (and beyond) - PowerPoint PPT Presentation

Parallelization of “stencil-based” methods: step-69 (and beyond) Martin Kronbichler 1 , Matthias Maier, Ignacio Tomas 2 Department of Mathematics Texas A&M University 1 Institute for Computational Mechanics, TUM, Munich 2 CSRI, Sandia National Laboratories 1 Kronbichler, Maier, Tomas | Parallelization of “stencil-based” methods

Motivation 2 Kronbichler, Maier, Tomas | Parallelization of “stencil-based” methods

Motivation Simulation of hyperbolic systems of conservation laws › robust (?!) › accurate (?!) › scalable (!) Goals: 3 Kronbichler, Maier, Tomas | Parallelization of “stencil-based” methods

Outline 1 Motivation 2 The compressible Euler equations 3 step-69 4 . . . and beyond 4 Kronbichler, Maier, Tomas | Parallelization of “stencil-based” methods

The compressible Euler equations 5 Kronbichler, Maier, Tomas | Parallelization of “stencil-based” methods

Compressible Euler equations Hyperbolic system of conservation laws u t + div f ( u ) = 0 ; where u ( x ; t ) : R d ˆ R ! R d + 2 , and f ( u ) : R d + 2 ! R ( d + 2 ) ˆ d , with space dim. d – 1. 2 3 m > 2 3  6 7  ` 1 m ˙ m + I p State: u = 6 m 7 Flux: f ( u ) = 6 5 ; 7 6 7 4 5 4 m > E  ( E + p ) Here  2 R + denotes the density, m 2 R d is the momentum, and E 2 R + is the total energy of the system. The pressure p 2 R + is determined by an equation of state, e. g., „ « E ` j m j 2 Polytropic gas equation: p ( u ) := ( ‚ ` 1 ) ; ‚ 2 ( 1 ; 5 = 3 ] : 2  6 Kronbichler, Maier, Tomas | Parallelization of “stencil-based” methods

Compressible Euler equations Variational principle viable ??? Solution theory Testing with u : We call u a viscosity solution if 1 d " ! 0 u " ; ‚ ‚ 2 “ ” u = lim ‚ ‚ u ‚ div f ( u ) ; u = 0 : L 2 Ω + ‚ L 2 (Ω) 2 d t | {z } with ??? u " t + div f ( u " ) = " ∆ u " : › no energy estimate available › no (quasi) best approximation ` ! global existence and uniqueness › etc. of viscosity solutions is an open problem. ` ! none of the classical finite element toolbox available. . . So what can we do? 7 Kronbichler, Maier, Tomas | Parallelization of “stencil-based” methods

Compressible Euler equations What can we meaningfully expect? Invariant set Conservation If u is a viscosity solution, then Conservation of mass, momentum and total energy: u ( x ; t ) 2 B 8 x 2 Ω ; 8 t – 0 ; d ˆ ˆ u d t + n ´ f ( u ) d o x = 0 : where u 2 B implies d t Ω @ Ω › positivity of density:  > 0 › positivity of internal energy: E ` j m j 2 > 0 2  › local minimum principle on specific ` ! Robustness: We want a numerical entropy: scheme that is conservative and maintains the invariant set. s ( u ) – min x 2 Ω s ( u 0 ( x )) 8 Kronbichler, Maier, Tomas | Parallelization of “stencil-based” methods

Compressible Euler equations How does (a subset of) the community judge our code? Picture norm comparison of flow characteristics to (known) experimental results. 9 Kronbichler, Maier, Tomas | Parallelization of “stencil-based” methods

Compressible Euler equations To summarize: › Robustness: conservative and invariant domain preserving Most importantly: This implies that our scheme will never crash (no matter what mesh was used, or (admissible) initial data was prescribed) › Accuracy (not discussed): second order convergent in case of known, smooth viscosity solutions › Scalability: Large-scale 3D computation 10 Kronbichler, Maier, Tomas | Parallelization of “stencil-based” methods

step-69 11 Kronbichler, Maier, Tomas | Parallelization of “stencil-based” methods

step-69 › Guermond & Popov, Invariant domains and first-order continuous finite element approximation for hyperbolic systems , SIAM J. Numer. Anal. 54(4):2466-2489. (Grossly oversimplified) formal derivation P Let u h ( x ; t ) = j ’ j ( x ) U j ( t ) be a finite element approximation: d “ ” “ ” ’ i ; u h + ’ i ; div f ( u h ) = 0 : d t Forward Euler: U n + 1 ` U n X X “ “ ”” j j U n ’ i ; div f = 0 : m ij + j fi n j 2I ( i ) j 2I ( i ) Lump mass matrix and approximate flux, f ( u n P j f ( U n h ) ı j ) ’ j : U n + 1 ` U n X ˆ ˆ “ ” i i U n m i + f ´ c ij = 0 ; m i = ’ i ; c ij = ’ i r ’ j : j fi n Ω Ω j 2I ( i ) 12 Kronbichler, Maier, Tomas | Parallelization of “stencil-based” methods

step-69 › Guermond & Popov, Invariant domains and first-order continuous finite element approximation for hyperbolic systems , SIAM J. Numer. Anal. 54(4):2466-2489. Scheme U n + 1 ` U n X X “ ” i i U n d n ij U n ´ c ij ` m i + f = 0 ; j j fi n j 2I ( i ) j 2I ( i ) where ˆ ˆ ’ i r ’ j ; m i = ’ i ; c ij = Ω Ω and where we have introduced a “graph viscosity” d n ij as stabilization term. The correct construction of d n ij is key for robustness 13 Kronbichler, Maier, Tomas | Parallelization of “stencil-based” methods

step-69 Compute two characteristic propagation speeds associated with U n i or U n j : v v " p ˜ ` ˜ " # # u p ˜ ` ˜ u p n p n t 1 + ‚ + 1 t 1 + ‚ + 1 u u j i – 1 ` ( U n i ; p ˜ ) := ˜ u n c n ; – 3 + ( U n j ; p ˜ ) := ˜ u n c n ` ˜ + ˜ ; i i p n j j p n 2 ‚ ˜ 2 ‚ ˜ i pos j pos and a two-rarefaction pressure: 2 ‚ 0 1 ‚ ` 1 ` ˜ ´ c j ` ‚ ` 1 ˜ c i + ˜ u j ` ˜ u i p ˜ ( U n 2 i ; U n B C ˜ j ) = ˜ ; p j @ A ` ˜ ´ ` ‚ ` 1 p i 2 ‚ + ˜ ˜ c i c j p j ˜ and a monotone increasing and concave down function p 8 2 ( p ` ˜ p ) > if p – ˜ p ; > > q > ˆ ˜ <  ˜ ( ‚ + 1 ) p + ( ‚ ` 1 ) ˜ p ( p ) := f ( U n i ; p ) + f ( U n j ; p ) + ˜ u j ` ˜ f ( U ; p ) := u i ; > » – > 2 ˜ ‚ ` 1 c > > 2 ‚ ` 1 ( p = ˜ p ) ‚ ` 1 ; otherwise : : Now:  ˜ p ˜ ( U n i ; U n “ˆ ” j ) if ( p max ) < 0 ; ˜ ˆ ˜ p ˜ := ˜ – 1 i ; p ˜ ) – 3 j ; p ˜ ) ` ( U n + ( U n – max = max neg ; ; p ˜ ( U n i ; U n min(˜ p max ; ˜ pos j )) otherwise : 14 Kronbichler, Maier, Tomas | Parallelization of “stencil-based” methods

step-69 Scheme „ « fi n X X “ ” U n + 1 = U n U n d n ij U n + f ´ c ij ` i i j j m i j 2I ( i ) j 2I ( i ) where ˆ ˆ “ ” d n – max ( n ij ; U n ˜ i ; U n j ) j c ij j ; ˜ – max ( n ji ; U n j ; U n ’ i ; c ij = ’ i r ’ j ; i ) j c ji j m i = ij = max Ω Ω I ( i ) — all column indices coupling to i . Stencil based: compare to matrix-vector multiplication X y i = a ij x j ; j 2I ( i ) . . . but highly nonlinear! 15 Kronbichler, Maier, Tomas | Parallelization of “stencil-based” methods

step-69 euler_step compute off-diagonal d n // Step 1: ij : for i = 1 , . . . , N do for j 2 I ( i ) , j 6 = i do “ ” d n – max ( n ij ; U n ˜ i ; U n j ) j c ij j ; ˜ – max ( n ji ; U n j ; U n max i ) j c ji j ij compute d n ii and fi n : // Step 2: fi n + 1 for i = 1 , . . . , N do “ ” ` P d n j 2I ( i ) ; j 6 = i d n m i ij , fi n min fi n ; ` c cfl 2 d n ii ii // Step 3: perform update for i = 1 , . . . , N do for j 2 I ( i ) do “ P ” ` ´ ´ c ij ` P U n + 1 fi n U n U n j 2I ( i ) d n ij U n + j 2I ( i ) f i i m i j j // Step 4: MPI synchronization U .update_ghost_values() 16 Kronbichler, Maier, Tomas | Parallelization of “stencil-based” methods

step-69 Parallelization approach ( step-69 ) › MPI parallelization via parallel::distributed::Triangulation LinearAlgebra::distributed::Vector › Precompute m i and c ij . (This is the only place where finite elements enter) › Distribute work onto threads with parallel::apply_to_subranges() › Avoids global to (MPI) local index translations. › (Asynchronous IO.) 17 Kronbichler, Maier, Tomas | Parallelization of “stencil-based” methods

MPI local index translations: Asynchronous IO: 18 Kronbichler, Maier, Tomas | Parallelization of “stencil-based” methods

19 Kronbichler, Maier, Tomas | Parallelization of “stencil-based” methods

. . . and beyond 20 Kronbichler, Maier, Tomas | Parallelization of “stencil-based” methods

. . . and beyond „ ”« X “ ” “ U L ; n + 1 j ) ´ c ij + d L ; n ` U n ` f ( U n U n j ` U n = fi n : m i i i ij i j 2I ( i ) ¸ n i + ¸ n „ ”« X X “ ” “ U H ; n + 1 j ) ´ c ij + d L ; n j ` U n ` f ( U n U n j ` U n m ij = fi n ; j j ij i 2 j 2I ( i ) j 2I ( i ) with a suitable indicator ¸ i . › The low-order U L ; n + 1 is robust, the high-order U H ; n + 1 is not. i i › Key idea: Invariant-domain preserving convex limiting 3 X U H ; n + 1 ` U L ; n + 1 P n P n = ij ; where ij = : : : i i j 2I ( i ) X ` U L ; n + 1 U n + 1 l n ij P n 0 » l n = ij ; ij » 1 : i i j 2I ( i ) 3 Guermond, Nazarov, Popov, Tomas, Second-order invariant domain preserving approximation of the Euler equations using convex limiting , SIAM J. Sci. Comput. 40 (2018) 21 Kronbichler, Maier, Tomas | Parallelization of “stencil-based” methods

Coming soon. . . https://github.com/conservation-laws/ 23 Kronbichler, Maier, Tomas | Parallelization of “stencil-based” methods