Transcript Level-Set Methods for Multimaterial Flows in Radiative

A Level-Set Method for Multimaterial
Radiative Shock Hydrodynamics
David Starinshak and Smadar Karni
Department of Mathematics, University of Michigan
Center for Radiative Shock Hydrodynamics (CRASH)
MULTIMAT 2011, September 5-9, Arcachon, France
Funding: DoE NNSA-ASC grant DE-FC52-08NA28616, NSF Award DMS 0609766 and
University of Michigan Rackham Travel Grant.
Goals and Outline
1. Develop a 2D multimaterial rad-hydro code which employs
state-of-the-art level set technology
2. Address concerns of species mass conservation and spurious
evaporation in the context of HEDP and rad-hydro
OUTLINE:
• Physics: Laser-driven HEDP experiments
• Mathematics: Model equations
• Numerics: Solver and Discretization Strategies
• Test Problems and Results
• Future Work and Extensions
The CRASH Experiment
Laser
Radiative
Shock
( ~ 3.8 kJ )
( ~ 200 km/s )
Xenon
Plastic
Beryllium
Key Physical Features:
1.
Strong radiative shock front
Ma ~ 102-103, radiation heats upstream media
2.
High energy density regime
ionization,
3.
Strong radiation-hydro coupling
4.
Multimaterial
> 1 g/cc, Tmat > 1 keV, Trad > 50 eV
prad ~ pmat , non-equilibrium
EOS and opacities vary across sharp interfaces
*X-ray
Radiograph from CRASH experiment (FW Doss, 2011)
Wall Ablation & Entrainment Problem
Wall Ablation
shocked
beryllium
unperturbed
xenon
• Irradiated wall heats up
and ablates
• Sends wall shock into Xe
• Thick-thin Interface
• Rad transport depends on
interface geometry
hydro
shock
plastic
( ~ 1 ns )
radiative precursor
entrained xenon
Entrainment
shocked
beryllium
plasma
• Shear wave instability
driven by shock front
• Interface dynamics drive
xenon entrainment
hydro
shock
dense shocked xenon
( ~ 10 ns )
Model Equations
Multi-Group Flux-Limited Diffusion with Ionization
Simplifications and Amendments
Monoenergetic “Gray” Radiation
Zero Explicit Ionization
Reduced Model:
Zero Electron Heat Flux
Single Opacity
Model Equations
Conservation of mass, momentum, and (material + radiation) energy
System Closure:
Opacity Model:
Material Parameterization:
Adiabatic index
Ideal gas constant
Opacity coefficient
Model Equations
Source
Blackbody Emission
Radiation Heat Flux
Absorption
•
Exchanges energy between matter and radiation
•
Attempts to equilibrate system:
The Level Set Model for Interfaces
Key Assumption
Materials are immiscible over timescales of interest
Sign of level curve determines material type
Level Curves advect at local flow velocity
Material-dependent quantities are discontinuous functions of the level curve
Numerical Solver Summary
Operator Splitting:
Explicit : 2nd order multimaterial HLL solver of MUSCL-Hancock type
•
•
Nonconservative products discretized to conserve total energy
Material designations determined from updated level set
Implicit : nonlinear root-finding + finite volume diffusion solver
•
•
Newton-Raphson iteration on T and ER
2nd order Crank-Nicolson discretization in time
Level Set Discretization Strategies
Initialize as signed distance function (satisfies Eikonal condition:
)
Many Attractive Methods Available
1.
2.
3.
4.
5.
6.
Black box hydro solver
HOUC
Particle LS
LS-CIR / Semi-Lagrangian
Fast Marching
ENO/WENO for H-J Equations
Advection Form
[HLL, upwind, PPM, etc.]
[Nourgaliev & Theofanous JCP 2006]
[Enright et al JCP 2002]
[Strain JCP 2000]
[Adalsteinsson & Sethian JCP 1999]
[Osher & Shu SIAM 1991]
Hamilton-Jacobi Form
Reinitialize:
Reinitialize:
where
where
Extension Velocity H-J
Extend normal velocity:
Efficient Implementation: Update only in a narrow band around interfaces
Interfaces and Material Terms
Sub-Zonal Interface reconstruction
• bi-linear interpolation from nodal level set values
• Interface normal obtained from interpolant
Sharp Volume Fractions
• Integrate Heaviside function of interpolant
Strategy for Material Terms
• Weighted by the computed volume fraction for each species
Diffusion Solver Discretization
Nonlinear, material-dependent, flux-limited heat conduction:
Coefficient
Flux limiter
Radiative Knudson number
Discretization: Crank-Nicolson in time; Finite volume, centered-differencing in space
Diffusion coefficient constructed cell faces
• Arithmetic average of ER
• Single-material faces:
-
Harmonic avg of opacity
(Arithmetic of DR)
• Two-material faces:
-
Arithmetic avg of opacity
(Harmonic of DR)
• “Mixed” cells treated as separate material
Diffusion Solver Discretization
Nonlinear, material-dependent, flux-limited heat conduction:
Coefficient
Flux limiter
Radiative Knudson number
Discretization: Crank-Nicolson in time; Finite volume, centered-differencing in space
Diffusion coefficient constructed cell faces
• Arithmetic average of ER
• Single-material faces:
-
Harmonic avg of opacity
(Arithmetic of DR)
• Two-material faces:
-
Arithmetic avg of opacity
(Harmonic of DR)
• “Mixed” cells treated as separate material
Spurious Pressure Oscillations
Pressure equilibrium not respected across diffused material fronts
•
•
•
Well-understood phenomenon for multi-fluid systems
Consequence of prescribing a mixed-cell EOS across (numerically) diffused interfaces
Does not occur in single-fluid systems
Oscillations
•
•
•
•
Occur at every time level
Occur in 1st order solvers
Not removed by mesh
refinement
Not removed by highorder solvers
Better control over pressure values is required near interfaces
Remedies:
Ghost Fluids
R Fedkiw et al, JCP, 1999
Pressure Evolution
S Karni, JCP, 1994
Single-fluid Solver
S Karni & R Abgrall, Proceedings, 2001
Single-Fluid Multimaterial Solver
Interface
Left State
Material 1
Right State
Material 2
Resolve waves at cell boundary so that WL and WR “see” the same fluid on the other side
Properties
Total mass and momentum perfectly conserved
(Material + radiation) energy essentially conserved
• Conserved almost everywhere
• Errors are small, do not accumulate, and tend to zero with mesh refinement
• NOTE: Material and radiation energies not individually conserved
S Karni & R Abgrall, Ghost-Fluids for the Poor: A Single Fluid Algorithm for Multifluids. Proceedings ICHP. 2001.
Pressure Defects in Pure Hydro
Multimaterial
Sod Shock Tube:
•
•
Pressure defect develops across interface
Single-fluid solver: defect removed with -0.17% error in energy conservation
Pressure Defects in Rad-Hydro
Source-Dominated Shock Tube:
Conservation Error
< 0.12%
Diffusion-Dominated Shock Tube:
Conservation Error
< 0.09%
1D Wall Ablation Problem
plastic
Hot ionized gas
Radiative precursor
Beryllium
plasma
Interface
xenon
Cold dense wall
1D shock tube initial conditions:
Radiation in equilibrium:
Temperature-dependent opacity:
Boundary conditions:
- Radiating left boundary (T = 100 eV)
- Zero gradient at right boundary
1D Wall Ablation Problem
Computed Solution
Semi-analytic scaling*
Dense
shock
Heat
front
Wall
shock
Interface
*Used by permission (Drake et al, 2010, preprint)
NOTE: System loses positivity if spurious pressure oscillations not addressed
1D Wall Ablation Problem
Species Mass Conservation
Interfaces are sharp, but material fronts diffuse numerically
•
Conservation of individual species masses not guaranteed
Complications to Rad-Hydro Models:
1. Density can vary by orders of magnitude across interfaces
– Small changes in interface position
large species mass errors
2. Opacity and EOS sensitive to density changes
– Radiative transfer rates
– Overflowing bounds of EOS / opacity table
3. Spurious evaporation of entrained fluid species
Species Mass Conservation Errors
Errors manifest as mass exchange between fluid species
+ 28.8%
+ 20.7%
+ 14.3%
Reinitialize Level Curves
+ 2.78%
+ 0.42%
+ 0.18%
Summary
Accomplishments
– Verified code for 1D, 2-material rad-hydro
– Implemented and begun testing 2D, N-material hydro with “modern” LS solvers
– Characterized species mass errors in 1D
Future Work
Implementation
– 2D extensions of source and diffusion solvers
– Generalize LS methods to 3+ materials
Sub-Grid Material Resolution
– Mixed-cell diffusion solver (tensor?)
– Multimaterial source (!)
– LS-informed adaptive grids
Species Mass Conservation Errors
– Consistency conditions: VOF and LS
– Modify reinitialization methods
Acknowledgements
Smadar Karni, Eric Myra, Bruce Fryxell, Ken Powell, and Paul Drake
The entire CRASH Team
University of Michigan, Texas A&M, and Simon Fraser contributors
U.S. Department of Energy’s NNSA-ASC Program, NSF, Rackham Graduate School
MULTIMAT 2011 organizing committee
References
1.
2.
3.
R Abgrall & S Karni, Computations of Compressible Multifluids, JCP, 169, 594 (2001).
R Drake et al, Behavior of Irradiated Low-Z Walls and Adjacent Plasma (preprint, 2010).
R Fedkiw, T Aslam, B Merriman, & S Osher, A Non-oscillatory Eulerian Approach to Interfaces in Multimaterial Flows
(the Ghost Fluid Method), JCP, 152, 457 (1999).
4. S Karni, Multicomponent Flow Calculations by a Consistent Primitive Algorithm, JCP, 112, 1 (1994).
5. S Karni & R Abgrall, Ghost-Fluids for the Poor: A Single Fluid Algorithm for Multifluids. Proceedings of the 10th
International Conference on Hyperbolic Problems, Theory and Numerics. 2001.
6. C Levermore & G Pomraning, A Flux-Limited Diffusion Theory, Astrophys. J., 248: 321-334 (1981).
7. R Lowrie, J Morel, & J Hittinger, The Coupling of Radiation and Hydrodynamics, Astrophys. J., 521: 432-450 (1999).
8. R Lowrie, R Rauenzahn, Radiative Shocks in the Equilibrium Diffusion Limit, Shock Waves, 16: 445-453 (2007).
9. R Lowrie, J Edwards, Radiative Shocks with Grey Nonequilibrium Diffusion, Shock Waves, 18: 129-143 (2008).
10. Mihalas & Mihalas, Foundations of Radiation Hydrodynamics, 1983.
11. W Mulder, S Osher, J Sethian, Computing Interface Motion in Compressible Gas Dynamics, JCP, 100, 2009 (1992).
12. B van der Holst, G Toth, I Sokolov, K Powell, J Holloway, E Myra, Q Stout, M Adams, J Morel, S Karni, R Drake, CRASH:
A Block-Adaptive-Mesh Code for Radiative Shock Hydrodynamics (preprint, 2011)
Implicit Solver
Solve:
Use EOS to transform E to T:
Kinetic terms do not vary in time, and material terms are frozen at * time level:
Implicit Backward Euler discretization in time:
Vectorize:
FV Spatial Discretization
- Tridiagonal matrix in 1D
- Banded tridiagonal in 2D
- Arithmetic averaging gives DR*
at cell boundaries
Implicit Solver (Cont)
Solve the nonlinear vector-operator equation:
Equivalently, perform nonlinear root-finding:
where
Newton-Raphson Iteration:
NOTE:
Inverting the Jacobian matrix
Initial iterate:
amounts to inverting the matrix from the
discretized diffusion operator. This is done
at each iteration.
Alternatives: preconditioned CG or GMRES
Pressure Evolution Solver
Away from Interfaces
Near Interfaces
Solve material energy equation:
Solve material pressure equation:
Form pressure using EOS:
Form material energy using EOS:
Solving for pressure directly ensures its continuity across material fronts
Solver Properties
•
•
Perfectly conserves total mass and momentum
Essentially conserves energy
– Conserved almost everywhere in computational domain
– Conservation errors are small and do not accumulate
– Errors decrease with mesh refinement
NOTE: Material and radiation energies not strictly conserved for this system; their sum is.
Spurious Pressure Oscillations
Pressure equilibrium not respected across material fronts
•
•
•
•
Pressure computed from diffused hydro quantities using the EOS
Two EOS exist across interfaces: a mixed-cell EOS is needed
As interface moves: mixed-fluid EOS becomes inconsistent with hydro variables
A pressure defect develops, sending signals into neighboring cells
Oscillations
•
•
•
•
Occur at every time level
Occur in 1st order solvers
Not removed by mesh
refinement
Not removed by highorder solvers
Better control over pressure values is required near interfaces
Remedies:
Ghost Fluids
R Fedkiw et al, JCP, 1999
Pressure Evolution
S Karni, JCP, 1994
Single-fluid Solvers
S Karni & R Abgrall, Proceedings, 2001
Single-Fluid Multimaterial Solver
Interface
Left State
Right State
Material 1
Material 2
Waves updating WL
Start with primitive variables
Form WL and WR using EOS from Mat. 1
Waves updating WR
Start with primitive variables:
Form WL and WR using EOS from Mat. 2
Properties
•
•
Perfectly conserves total mass and momentum
Essentially conserves energy
– Conserved almost everywhere in computational domain
– Conservation errors are small and do not accumulate
– Errors decrease with mesh refinement
NOTE: Material and radiation energies not strictly conserved for this system; their sum is.