1 / 42

IMPLICIT SCHEMES FOR THE SIMULATION OF UNDERWATER IMPLOSION AND EXPLOSION PROBLEMS?

IMPLICIT SCHEMES FOR THE SIMULATION OF UNDERWATER IMPLOSION AND EXPLOSION PROBLEMS? Charbel Farhat, Alex Main and Kevin Wang Department of Aeronautics and Astronautics Department of Mechanical Engineering Institute for Computational and Mathematical Engineering Stanford University

saman
Download Presentation

IMPLICIT SCHEMES FOR THE SIMULATION OF UNDERWATER IMPLOSION AND EXPLOSION PROBLEMS?

An Image/Link below is provided (as is) to download presentation Download Policy: Content on the Website is provided to you AS IS for your information and personal use and may not be sold / licensed / shared on other websites without getting consent from its author. Content is provided to you AS IS for your information and personal use only. Download presentation by click this link. While downloading, if for some reason you are not able to download a presentation, the publisher may have deleted the file from their server. During download, if you can't get a presentation, the file might be deleted by the publisher.

E N D

Presentation Transcript

  1. IMPLICIT SCHEMES FOR THE SIMULATION OF UNDERWATER IMPLOSION AND EXPLOSION PROBLEMS? Charbel Farhat, Alex Main and Kevin Wang Department of Aeronautics and Astronautics Department of Mechanical Engineering Institute for Computational and Mathematical Engineering Stanford University Stanford, CA 94305

  2. MOTIVATION Are explicit schemes really required to perform such simulations?

  3. SHOCK TUBE (GAS-GAS) • MURI Review, April 2007 0.5 1.0 P = 1 r = 1 M = 0 g = 1.6 P = 0.1 r = 0.125 M = 0 g = 1.2

  4. POTENTIAL OF IMPLICIT SCHEMES • MURI Review, April 2007 t1 = 1e-1 s 5.0e-01 Initial CFL = 1 4.0e-01 CFL= 20 3.0e-01 CFL = 40 Density (kg/m3) 2.0e-01 1.0e-01 0.0e+00 5.0e-01 6.0e-01 7.0e-01 X (m)

  5. OBJECTIVES 2010 • Develop fully-implicit schemes for the solution of: - Multi-fluid and multi-phase flow problems - Multi-fluid and multi-phase fluid-structure interaction problems • Implement them in AERO-F and AERO-F/AERO-S • Benchmark them using representative 3D underwater • implosion problems

  6. 1 1 Fj,j+1 = Fj+1/2 (nj,j+1) = (Fj + Fj+1 )- | F’ |j+1/2 (Wj+1 – Wj) = Roe (Wj, Wj+1, gs, ps) (stiffened gas) 2 2 @(rf) (ruf) @ + = 0 @t @x COMPUTATIONAL FRAMEWORK • MUSCL-based solver with Roe Flux j + 1/2 j j + 1 • Interface capturing via the level-set equation (conservation form)

  7. FVM-ERS • FVM with exact local Riemann solver for multi-phase flows Wjn W*n W*n Wj+1n j - 1 j - 1/2 j j + 1/2 j + 1 - Fj,j+1 = Roe (Wjn, W*n, EOSj) Fj+1,j = Roe (Wj+1n, W*n, EOSj+1) - W*n and W*n determined from the exact solution of local two-phase Riemann problems C. Farhat, A. Rallu and S. Shankaran, "A Higher-Order Generalized Ghost Fluid Method for the Poor for the Three-Dimensional Two-Phase Flow Computation of Underwater Implosions", Journal of Computational Physics, Vol. 227, pp. 7674-7700 (2008)

  8. Wnpj+1 Wnpj • Exact solution of the analytical problem 1 1 (RR(pI; pR,rR) - RL(pI; pL,rL)) uI = (uL + uR) + 2 2 RL(pI; pL,rL) + RR(pI; pR,rR) + uR – uL = 0 pI, rIL, rIR, uI - Newton’s method LOCAL RIEMANN SOLVER • Wave structure and Riemann problem rIL,pI,uI ,rIR contact discontinuity rarefaction shock t gas water x j j + 1/2 j + 1 rLuL pL rRuR pR

  9. FVM-ERS (EXPLICIT) • GFMP with exact local Riemann solver j - 1/2 j + 1/2 j - 1 j j + 1 • If fjnfj+1n> 0 then Fj,j+1 = Fj+1,j = Roe (Wjn, Wj+1n, EOSj= EOSj+1) If fjnfj+1n< 0 then Fj,j+1 = Roe (Wjn, WjRn(rIL, pI, uI), EOSj) Fj+1,j = Roe (Wj+1n, W(j+1)Rn(rIR, pI, uI), EOSj+1) Dt ~ - Wjn+1 = Wjn - (Fj,j+1 - Fj,j-1) (forward Euler) Dx ~ - Unpack Wn+1 using fn and solve the level-set equation to get fn+1 - Pack Wpn+1 using fn+1 to get the updated solution Wn+1

  10. FVM-ERS (IMPLICIT) • Implicit Extension of FVM-ERS method j - 1/2 j + 1/2 j - 1 j j + 1 • If fjnfj+1n> 0 then Fj,j+1 = Fj+1,j = Roe (Wjn+1, Wj+1n+1,EOSj= EOSj+1) If fjnfj+1n< 0 then Fj,j+1 = Roe (Wjn+1, WjRn+1,EOSj) Fj+1,j = Roe (Wj+1n+1, W(j+1)Rn+1, EOSj+1) Dt ~ - Wjn+1 = Wjn - (Fj,j+1 - Fj,j-1) (backward Euler) Dx ~ - Unpack Wn+1 using fn and solve the level-set equation to get fn+1 - Pack Wpn+1 using fn+1 to get the updated solution Wn+1

  11. IMPLICIT FLUID-FLUID • Implicit scheme  nonlinear equation •  Jacobians of the two-phase • Riemann problem dFj,j+1 @Fj,j+1 @Fj,j+1@WjRn+1 + = dWjn+1@Wjn+1@WjRn+1 @Wjn+1 dFj,j+1 @Fj,j+1@WjRn+1 = dWj+1n+1@WjRn+1 @Wj+1n+1 dFj+1,j @Fj+1,j @Fj+1,j@W(j+1)Rn+1 + = dWj+1n+1@Wj+1n+1 @W(j+1)Rn+1 @Wj+1n+1 dFj+1,j @Fj+1,j@W(j+1)Rn+1 = dWjn+1@W(j+1)Rn+1 @Wjn+1

  12. STIFFENED GAS (SG) • Local two-phase Riemann solver for SG-SG requires the • solution of the equation uL + FL(rL, pL;pI) = uIL = uIR = uR + FR(rR, pR; pI) • Differentiating the above contact equation gives @FL @FL @FL duL + drL dpL dpI + + @rL @pL@pI @FR @FR @FR = duR + drR dpR dpI + + @rR @pR@pI dpI dpI dpI dpI dpI dpI , , , , , dpL drL dpR drR duL duR

  13. STIFFENED GAS • Since , the derivatives • can be obtained by straightforward differentiation uI = uL + FL(rL, pL;pI) = uR + FR(rR, pR;pI) duI duI duI duI duI duI , , , , , dpL drL dpR drR duL duR rIR, rIL • Finally, the derivatives of the interfacial densities • are obtained from the straightforward differentiation of rIR = RR(rR,uI,pI) rIL = RL(rL,uI,pI)

  14. OTHER EOSs • Implicitization was also performed for - Tait-Tait - Tait-SG (which also implies Tait-PG)

  15. p = A(1 - )e-R1+ B(1 - )e-R2 + wre wr wr R1r0 R2r0 r0 r0 r r JWL EOS • Jones-Wilkins-Lee (JWL) equation of state for modeling explosive products of combustion (and in particular Trinitrotoluene — a.k.a. TNT) where A, B, R1, R2, w and r0 are material constants - Highly nonlinear function p(r,e) - Presence of exponentials

  16. uL + FL(rL, pL;rIL) = uIL = uIR = uR + FR(rR, pR; rIR) GL(rL, pL; rIL) = pIL pIR= GR(rR, pR; rIR) = JWL EOS • Solution of exact Riemann problem involves a • system of two nonlinear equations (1) (2) • FL and GL depend on the nature of the interaction in the • phase modeled by the JWL EOS • shock algebraic equation • rarefaction differential equation

  17. rIR,uIR ,pIR t rarefaction c(r,p) rR,uR ,pR = r x du + _ dr = s rw+1 p - Ae-R1+ Be-R2 r0 r0 r r SG-JWL RIEMANN SOLVER • Rarefaction wave in a JWL medium (k) • The isentropic evolution in the • rarefaction fan between two • constant states is given by (1) (2) complex Riemann problem • Algebraic entropy (s) formula for the JWL EOS • No obvious algebraic Riemann invariants for the JWL EOS • No analytical Jacobians of the invariants either

  18. JWL EOS • Riemann invariants are tabulated for the explicit • time-stepping scheme • For implicit time-stepping where Jacobians are • required, they are currently computed on-line by • solving an ODE • Implicitization was performed for - JWL-JWL - JWL-SG (which also implies JWL-PG)

  19. TIME-INTEGRATORS • AERO-F is equipped with two implicit time-integrators - Backward-Euler - Three-point Backward Differencing Formula (3PBDF) • Backward Euler estimates the time-derivative as follows ~ dWi Win+1 - Win = dt Dt • Node i has the same fluid ID at both tn and tn+1 • and therefore the above estimation is meaningful

  20. 3PBDF • 3PBDF approximates the time-derivative as follows ~ dWi a0Win+1 -a1Win + a2Win-1 = dt Dt where are constants a0, a1, a2 • Problem: at tn-1, node i may have a different fluid ID • than at tn and tn+1: since the density can be • discontinuous across a fluid-fluid interface (contact • surface), Win-1 and Win are not necessarily related in • this case and therefore the above approximation is • invalid

  21. 3PBDF • Solution - When node i has changed phase between tn-1 and tn, replace Win-1 by W(i-1)Rn-1, the exact solution of the two-phase Riemann problem on the upstream side of the interface at node i-1 and time tn-1 n-1 n W(i-1)Rn-1 i+1 i-2 i i-1

  22. 2fin+1 - 2fin = Dt d fi dt 3PBDF FOR LEVEL SET • Similar issue arises when 3PBDF is applied to the level • set equation: after re-initialization, fn-1no longer exists! • Solution: 1 dfin - 2 dt where the last term can be estimated from the fluxes at tn

  23. LIMITATION • The fluid-fluid interface may cross no more than one cell • per time-step - Required to address phase change • AERO-F enforces this condition by adjusting the time-step as necessary

  24. r = 50 (kg/m3) r = 1000.0 (kg/m3) u= 0.0 (m/s) u= 0.0 (m/s) p= 105 (Pa) p= 109 (Pa) SHOCK TUBE PROBLEM • Shock tube: air to the left, water to the right • Air modeled as a perfect gas (g = 1.4); water modeled • as a stiffened gas (g = 4.4, p = 6.0x 108) Air Water • Simulation up to t = 1 x 10-5 s

  25. SHOCK TUBE PROBLEM

  26. SHOCK TUBE PROBLEM

  27. TURNER IMPLOSION PROBLEM • Implosion of a glass sphere (D = 0.0762 m) (0.5m, 0.5m) Air (P = 105 Pa) (0, 0) Sensor Water (P = 6.996 MPa) z (0.5m, -0.5m) x

  28. TURNER IMPLOSION PROBLEM • Air modeled as a perfect gas (g = 1.4); water modeled • as a stiffened gas (g =7.15, p = 2.89 x 108 Pa) • 780,000 grid points • Simulation up to t = 0.5 ms

  29. TURNER IMPLOSION PROBLEM Explicit (RK2), CFL = 0.5 Implicit (3BDF), CFL = 100

  30. TURNER IMPLOSION PROBLEM

  31. TURNER IMPLOSION PROBLEM

  32. TURNER IMPLOSION PROBLEM • CPU performance on 168 cores of a Linux cluster - All runs with I/O on were performed with equal amount of I/O speedup factor ~ 8.74 (FE/BE) ~ 25 (3PBDF/RK2)

  33. pI, rIRus x = x(t) contact discontinuity not involved rarefaction* t structure fluid x i j Mij Wnj rRuR pR EMBEDDED FSI FRAMEWORK • Fluid-structure Riemann problem (1/2 problem) n w(x,0) =W , if x ≥ 0 j w F (w) = 0 + x t * could also be a shock u(x(t), t) = u (Mij)∙ nG(Mij) s

  34. pI, rIRus x = x(t) contact discontinuity not involved rarefaction* t structure fluid x i j Mij Wnj rRuR pR ONE-SIDED RIEMANN PROBLEM • Fluid-shell problem • Closed form algebraic solution of the problem exists (SG, Tait) us = uR + R2(pI(2); pR , rR) - Closed form Jacobians exist as well (SG, Tait)

  35. FLUX COMPUTATION • The flux across the face at Mij is given by G Mij i j fluid 1 fluid 2 Fij= Roe (us, pI(1), Wni , EOS(1), uij) Fji= Roe (us, pI(2), Wnj , EOS(2), uji)

  36. EMBEDDED FSI (IMPLICIT) • Implicitization was performed for - SG-structure (which also implies PG-structure) - Tait-structure remains to be done • 3PBDF n-1 n W(i-1)Rn-1 structure i+1 i-2 i i-1

  37. 2D IMP45 • Simplified IMP45 using a thin slice of the aluminum tube air ( p = 14.5 psi ) water ( p = 1500 psi) X 400 AERO-F/DYNA3D • Explicit = Explicit/Explicit simulation uses dt = 0.75 x 10-8 • Implicit = Implicit/Explicit simulation uses dt = 3.0 x 10-6

  38. 2D IMP45 • Pressure at t=0.4 ms • Solution - Explicit (RK2) - Implicit (3PBDF) - Implicit (BE)

  39. 2D IMP45 • Pressure at a sensing node

  40. 2D IMP45 Explicit (RK2), CFL=0.5 Implicit (3BDF), CFL=100

  41. 2D IMP45 • CPU performance on 64 cores of a Linux cluster - AERO-S (nonlinear): implicit midpoint rule time-integrator - Same amount of I/O performed in all runs with I/O on speedup factor = 35.5 (FE/BE) = 43.0 (RK2/3PBDF)

  42. SUMMARY • Implicitization of fluid-fluid schemes in AERO-F • Development of a new multi-phase flow solver based on the • three point backward difference formula - Verification and validation using shock tube and Turner’s implosion problems - Achievement of speedup factors of 4-5 • Implicitization of fluid-structure schemes using AERO-F - Verification on a 2D implosion problem - Achievement of speedup factor of 43 using the midpoint rule in AERO-S

More Related