Skip to main content

Module spherical

Module spherical 

Source
Expand description

§Spherical coordinates

The Euler equations in spherical symmetry: a one-dimensional radial problem (geometry = "spherical"). The geometric terms are those of the cylindrical case with face areas $4\pi r^2$ and volumes $\tfrac43\pi r^3$ per unit solid angle, and the validation is the same set of questions. Problem files are in problems/spherical/. Test paths are relative to crates/merlin-finite-volume/tests/ unless they start with crates/.

CaseReferenceHeadline result
Sound waveexact linear solution $f(r-ct)/r$order 1.72, 1.67, 1.88 (MUSCL); tends to 1 without reconstruction
Noh’s implosionexactshock at 0.19875 (exact 0.2); plateau 63.12 (exact 64)
Sedov-Taylor blastexact self-similar solution (in Python)shock radius 0.69375 (exact 0.69096); error 2.6e-3 at 800 cells
Blasta fine run in the same geometryorder 1.09 on the finest pair
Conservationinvariantsdrift 1.1e-13 / 2.4e-13

§Outgoing sound wave

A spherical wave of small amplitude has the exact solution $p'=\rho_0 c\,f'(r-ct)/r$ from the potential $\varphi=f(r-ct)/r$ (Landau and Lifshitz, Fluid Mechanics, Section 69; Equations). It solves the nonlinear equations up to terms of second order in the amplitude, which the relative amplitude $10^{-7}$ makes negligible. It is the Rust test of the radial source term: an inconsistent geometric term would show up as an error that does not shrink with the grid. A Gaussian of width 0.06, $t=0.3$, HLLC, SSP-RK3; errors in units of the amplitude:

SchemeCells$L_1$ errorsObserved orders
MUSCL (van Leer)50, 100, 200, 4001.246, 0.378, 0.119, 0.03231.72, 1.67, 1.88
first order200, 400, 800, 16001.422, 0.875, 0.495, 0.2650.70, 0.82, 0.90, tending to 1

The same wave computed through the Python reference (pymerlin convergence spherical/euler/acoustic, errors in units of the pulse) gives 5.24e-3, 1.65e-3, 4.46e-4, 1.09e-4 and orders 1.67, 1.88, 2.03, fitted 1.86, in agreement with the Rust test.

Density, spherical sound wave

Tests that enforce this

Enforced by geometry/acoustic.rs: a_spherical_sound_wave_converges_at_second_order_with_muscl, a_spherical_sound_wave_converges_at_first_order_without_reconstruction; and in continuous integration pymerlin check --self-test acoustic and pymerlin convergence spherical/euler/acoustic --expect-order=2 --tolerance=0.4.

§Noh’s implosion

The same exact problem as in the cylinder, with $\nu=3$. For $\gamma=5/3$ the post-shock density is $4^3=64$ and the shock reaches $r=0.2$ at $t=0.6$ (spherical/euler/noh.toml, 400 cells).

Noh, sphere, t = 0.6

  • Shock at 0.19875 against 0.2.
  • Mean density behind the shock 63.12 against 64.
  • $L_1$ density error over $0.05<r<0.35$ at 100, 200, 400 cells: 0.649, 0.339, 0.178; observed orders 0.94, 0.93, a shock and the wall-heated centre limiting the order to about one.
Tests that enforce this

Enforced by geometry/noh.rs (parametrised over both geometries): the_density_behind_the_shock_is_the_exact_plateau, the_shock_is_at_the_exact_radius, the_error_falls_as_the_grid_is_refined, muscl_is_more_accurate_than_first_order_on_the_implosion; crates/merlin-cli/tests/run/curvilinear_problems.rs: noh_files_show_the_exact_gas_behind_the_shock.

§Sedov-Taylor blast

A hot ball of radius 0.02 ($p=10^4$ against $10^{-3}$, energy $E=0.838$) in a cold gas at rest drives a strong spherical shock. Once the shock is several radii of the ball away, the flow is the self-similar solution of Sedov (1959), $R(t)=(Et^2/\alpha\rho_0)^{1/5}$ with $\alpha=0.85107$ for $\gamma=1.4$, which is Sedov’s published 0.851; it is computed in Python (src/pymerlin/exact_solutions/sedov.py, derived on the Equations page). 400 cells, $t=0.4$ (spherical/euler/sedov.toml):

Sedov blast, sphere, t = 0.4

  • Shock radius 0.69375 against 0.69096: 0.4 %, a little more than a cell.
  • The density peak is 5.00 against the exact 6, from the smearing of the shock.
  • Mass and energy drift 3.3e-13 and 2.3e-12.
  • $L_1$ density error against the exact solution (units of the density range) at 100, 200, 400, 800 cells: 1.34e-2, 8.55e-3, 5.39e-3, 2.59e-3, observed orders 0.65, 0.67, 1.05 (fit 0.78).

As in the cylindrical case the check is a bound, not an order: the simulation starts from a ball of finite size, so the error does not go to zero as the grid is refined. It is checked from Python, in continuous integration: uv run pymerlin check --self-test sedov and uv run pymerlin convergence spherical/euler/sedov --max-error=0.005.

§Blast against a fine run

The hot core $r<0.25$ in a closed tube, $t=0.2$ (spherical/euler/blast.toml); the reference is a 1600-cell run in the same geometry.

Cells50100200400
$L_1$ density error6.9e-26.0e-23.4e-21.6e-2

Order 1.09 on the finest pair.

Spherical blast

Tests that enforce this

Enforced by geometry/blast.rs: a_blast_converges_to_the_fine_run_in_the_same_geometry.

§Orders at a glance

1 2 4 8 1 1e−1 1e−2 slope 1 slope 2 sound wave (MUSCL) sound wave (first order) Noh, away from the centre Sedov, against the exact solution blast, against a fine run refinement factor (cells ÷ coarsest grid) error ÷ coarsest error

Error against refinement for the spherical cases, each series divided by its own coarsest error. The exact profile of Noh’s implosion is on the cylindrical page.

§Conservation and equilibrium

  • A gas at rest has exactly zero right-hand side, as in the cylinder (geometry/operator.rs).
  • A closed ball conserves mass and energy: relative drift 1.1e-13 and 2.4e-13 over $t=0.5$ with 400 cells (5478 steps), with the total mass $4\pi/3$ required by the volumes (geometry/conservation.rs: a_closed_cylinder_or_ball_conserves_mass_and_energy; crates/merlin-cli/tests/run/curvilinear_problems.rs: closed_problems_conserve_the_volume_integrals_of_mass_and_energy).
  • Results do not depend on the thread count (1 to 8) or the block layout.
  • A gas with swirl is refused on a spherical grid, because the swirl needs an axis to rotate about (geometry/swirl.rs: swirl_needs_an_axis_to_rotate_about).
  • A strong outward flow at the centre keeps a positive density, because the time step accounts for the cell next to the centre, whose face area over volume is $3/\Delta r$ (geometry/conservation.rs: the_time_step_accounts_for_the_small_cells_next_to_the_axis, a_strong_outward_flow_at_the_axis_keeps_a_positive_density).

Integrals of the spherical blast

§References

  • Noh, W. F. (1987), Errors for calculations of strong shocks using an artificial viscosity and an artificial heat flux, J. Comput. Phys. 72, 78-120.
  • Landau, L. D. and Lifshitz, E. M., Fluid Mechanics, Sections 69 and 106.
  • Sedov, L. I. (1959), Similarity and Dimensional Methods in Mechanics, Academic Press, Chapter IV.
  • Toro, E. F. (2009), Riemann Solvers and Numerical Methods for Fluid Dynamics, Chapter 1.