STAGGERED MESH GODUNOV (SMG) SCHEMES FOR ALE HYDRODYNAMICS

1. STAGGERED MESH GODUNOV (SMG) SCHEMES FOR ALE HYDRODYNAMICS Gabi Luttwak1 and Joseph Falcovitz2 1Rafael, P.O. Box 2250, Haifa 31021, Israel 2Institute of Mathematics, The Hebrew University of Jerusalem, Israel

2. Hydro-Schemes • Classical Lagrange (CL) & ALE • staggered in space & time • even schemes are staggered in space with predictor-corrector time integration • pseudo-viscosity used for shock capturing • Second Order Godunov (Eulerian) • zone centered variables & gradients • gradients limited to preserve monotonic profiles • Riemann Problems (RP) solved at zone faces • We seek the best approach for an MM-ALE code • It is required to treat the flow in both Lagrangian and (Multi-material) Eulerian framework

3. Staggered Mesh Godunov - SMG SMG schemes: • staggered in space & time • or staggered in space and even-time integration • use RP to capture shocks • bridge the Lagrange-to-Godunov “conceptual gap”

4. (MM)ALE code • Geometry Update (not necessary for pure Euler) • compute areas, volumes etc. • Lagrangian Phase • CL/SMG Scheme • Definition of pressure and stress in a MM zone • Advection Phase (not present for pure Lagrange) • Second order advection for cell-centered variables • Staggered Mesh Momentum advection • SALE • HIS (Benson) • Single step /Partially-split advected-volume integration • Interface tracking: VOF (Luttwak 1984,2002) similar to Youngs scheme • Grid Motion Algorithm

5. Classical Lagrange (CL) Schemes • simple, efficient and robust • ~ second order accurate for smooth flow • on a uniform mesh with constant time steps • until grid gets distorted • pseudo-viscosity is used to capture shocks • over 3-4 zones • + linear terms required for weak shocks

6. Pseudo-viscosity • form and coefficients somewhat arbitrary • not from first principles • unphysical heating in regions of isentropic flow • especially in Multi-D • second order Godunov-like methods help overcome deficiencies of classical Lagrange Q • e.g. Christensen, Benson, Caramana et al.

7. Second-order zone-centered Godunov sharp shock resolution < 2 zone thick second order accurate physically intuitive BUT hard to add more physics Elastic-plastic flow, chemical reaction etc. Multi-D is well formulated only in Eulerian coordinates (operator splitting)

8. Zone-Centered Godunov with Lagrange/ALE mesh • Vertex velocities obtained by interpolating either • zone-centered velocities (from momentum equation), or • face velocities (from RP solution) • Result: Weak coupling between node motion and zone deformation •  Spurious grid deformation • Remedy: SMG scheme

9. SMG schemes • SMG/L scheme • RP solved only at zone faces • Zone-centered velocity obtained by interpolation of vertex velocities • Face-centered data extrapolated from the zone-centered data using the (limited) zone–centered gradients • Additional impulse due to RP solution distributed to neighboring vertices • SMG/G scheme • “collision” RP solved at corner-zone faces is used to update the momentum equation • “contact discontinuity” RP solved at zone faces serves to update the (total) energy equation • SMG/Q scheme (closest to CL)

10. Staggered Mesh • Zone-centered but node-centered velocity • velocity jump on corner-zone faces • Jump in only at the zone faces • zones k, edge-neighbor nodes i,i+1 • Corner-zones (i,k). • The faces separate corner-zones (i,k) and (i+1,k)

11. SMG/Q scheme • Data for “collision” RP at corner zone faces: • velocity jump at face : • Continuous in zone k • depends on velocity gradient in the zone and on velocity in neighboring nodes. for first order accuracy (sharp velocity jump). in regions of smooth flow.  solution

12. “collision” Riemann problem • 1D symmetric collision RP (data jump only in velocity). • Solution: either 2 shock waves or 2 centered rarefaction waves

13. SMG/Q scheme • momentum equation: • sum up “corner face forces” on all faces surrounding node i • internal energy equation: • assume all work performed by against the deformation of zone k is deposited in the zone

14. SMG/Q characteristics • As each depends only on the velocity difference the Lagrange scheme is Galilei-Invariant (like CL). • As in “traditional” Godunov schemes no pseudo-viscosity is required. • Needs solution to a simple “collision” RP. • Formulated with internal energy. • see also Christensen, Caramana et al.

15. SMG/Q and pseudo-viscosity • We can regard as a kind of pseudo-viscosity. • as at a shock • except at phase changes or for reacting HE. • define it as a uni-axial tensor acting along the velocity difference (Not along ).

16. Uni-axial tensor Q • Edge-defined Q is better if related lumped forces are aligned with the velocity difference, instead of the edge (Margolin, Caramana et al.). • a moment opposing edge rotation is thereby generated (a grid stabilizing effect). • Summing up “corner face forces” over all faces separating edge-neighbor nodes, produces normal-to-face forces. • Alternately define as uni-axial tensor along the direction of :

17. The Riemann Problem solutions • The RP solution is aimed to provide the necessary dissipation at shocks. • An approximate Riemann solution may be adequate for moderate shocks. • SMG/Q can be used for elastic-plastic and reactive flow even without solving the exact RP (with results at least similar to Q-based CL). • Using the specific RP solution in each case might further improve the results.

18. “Collision” RP solution • The Hugoniot is and the RP can be solved assuming a linear shock-to-particle-velocity relationship • The exact RP can be solved numerically (or, for an ideal gas, analytically) at strong shocks or strong rarefactions. • Note: two waves, each with half the velocity jump !! • Kuropatenko’s Q: single shock with full velocity jump!!

19. Second-Order SMG/Q Velocity jump at faces computed using the zone-centered velocity gradient limited to preserve a smooth velocity profile: with the largest such that the velocities extrapolated to the center of every single-vertex-neighbor zone stay in the min/max range of nodes in surrounding zones. In 1D: this algorithm (QG ) is similar to Christensen’s, which we denote QT

20. The slope limiter • This limiter extends van Leer 1D scheme to a multi-D staggered mesh • Uses average zone-centered velocity gradient • The limiter is defined on either structured or unstructured mesh • In 1D it reduces to , while Christensen’s TVD-like limiter is:

21. Planar Noh Problem • with the data: • and analytic solution: • The 1D calculation had 100 unit length zones • An exact Rieman solver was used in the 1D simulations • Sharp and accurate shock capture and almost identical results for QG and QT

22. Fig.1 Planar Noh problem.

23. Centered Compression Wave (CCW) • The initial (t=0) smooth profile has steepened up to a sharp discontinuity at t=20 • The t=0 profile is the “inverse” of a centered rarefaction with a linear velocity distribution • Results show an almost isentropic compression

24. Fig.2 Centered Compression Wave (CCW).

25. The Sod Problem • The initial (shock tube) data: • It was run with 100 unit-size zones. • The results at t=20 compare favorably to those of the GRP (a genuine second-order Godunov scheme see Ben Artzi & Falcovitz) • SMG with even-time integration gave virtually identical results for this problem

26. Sod Problem. SMG/QG and GRP

27. Sedov-Taylor blast wave problem • The initial data is: cold ideal gas with a “point source” represented by setting in a single zone at the origin (see Caramana et al.) • The physical space is a 253 boxdivided into a uniform 903 mesh • 3D computation conducted in the octant with the symmetry planes • The exact solution has a post-shock density and the front radius is at .

28. Fig.4 Sedov-Taylor blast wave. The 453 mesh spans a 1.253 box. The velocity plot at t=1. (a) xy plane (b) 3D view

29. Fig.5 The Sedov-Taylor blast wave. The xy plane mesh plot at t=1, (left) without NC, (right) with NC.

30. Fig. 6 The Sedov-Taylor blast wave. Isobar plot at t=1.

31. Fig. 7 The Sedov-Taylor blast wave. Isodensity plot at t=1.

32. Fig. 8 The Sedov-Taylor blast wave. Density profiles at t=1.

33. The Saltzman Problem • A piston moves with a velocity of u=10 into an almost cold ideal gas • This 1D problem is solved on a 2D mesh with planar symmetry • The 100x10 mesh shown above tests the propagation of a planar shock into an obliquely-inclined non-uniform mesh

34. Fig.9 The Saltzman Problem. 100x10 Lagrange mesh at t=0.7 (a) Scalar Q (b) SMG with uni-axial tensor q and NC (a) (b)

35. Fig.9 Saltzman Problem. 100x10 Lagrange mesh at t=0.7 Improved SMG/Q (RP solved along velocity vector difference of edge-neighbor pairs). It seems no NC is required!!!

36. Cylindrical Noh Problem • with the data: • and analytic solution:

37. Cylindrical Noh problem 1D Lagrangian mesh in 3D code • 3D code can handle 1D,2D mesh with axial symmetry. • The B.C. and the stability condition must be adjusted • A 100 equal size grid was used. • The over-shoot at the front might be caused by the use of the approximate RP solver in the 3D code. An infinite strength shock in an ideal gas does not obey the linear shock to particle relationship. t=0.6

38. Cylindrical Noh problem 2D ALE mesh in 3D code • We consider a cylinder with a body-fitted grid which has 32x32x1 zones • There are symmetry planes at z=0,z=0.05 • The ALE grid motion keeps the mesh smooth while the outer radius shrinks

39. Cylindrical Noh problem 2D ALE mesh in 3D code at T=0.6 (a) velocity plot (b) isobars

40. Cylindrical Noh problem 2D ALE mesh in 3D code at T=0.6. isodensity plot

41. Cylindrical Noh problem 2D ALE mesh in 3D code at T=0.6 (a) velocity plot (b) isobars

42. Spherical Noh Problem • with the data: • and analytic solution: see Noh, Birman et. al.

43. 3D ALE mesh for Spherical Noh problem • The mesh is Cartesian. The ALE grid is uniformly contracting at a rate given by the grid velocity: • where is the velocity of a grid point at . is the spherical boundary at time t. • The 3D mesh has 503 zones. The material boundary cuts through mesh zones.

44. Spherical Noh problem 3D (Eulerian-like) ALE mesh at T=0

45. Spherical Noh problem 3D (Eulerian-like) ALE mesh at T=0.6. isobar plot

46. Spherical Noh problem 3D (Eulerian-like) ALE mesh at T=0.6. isodensity plot

47. Concluding Remarks • SMG is a ‘dual scheme’. It is both a Godunov method and a classical-Lagrange method. • SMG/Q bridges the Lagrange-to-Godunov “conceptual gap”. • The second-order extension produces sharp shock capturing while minimizing entropy errors in regions of smooth compression. • The method can be applied to structured mesh, as well as to unstructured mesh.

48. Future Work • Improve the velocity limiter for multi-D structured and unstructured mesh • Explore alternate Q formulations: • Zone-centered “directional scalar Q” • Zone-centered “directional uni-axial tensor Q” • Elastic-plastic flow

49. References • Wilkins M. L., "Calculation of Elastic-Plastic Flow”, Meth. Comp. Phys., 3, p211, B. Alder et al. eds., Academic Press (1964). • van Leer B., J. Comp. Phys., 32, p101, (1979). • Ben Artzi M., Falcovitz J., "Generalized Riemann Problems in Computational Fluid Dynamics", Cambridge Univ. Press, London, (2003). • Luttwak, G., "Comparing Lagrangian Godunov and Pseudo-Viscosity Schemes for Multi-Dimensional Impact Simulations”, p255-258, Shock Compression of Condensed Matter-2001, Furnish M. D. et al Eds, (2002), AIP, CP620. • Luttwak, G. ”Staggered Mesh Godunov (SMG) Schemes for Lagrangian Hydrodynamics”, presented at the APS Topical Conference on Shock Compression of Condensed matter,Baltimore,MD,July31-Aug5,2005

50. References • Christensen R.B., "Godunov Methods on a Staggered Mesh. An Improved Artificial Viscosity", L.L.N.L report UCRL-JC-105269, (1990). • Benson D. J., Schoenfeld S., Comp. Mech., 11, p107-121, (1993). • Caramana E.J., Shaskov M.J., Whalen P.P., J. Comp. Phys. 144, p70, (1998). • Margolin L.G., LLNL report UCRL-5382,(1988) • Noh F. W.,J. Comp. Phys., 72,p78,(1987) • Sod G.,A., J. Comp. Phys., 27,p1,(1978) • Birman A., Har’el N. Y. ,Falcovitz J. , Ben-Artzi M., Feldman U., 22nd Int. Symp. on Shock Waves, Imperial College, London, UK, July 18-23, (1999)