%20-%20https://sketch.com%20--%3e%3ctitle%3e3ECBBF50-FD08-412A-943D-2CFC5E025A16%3c/title%3e%3cdesc%3eCreated%20with%20sketchtool.%3c/desc%3e%3cg%20id='Stylesheet'%20stroke='none'%20stroke-width='1'%20fill='none'%20fill-rule='evenodd'%3e%3cg%20transform='translate(-926.000000,%20-668.000000)'%20id='BUTTONS'%3e%3cg%20transform='translate(100.000000,%20597.000000)'%3e%3cg%20id='BTN-/-RIGHTS'%20transform='translate(816.000000,%2065.000000)'%3e%3cg%20id='ICON-/-Rights'%20transform='translate(10.000000,%206.000000)'%3e%3cpath%20d='M13.5752841,0%20L1.34943182,0%20C0.59375,0%200,0.586666667%200,1.33333333%20L0,22.6666667%20C0,23.4133333%200.620738636,24%201.34943182,24%20L17.6505682,24%20C18.40625,24%2019,23.3866667%2019,22.6666667%20L19,5.36%20L13.5752841,0%20Z%20M13.5752841,1.89333333%20L17.0838068,5.36%20L13.5752841,5.36%20L13.5752841,1.89333333%20Z%20M1.34943182,22.6666667%20L1.34943182,1.33333333%20L12.1988636,1.33333333%20L12.1988636,5.36%20C12.1988636,6.10666667%2012.7926136,6.69333333%2013.5482955,6.69333333%20L17.6235795,6.69333333%20L17.6505682,22.6666667%20L1.34943182,22.6666667%20Z'%20id='Shape'%20fill='%23595959'%20fill-rule='nonzero'%3e%3c/path%3e%3ccircle%20id='Oval'%20stroke='%23595959'%20stroke-width='1.5'%20fill='%23595959'%20cx='9.5'%20cy='13.5355339'%20r='5.5'%3e%3c/circle%3e%3cpolyline%20id='Path-2'%20stroke='%23FFFFFF'%20stroke-width='1.5'%20transform='translate(9.500000,%2013.035534)%20rotate(45.000000)%20translate(-9.500000,%20-13.035534)%20'%20points='8%2016.0355339%2011%2016.0355339%2011%2010.0355339'%3e%3c/polyline%3e%3c/g%3e%3c/g%3e%3c/g%3e%3c/g%3e%3c/g%3e%3c/svg%3e)
1 Introduction
A large number of authors have studied jets of various kinds in the atmospheres and from volcanic eruptions of celestial bodies. Volcanic eruptions are the most dominant motivation for these, but a significant group concerns jets that reach very large distances or even escape gravity.
An experimental study of volcanic eruptions, Kieffer & Sturtevant (Reference Kieffer and Sturtevant1984), simulated possible eruption scenarios by using various gases and conditions, ranging from high molecular weight refrigerants to helium and from subsonic to supersonic flow. The barrel shocks and Mach disks characteristic of supersonic jets were demonstrated in these experiments and the importance of particulate loading in reducing the sound speed was recognised. Another analysis of the debris landscape of the eruption of Mt St Helens by Kieffer & Sturtevant (Reference Kieffer and Sturtevant1988) also demonstrated features of supersonic jets in volcanic scenarios.
Also concerned with volcanism on the Earth, a computational study of jets that issue vertically from the surface and are sonic at the surface and then expand to become supersonic was published in Ogden (Reference Ogden2011). The interest in this study was in jets within the Earth’s atmosphere, so that a constant gravitational acceleration and a flat planet surface could be assumed. Also, the buoyancy provided by the atmosphere is a significant factor in this regime. Notably, previous one-dimensional results were extended to two-dimensional investigations and various geometrical parameters were examined.
In an earlier paper, Ogden, Glatzmaier & Wohletz (Reference Ogden, Glatzmaier and Wohletz2008), another computational study examined the effect of pressure ratio between the reservoir supplying the jet and the atmosphere, a parameter that also plays a significant part in this study. Stationary shock waves and in some cases oscillatory behaviour were observed. Again, the parameters were such that uniform gravity and a flat surface sufficed.
Another examination of the phenomena exhibited by the Mt St Helens eruption was provided by Orescanin et al. (Reference Orescanin, Prisco, Austin and Kieffer2008). This was an experimental simulation of a supersonic jet impinging on a plane surface. The eruption from Mt St Helens was virtually horizontal, so that the supersonic jet impinged on the surrounding surface. This extensive experimental investigation gives details of how the jet structure is modified by the impingement and how its shock waves interact with the boundary layer created by the eruption flow.
More closely relevant to the interest of the present study of flows that extend to distances comparable to or greater than the radius of the moon concerned, Ingersoll & Ewald (Reference Ingersoll and Ewald2011) presented an analysis of the data from the Cassini mission, in which the plume from a jet issuing from Saturn’s moon Enceladus was observed. Detailed estimates of the particle sizes and mass fluxes led to the conclusion that the jet, which escapes Enceladus’ gravity, feeds the E-ring of Saturn with ice.
A very detailed computational study of the jets on Enceladus was given by Yeoh et al. (Reference Yeoh, Chapman, Goldstein, Varghese and Trafton2015). In this study, the direct simulation Monte Carlo (DSMC) method was used to compute the flow of water vapour laden with tiny ice grains. The DSMC method permits accurate simulation of flows ranging from continuum to rarefied and collisionless conditions. The effects of vent size on the altitude at which the transition from continuum to rarefied conditions occurs and on the escape condition relating to particle size were among the important results.
From the same group, McDoniel et al. (Reference McDoniel, Goldstein, Varghese and Trafton2015) presented DSMC simulation results of the giant Pele plume of Jupiter’s moon Io. This sulphur-dioxide and silicate ash plume rises high above the surface and spills material from a mushrooming plume onto the surface in a characteristic pattern surrounded by a red surface ring. This elegant study showed that the asymmetry of the surface pattern and the geometry of the ring are caused by gas-dynamic interactions between the jet and the lava lake from which it erupts. By varying the ash particle size in the simulation they were able to estimate the particle size from the way it affects the observed features on the moon’s surface. See also McDoniel (Reference McDoniel2015).
The references cited are, of course, only a small excerpt of the extensive literature on this subject. They have been selected to show that very detailed studies of the kind of flows of interest here have been made. It is not the aim here to provide comparable studies that include details of the physical and chemical processes that may be, and often certainly are, active in such flows. Instead, the intent is to consider a model situation, in which the medium of the jet is considered to be an inviscid non-heat-conducting perfect gas. Moreover, the flow is restricted to being symmetrical about the axis of the jet which is taken to issue radially from a spherical, non-rotating moon with inverse-square gravity. Even within these narrow confines, the problem is determined by five independent dimensionless parameters. The exploration of this parameter space is the purpose of this work.
2 Theory
2.1 Equations of motion
The continuity, momentum and energy equations for inviscid, non-heat-conducting flow of a perfect gas, pressure
$p$
, density
${\it\rho}$
, velocity
$\boldsymbol{u}$
, specific reservoir enthalpy
$h$
and specific heat ratio
${\it\gamma}$
acted on by a body force per unit mass
$\boldsymbol{f}$
in an inertial frame and time
$t$
are:
$$\begin{eqnarray}\displaystyle & \displaystyle \frac{\text{D}{\it\rho}}{\text{D}t}+{\it\rho}\boldsymbol{{\rm\nabla}}\boldsymbol{\cdot }\boldsymbol{u}=0, & \displaystyle\end{eqnarray}$$
$$\begin{eqnarray}\displaystyle & \displaystyle {\it\rho}\,\frac{\text{D}\boldsymbol{u}}{\text{D}t}=-\boldsymbol{{\rm\nabla}}p+{\it\rho}\,\boldsymbol{f}, & \displaystyle\end{eqnarray}$$
$$\begin{eqnarray}\displaystyle & \displaystyle {\it\rho}\,\frac{\text{D}h}{\text{D}t}=\frac{\partial p}{\partial t}+{\it\rho}\,\boldsymbol{u}\boldsymbol{\cdot }\boldsymbol{f}, & \displaystyle\end{eqnarray}$$
with
$$\begin{eqnarray}h=\frac{{\it\gamma}}{{\it\gamma}-1}\,\frac{p}{{\it\rho}}+\frac{|\boldsymbol{u}|^{2}}{2}.\end{eqnarray}$$
The specific reservoir enthalpy, denoted here by
$h$
, is sometimes called the specific total enthalpy and given the subscript
$t$
or
$0$
. However, since specific enthalpy only occurs as its reservoir value throughout this study, the subscript is omitted.
In accordance with the assumption of axial symmetry, introduce cylindrical polar coordinates
$(x,y)$
, velocity components
$(u,v)$
, and choose the axial coordinate
$x$
to lie on the axis of the jet issuing radially from the spherical moon’s surface.
The equations become:
$$\begin{eqnarray}\displaystyle & \displaystyle \frac{\partial {\it\rho}}{\partial t}+\frac{\partial {\it\rho}u}{\partial x}+\frac{\partial {\it\rho}v}{\partial y}+\frac{{\it\rho}v}{y}=0, & \displaystyle\end{eqnarray}$$
$$\begin{eqnarray}\displaystyle & \displaystyle \frac{\partial u}{\partial t}+u\,\frac{\partial u}{\partial x}+v\,\frac{\partial u}{\partial y}=-\frac{1}{{\it\rho}}\,\frac{\partial p}{\partial x}+f_{x}, & \displaystyle\end{eqnarray}$$
$$\begin{eqnarray}\displaystyle & \displaystyle \frac{\partial v}{\partial t}+u\,\frac{\partial v}{\partial x}+v\,\frac{\partial v}{\partial y}=-\frac{1}{{\it\rho}}\,\frac{\partial p}{\partial y}+f_{y}, & \displaystyle\end{eqnarray}$$
$$\begin{eqnarray}\displaystyle & \displaystyle \frac{\partial h}{\partial t}+u\,\frac{\partial h}{\partial x}+v\,\frac{\partial h}{\partial y}=\frac{1}{{\it\rho}}\,\frac{\partial p}{\partial t}+u\,f_{x}+v\,f_{y}. & \displaystyle\end{eqnarray}$$
With an inverse-square gravity field, and acceleration due to gravity
$g$
at the moon’s surface (
$r=r_{0}$
),
$$\begin{eqnarray}f_{x}=-g\,\frac{r_{0}^{2}}{r^{2}}\,\frac{x}{r},\quad f_{y}=-g\,\frac{r_{0}^{2}}{r^{2}}\,\frac{y}{r},\quad r^{2}=x^{2}+y^{2}.\end{eqnarray}$$
Although the computational part deals with unsteady flow, it is instructive to consider steady flow first. A possible scenario is that an erupting jet flow eventually reaches a steady state. The appropriate equations are (2.5)–(2.8) with the five time derivatives set to zero. Along the
$x$
-axis
$v=0$
,
$y=0$
and
$x=r$
, so that, along the jet axis, which by assumption is also a streamline, the steady-flow energy equation becomes (
$u\neq 0$
):
$$\begin{eqnarray}\frac{\text{d}h}{\text{d}r}=-g\,\frac{r_{0}^{2}}{r^{2}},\end{eqnarray}$$
stating that the only agency that reduces
$h$
along
$x=r$
is the work done against gravity. This may be integrated to give
$$\begin{eqnarray}h(r)=h(r_{0})-g\,r_{0}\left(1-\frac{r_{0}}{r}\right).\end{eqnarray}$$
2.2 Escape velocity
Assume that there is a constriction at the surface of the moon, so that the flow from a reservoir inside the moon to the outside passes through a minimum cross-sectional flow area. Assume also that the pressure ratio between the reservoir and the moon’s surface (
$p_{r}/p_{0}$
) is sufficiently large to make the constriction a sonic throat, so that (
$a=$
speed of sound)
$$\begin{eqnarray}u_{0}^{2}=a_{\ast }^{2}=\frac{{\it\gamma}p_{\ast }}{{\it\rho}_{\ast }}.\end{eqnarray}$$
Here the asterisk denotes the state at the sonic throat. Equation (2.12) identifies the jet speed at the surface (
$u_{0}$
) with the sonic speed within the jet at the throat (
$a_{\ast }$
). Note that the state at the throat (
$p_{\ast },{\it\rho}_{\ast }$
) is determined by the isentropic expansion from the reservoir state (
$p_{r},{\it\rho}_{r}$
) and is independent of the surface state (
$p_{0},{\it\rho}_{0}$
).
Substituting in (2.11) for
$h(r_{0})=h_{\ast }$
from (2.4) and setting
$h(r)=0$
at
$r=\infty$
then gives
$$\begin{eqnarray}0=\frac{{\it\gamma}}{{\it\gamma}-1}\,\frac{p_{\ast }}{{\it\rho}_{\ast }}+\frac{u_{0}^{2}}{2}-g\,r_{0}=\frac{a_{\ast }^{2}}{{\it\gamma}-1}+\frac{u_{0}^{2}}{2}-g\,r_{0}=\frac{{\it\gamma}+1}{{\it\gamma}-1}\,\frac{u_{0}^{2}}{2}-g\,r_{0}.\end{eqnarray}$$
This represents the escape condition at which the jet speed at the surface is just large enough to overcome gravity. Solving for this special speed,
$$\begin{eqnarray}u_{0}^{2}=u_{eg}^{2}=2\,g\,r_{0}\,\frac{{\it\gamma}-1}{{\it\gamma}+1}.\end{eqnarray}$$
The ‘
$g$
’ in the subscript ‘
$eg$
’ distinguishes
$u_{eg}$
as the escape velocity of a gas, because it is significantly smaller than (
${\it\gamma}>1$
) the escape velocity
$u_{es}=\sqrt{2gr_{0}}$
of a solid body. The reason is that the thermal energy of the gas at
$r=r_{0}$
is converted to ordered kinetic energy in the isentropic expansion of the gas as it flows to increasing
$r$
; isentropic, because dissipative processes have been excluded by assumption.
This result is not to be confused with the so-called Jeans escape. The latter describes the escape of molecules from the collision-rare condition at high altitude by virtue of their thermal speed. The result in (2.14) represents a (collision-rich) continuum phenomenon. Some particular values of the escape velocities
$u_{eg}$
and
$u_{es}$
on moons in our solar system are given in the Appendix.
2.3 Maximum penetration radius
Still within the confines of steady flow, consider the case when the jet velocity at the surface is smaller than the escape velocity,
$u_{0}<u_{eg}$
. In that case,
$h\rightarrow 0$
at a finite value of
$r$
,
$r_{p}$
say. Substituting
$h(r)=0$
and
$r=r_{p}$
in (2.11) and using (2.12) and (2.14), the maximum steady-state penetration radius is obtained:
$$\begin{eqnarray}\frac{r_{p}}{r_{0}}=\left(1-\frac{u_{0}^{2}}{u_{eg}^{2}}\right)^{-1}.\end{eqnarray}$$
The surface to escape velocity ratio may be related to the reservoir conditions by using the relations for isentropic flow from the reservoir (at zero velocity) to the throat (sonic condition):
$$\begin{eqnarray}\frac{u_{0}^{2}}{u_{eg}^{2}}=\frac{{\it\gamma}}{{\it\gamma}-1}\,\frac{p_{r}}{{\it\rho}_{r}}\,\frac{1}{gr_{0}}.\end{eqnarray}$$
2.4 Atmosphere
The earlier mention of the pressure ratio
$p_{r}/p_{0}$
implies that the moon has an atmosphere. Two idealised cases of the structure of the atmosphere are considered: isothermal and adiabatic. In an isothermal atmosphere,
$p/{\it\rho}=p_{0}/{\it\rho}_{0}$
. Using this in the static equilibrium of the atmosphere gives
$$\begin{eqnarray}\frac{\text{d}p}{\text{d}r}=-{\it\rho}g\,\frac{r_{0}^{2}}{r^{2}}=-\frac{{\it\rho}_{0}gr_{0}^{2}}{p_{0}}\,\frac{p}{r^{2}},\end{eqnarray}$$
and integrating,
$$\begin{eqnarray}\frac{p}{p_{0}}=\frac{{\it\rho}}{{\it\rho}_{0}}=\exp \left[-\frac{{\it\rho}_{0}gr_{0}}{p_{0}}\left(1-\frac{r_{0}}{r}\right)\right].\end{eqnarray}$$
The corresponding result for an adiabatic atmosphere, for which
$p/{\it\rho}^{{\it\gamma}}=p_{0}/{\it\rho}_{0}^{{\it\gamma}}$
, is
$$\begin{eqnarray}\frac{p}{p_{0}}=\left(\frac{{\it\rho}}{{\it\rho}_{0}}\right)^{{\it\gamma}}=\left[1-\frac{{\it\gamma}-1}{{\it\gamma}}\,\,\frac{{\it\rho}_{0}gr_{0}}{p_{0}}\left(1-\frac{r_{0}}{r}\right)\right]^{{\it\gamma}/({\it\gamma}-1)}.\end{eqnarray}$$
Equation (2.18) is, of course, the limiting case of (2.19) as
${\it\gamma}\rightarrow 1$
. In both cases, the scale height
$\ell _{s}$
, defined as the height above the surface where the pressure would go to zero if it were to decrease from
$p_{0}$
with its surface gradient, is
$$\begin{eqnarray}\ell _{s}=\frac{p_{0}}{{\it\rho}_{0}g}.\end{eqnarray}$$
Figure 1 shows examples of the atmospheric structure.
Examples of the structure of isothermal (full line) and adiabatic (chain-dotted line,
${\it\gamma}=1.4$
) atmospheres with
$\ell _{s}/r_{0}=0.1$
.
2.5 Parameters of the problem
Within the scheme of the assumptions made, the problem reduces to one in which any dimensionless quantity
$Q$
depends on five independent parameters:
$$\begin{eqnarray}Q=Q\left(\frac{p_{r}}{p_{0}},\,\frac{{\it\rho}_{r}}{{\it\rho}_{0}},\,\frac{p_{0}}{{\it\rho}_{0}\,g\,r_{0}},\,{\it\gamma},\,\frac{r_{\ast }}{r_{0}}\right),\end{eqnarray}$$
where
$r_{\ast }$
is the radius of the opening in the moon’s surface. The third of these is
$\ell _{s}/r_{0}$
, which is held constant at 0.1 throughout this study. The remaining four independent parameters represent the space that is explored by the computational study. It turns out that the variable that influences the results most strongly is the velocity ratio
$u_{0}/u_{eg}$
which depends on all of the first four parameters in (2.21), see (2.16). When describing the results this variable is therefore chosen in place of
${\it\rho}_{r}/{\it\rho}_{0}$
. It is also convenient for comparison with the analytical result of (2.15). Most of the results presented in the following have been obtained with
$r_{\ast }/r_{0}=0.1$
. However, some results showing the effect of the opening radius are presented in § 4.6.
The space and time scaling may be expressed by
$$\begin{eqnarray}\left(\frac{r}{r_{0}},\sqrt{\frac{{\it\rho}_{r}}{p_{r}}}gt\right).\end{eqnarray}$$
3 Computational setup
3.1 Equations
For numerical solutions of the Euler equations (2.5)–(2.8) they are usually written in the equivalent form
$$\begin{eqnarray}\frac{\partial \boldsymbol{W}}{\partial t}+\frac{\partial \boldsymbol{F}}{\partial x}+\frac{\partial \boldsymbol{G}}{\partial y}={\it\rho}\boldsymbol{B}-\frac{{\it\rho}v}{y}\boldsymbol{S},\end{eqnarray}$$
in which the vector of dependent variables,
$\boldsymbol{W}$
, the flux vectors
$\boldsymbol{F}$
and
$\boldsymbol{G}$
, specific body-force vector
$\boldsymbol{B}$
and axial symmetry vector
$({\it\rho}v/y)\boldsymbol{S}$
are determined by
$$\begin{eqnarray}\boldsymbol{W}=\left[\begin{array}{@{}c@{}}{\it\rho}\\ {\it\rho}u\\ {\it\rho}v\\ {\it\rho}h-p\end{array}\right],\quad \boldsymbol{F}=\left[\begin{array}{@{}c@{}}{\it\rho}u\\ {\it\rho}u^{2}+p\\ {\it\rho}uv\\ {\it\rho}uh\end{array}\right],\quad \boldsymbol{G}=\left[\begin{array}{@{}c@{}}{\it\rho}v\\ {\it\rho}vu\\ {\it\rho}v^{2}+p\\ {\it\rho}vh\end{array}\right],\quad \boldsymbol{B}=\left[\,\begin{array}{@{}c@{}}0\\ f_{x}\\ f_{y}\\ f_{x}u+f_{y}v\end{array}\right],\quad \boldsymbol{S}=\left[\begin{array}{@{}c@{}}1\\ u\\ v\\ h\end{array}\right],\end{eqnarray}$$
with
$f_{x}$
and
$f_{y}$
as given by (2.9).
3.2 Software
The software system Amrita, constructed by James Quirk, see Quirk (Reference Quirk1998), was used. Amrita is a system that automates and packages computational tasks in such a way that the packages can be combined (dynamically linked) according to instructions written in a high-level scripting language. The present application uses features of Amrita that include the automatic construction of the Euler solver, documentation of the code, adaptive mesh refinement according to simply chosen criteria, and scripting-language-driven computation, archiving and post-processing of the results. The automation of the assembly and sequencing of the tasks makes for dramatically reduced possibility of hidden errors. It also makes computational investigations transparent and testable by others. The ability to change one package at a time, without changing the rest of the scheme, facilitates detection of sources of error. The Euler solver generated for the present computations was an operator-split scheme with HLLE flux and kappa-MUSCL reconstruction. The bare bones Amrita system did not include a body force, so that the necessary code to include
$\boldsymbol{B}$
was constructed and extensively debugged.
3.3 Initial and boundary conditions, discretisation
Figure 2 shows a pseudo-schlieren image of an example of the initial condition of a computation. The (
$x,y$
) plane is discretised by a Cartesian grid of
$600\times 600$
coarse-grid cells that are adaptively refined by a factor of 3 to make an effective grid of
$1800\times 1800$
cells. The criterion for adaptation is a chosen threshold of the magnitude of the fractional density gradient (density gradient divided by local density). The grey-shading of the visualisation is a monotonic function of the magnitude of the fractional density gradient; fractional in order to maintain sensitivity at small densities. In the initial state shown in figure 2 the spatial variation of density is solely due to the atmospheric structure, in this instance isothermal. The graph above shows the distribution of pressure ratio
$p/p_{0}$
(solid line) and density ratio
${\it\rho}/{\it\rho}_{0}$
(dashed line) along the symmetry axis, which is the bottom edge of the figure. The strength of the gravity field is chosen to give
$\ell _{s}=0.1r_{0}$
, as in figure 1.
Meridional section through one quadrant of the flow field around a hollow spherical moon with an opening at a pole, showing the initial conditions for the computation (see text).
$p_{r}/p_{0}=15$
,
${\it\gamma}=1.4$
,
${\it\rho}_{r}/{\it\rho}_{0}=7.500$
,
$u_{0}/u_{e}=0.8367$
.
One quadrant of the moon is shown with a radius of 100 coarse grid cells. It is hollow with the reservoir inside, in which the pressure and density are initially uniform (see the graph). The solid moon shell is modelled by a level set that specifies the smallest distance of a field point from any solid boundary. It makes the solid surface a no-through-flow boundary. The opening on the axis is initially closed. The conditions on most of the domain boundary are open. The exceptions are the symmetry axis and the portion of the left boundary that lies inside the reservoir. The latter is reflecting, i.e. like a solid boundary.
4 Computational results
In all the computational results presented, the value of
$p_{0}/({\it\rho}_{0}r_{0}g)$
is held constant at 0.1 thus fixing the dimensionless scale height of the atmosphere. The effect of the atmosphere, both isothermal and adiabatic, can then be studied by varying
$p_{r}/p_{0}$
. Although most of the results that will be shown are those with
${\it\gamma}=1.4$
, the effect of
${\it\gamma}$
is discussed at the end of § 4.5.
4.1 Some features of the flow
To set the scene, it is useful to describe some of the features of the jet. Figure 3 shows an image of a late phase of a jet with
$p_{r}/p_{0}=60$
,
$u_{0}/u_{eg}=0.8367$
,
${\it\gamma}=1.4$
,
$r_{\ast }/r_{0}=0.1$
. The expansion of the gas from the hot, high-pressure reservoir through the sonic constriction to a supersonic flow encounters an opposing pressure from the atmosphere that causes the flow to be deflected back closer to the
$x$
-direction by the ‘barrel’ shock, and parallel to a shear layer that issues from the lip of the constriction. The resistance provided by gravity causes a nearly normal shock to appear in the flow close to the axis. This is sometimes referred to as the Mach disk in jet parlance.
Example of result of a computation showing features of a jet that has reached near-steady state.
$p_{r}/p_{0}=60$
,
$u_{0}/u_{eg}=0.8367$
,
${\it\gamma}=1.4$
,
$r_{\ast }/r_{0}=0.1$
.
The normal shock and the barrel shock intersect at a triple shock point, as a third shock results from the clash of the flows that have passed through the normal and barrel shocks. These two flows are separated by a shear layer that issues from the triple point on the downstream side of the third shock. The jet reaches a maximum penetration radius on its axis, and the material that has passed through the normal shock – as well as the material that has come through the barrel and third shocks – spills to lower latitudes (pole on
$x$
-axis) and falls down onto the moon surface at an impingement area. The outside field shows some weak shock waves that have been emitted into the tenuous atmosphere by the unsteadiness of the jet.
The triple ‘point’ is, of course a triple ‘line’. In this axisymmetric flow the shock waves are axisymmetric surfaces and the intersection of these three shock surfaces is the triple line, a perfect circle around the symmetry axis. Its trace in the meridional plane is the triple point.
4.2 Development of the jet with time
Figure 4 shows a time sequence of the development of the flow with the same parameters as in the case of figure 3. After the constriction is opened to start the jet erupting, figure 4(a) shows a nearly spherical shock wave that is generated by the initial eruption and that has propagated into the atmosphere to
$x\simeq 3r_{0}$
on the axis. This is followed by the mushrooming jet that already contains some of the features illustrated in figure 3. In figure 4(a–c) this develops further and the rim of the mushroom begins to fall. In the second row, the jet height on the axis decreases somewhat as impingement wanders toward the equator, and the barrel shock exhibits some unsteadiness. In the bottom row, the supersonic part of the jet becomes virtually steady, but the rest of the jet remains unsteady, so that the ‘maximum steady-state penetration radius’ necessarily has to be quoted with an error bar that characterises the fluctuation around a mean value. Note that, during the unsteady phase of the flow, the mushroom overshoots the eventual steady-state penetration radius.
Time sequence of pseudo-schlieren images.
$p_{r}/p_{0}=60$
,
$u_{0}/u_{eg}=0.8367$
,
${\it\gamma}=1.4$
,
$r_{\ast }/r_{0}=0.1$
. The time intervals between images are equal. See supplementary movie 1 available at http://dx.doi.org/10.1017/jfm.2016.184.
It should be pointed out that during the whole of the development of the jet the pressure and density changes in the reservoir are extremely small in all cases, indicating that the drainage of the reservoir is insignificant. The reservoir pressure and density also remain uniform except in the immediate vicinity of the throat where they fall in accordance with steady isentropic expansion.
Figure 5 shows another, similar, time sequence, with
$p_{r}/p_{0}=250$
. Particularly evident at this high pressure ratio, the initial eruption sends a narrow axial jet far ahead of the main mushroom. As time evolves, this axial jet falls back into the mushrooming main jet. Also, the transverse expansion of the jet is greater than at lower pressure ratio, so that the shear layer leaves the throat lip at a larger angle from the axis and the barrel shock is conical rather than barrel-shaped. Furthermore, the impingement area moves beyond the equator at later times. Because this is outside the computational domain and flow across the boundary can occur in both directions, it has to be considered with some care. The upstream influence from the left domain boundary is less important if the impingement area is inside the domain.
Time sequence of pseudo-schlieren images.
$p_{r}/p_{0}=250$
,
$u_{0}/u_{eg}=0.8367$
,
${\it\gamma}=1.4$
,
$r_{\ast }/r_{0}=0.1$
. The time intervals between images are the same as in figure 4.
4.3 Effect of the velocity ratio
Figure 6 shows the effect of changing the velocity ratio
$u_{0}/u_{eg}$
in late-phase images at a pressure ratio of 250. As the velocity ratio is increased, the maximum steady-state penetration radius of the jet increases dramatically. The conical shape of the barrel shock that is characteristic for this high pressure ratio is evident in all frames and its angle is approximately constant. The impingement area moves to lower latitudes with increasing
$u_{0}/u_{eg}$
at low values, but disappears entirely as
$u_{0}/u_{eg}\rightarrow 1$
. In the last image, at
$u_{0}/u_{eg}=0.9354$
, the jet continues to expand and the normal shock is outside the region shown. In fact, at this condition it even leaves the computational domain completely.
Effect of the velocity ratio.
$p_{r}/p_{0}=250$
,
${\it\gamma}=1.4$
,
$r_{\ast }/r_{0}=0.1$
.
$u_{0}/u_{eg}=0.5401$
(a), 0.6236 (b), 0.7071 (c), 0.7638 (d), 0.8367 (e), 0.9354 (f).
Figure 7 is the same as figure 6, but for a pressure ratio of 60. The features are qualitatively similar, except that the barrel shock exhibits a barrel shape and the angle at which it leaves the throat lip is smaller.
Effect of the velocity ratio.
$p_{r}/p_{0}=60$
,
${\it\gamma}=1.4$
,
$r_{\ast }/r_{0}=0.1$
.
$u_{0}/u_{eg}=0.5401$
(a), 0.6614 (b), 0.7071 (c), 0.7638 (d), 0.8367 (e), 0.9354 (f).
4.4 Effect of the pressure ratio
Next, figure 8 shows the effect of changing the pressure ratio in late-phase images at a constant velocity ratio of 0.7071. Figure 8(a) with
$p_{r}/p_{0}=4$
exhibits a tiny barrel shock and normal shock. The shock pattern increases dramatically in size, the shear layer angle increases and the impingement area moves to lower latitudes with increasing pressure ratio. These features are qualitatively the same in figure 9, showing late-phase images at a larger velocity ratio.
Effect of the pressure ratio.
$u_{0}/u_{eg}=0.7071$
,
${\it\gamma}=1.4$
,
$r_{\ast }/r_{0}=0.1$
.
$p_{r}/p_{0}=4$
(a), 15 (b), 60 (c), 250 (d).
Effect of the pressure ratio.
$u_{0}/u_{eg}=0.8367$
,
${\it\gamma}=1.4$
,
$r_{\ast }/r_{0}=0.1$
.
$p_{r}/p_{0}=4$
(a), 15 (b), 60 (c), 250 (d).
Figure 10 shows two images of jets that escape gravity. A feature of such cases is that the jet spills little material back onto the moon’s surface. Thus, the decrease of the impingement latitude that accompanies increase of
$u_{0}/u_{eg}$
at small values is followed by an impingement-latitude increase at higher velocity ratios.
Two cases where
$u_{0}/u_{eg}>1$
,
${\it\gamma}=1.4$
,
$r_{\ast }/r_{0}=0.1$
,
$u_{0}/u_{e}=1.0801$
. (a)
$p_{r}/p_{0}=15$
,
${\it\rho}_{r}/{\it\rho}_{0}=4.500$
; (b)
$p_{r}/p_{0}=60$
,
${\it\rho}_{r}/{\it\rho}_{0}=18$
. Note how, when the jet escapes gravity, very little of it spills back onto the moon. The ripples in the jet, as in other cases at later times, are caused by the reverberation of waves in the reservoir.
4.5 Maximum steady-state penetration radius
While figures 6 and 7 show a substantial increase in
$r_{p}$
as the velocity ratio is increased, it is interesting to compare this growth quantitatively with (2.15). This is done in figure 11 for pressure ratios of 15 and 3. Each of the square points in the figure represents an average of measurements from several late-phase computed images with an error bar representative of the fluctuation amplitude. The dashed line is arbitrarily drawn at 70 % of the altitude above the moon surface of the theoretical line because it passes approximately through the points of the case
$p_{r}/p_{0}=15$
. This also serves as a reference in the case
$p_{r}/p_{0}=3$
in figure 11(b) and in figure 12.
Plot of dimensionless steady-state jet penetration radius versus
$u_{0}/u_{eg}$
for
${\it\gamma}=1.4$
,
$r_{\ast }/r_{0}=0.1$
. (a)
$p_{r}/p_{0}=15$
; (b)
$p_{r}/p_{0}=3$
. For notation, see text.
Plot of dimensionless steady-state jet penetration radius versus
$u_{0}/u_{eg}$
for
${\it\gamma}=1.4$
,
$r_{\ast }/r_{0}=0.1$
. (a)
$p_{r}/p_{0}=50$
; (b)
$p_{r}/p_{0}=250$
. For notation, see text.
Also shown in these figures as a chain-dotted line is the dimensionless penetration radius that would be reached by a solid object leaving the moon’s surface at the same speed. This is very much smaller than that of the gas because the escape velocity of the solid is larger by a factor of
$\sqrt{({\it\gamma}+1)/({\it\gamma}-1)}$
than that of the gas.
Figure 12 repeats figure 11 for the case of pressure ratios of 50 and 250. The curves are exactly the same as those in figure 11. Note that, as the pressure ratio is increased from 3 to 250, the points approach the theoretical isentropic curve more closely. As
$p_{r}/p_{0}$
is increased, the impediment to the jet that is provided by the presence of the atmosphere is reduced. This is shown in another form in figure 13(a).
Plot of dimensionless steady-state jet penetration radius. (a) Effect of
$p_{r}/p_{0}$
for
$u_{0}/u_{eg}=0.71$
,
${\it\gamma}=1.4$
,
$r_{\ast }/r_{0}=0.1$
. (b) Effect of
${\it\gamma}$
for
$u_{0}/u_{eg}=0.71$
and
$p_{r}/p_{0}=15$
,
$r_{\ast }/r_{0}=0.1$
. The heavy black line represents the value given by (2.15) in both graphs.
It is interesting to observe how strongly
${\it\gamma}$
affects
$r_{p}/r_{0}$
(see figure 13
b). One might speculate that the discrepancy between the computed penetration radius and the theoretical isentropic value could be related in part to the entropy increase across the normal shock which is given by
$$\begin{eqnarray}\frac{{\rm\Delta}s}{R}=\frac{1}{{\it\gamma}-1}\,f(M,{\it\gamma}),\end{eqnarray}$$
where
$s$
is specific entropy,
$R$
is the specific gas constant,
$M$
is the Mach number upstream of the shock and the function
$f$
depends relatively weakly on
${\it\gamma}$
. The behaviour exhibited by figure 13(b), showing the dependence of
$r_{p}/r_{0}$
on
${\it\gamma}$
with
$p_{r}/p_{0}$
,
$u_{0}/u_{eg}$
and
$r_{\ast }/r_{0}$
held constant, seems to mirror the prefactor of
$f$
. However, since
$M$
is not quite the same in the five cases plotted, a quantitative comparison is not warranted.
4.6 Effect of throat radius
In all of the results presented so far,
$r_{\ast }/r_{0}=0.1$
. To examine the effect of this parameter, several computations were performed with the reduced value of 0.03. At this value the radius of the opening amounts to only nine fine grid cells and is near the limit of reasonable resolution. The throat radius does not enter the derivation of (2.15), so that it is not expected to influence the steady-state penetration radius significantly.
Figure 14 shows four late-phase images from computations with
$r_{\ast }/r_{0}=0.1$
in (a,c) and 0.03 in (b,d). Figure 14(a,b) is at
$p_{r}/p_{0}=60$
, (c,d) at 250. The high velocity ratio has been chosen for this comparison because the effect of throat radius is largest there. However, note that, though the throat radius is changed by a factor of more than 3, the scale of the jet geometry is clearly not determined by
$r_{\ast }$
. The throat radius does influence the maximum steady-state penetration radius to some extent, as can be seen in figure 15. Reduction of
$r_{\ast }$
also reduces the transverse extent of the jet somewhat.
(a,c) From computations with a throat diameter
$r_{\ast }/r_{0}=0.1$
and (b,d) with
$r_{\ast }/r_{0}=0.03$
. (a,b) From computations with
$p_{r}/p_{0}=60$
,
${\it\rho}_{r}/{\it\rho}_{0}=27$
; (c,d) with
$p_{r}/p_{0}=250$
,
${\it\rho}_{r}/{\it\rho}_{0}=112.500$
. Velocity ratio
$u_{0}/u_{eg}=0.8819$
,
${\it\gamma}=1.4$
.
Maximum steady-state penetration radius obtained with the reduced throat diameter (
$r_{\ast }/r_{0}=0.03$
). Because the curves in this figure are the same as those in figures 11 and 12, it can be seen that the smaller throat gives slightly reduced penetration relative to the case
$r_{\ast }/r_{0}=0.1$
. Since error bars would be confusing here, they have been omitted.
4.7 Atmosphere structure
All the computational results shown are for an isothermal atmosphere. The same parameter space has been covered by computations with an adiabatic atmosphere. However, the results differ so little from the cases with isothermal atmosphere that no additional information would be conveyed by presenting them. The insensitivity of the results to atmospheric structure indicates that the effect of the atmosphere comes about through the atmospheric pressure and density in a region close to the moon’s surface, where the two structures are similar, see figure 1.
5 Experimental analogy
The problem of simulating a spherical inverse-square gravity field in the laboratory poses enormous difficulties. The problem at hand does, however, suggest that one might make use of the shallow water wave analogy for compressible flow in a very approximate way. To this end, consider a solid surface that is symmetrical about a vertical axis and has the shape
$$\begin{eqnarray}\frac{z}{b}=\frac{z_{0}}{b}-\frac{b}{r},\end{eqnarray}$$
where
$r$
is the horizontal radius from the axis,
$z$
is the vertical elevation of the surface with asymptote
$z_{0}$
at
$r=\infty$
and
$b$
is a constant. In a two-dimensional fashion
$z(r)$
may be considered as the gravitational potential well of an inverse-square gravity field.
If, at some finite radius, a water jet issues tangentially up onto this surface in a radial direction, i.e. with no initial circumferential velocity component, then its kinetic energy is converted to potential energy as it flows up along the surface; it slows down and, if its initial speed is insufficient to escape the gravitational potential well, reaches only a finite
$r$
. If its speed at small
$r$
is greater than the shallow water wave speed (‘super-undal flow’ corresponding to supersonic flow in the analogy) this means that it will encounter a hydraulic jump (corresponding to a normal shock in the analogy) to slow it down abruptly to sub-undal flow, allowing it to spill transversely to fall back down into the potential well on both sides.
In order to examine the analogy, an apparatus was constructed that consists of a cylindrical tube of 115 mm outside diameter standing in a water bath and fitted with an axisymmetric surface of the shape of (5.1) that was made by three-dimensional printing (see figure 16). In the experimental apparatus the value of
$b$
in (5.1) is 58 mm and the outer radius of the surface is 200 mm. The tube, representing the moon, was fitted with an orifice directed tangentially to the inside slope (
$\simeq 45^{\circ }$
) of the surface. The cross-section of the orifice blends smoothly from the inside wall of the tube to a circular and then to a 3 : 1 elliptical shape, with long axis horizontal, at the outside wall. A small pump feeds the inside of the tube with water and the level inside the tube is kept constant by bleeding off at the appropriate level. There is a small gap between the inside of the potential well surface and the outside of the tube.
Experimental apparatus. The vertical tube represents the moon, the partial surrounding surface acts as the potential well. The little black pump tops up the water level in the tube.
In figure 17(a,b) an experimental jet flow is shown in a case with a moderate velocity ratio in two views, one from the side and one from the top. Figure 17(c,d) shows top views of two cases with larger velocity ratio. In (d), the jet just escapes ‘gravity’. It is clear in these cases that the normal shock is there in the form of a normal hydraulic jump that terminates the super-undal flow. However, the barrel shock, though present, does not come from the lip of the jet, but touches the ‘moon’ at a finite distance from the jet lip. This is because of the absence of an ‘atmosphere’, permitting the flow to expand up to a right angle from the symmetry plane.
(a) Side view of experimental jet with the water level just level with the jet exit. (b) A similar case in top view. Note the super-undal flow and hydraulic jump. The water has been dyed with red food dye for contrast. (c) Top view of a case with larger velocity ratio. (d) Top view of a case where the jet just escapes the gravitational potential well. See supplementary movie 2 available at http://dx.doi.org/10.1017/jfm.2016.184.
The qualitative features of these jets in the crude, quasi-two-dimensional analogy are remarkably similar to the gas-dynamic axisymmetric moon jets. The size of the apparatus is too small to show the triple point and third shock. The scale of surface tension effects masks these features. In the analogy, there is, of course, no equivalent of the reservoir temperature, so that one less parameter is available. However, it is possible to give the jet an equivalent of the impediment of an atmosphere by raising the water level to a point at which the jet exit is submerged.
Results showing such cases are shown in figure 18 with different velocity ratios. The barrel shocks in these clearly show the effect of the ‘atmosphere’, though its position is influenced by the entrainment of two air bubbles near the ‘moon’ surface in some of them.
Four top views of flows in which the jet exit is submerged. Increasing velocity ratio from (a–d).
6 Discussion
In the large number of computations that were performed to explore the parameter space of the model problem, there are, of course, many features that could be studied in detail. Among these are the location of the normal shock on the axis, the angles of the barrel shock and shear layer at the jet lip, the maximum width of the barrel shock and the impingement location. Relations of some of these (such as e.g. the shear layer angle) to the independent parameters might even be determined theoretically. However, the theoretical and computational results bring out the maximum steady-state penetration radius of the jet on its axis as the quantity of most prominent interest, so that it dominates the scope of this study.
6.1 Penetration radius
The results show that, even within the assumptions of the computational scheme, there are four effects that cause the computed penetration radius to depart from the theoretical isentropic value given by (2.15).
The first is the effect of pressure unsteadiness. The term
$\partial p/\partial t$
in (2.8) is positive in almost the whole flow field during the early phase of the flow after first eruption. Consequently, it contributes positively to the change in
$h$
, and a part of the jet penetrates beyond the steady value before falling back into the mushroom. It even causes the mushroom itself to overshoot its late-phase penetration radius before settling down.
The second is the effect of the atmosphere. It manifests itself through the parameter
$p_{r}/p_{0}$
. The impediment of the jet by the atmosphere reduces the steady penetration radius. Thus, as
$p_{r}/p_{0}$
is increased, the computed penetration radius approaches the isentropic value more closely (see figure 13
a).
The third is a consequence of the assumption of zero dissipation. As discussed in § 4.5, some dissipative effects are related to the entropy rise through shock waves, a feature that is correctly reproduced in the computation. It manifests itself most clearly through the effect of the parameter
${\it\gamma}$
(see figure 13
b).
Finally, reducing the diameter of the opening in the moon’s surface causes the penetration radius to be reduced somewhat, see figure 15. It also causes the jet to be a little narrower.
Though all four of these effects influence the geometry of the jet to some extent, by far the dominant and scale-determining dimension of the flow is
$r_{p}$
as given by (2.15).
6.2 Relation of the model to real flows
Some effects that have been neglected were mentioned earlier; among them gas rarefaction, radiation, multi-phase flow. Three particular assumptions are discussed here.
The assumption of axial symmetry constrains the results. It means that in the computed results all vortex lines (integral curves of the vorticity field) are perfect circles around the symmetry axis. In three-dimensional flow such circular vortex lines would be broken up into more irregular forms. Since vortex lines resist stretching, the constraint of axial symmetry may make the jet narrower than a three-dimensional one would be.
It is appropriate also to note a consequence of the assumption that the medium is a perfect gas. Even if there are no liquid or solid ejecta accompanying the gas from within the reservoir, the gas cools down during the expansion and eventually condenses. Just like the conversion of thermal energy into ordered kinetic energy that is included in this study, condensation would also free latent heat and give additional impetus to the jet. However, the condensate is very much heavier than the gas and will rain down relative to the gas jet, thus impeding it. Including this effect in the problem would be a major complication that would likely necessitate studies of special cases; studies nevertheless well worth undertaking.
The assumption that the flow is inviscid and non-heat-conducting is probably satisfactory for celestial bodies in most of the flow field. Obvious regions where viscous effects enter the problem are the shear layers and the turbulent regions in the subsonic parts of the flow. These are considered to be of minor importance to quantities such as
$r_{p}/r_{0}$
. The boundary layer on the walls in the throat region merely causes the effective throat radius to be slightly smaller than the actual one. In regions at high altitude gas rarefaction effects do become important. Though heat conduction, like viscosity, is of minor importance, heat loss by radiation may become important. Its effect would be to reduce
$r_{p}/r_{0}$
.
7 Conclusions
A simplified form of the complex flow that occurs when supersonic jets are emitted from the surface of a celestial body (moon) is considered by treating an abstracted model problem. In the model problem the medium is an inviscid non-heat-conducting perfect gas that leaves the surface of the spherical, non-rotating moon radially at sonic speed. The flow is assumed to be symmetrical about the axis of the jet. Theoretical considerations lead to an analytical expression for the escape velocity of a gas that is much smaller than the escape velocity of a solid body. An analytical expression is also obtained for the maximum steady-state penetration radius of such a jet. The model problem is determined by five independent dimensionless parameters that include the effect of an atmosphere. An extensive computational study of this parameter space provides a rich picture of the features of such jets. The computations broadly confirm the theoretical results, but show that minor deviations from them result from the impediment by the atmosphere, by dissipation in shock waves and by effects of the opening size in the moon’s surface. An experiment that uses the shallow water analogy and an axisymmetric surface to simulate a gravitational potential well is performed. The surface surrounds a vertical circular cylinder that simulates the moon. It is filled with water and fitted with a nozzle that provides a jet issuing tangentially to the surface. The experimental flows provide a good mirror of the qualitative features of the gas-dynamical jets.
Acknowledgements
My interest in this problem was sparked by discussions with D. Goldstein during his visit to Caltech in the fall of 2014. B. Valiferdowsi designed and built the experimental apparatus and also helped me with the performance of the experiments. I thank them both.
Supplementary movies
Supplementary movies are available at http://dx.doi.org/10.1017/jfm.2016.184.
Appendix
Table 1 gives examples of moons in our solar system together with the gases in their atmospheres and the solid and gas escape velocities.
Examples of parameters of moons in our solar system.
'/%3e%3cpath%20d='M17%207.971L15.711%206.677%2010.912%2011.495%2010.912%200%209.088%200%209.088%2011.495%204.289%206.677%203%207.971%2010%2015z'%20transform='translate(8%208)'/%3e%3c/g%3e%3c/g%3e%3c/svg%3e)
