Large Scale Simulation of Cloud Cavitation Collapse Ursula - - PowerPoint PPT Presentation

CSElab

Computational Science & Engineering Laboratory http://www.cse-lab.ethz.ch

Large Scale Simulation of Cloud Cavitation Collapse

Ursula Rasthofer

4th Cavitation Workshop Crete, 30 May - 1 June 2016

with: Fabian Wermelinger, Jonas Sukys, Diego Rossinelli, Panagiotis Hadjidoukas, Petros Koumoutsakos

Outline

• 1. Motivation
• 2. Cubism-MPCF
• 3. Modeling Approach and Governing Equations
• 4. Simulations
• 5. Conclusions

2

Outline

• 1. Motivation
• 2. Cubism-MPCF
• 3. Modeling Approach and Governing Equations
• 4. Simulations
• 5. Conclusions

3

Power of Cavitation

4

Engineering

https://de.wikipedia.org/wiki/Kavitation extracted from: Brennen, Hydrodynamics of pumps, Oxford University Press, 1994 http://www.forddoctorsdts.com http://www.forddoctorsdts.com extracted from: Bazan-Peregrino et al., Cavitation-enhanced delivery of a replicating oncolytic adenovirus to tumors using focused ultrasound , Journal of Controlled Release Volume 169, Issues 1–2, 2013, 40-47

Biomedical Applications

https://en.wikipedia.org/wiki/ Extracorporeal_shock_wave_lithotripsy

Nature

https://de.wikipedia.org/wiki/ Knallkrebse herve.cochard.free.fr

• Thousands of bubbles
• Cloud dynamics dominated by bubble

interactions

• Experiments
• Theory
• based on Rayleigh-Plesset equation or similar

relations

• usually assume spherical bubbles
• Simulations
• Lagrangian approaches for bubbles
• fully resolved bubbles usually restricted to small

clouds (~100 bubbles)

Cloud Cavitation Collapse

5

esource.alstrbology.com www.jotformeu.com www.amc.edu.au

120 Collapsing Bubbles (TUM)

6

Bubbles = 120 Np = 1.2E+08 Tw = 29.9 CSElab: 15E+03 Bubbles, Np = 1.3E+13, Tw = 1.8

Outline

• 1. Motivation
• 2. Cubism-MPCF
• 3. Modeling Approach and Governing Equations
• 4. Simulations
• 5. Conclusions

7

Cubism-MPCF

• compressible multicomponent flow solver

[Rossinelli et al. 2013]

• 3D structured-grid finite volume solver
• wavelet-based compression of simulation data

8

C u b i s m

• MPCF

C S E l a b

• conservation-law form
• cell-average values
• semi-discrete form

Finite Volume Method

9

Uijk(t) = 1 h3 Z

Ωc

U(x, t)dΩc

dUijk(t) dt + 1 h ⇣ Fx

i+ 1

2 − Fx

i− 1

2 + Fy

j+ 1

2 − Fy

j− 1

2 + Fz

k+ 1

2 − Fz

k− 1

2

⌘ = 0

ci ci-1 ci+1 xi-1/2 xi+1/2 Ui-1 Ui Ui+1 Fi+1/2 Fi-1/2 h

∂U(x, t) ∂t + r · F(U) = 0 U = (α1ρ1, α2ρ2, ρu, ρv, ρw, E)T

Godunov-Type Finite Volume Method

• reformulation of gas-volume-

fraction equation

[Johnsen & Colonius 2006]

• high-order reconstruction of

face values of primitive variables using WENO3/5

[Liu et al. 1994, Jiang & Shu 1998]

• approximate HLLC Riemann

solver for flux reconstruction

[Toro et al. 1994]

• low-storage 3rd-order Runge-

Kutta scheme for time discretization

[Gottlieb et al. 2001]

10

∂α2 ∂t + r · (α2u) = α2ρ1c2

1

α1ρ2c2

2 + α2ρ1c2 1

r · u

ci ci-1 ci+1 xi-1/2 xi+1/2 Ui-1 Ui Ui+1 Fi+1/2 Fi-1/2 h WL WR

s r c x y U- U1 U2 U+ UL* UR*

Software Layout

• cluster (MPI)
• node (OpenMP)
• core (SIMD): WENO, HLLC

11

Core Node Cluster

grid cell computational domain subdomain block

www.cscs.ch www.alcf.anl.gov

Performance and Comparison

12

PFLOPS (% Peak) 14.4 PFLOPS (72%) 0.1 - 3% (TUM) 1.3 - 6.4% (Stanford) SIZE 1.3 E13 - 15K bubbles 1.2 E08 - 0.15K bubbles (TUM) 0.4 E13 - turbulence (Stanford) I/O COMPRESSION 10-100x

• (TUM)
• (Stanford)

TIME TO SOLUTION (no I/O) Tw = 1.8 Tw= 29.7 (TUM) Tw= 16.3 - 39.0 (Stanford)

Outline

• 1. Motivation
• 2. Cubism-MPCF
• 3. Modeling Approach and Governing Equations
• 4. Simulations
• 5. Conclusions

13

Numerical Challenges

• O(100) to O(104) bubbles
• broad range of length and time scales
• complex interface evolution due to

bubble-bubble interactions

• high density ratios O(103)
• surface-tension effects
• coupling of gas and nearly-

incompressible liquid

• secondary cavitation

14

Modeling Strategies

15

complexity/fidelity computational efficiency

two-fluid models

Lagrangian/Eulerian methods, e.g., based

• n level-set description

single-fluid models

e.g., diffuse interface methods

mixture approches

e.g., thermodynamic equilibrium cavitation model

Diffuse Interface Method

16

• based on two-fluid model with 7 equations

[Baer & Nunziato 1986, Saurel & Abgrall 1999]

• mass, momentum, energy conservation for each phase
• additional topological equation for interface
• reduced to 5-equation model by asymptotic analysis in limit
• f zero relaxation time (velocity and pressure equilibrium at

interfaces)

[Kapila et al. 2001, Allaire et al. 2002, Murrone & Guillard 2005, ….]

Governing Equations

where and

17

∂α1ρ1 ∂t + r · (α1ρ1u) = 0 ∂α2ρ2 ∂t + r · (α2ρ2u) = 0 ∂ (ρu) ∂t + r · (ρu ⌦ u + pI) = 0 ∂E ∂t + r · ((E + p) u) = 0 ∂α2 ∂t + u · rα2 = K r · u

K = α1α2(ρ1c2

1 − ρ2c2 2)

α1ρ2c2

2 + α2ρ1c2 1

[Kapila et al. 2001, Murrone & Guillard 2005, Saurel et al. 2009, …]

O(h)

α1 + α2 = 1

• stiffened equation-of-

state (EoS) for pure fluids

Thermodynamic Closure

18

1 ρc2 = α1 ρ1c2

1

+ α2 ρ2c2

2

[Wood 1930] [Menikoff & Plohr 1989, …]

pk = (γk − 1) ρkek − γkpc,k

ρe = α1ρ1e1 + α2ρ2e2

ρ = α1ρ1 + α2ρ2

• pressure equilibrium p=pk
• mixture speed of sound
• mixture relations

• source term of gas-volume-fraction equation neglected
• mixture mass conservation
• equation-of-state must be valid in all flow regions (pure fluids

and artificial mixing zone) and

Further Variants

19

Γ = 1 γ − 1 and Π = γpc γ − 1.

[Shyue 1998, Johnsen & Colonius 2006,…]

∂ρ ∂t + r · (ρu) = 0 ∂ (ρu) ∂t + r · (ρu ⌦ u + pI) = 0 ∂E ∂t + r · ((E + p) u) = 0 ∂φ ∂t + u · rφ = 0 where φ 2 {Γ, Π}

ρc2 =

α1ρ1c2

1

γ1−1 + α2ρ2c2

2

γ2−1 α1 γ1−1 + α2 γ2−1

Model Properties

• non-monotonic speed of sound with

respect to gas volume fraction

• right-hand-side term of gas volume

fraction equation

• increase of gas volume fraction across

rarefaction waves

• decrease of gas volume fraction across

compression waves

20

• α1α2(ρ1c2

1 ρ2c2 2)

α1ρ2c2

2 + α2ρ1c2 1

r · u

• Key Issue I: Negative Pressure
• High Pressure Ratios
• Large Clouds
• amplification of waves due to

many bubble layers

• high peak pressures in cloud

resulting in strong rarefaction in cloud interior

21

Key Issue I: Negative Pressure

• effects of choice for EoS for pure fluids
• negative pressure in liquid due to stiffened EoS
• tensile stresses appearing as negative pressure regions

due to strong rarefaction waves

• indications of secondary cavitation?
• ther?

22

pk = (γk − 1) ρkek − γkpc,k

Note: Negative pressures also observed by other groups (private communication)

Negative Pressure - Remedies?

23

• right-hand-side term of gas volume fraction equation
• dynamic creation of interfaces by mechanical relaxation if small

portion of gas available [Saurel et al. 2008]

• integration of mass transfer directly enabled by separate

continuity equation for each fluid [Saurel et al. 2008]

∂α2 ∂t + u · rα2 = α1α2(ρ1c2

1 ρ2c2 2)

α1ρ2c2

2 + α2ρ1c2 1

r · u

• two diverging rarefaction waves in water with small amount of

air volume fraction content

• initial conditions

Cavitation Tube - Setup

24

rarefaction wave rarefaction wave stationary contact wave

x

t

WL

WR

(ρ, u, p, α2)0 = ( (1150.0 kg/m3, −100.0 m/s, 1.0 · 105 Pa, 1.0 · 10−2) for 0 m ≤ x ≤ 0.5 m (1150.0 kg/m3, 100.0 m/s, 1.0 · 105 Pa, 1.0 · 10−2) for 0.5 m < x ≤ 1.0 m

Cavitation Tube - Results

25

Key Issue II: Numerically Diffusing Interface

• numerical thickening of interface by

inherent diffusion of shock capturing schemes

• maintaining approximately constant

interface thickness

• compatibility with thermodynamic

mixture model in interface zone

26

Key Issue II: Numerically Diffusing Interface

where

27

∂α1ρ1 ∂t + r · (α1ρ1u) = R1 ∂α2ρ2 ∂t + r · (α2ρ2u) = R2 ∂ (ρu) ∂t + r · (ρu ⌦ u + pI) = u (R1 + R2) ∂E ∂t + r · ((E + p) u) = 1 2kuk2 (R1 + R2) + (ρ2e2 ρ1e1) R ∂α2 ∂t + u · rα2 = K r · u + R

R(↵2) = L(↵2)U0n · r (✏kr↵2k ↵2 (1 ↵2)) R1(↵1⇢1, ↵2) = L(↵2)U0n · (r (✏n · r(↵1⇢1)) (1 2↵2)r(↵1⇢1)) R2(↵2⇢2, ↵2) = L(↵2)U0n · (r (✏n · r(↵2⇢2)) (1 2↵2)r(↵2⇢2))

[Tiwari et al. 2013]

Interface Sharpening

Outline

• 1. Motivation
• 2. Cubism-MPCF
• 3. Modeling Approach and Governing Equations
• 4. Simulations
• 5. Conclusions

28

Cloud Cavitation Collapse

29

Cloud Cavitation Collapse

30

Water-Air Shock Tube

31

0.1 1 10 100 100 1000 10000 error ncell L2 wo Kdiv L1 wo Kdiv L2 w Kdiv L1 w Kdiv L2 =0.75 wo Kdiv L1 =0.75 wo Kdiv O(h)

1000 2000 3000 4000 5000 6000 7000 8000 9000 10000 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 p [105 Pa] x [m] exact Cubism-MPCF 100 200 300 400 500 600 700 800 900 1000 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 [kg/m3] x [m] exact Cubism-MPCF

50 100 150 200 250 300 350 400 450 500 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 u [m/s] x [m] exact Cubism-MPCF

density

(ρ, u, p, α2)0 = ( (1000.0 kg/m3, 0 m/s, 1.0 · 109 Pa, 0) for 0 m ≤ x ≤ 0.75 m (10.0 kg/m3, 0 m/s, 1.0 · 105 Pa, 1) for 0.75 m < x ≤ 1.0 m

Spherical Bubble Collapse - Setup

• pressure ratios 10 and 50
• bubble radius and pressure
• Keller-Miksis equation

32

¯ p = R

Ω

α2pdΩ R

Ω

α2dΩ Req = @ 3 4π Z

Ω

α2dΩ 1 A

1 3

[Keller & Miksis 1980]

p0 = ( pB,0 if 0 ≤ r ≤ RB,0 p∞ + RB,0

r

(pB,0 − p∞)

• therwise

1 − ˙ RB c1 ! RB ¨ RB + 3 2 − ˙ RB 2c1 ! ˙ R2

B = 1

ρ1 1 + ˙ RB c1 ! (pB − p∞) + ˙ RB ρ1c1 d dt (pB − p∞) pB = pG,0 ✓RB,0 RB ◆3γ

credit: DynaFlow Inc.

Spherical Bubble Collapse - Results

33

pressure ratio 10 pressure ratio 50

2 4 6 8 10 12 14 16 5 10 15 20 pB/p t [s] Keller-Miksis L=15R R0/h=12.8 L=15R R0/h=12.8 wo Kdiv L=15R R0/h=25.6 L=15R R0/h=25.6 wo Kdiv 0.2 0.4 0.6 0.8 1 5 10 15 20 RB [mm] t [s] Keller-Miksis L=15R R0/h=12.8 L=15R R0/h=12.8 wo Kdiv L=15R R0/h=25.6 L=15R R0/h=25.6 wo Kdiv

reduced equation system without “K div u”-term reduced equation system without “K div u”-term equation system with “K div u”-term equation system with “K div u”-term

Spherical Bubble Collapse - Interface

34

• 0.5

0.5 1 1.5 2 5 10 15 20 (dI-dI,0)/dI,0 t [s] L=15R R0/h=12.8 L=15R R0/h=12.8 wo Kdiv L=15R R0/h=25.6 L=15R R0/h=25.6 wo Kdiv

diso

I

= R0.1 − R0.9

without sharpening with sharpening

• 0.5

0.5 1 1.5 2 5 10 15 20 (dI-dI,0)/dI,0 t [s] L=15R R0/h=12.8 L=15R R0/h=12.8 w R L=15R R0/h=25.6 L=15R R0/h=25.6 w R

• 4
• 3
• 2
• 1

1 2 3 4 5 10 15 20 K / div u / K div u t [s] L=15R R0/h=12.8 K L=15R R0/h=12.8 div u L=15R R0/h=12.8 K div u

decrease of gas fraction in interface during compression increase of gas fraction in interface during expansion reduced equation system without “K div u”-term equation system with “K div u”-term without sharpening with sharpening

• average radius
• equivalent radius
• gas volume fraction
• cloud interaction parameter
• cloud generation
• locations: uniform distribution
• radii: log-Gaussian

Collapsing Bubble Clouds

35

water bubbles cloud radius bubble radius

Req =

3

v u u t

N

X

i=1

R3

i

Ravg = 1 N

N

X

i=1

Ri α = 1 Vc

N

X

i=1

4 3πR3

i

β = α ✓ Req Ravg ◆2

50K Bubbles, 64 Billion Cells

36

25K time steps, 72h x 32K cores

• 2500 Bubbles - Collapse Process

37

• bubble radius between 0.1 and 0.15 mm
• cloud interaction parameter of 15
• 100 bar ambient pressure
• 45 bar bubble pressure
• bubble radius resolved by 10 to 15 cells

• 2500 Bubbles - Pressure Evolution

38

2500 Bubbles - Averaged Fields

• averaging of pressure and

velocity over shells of equal thickness

• thickness approximately

average bubble diameter

39

spherical cloud r

Evolution of Bubble Shapes

40

S - surface area of the cavity V - volume of the cavity Vh - volume of the convex hull

Sphericity Ψ = π

1 3 (6V ) 2 3

S

Porosity

φ = V Vh

Non-Sphericity and Non-Porosity

41

Shock-Induced Collapse

42

Outline

• 1. Motivation
• 2. Cubism-MPCF
• 3. Modeling Approach and Governing Equations
• 4. Simulations
• 5. Conclusions

43

Summary

• Cubism-MPCF : finite volume solver capable of ~13T cells for

compressible, multicomponent flow simulations

• key issues on models for cloud cavitation collapse

✓ important to account for gas expansion and compression in interface ✓ presently limited to weak-to-moderate collapse processes

• HPC for clouds of collapsing bubbles: 2K to 50K

✓ pressure and velocity fields as well as bubble shape evolution

44

Outlook

• data analysis of large scale clouds
• coarse grained models of cavitation
• development of cavitation models
• turbulent flow with cavitation

45

46

Thank you for your attention!