top of page

Test 23

The Wind-Bubble Shock

(R. Cabézon, L. Schmidt)

Introduction

The wind-cloud (or “blob”) test, popularized by Agertz et al. (2007), is a benchmark problem that evaluates a code’s ability to capture strong shocks and the ensuing fluid mixing driven by Kelvin-Helmholtz (KH) instabilities across a large density contrast. The standard setup places a cold, dense, spherical cloud in pressure equilibrium with a hotter, more diffuse ambient medium. A supersonic wind in the background frame impacts the initially stationary cloud. Traditional density-entropy SPH formulations have historically struggled with this problem: spurious surface-tension-like forces originating from zeroth-order (E0) errors at contact discontinuities suppress KH growth and artificially delay the cloud’s erosion and mixing. This issue, as well as various strategies to overcome it, such as improved gradient estimators, pairing-resistant kernels, better volume elements, artificial conductivity/heat flux, and modern viscosity treatments, have been extensively documented in the SPH literature.

We quantify cloud disruption and mixing via the surviving dense-gas mass fraction.

Unless otherwise specified, all quantities are given in normalized code units.



Flow phenomena

  • three-dimensional

  • mixing

  • instabilities

  • supersonic

  • compressible flow



Geometry

The global domain is a 3D box with extents [Lx, Ly, Lz] = [1, 0.25, 0.25] filled with a uniform wind. A denser spherical cloud of radius 𝑟_𝑐 = 0.025 is centered at [0.125, 0.125, 0.125]



Boundary conditions

Periodic boundary conditions in all three spatial directions.



Initial conditions

The density contrast is χ = ϱ_𝑐/ϱ_𝑤 = 10 between the wind and the cloud. For simplicity, ϱ_𝑤 = 1 and ϱ_𝑐 = 10 are used. To enforce an initial pressure equilibrium with an ideal gas (γ = 5/3), the internal energies are 𝑢_𝑤 = 3/2 and 𝑢_𝑐 = 3/20, so that 𝑃 = ϱ_𝑢 (γ−1) = 1 in both media.

The cloud is initially at rest. The wind moves along +X with 𝑣_𝑥/ = 𝑀 𝑐_𝑠, with a Mach number 𝑀 = 2.7.

To reduce discretization noise, the smoothing length is tapered smoothly around the cloud:


ℎ = ℎ_𝑐 + 12 (ℎ_𝑤 − ℎ_𝑐) [1 + 𝑡𝑎𝑛ℎ (𝑘 (𝑟 − 𝑟_𝑐 − ℎ_𝑤))],


applied to particles with 𝑟 <= 𝑟_𝑐 + 2 ℎ_𝑤, where ℎ_𝑐 and ℎ_𝑤 are the target smoothing lengths inside the cloud and in the wind, respectively; r is the distance from the cloud center, 𝑟_𝑐 is the radius of the cloud, and 𝑘 = 𝑁_𝑛𝑔𝑏,𝑚𝑎𝑥/0.125 with 𝑁_𝑛𝑔𝑏,𝑚𝑎𝑥 = 150 being the maximum number of neighbors. Particle masses are equal everywhere, and the targeted number of neighbors is 100 for all particles.

Simulations were performed with SPH-EXA and can be reproduced with:


./sphexa-cuda --init wind-shock --glass 50c.h5 -n 150 -s 0.5 -w <snapshot-frequency> --avclean


where <snapshot-frequency> is an integer (iterations between snapshots) or a float (physical time between snapshots). The -n parameter controls the resolution (default 150 corresponds to the particle counts specified below). The -s parameter controls the simulation duration (integer = iterations, float = physical time), analogous to -w. For CPU runs, use ./sphexa instead of ./sphexa-cuda. The --glass argument supplies a pre-relaxed glass-like cube template to construct the initial conditions (see Discretisation for more details). --avclean activates the slope linear reconstruction method to have a better control over numerical dissipation.



Discretization

The initial geometry is assembled from a pre-relaxed glass-like cubic template (stacked four times along X for the wind box). A spherical cavity is cut at the cloud location and filled with a sphere carved from the same template and stretched to achieve χ = 10 with equal particle masses. The glass template can be found on the SPH-EXA Github repository https://github.com/sphexa-org/sphexa.


Results specification

Quantitative comparison focuses on cloud fragmentation and mixing via the surviving mass fraction of the cloud, measured on a characteristic instability timescale. The characteristic time scale in our setup of the emerging Kelvin-Helmholtz instabilities for χ = 10 is 𝑡_𝑘ℎ = 0.0726 (Agertz et al. 2007).

A particle is considered to belong to the cloud if it has a density ϱ >= 0.64 ϱ_𝑐 and an internal energy 𝑢 <= 0.9 𝑢_𝑤, as defined by García-Senz et al. (2022).



Results format

The results are provided in an ASCII file. First column: time, normalized to the Kelvin-Helmholtz characteristic time (𝑡_𝑘ℎ); second to sixth columns: surviving mass fraction of the cloud at different resolutions (in millions of particles): 0.5, 4, 13.6, 32.3, and 258.4. These are obtained with the following values for the argument -n: 50, 100, 150, 200, 400.



Benchmark results

In agreement with grid-based codes and results of Agertz et al. (2007), the cloud is effectively destroyed by 𝑡 ∼ 2.5 𝑡_𝑘ℎ (i.e. the surviving mass fraction tends to zero by that time). Data is provided in an ASCII file.


ree

Figure 1. The surviving mass fraction of the cloud as a function of the time normalized to 𝑡_𝑘ℎ for all resolutions, serving also as a convergence test.



ree

Figure 2. Converged results are compared to those obtained with the mesh code ENZO, extracted from Agertz et al. (2007), and those obtained with the SPH code CRKSPH, extracted from Frontiere et al. (2017).



ree

Figure 3. A thin slice of the box around the plane Z = 0.125 (i.e. |𝑧_𝑎 − 0.125|/ℎ_𝑎 <2) at different times. Color shows density.


ree

Figure 4. A thin slice of the box around the plane Z = 0.125 at t = 0.



Download

You can download the full test case below:



References

bottom of page