# Two- and three-dimensional explosions

We solve the Euler equations for a compressible gas. We also need to compute volume fractions for the initial condition.

``````#include "compressible.h"
#include "fractions.h"

#if dimension == 2
# define LEVEL 7
#else // 3D
# define LEVEL 6
#endif

int main() {``````

We make boundary conditions free outflow.

``````  foreach_dimension() {
w.n[right] = neumann(0);
w.n[left]  = neumann(0);
}``````

The domain spans $\left[-1:1\right]×\left[-1:1\right]×\left[-1:1\right]$.

``````  origin (-1, -1, -1);
size (2.);
init_grid (1 << LEVEL);
run();
}``````

Initial conditions come from Toro’s book (Riemann Solvers and Numerical Methods for Fluid Dynamics, 3rd Edition, Springer Ed.) Chapter 17 section 17.1.1 are given in terms of density ($\rho$), pression ($p$), velocity ($u$) both at the left and right side of the discontinuity placed at $R=0.4$.

``````event init (t = 0)
{
double R = 0.4 ;
double rhoL = 1., rhoR = 0.125 ;
double pL = 1.0,  pR = 0.1 ;``````

To define an accurate (i.e. symmetrical) initial sphere of rayon $R$, we compute the volume fraction corresponding to the spherical interface.

``````  scalar f[];
fraction (f, sq(x) + sq(y) + sq(z) - sq(R));``````

Left and right initial states for $\rho$, $\mathbf{\text{w}}$ and energy $E=\rho {\mathbf{\text{u}}}^{2}/2+p/\left(\gamma -1\right)$.

``````  foreach() {
ρ[] = rhoR*f[] + rhoL*(1. - f[]);
foreach_dimension()
w.x[] = 0.;
E[] = (pR*f[] + pL*(1. - f[]))/(gammao - 1.);
}

θ = 1.3; // tune limiting from the default minmod
}

event print (t = 0.25)
{``````

At $t=0.25$ we output the values of $\rho$ and the normal velocity ${\mathbf{\text{u}}}_{n}$ as functions of the radial coordinate.

``````  foreach() {
double r = sqrt(sq(x) + sq(y) + sq(z));
double wn = (w.x[]*x + w.y[]*y + w.z[]*z)/r;
printf ("%g %g %g\n", r, ρ[], wn/ρ[]);
}``````

For reference we also output a cross-section at $y=0$.

``````  for (double x = 0; x <= 1; x += 1e-2)
fprintf (stderr, "%g %.4f %.4f\n", x,
interpolate (ρ, x, 0.),
interpolate (w.x, x, 0.));
}``````

On trees, we adapt the mesh by controlling the error on the density field.

``````#if TREE
Results are presented in terms of $\rho$ and normal velocity ${u}_{n}$ for Cartesian (7 levels in 2D and 6 levels in 3D) and adaptive (8 levels in 2D and 7 levels in 3D) computations. The numerical results compare very well with Toro’s numerical experiments.