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
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
DG finite element solution of steady state linear neutron transport equations, 1974.
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.
porous media equations with explicit Euler time stepping, 1982.
approximating shocks for scalar conservation laws, INRIA 1989.
with DG space discretization for scalar conservation laws, 1991.
Programming + Lapidus artificial viscosity.
Engineering 158, 81-116 (1998).
parallel platforms with distributed and shared memory.
tracing for 3D ICF simulations: Computer Methods in Applied Mechanics and Engineering 187, 181-200 (2000).
Tetrahedron Pyramid Prism Hexahedron
( ) ( ) ( ) ( ) ( )
Mass / Force Body External (EOS) Mass Energy / Internal = ) , ( Mass Energy / Total = ) , ( 2 1 = Pressure = P Velocity, = Density,
3 2 1 3 3 2 2 1 1 3 2 1
= + = = − + − + − + − + − = = = ∂ ∂ + ∂ ∂
i i i i i g j j i g j j j g j j j g j j j g j j ji i i i
G P I P I v v E v G v G G G J S v v E v P v v v P v v v P v v v P v v F E v v v J A S x F t A ρ ρ ρ ρ ρ ρ ρ ρ δ ρ δ ρ δ ρ η ρ ρ ρ ρ ρ
α α α α α α
r
i g i g i ij ji i j i i i j i i i
v v v J J JJ x x J t x v t x x t x x x x ≈ = = = ≡ ∂ ∂ = = =
−
: grid Lagrangian , : grid Eulerian | | , , , : velocity Grid ) ( , ) , ( = : mesh moving ALE t, : s t variable Independen
1 ji ij g i
η ∂ ∂
faces element across continuous be to required is faces element across
discontinu general in are P , , element in nodes
number =
all at and node at 1 fns basis Linear = , node at ), ( = element each in tion representa ) ( eric Isoparamet P, , , : Variables Dependent Linear Piecewise
1 =
x v n Q Q Q Q x v Q
n
r r l l r r r r
l l l l l
ρ φ ξ φ ξ ρ = ≡
+ − ≠
Γ Γ ∇ ∇ Γ
A A N F N F A
i x here upwind flu use Roe's i i i i K t
r r 4 4 3 4 4 2 1 r r r
l l l l
) elements finite uous (discontin general In side. (+) the to side (-) the from surface to normal pointing
is where ) (d
x (d F = ) x (d
) x (d : surface K with element each for Equations Galerkin
K K K K i i K K α ∂ α α α
φ φ φ φ ∂
) , , max( , |) | , max( | | | | : fix Entropy Hyman
), ( : Average Roe by the defined matrix ation diagonaliz the speed, sound the is ) ( ), ( , matrix eigenvalue diagonal a is where ) ( 2 1
* * * * * 1 * * * * * * * i i * * * * i * i * * i i 1 * * * i i i
λ λ λ λ ε λ ε λ λ ∂ λ
β α αβ β β αβ α α α α α α α
− − = − + → Λ ≡ ∂ ≡ − ≡ − ± − − = Λ − Λ + + ≡
+ − − + − + − + − − + −
R R A F N W A A W F N F N R c c v v N v v N F N F N R sign R F N F N F N
i i i i g i i g i i i i i i Roe i
N
g i i bndry g i i i i g i i i g j j j bndry bndry g i i bndry bndry bndry Roe i i g i i i i g j j j i i ghost i ghost ghost
v N P v N v N v v N c v v N P P P v N P N P N P N F N v N v N v v N N v v P P for solve given, If for it use bndry,
given If )] ( )][ ( [ , : boundary
flux Roe : face boundary
velocity normal Average Roe )] ( [ 2 , : faces boundary
states Ghost
int * int int int 3 2 1 * int int int int
− + − + = = = − − = = = ρ ρ ρ
α
must be continuous while fluid velocity is generally discontinuous across elements.
estimate to determine the grid velocity from the fluid velocity.
constraints on the grid velocity.
Lagrangian code. In general we get an “almost” Lagrange code.
estimate. squares least the in
s constraint linear 3 , 2 , 1 impose also can We from determined
BC velocity normal by given is : faces boundary For for ] [ solve : faces interior For :
across flux mass no requiring by determined is node share that faces all : where , ] [ Minimize : for node at velocity fluid
estimate Squares Least
1 2 } { g in c bndry f n f f Roe i i n n f n g if if f f f g in if g in
v n P Y f Y Y F N f f Y n f v N Y Y v N v n
n n n n n n n n n n
= = ≡ ≡ −
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
|) ) ( | |, ) ( | |, ) ( max(| | | | | Min 3 . ,
* * * * * max * max * 3 K
c v v N c v v N v v N d x d t CFL t CFL t
g i i i g i i i g i i i K K Courant Courant
+ − − − − = Γ = ≈ ≤ ∆
∫ ∫
∂
λ λ
each element K, reduce to a system of n ODE s for the moments:
Kutta (RK) time stepping.
with CFL number of about 1/3.
, ~ , ) ~ ~ ( , ~ ~ : ve Conservati , ) P , , ( : Primitive ), ( = Variables Linear Piecewise
3 3 1 3 1 n 1 =
∫ ∫ ∫ ∑
>= < ∂ ∂ = > < − > < ≈ > < − = Φ Φ = =
− − K K m K l lm m lm l
x d x d Q Q B A T A A T B B x d M A A v B Q Q
β α αβ β β αβ α α
ϕ ϕ ρ ξ φ r r
l l l
x d A M
K l l 3
= ϕ
limiter slope MINMOD s VanLeer' to reduces this 2 = 1D, For and : subject to ) ( 2 1 Minimizing by Find : problem g Programmin Quadratic e it thru th modify we s constraint s VanLeers' violates
1 VL max VL min 1 2 VL VL
n Q Q w Q Q Q Q Q w Q Q If
n n
〉 〈 = 〉 〈 ≤ ≤ 〉 〈 −
∑ ∑
= = l l l l l l l l l l
boundaries
nodes for boundaries general
nodes for 1 ) 1 ( : Nodes For node sharing elements all
max) min, ( : limits ion stabilizat VanLeer node at
around n fluctuatio ) 2 1 = ( let K, element each In , , : t requiremen Physical
max min, max min
symmetry Q Q Q Boundary l Q Q Q Q Q Q Q Q Q ,...n , Q Q Q E P
face bndry l l l
= = 〉 〈 + 〉 〈 − → 〉 〈 〉 〈 = 〉 〈 〉 〈 ≤ ≤ 〉 〈 〉 〈 ≡ + 〉 〈 = ≥ ω ω ω ω δ δ δ ρ
l l
l l
pressure must be enforced.
stabilization limits the fluctuations of variables around their average values within an element.
mesh generalization of VanLeer’s stabilization has a unique and simple solution.
diffusive artificial viscosity.
viscosity.
shocks to keep 2nd order accuracy in smooth regions.
Correction to the Roe Flux for interior faces.
∫ ∫ ∫ ∫ ∫
∂ ∂ ∂ + − ∂ ∂ +
Γ Γ ± − − = ≈ = ± Γ > < − > < − Γ → Γ → ∆ ∆ → ∇ ⋅ ∇ + = ⋅ ∇ + ∂
K g i i g i i K i L L Roe i Roe i L ity Vis Lapidus L t
c v v v v N l D A A l D F N F N t x D A D S F A
K K * * * * * K K K K K K cos
d d ) | ) ( | |, ) ( | ( max .3 , face the
side either
elements 2 the to refers ) (d ] [ ) (d ) (d : shocks strong near faces interior For : tion implementa DGM Our ) , ( as ) ( λ κ λ κ φ φ φ
α α α α α α α α l l l
r r r 43 42 1 r r r r
2 . at time viscosity Lapidus With 13 . at time viscosity pidus Without La : ) at vs. ( mesh
Plot ) 26 , ( at 6) (0,-70.746 velocity Normal : BC 1 , 619 . , 26 hexes
mesh 1 x 6 x 208 Limit Godunov mode Lagrangian Geometry ) , , ( Cartesian 3D )] 1 ( 10000 , 7466 . 70 , 4 [ ] , , [ ) 1 , , 1 ( ] , , [ 26 13 : , 13 : , 5/3 with gas Law Ideal : State Initial : contact no and shock strong Very : Problem Tube Shock = = = = = ≤ ≤ ≤ ≤ ≤ ≤ − − = = ≤ ≤ ≤ ≤ = = = t t z x y x z y x z y x P v P v x Right x Left v v
Right x Left x z y
γ ρ ρ γ γ
To avoid tampering with the solution in situations where accuracy is required, we implemented a hybridization scheme that smoothly interpolates from adiabatic (Rayleigh-Taylor) flows to strong shock regions
flows adiabatic s s shocks strong r r s r B B sB B r rB s B Hybrid
Godunov RK Godunov VL hybrid
near 1 with rate production entropy measures near with strength shock measures 1 , , : where ] ) 1 ( [ ) 1 ( : variables primitive
blend a Construct → → ≤ ≤ > < ≡ + − + − =
α α α α α α
various classes
pointers to access data across classes
geomtries ) ( for ) 2 , 1 , ( where ) / 1 ( , , 1 : h region wit unshocked Outer , , : h region wit shock
Central : shock expanding an by separated regions 2 : solution analytical Similar
boundary. symmetry a into ) 1 ( ity unit veloc with ) 1, ( gas cold initially
Impact pherical indrical,s planar,cyl g r t v P v P P v v v P
g s s s
= − = = − = = = = = − = = = ρ ρ ρ ρ ρ
at vs. and
Plots 3 / 64 , 64 3 / 4 sphere
Velocity Normal : BC 4 / , 2 / , 1 hexes
mesh 1 x 1 x 128 , 1 , 1 , 5/3 : with gas Law Ideal : Condition Initial gas imploding y sphericall a
stagnation simulates ) , , ( Geometry Spherical 1D = = = = = = = ≤ ≤ ≤ ≤ ≤ ≤ = − = = = φ θ ρ ρ π φ π θ ρ γ γ φ θ r P P U d Shock Spee r P v r
s s s r
− − −
→ + = → → + = → → − + = → = = = = = = = = 1 / , ) 1 /( 2 1 / , ) 1 /( 2 1 / , ) 1 /( ) 1 ( 151667 . 1 : 3 / 5 For ) / ( ) 5 / 2 ( : ) / ( : : Solution Analytical Similar
Shock. Strength Infinite ) ( , , , : Gas Cold Law Ideal ensity Constant D : Condition Initial
2 5 / 1 3 5 / 1 2 3 s s s s s s s s s s
r r U P P r r U u u r r k t E k U Velocity Shock t E k r Position Shock x E E v P γ ρ γ γ γ ρ ρ ρ γ ρ ρ δ ρ ρ γ r r 1 at time for vs. / , / , ) vs. ( mesh
Plot 3 / 5 , 4935889 . , 1 Velocity Normal : BC 2 , 125 . 1 , Hexes
Mesh 45 x 1 x 45 ) , , ( Geometry l Cylindrica 3D Mode Lagrangian = = = = = = ≤ ≤ ≤ ≤ z r P P z r E z r z r
s s
ρ ρ γ ρ π θ θ
10 Interface at Ratio Density R 1) 1)/(R
, ) / g (2 Rate Growth Linear z .00045 .001, , / 2 ) ( ) cos( ) cos(k .5 z :
Perturbati : Interface .33671717 g , 1 . 1 , d g P P 2 for , ) / ( 10 for , ) / ( 100 10 : Gas Ideal
1/2 2 / 1 2 2 x 9 / 1 9 / 1
= = + = = ∆ = = = + = = = = = ∫ + = ≤ ≤ = ≤ ≤ = = α λ α π χ λ λ π χ δ λ λ ρ λ λ ρ λ ρ γ
y x y
k k k y k x z l z z l z z l z
modes r rectangula & square 3D 2D, : for time vs. ude) Log(amplit
Plot .1 amplitude when interface
Plot : 2D for 3 : mode r Rectangula : mode Square Velocity Normal : BC Mode Lagrange 2 , 5 . , hexes 40 x 20 x 20 : Grid Geometry ) ( Cartesian 3D λ λ λ ≈ = = = = ≤ ≤ ≤ ≤
y y x y x
k k k k k z y x x,y,z
(www.t12.lanl.gov/~lagrit)
(www-users.cs.umn.edu/~karypis/metis)
) 5 / , 5 1 (cos ) , ( 5 / : 2 1 :
1
π φ θ π ϕ ± = ± = =
−
and
through plane planes azimuthal r sphere
27.46 ) max( run 1D 26, ) max( run 3D .6 .58, , 56 . at s. Density v
Plots 4 / /2, , 1 hexes, 1 x 1 x 200 run comparison ) , , ( Geometry D 1 563 . ) max( with 6 . at Density
view
:
Plot 3 / 4 : BC motion Grid ALE mode Lagrange 1, : IC Gas Law 5/3 Ideal = ≈ = ≤ ≤ ≤ ≤ ≤ ≤ = = = = = = = ρ ρ π ϕ π θ ϕ θ ρ γ t r r r Spherical z t P d Constraine v P
bndry
r
v R P
e
field radiation the
t coefficien diffusion derivative Lagrangian coupling matter to radiation laser) a to due (e.g. source external flux, heat
divergence Force, gravity ) ( ) ( ) ( ( = = = = = = − ∇ ⋅ ∇ = + + + ⋅ = ⋅ ∇ + ∂ = ⋅ ∇ + ∂ = ⋅ ∇ + ∂
ν νε ε ε ν ν ν νε ε ε ρ ρ
ρ ν ρ ρ ρ ρ ρ ρ D dt d K S H K u D u dt d d K S H E
νε E t t t
g v g F g F v v)
v
2 2 2 2 2 2 2
c e c e pe
4
r P t v r P r r r t
r r r r R r
physics packages (laser, heat conduction, hydro)
shock wave arises
KeV) ergs/(g 10 c , 4 . 1
15 v =
= γ
19 2 / 5
10 , , = = ∇ − χ χ χ χ T T
3
g/cm 001668 . =
c
ρ
c
r r r r ρ ρ ρ ρ ρ ρ 01 . ) ( , 2 ) 16 15 ( , 10 ) 8 7 (
c c
= = = ≤
3 7
erg/cm 10 67 . 6 ) , 10 / ( ⋅ = = T p p
c b
ρ
2 13 W/cm
10 875 . 2 ⋅
cm 2 .
0 =
r
ablator (5104 tets, 1246 nodes, 10,915 faces)
– Hydro: Pbdry=58. GPa (corresponds to Be at ρ=1.85 g/cc and T) – Heat Conduction: F=0 – Radiation: Energy source (correspond to Tr = 0.16 keV)
with METIS
combination.
processors.
meshes