This work was performed under the auspices of the U.S. Department of Energy

1. NUMERICAL MODELLING OF EXPLOSIONS IN UNDERGROUND CHAMBERS USING INTERFACE TRACKING AND MATERIAL MIXING Numerical Methods for Multi-Material Fluid Flows September 5th-8th, 2005 St. Catherine’s College, Oxford, UK Benjamin T. Liu and Ilya Lomov Energy and Environment Directorate Lawrence Livermore National Laboratory This work was performed under the auspices of the U.S. Department of Energy by University of California, Lawrence Livermore National Laboratory under contract No. W-7405-Eng-48.

2. Gas-phase mixing Droplets or bubbles Introduction: Sharp and Diffusive Interfaces Sharp Interfaces Diffusive “Interface”

3. Outline • Treatment of sharp interfaces • Treatment of diffusive interfaces • Simulations combining sharp and diffusive interfaces

4. GEODYN • High-order Godunov Eulerian code • Able to model large deformations • Able to capture shocks • Treatment of interfaces is important • Structured rectangular grids with adaptive mesh refinement • Multi-material with a fully integrated stress tensor • Characteristic tracing of stress tensor • Acoustic approximation for shear waves • Flexible material library • Analytic and tabular EOS • Wide range of constitutive models • Especially designed to model the response of geophysical media • Includes a variety of yield strength models

5. Outline: Treatment of Sharp Interfaces • Treatment of sharp interfaces • Standard treatment • Hybrid energy update • Stress equilibration • Treatment of diffusive interfaces • Simulations combining sharp and diffusive interfaces

6. Standard Treatment of Sharp Interfaces • Volume-of-fluid approach • High order interface reconstruction • used to calculate transport volumes • preserves linear interface during translation • Thermodynamics based equations for the mixed cell update

7. Standard Pressure Relaxation Scheme • Iterative adjustment of volume fractions • - bulk modulus • - numerically or physically based limiter

8. Problems with the Standard Treatment • Conservative energy update is not robust for mixed materials • Materials with drastically different properties (r, K, etc) • Most severe when kinetic energy is large relative to internal energy • Pressure relaxation unsuitable when strength in mixed cells is important • Relaxation scheme ignores strength • Effective strength for material in mixed cells with fluid is zero

9. Hybrid Energy Update • Conservative equation • Non-conservative equation • Hybrid (conserves energy)

10. Hybrid Energy Update Test: Aluminum Flyer Plate (3 km/s) in Air GPa GPa ConservativeNon-conservativeHybrid mm mm Position of flyer plate

11. no strength no strength no strength no strength Pressure Equilibration in Mixed Cells with Strength • Pressure relaxation ignores strength • Problem in mixed cells with solid and fluid • Solid w/strength and fluid w/o strength • Pressures in solid and fluid are equal • Mixed cells containing fluid have no strength • Material is weaker near interfaces • Introduces strong mesh dependence • Results in cells containing differing solids w/strength are also wrong fluid no strength P solid w/ strength

12. Normal Stress Equilibration in Mixed Cells with Strength • Equilibrate normal stress instead of pressure • Information within mixed cell insufficient • Need to calculate T’nn • Requires elastic hoop strain (ett) • Solution: Use properties from single-material cells in the vicinity of the mixed cell • Consistency conditions: • Stress normal to interface is continuous • Elastic strains in transverse direction taken from single-material cells • Interfacial shear stress can be calculated using a friction law • Fall back to pressure relaxation scheme when: • No single-material cells in the direction of normal • More then 2 materials in the cells ett2 Tnn ett1

13. Stress Relaxation • Elastic hoop strain in the single material cell: • Normal component of the stress deviator in the mixed cell: Relax total normal stress in each material to the average across the cell: Constraint modulus

14. Stress Relaxation Test: Aluminum Flyer Plate (3 km/s) in Air Pressure Relaxation Stress Relaxation Pressure Normal Stress (-Tnn) Pressure Normal Stress (-Tnn) For elastic 1D strain: for the aluminum plate

15. Test Problem - Cylindrical Cavity Expansion Pressure Relaxation Results Vacuum Aluminum + 1 bar 0 -1 bar Air 1 bar Radial Stress

16. Stress Relaxation Pressure Relaxation Radial Stress + 1 bar 0 -1 bar Hoop Stress

17. Problems Requiring Stress Equilibration • Quasi-static solution after initial waves have passed • Cavity expansion • Blast or impact loading of deeply buried structures • Overall response driven by deformation in the mixed zones • Fast moving solids undergoing “slow” deformation • Void nucleation and growth under positive pressure • Pressure relaxation will cause voids to immediately close • Strength in the material surrounding voids is important

18. Outline: Treatment of Diffusive Interfaces • Treatment of sharp interfaces • Standard treatment • Hybrid energy update • Stress equilibration • Treatment of diffusive interfaces • Track mass fractions of components • Use effective mixture gamma • Iterate for real materials • Simulations combining sharp and diffusive interfaces

19. Diffusive Material Interface Treatment • Consider materials that diffuse into one another • Separate components within a single computational “material” • Mass fractions (with total r, e) sufficient to reconstruct mixture state variables • Should enforce pressure and temperature equilibrium between components

20. Ideal Gas Mixing Ideal Gas Mixture • Internal energy • Effective molecular weight • Effective gamma

21. Ideal Gas Pressure Calculation Molecular mixture Droplets or bubbles fi: fraction of mixture volume occupied by component i Ideal Gas Pressure Applying Dalton’s Law: For an ideal gas: Enforcing pressure equilibrium: Pressure for ideal gas mixture independent of spatial component distribution

22. Non-Ideal Equations of State Non-Ideal Equations of State • Define an effective (component) gamma: • a constant for ideal gases • a relatively slowly varying parameter for a wide range of densities and temperatures for many real materials • Calculate pressure based on mixture gamma: • Similarly calculate temperature: • Zeroth order approximation: ei = e, ri = mir • Yields correct averages for ideal gases

23. Iterative Refinement for Non-Ideal Gases Non-Ideal Equations of State • Initial guess: ei = e, ri = mir • Iterate on component densities and energies • Iterative estimate for energy • Pressure relaxation scheme for density • Two-phase region may be singular and non-convergent • Solution has oscillations • Saurel & Abgrall (1999), Karni (1994), et al • Zeroth order approximation good when gamma is changing slowly

24. Outline: Simulations • Treatment of sharp interfaces • Standard treatment • Hybrid energy update • Stress equilibration • Treatment of diffusive interfaces • Track mass fractions of components • Use effective mixture gamma • Iterate for real materials • Simulations combining sharp and diffusive interfaces • Mixing and heating in underground chambers • 2D simulation • Large-scale 3D simulation

25. Explosions in Underground Chambers • Fundamental study of multi-material mixing and heating • Demonstrate combination of diffuse and sharp interfaces • No explicit subgrid model • Turbulence implicitly modeled by truncation errors • Monotone Integrated Large Eddy Simulation (MILES) [J. Boris, 1992] • Physical rationale by L. Margolin and W. Rider in 2002 • Examine heating of water contained in underground chambers • Consider different modes of heating after an explosion • Shock heating (PdV work) • Convective mixing • Measure degree of heating by fraction of water above 650K • Critical point for water • Vapor and liquid indistinguishable

26. 2D Problem Setup 167 GJ source 4 tons water 1.5mm steel liner

27. Density/Temperature Profiles

28. Convective mixing dominates heat transfer Expansion and cooling Shock heating Temperature Distribution T < 650 K 650 K  T < 2600 K T  2600 K

29. 3D Calculation • Run on LLNL’s Thunder supercomputer • Utilized 960 nodes (3840 Itanium CPU’s) • Used almost 1 TB of total memory • Largest problem of its kind to date • Two levels of refinement • 16.8 million zones (6 cm resolution) on the coarse level • ~160 million zones (1.5 cm resolution) on the fine level 0.5 m DOB 60 m x 10 m x 10 m chamber

30. 0.5 m roof 8.4 TJ 60 m x 10 m x 10 m chamber 200 tons water

31. Conclusions • Improved treatment of sharp interfaces • Hybrid energy update robustly captures shocks while conserving energy • Stress equilibration improves modelling of material with strength • Implemented simple treatment of diffusive interfaces • Store mass fractions and calculate an effective gamma • Zeroth order approximation sufficient for many applications • Successfully simulated problems including sharp and diffusive interfaces • Performed both 2D and 3D simulations • Examined mixing and heating of explosions in bunkers