Implementing ALE Motion in a Discontinuous Finite Element Hydro - PowerPoint PPT Presentation

Implementing ALE Motion in a Discontinuous Finite Element Hydro Code* Manoj K. Prasad, Jose L. Milovich, Aleksei I. Shestakov, David S. Kershaw, and Michael J. Shaw Lawrence Livermore National Laboratory, Livermore, CA 94550, USA *Work performed under the auspices of the U.S. Department of Energy by the Lawrence Livermore National Laboratory under contract number W-7405-ENG-48

Motivation ALE Hydro code that combines the accuracy of a higher order Godunov scheme with the unstructured mesh capabilities of finite elements which can be explicitly evolved in time.

Brief History of Discontinuous Galerkin (DG) Finite Element Methods • P. LeSaint & P. A. Raviart: rigorous error analysis & convergence rates for DG finite element solution of steady state linear neutron transport equations, 1974 . D. S. Kershaw & J. A. Harte , implemented a fully implicit 2D time dependent linear neutron transport on triangular mesh – solving (3 x 3) block lower triangular linear system at each time step – using the partial ordering scheme of LeSaint & Raviart with no cycles, LLNL 1993. • G. Chavent & G. Salazano: applied DG methods to time dependent nonlinear porous media equations with explicit Euler time stepping, 1982 . • G. Chavent & B. Cockburn: linear error analysis and slope limiters for approximating shocks for scalar conservation laws, INRIA 1989. • B. Cockburn & C.-W. Shu: combined Runge-Kutta explicit time discretization with DG space discretization for scalar conservation laws, 1991.

Model • Godunov scheme: 2nd order. • Piecewise Linear Finite Element discretization. • Roe upwind surface flux & Harten-Hyman Entropy fix • Explicit time stepping : 2nd order Runge-Kutta. • 3D shock stabilization : VanLeer “minmod” by Quadratic Programming + Lapidus artificial viscosity. • Algorithm appears in : Computer Methods in Applied Mechanics and Engineering 158 , 81-116 (1998).

Implementation • 3D Unstructured Mesh : tets, pyramids, prisms & hexes. • 1-Step Arbitrary Lagrangian-Eulerian (ALE) moving mesh . • 3D Geometries : Cartesian, Cylindrical & Spherical. • Object oriented C++ design untangles mesh from physics. • Parallelized using Domain Decomposition & MPI message passing. • Portability: code runs on uniprocessor workstations and massively parallel platforms with distributed and shared memory. • Integrated with electron, radiation diffusion transport & laser ray tracing for 3D ICF simulations: Computer Methods in Applied Mechanics and Engineering 187 , 181-200 (2000).

Allowable Cell Types Pyramid Prism Tetrahedron Hexahedron Mesh connectivity requires that cells share like-kind faces. No slide lines allowed.

3D ALE Moving Mesh Hydro Equations Conservative form: Grid Motion: ∂ ∂ A F + = α α i S α ∂ ∂ 0 t x 0 Independen t variable s : t, x i i ( )   ρ ρ − g   v v  0  0 ALE moving mesh : x = x ( x , t ) , j j  ( )      i i j ρ δ + ρ − g ρ v P v v v G       1 j 1 j j 1 1 ( ) = = 0 x x ( t 0 )       = ρ = η δ + ρ − g = ρ A J v F P v v v S J G i i α α α i ji 2 j 2 j j 2 2     ( )   ∂ x ρ δ + ρ − g ρ v P v v v G       = g Grid velocity : v i , 3 j 3 j j 3 3 ( )       i ∂ ρ + − g ρ t E P v E v v v G       i j j i i r ∂ x ρ = Density, v = Velocity, P = Pressure ≡ j η = − = 1 J , JJ , J | J | ij ji ji ij ∂ 0 1 x + ρ E = v v I ( , P ) = Total Energy / Mass i i i 2 g = Eulerian grid : v 0 , ρ I ( , P ) = Internal Energy / Mass (EOS) i ≈ g Lagrangian grid : v v = G External Body Force / Mass i i i

Finite Element Discretization r r ρ Piecewise Linear Dependent Variables Q : , v , P, x r ξ Isoparamet eric ( ) representa tion in each element r n ∑ φ ξ ≡ Q = Q ( ), Q Q at node , l l l l = 1 l φ = = Linear basis fns 1 at node and 0 at all others l l n = number of nodes in element r ρ , v , P are in general discontinu ous across element faces r x is required to be continuous across element faces

Roe Upwind Surface Flux Γ Galerkin Equations for each element K with surface : K r r ∫ ∫ ∂ φ φ ∇ A (d x ) = - F (d x ) α α t i i l l K K r ∫ ∫ ∇ φ φ Γ = F (d x ) - N F (d ) α α i i i i K l l K ∂ K 1 4 4 2 4 4 3 use Roe's upwind flu x here Γ where N is outward pointing normal to surface from the (-) side i K r r − ≠ + to the (+) side. In general (discontin uous finite elements ) A A [ ] { } 1 ≡ − + + + Λ − − − + Roe * * * 1 N F N F N F R sign ( ) R N F N F α α α α α i i i i i i i i i i 2 Λ λ = − − ± * * * g * g * where is a diagonal eigenvalue matrix , N ( v v ), N ( v v ) c i i i i i i * * c is the sound speed, R the diagonaliz ation matrix defined by the Roe Average : ∂ * N F − + − − + ≡ * − * ≡ ≡ Λ − α * * * 1 N F N F W ( A A ), W i i R R α α αβ β β αβ i i i i ∂ * A β Harten - Hyman Entropy fix : − + λ * → λ * + ε − λ * ε = λ * − λ λ − λ * | | | | max( 0 , | |) , max( 0 , , )

+ Roe Averages: N * _ Roe Flux: =

Boundary Conditions Ghost states on outside of boundary faces : ρ ghost = ρ int ghost = int , P P • Lagrange BC: specified ghost = int − int − g v v 2 N [ N ( v v )] i i i j j j pressure or normal Roe Average normal velocity on boundary face : velocity with no mass flux * = g N v N v i i i i • We use “ghost” states on Roe flux on boundary : the outside of boundary   0 faces   bndry N P   1   = Roe bndry • Roe flux from the “ghost” N F N P ,   α i i 2   bndry N P and interior states is the 3    g bndry  N v P   i i required boundary flux = + ρ − + − bndry int int int g * int g P P [ N ( v v )][ c N ( v v )] j j j i i i g If N v given on bndry, use it for N v i i i i bndry g If P given, solve for N v i i

Almost Lagrangian Grid Velocity • Grid Position and Velocity g Least Squares estimate of fluid velocity at node n for v : in must be continuous while fluid velocity is generally ∑ − ≡ g 2 g Minimize [ N v Y ] , Y N v if in f f if if n n n n n discontinuous across elements. { f } n ≡ where : f all faces that share node n n • We use a “least squares” estimate to determine the grid Y is determined by requiring no mass flux across f or BC : f n n velocity from the fluid velocity. = Roe For interior faces f : solve N F [ Y ] 0 for Y n i 1 i f f n n • We can impose 1,2, or 3 linear For boundary faces f : Y is given by normal velocity BC n f n constraints on the grid velocity. bndry or determined from P • For 1D this gives an exact Lagrangian code. In general we = g We can also impose n 1 , 2 , 3 linear constraint s on v in the c in get an “almost” Lagrange code. least squares estimate.

Roe Mass Flux Normal Grid Velocity giving Zero Mass Flux on Face f Least squares estimate of velocity at vertex p from all faces f around it

Explicit 2nd order Runge-Kutta Time Advancement • The DGM equations, for ∆ ≤ ≈ t CFL t , CFL . 3 Courant each element K , reduce to a ∫ 3 d x system of n ODE s for the = t Min K Courant ∫ λ * Γ K | | d moments: max ∂ K λ * = * − g * − g − * * − g + * | | max(| N ( v v ) |, | N ( v v ) c |, | N ( v v ) c |) ∫ = ϕ 3 M A d x max i i i i i i i i i l l K r n ∑ φ ξ Piecewise Linear Variables Q = Q ( ), • We use 2nd order Runge- l l = 1 l ~ r Kutta (RK) time stepping. = ρ Primitive : B ( , v , P ) , Conservati ve : A ~ ∫ − = Φ 1 Φ = ϕ ϕ 3 A M , d x • Courant time step control l lm m lm l m K with CFL number of about ~ ~ − < > ≈ < − > − < > 1 B B T ( A A ) , α α αβ β β 1/3. ∫ 3 Q d x ~ ∂ A • Cockburn & Shu, Math. = α < >= K T , Q αβ ∫ ∂ B 3 d x Comput. 52 (1989) 411. β K