Aps/MaclaurinSpheroidFreeFall

Free-Fall Collapse of an Homogeneous Spheroid

Free-Fall Collapse of an Homogeneous Spheroid

"What is the form of the collapse under gravitational forces of a uniformly rotating spheroidal gas cloud? In the special case where initially the gas is absolutely cold and of uniform density within the spheroid, we show that the collapse proceeds through a series of uniform, uniformly rotating spheroids until a disk is formed."

— Drawn from the first paragraph of 📚 Lynden-Bell (1962)

Simplified Governing Relations

When studying the dynamical evolution of strictly axisymmetric configurations, it proves useful to write the spatial operators in our overarching set of principal governing equations in terms of cylindrical coordinates, $(\varpi ,\varphi ,z)$, and to simplify the individual equations as described in our accompanying discussion. The resulting set of simplified governing relations is …

Equation of Continuity

$\frac{d\rho }{dt}+\frac{\rho }{\varpi }\frac{\partial }{\partial \varpi }\left[\varpi \dot{\varpi }\right]+\rho \frac{\partial }{\partial z}\left[\rho \dot{z}\right]=0$

Euler Equation

${\hat{e}}_{\varpi }[\frac{d\dot{\varpi }}{dt}-\frac{j^2}{{\varpi }^3}]+{\hat{e}}_z\left[\frac{d\dot{z}}{dt}\right]=-{\hat{e}}_{\varpi }[\frac{1}{\rho }\frac{\partial P}{\partial \varpi }+\frac{\partial \mathrm{\Phi }}{\partial \varpi }]-{\hat{e}}_z[\frac{1}{\rho }\frac{\partial P}{\partial z}+\frac{\partial \mathrm{\Phi }}{\partial z}]$
where, the specific angular momentum, $j(\varpi ,z)\equiv {\varpi }^2\dot{\varphi }=\mathrm{c}\mathrm{o}\mathrm{n}\mathrm{s}\mathrm{t}\mathrm{a}\mathrm{n}\mathrm{t}(\mathrm{i}\mathrm{.}\mathrm{e}\mathrm{.},\mathrm{i}\mathrm{n}\mathrm{d}\mathrm{e}\mathrm{p}\mathrm{e}\mathrm{n}\mathrm{d}\mathrm{e}\mathrm{n}\mathrm{t}\mathrm{o}\mathrm{f}\mathrm{t}\mathrm{i}\mathrm{m}\mathrm{e})$

Adiabatic Form of the
First Law of Thermodynamics

$\frac{dϵ}{dt}+P\frac{d}{dt}\left(\frac{1}{\rho }\right)=0$

Poisson Equation

$\frac{1}{\varpi }\frac{\partial }{\partial \varpi }\left[\varpi \frac{\partial \mathrm{\Phi }}{\partial \varpi }\right]+\frac{{\partial }^2\mathrm{\Phi }}{\partial z^2}=4\pi G\rho .$

This study is closely tied to our separate discussion of the free-fall collapse of uniform-density spheres. For example, by definition, an element of fluid is in "free fall" if its motion in a gravitational field is unimpeded by pressure gradients. The most straightforward way to illustrate how such a system evolves is to set $P=0$ in all of the governing equations. In doing this, the continuity equation and the Poisson equation remain unchanged; the equation formulated by the first law of thermodynamics becomes irrelevant; and the two components of the Euler equation become,

${\widehat{\mathbf{𝐞}}}_{\varpi }$:	$\frac{d\dot{\varpi }}{dt}-\frac{j^2}{{\varpi }^3}$	$=$	$-\frac{\partial \mathrm{\Phi }}{\partial \varpi },$
${\widehat{\mathbf{𝐞}}}_z$:	$\frac{d\dot{z}}{dt}$	$=$	$-\frac{\partial \mathrm{\Phi }}{\partial z}.$

Gravitational Potential

Here, our specific interest is in modeling the free-fall collapse of a uniform-density spheroid. Ignoring, for the moment, the time-dependent nature of this problem, we appreciate from a separate, detailed derivation that the gravitational potential inside (or on the surface) of an homogeneous, triaxial ellipsoid with semi-axes $(x,y,z)=(a_1,a_2,a_3)$ is given, to within an arbitrary additive constant, by the expression,

$\mathrm{\Phi }(\overrightarrow{x})=\pi G\rho [A_1{x}^2+A_2{y}^2+A_3{z}^2],$

where the three, spatially independent coefficients, $A_1,A_2,$ and $A_3$ are functions of the chosen lengths of the three semi-axes. When deriving mathematical expressions for each of the three $A_i$ coefficients, in our accompanying discussion we have found it useful to initially attach a subscript, $(\ell ,m,\mathrm{o}\mathrm{r}s)$ — indicating whether the coefficient is associated with the (largest, medium-length, or smallest) semi-axis — before specifying how, for a given problem, $(\ell ,m,s)$ are appropriately associated with the three $(x,y,z)$ coordinate axes.

Oblate Spheroids

For example, for an oblate-spheroidal mass distribution, by definition the "largest" and "medium-length" semi-axes are equal to one another. Hence, $a_{\ell }=a_m$ and, according to our associated derivation,

$A_{\ell }=A_m$

$=$

$-\frac{(1-e^2)}{e^2}+\frac{(1-e^2{)}^{1/2}}{e^3}\left[{\sin }^{-1}e\right],$

where, $e\equiv (1-a_s^2/a_{\ell }^2{)}^{1/2}$; and,

$A_s$

$=$

$\frac{2}{e^2}-\frac{2(1-e^2{)}^{1/2}}{e^3}\left[{\sin }^{-1}e\right].$

Conventionally, the $z$-axis is aligned with the symmetry (in this case, shortest) axis of the mass distribution, so we set $A_1=A_2=A_{\ell }$, and $A_3=A_s$. Therefore — see also our parallel discussion — we appreciate that, for oblate-spheroidal mass distributions,

$\mathrm{\Phi }(\overrightarrow{x})=\pi G\rho [A_{\ell }{\varpi }^2+A_s{z}^2].$

These same coefficient expressions may also be found in, for example:  Chapter 3, Eq. (36) of [EFE]; §4.5, Eqs. (48) & (49) of [T78]; and the first column of Table 2-1 (p. 57) of [BT87]. Note that, as we have pointed out in a separate discussion of Maclaurin Spheroids, the expressions for $A_{\ell }$ and $A_s$ have the following values in the limit of a sphere $(e=0)$ or in the limit of an infinitesimally thin disk $(e=1)$:

Table 1:  Limiting Values (for oblate spheroids)
	$e\to 0$	$\frac{a_3}{a_1}\to 0$
$\frac{{\sin }^{-1}e}{e}$	$1+\frac{e^2}{6}+{𝒪}\left(e^4\right)$	$\frac{\pi }{2}-\left(\frac{a_3}{a_1}\right)+\frac{\pi }{4}(\frac{a_3}{a_1}{)}^2-{𝒪}(\frac{a_3^3}{a_1^3})$
$A_{\ell }$	$\frac{2}{3}[1-\frac{e^2}{5}-{𝒪}(e^4\left)\right]$	$\frac{\pi }{2}\left(\frac{a_3}{a_1}\right)-2(\frac{a_3}{a_1}{)}^2+{𝒪}(\frac{a_3^3}{a_1^3})$
$A_s$	$\frac{2}{3}[1+\frac{2e^2}{5}+{𝒪}(e^4\left)\right]$	$2-\pi \left(\frac{a_3}{a_1}\right)+4(\frac{a_3}{a_1}{)}^2-{𝒪}(\frac{a_3^3}{a_1^3})$

Prolate Spheroids

For a prolate-spheroidal mass distribution, by definition the "smallest" and "medium-length" semi-axes are equal to one another. Hence, $a_m=a_s$ and, according to our associated derivation,

$A_{\ell }$

$=$

$\frac{(1-e^2)}{e^3}\cdot \ln \left[\frac{1+e}{1-e}\right]-\frac{2(1-e^2)}{e^2},$

where, as above, $e\equiv (1-a_s^2/a_{\ell }^2{)}^{1/2}$; and,

$A_s=A_m$

$=$

$\frac{1}{e^2}-\frac{(1-e^2)}{2e^3}\cdot \ln \left(\frac{1+e}{1-e}\right).$

If the symmetry (in this case, longest) axis of this prolate-spheroidal mass distribution is aligned with the $z$-axis of the coordinate system, then we should set $A_1=A_2=A_s$, and $A_3=A_{\ell }$. This means that the expression for the gravitational potential is,

$\mathrm{\Phi }(\overrightarrow{x})=\pi G\rho [A_s{\varpi }^2+A_{\ell }{z}^2].$

These same coefficient expressions may also be found in, for example, the second column of Table 2-1 (p. 57) of [BT87].

NOTE:  If, following [EFE], we instead align the longest (and, in this case, symmetry) axis of the prolate mass distribution with the $x$-axis, then $A_2=A_3=A_s$ and $A_1=A_{\ell }$. This matches the coefficient expressions presented in our parallel discussion of the potential inside and on the surface of a prolate-spheroidal mass distribution.

Consider a Time-Varying Eccentricity

If the eccentricity of an homogeneous spheroid varies with time — that is, if $e\to e(t)$ — while it remains homogeneous, the result will be a potential of the form,

$\mathrm{\Phi }(\overrightarrow{x},t)=A(t){\varpi }^2+C(t)z^2,$

whether the spheroid is oblate or the spheroid is prolate.

Lynden-Bell's (1962) Insight

Let's examine the analysis by 📚 D. Lynden-Bell (1962, Math. Proc. Cambridge Phil. Soc., Vol. 58, Issue 4, pp. 709 - 711) — hereafter, LB62 — of the "Gravitational Collapse of a Cold Rotating Gas Cloud."

Motion of a Single Particle

Consider a particle that, at time $t=0$, is at position $({\varpi }_0,{\varphi }_0,z_0)$ and is moving about the $z$-axis with velocity, ${\varpi }_0\mathrm{\Omega }\Rightarrow j_0={\varpi }_0^2\mathrm{\Omega }$. Consider furthermore that its acceleration is subject to the force arising from an axisymmetric gravitational potential of the form,

$\mathrm{\Phi }(\overrightarrow{x},t)=A(t){\varpi }^2+C(t)z^2.$

[This is the gravitational potential adopted by LB62 — see his equation (1) — except he adopted a different sign convention to ours. He would therefore have also attached a sign to the gradient of the potential that is the opposite of the sign that appears on the right-hand side of our Euler equation expression.] In this case, the two components of the Euler equation that govern the particle's motion are,

${\widehat{\mathbf{𝐞}}}_{\varpi }$:	$\ddot{\varpi }-\frac{j_0^2}{{\varpi }^3}$	$=$	$-2A\varpi ,$
${\widehat{\mathbf{𝐞}}}_z$:	$\ddot{z}$	$=$	$-2Cz,$

where we have adopted the familiar shorthand notation, $d\dot{\varpi }/dt\to \ddot{\varpi }$ and $d\dot{z}/dt\to \ddot{z}$. If we divide the first of these relations by ${\varpi }_0$ and the second by $z_0$, then adopt the dimensionless variables, $R\equiv \varpi /{\varpi }_0$ and $Z\equiv z/z_0$, we can write,

$\frac{\ddot{\varpi }}{{\varpi }_0}-\frac{j_0^2}{{\varpi }_0{\varpi }^3}$	$=$	$-2A\left(\frac{\varpi }{{\varpi }_0}\right)$
$\Rightarrow \ddot{R}-\frac{{\mathrm{\Omega }}^2}{R^3}$	$=$	$-2AR,$
LB62, p. 710, Eq. (10)

and,

$\frac{\ddot{z}}{z_0}$	$=$	$-2C\left(\frac{z}{z_0}\right)$
$\Rightarrow \ddot{Z}$	$=$	$-2CZ.$
LB62, p. 710, Eq. (11)

Finally, we use $\phi (t)$ to represent the particle's time-varying angular-coordinate position relative to its initial position — that is, we adopt the definition, $\phi (t)\equiv \varphi (t)-{\varphi }_0$. Then, conservation of angular momentum implies that, at any moment, the particle's rotation frequency about the symmetry axis will be,

$\dot{\phi }=\dot{\varphi }$	$=$	$\frac{\mathrm{\Omega }}{R^2}.$
LB62, p. 710, Eq. (9)

This governing set of evolutionary equations has been set up such that at time, $t=0$: $R=1$, $Z=1$, $\phi =0$, $\dot{R}=0$, $\dot{Z}=0$, and $\dot{\phi }=\mathrm{\Omega }$. With this set of initial conditions in hand, along with an appropriate specification of the two time-dependent coefficients, $A(t)$ and $C(t)$, the set of governing relations can be integrated (numerically) to give $R(t),\phi (t)$, and $Z(t)$. This is the result that 📚 Lynden-Bell (1962) established for the motion of one particle.

Evolution of the Spheroid

Following his examination of the motion of an individual particle, LB62 recognized that, "$R(t),\phi (t)$, and $Z(t)$ are all independent of ${\varpi }_0,{\varphi }_0$, and $z_0$ because [none of the three evolutionary equations] nor the above initial conditions mention them." You only need to integrate the coupled set of governing relations once then — assuming that the functions, $A(t)$ and $C(t)$, are the same in all cases — the time-dependent coordinates of any particle are given by $({\varpi }_{0}R,{\varphi }_0+\phi ,z_{0}Z)$, where $({\varpi }_0,{\varphi }_0,z_0)$ are the initial coordinates of that particle. "Thus the result of the motion is merely a change of scales."

Consider then, as did LB62, the evolution of a spheroid that is initially uniformly filled with free particles and whose only motion, initially, is uniform rotation, $\mathrm{\Omega }$, about the z-axis. As LB62 puts it, since the motion of each particle can be described merely via a change of scales: "… the distribution of the particles remains uniform, and the boundary remains spheroidal"; and, while the angular frequency of each particle, $\dot{\varphi }$, varies with time, "… since $\dot{\varphi }=\dot{\phi }$ the rotation remains uniform in space."

It should be clear, as well, that the eccentricity of the evolving spheroid will vary with time. Specifically in the case of an oblate spheroid, the time-dependent semi-axes are $(a_{\ell }R(t),a_{\ell }R(t),a_{s}Z(t))$; hence,

$e(t)$

$=$

$[1-(\frac{a_{s}Z}{a_{\ell }R}{)}^2{]}^{1/2}=[1-(1-e_0^2)(\frac{Z}{R}{)}^2{]}^{1/2},$

where, $e_0=(1-a_s^2/a_{\ell }^2{)}^{1/2}$ is the eccentricity of the spheroid initially, and the time-variation enters via the pair of functions, $Z(t)$ and $R(t)$. In the case of an prolate spheroid, the time-dependent semi-axes are $(a_{s}R(t),a_{s}R(t),a_{\ell }Z(t))$; hence,

$e(t)$

$=$

$[1-(\frac{a_{s}R}{a_{\ell }Z}{)}^2{]}^{1/2}=[1-(1-e_0^2)(\frac{R}{Z}{)}^2{]}^{1/2}.$

In turn, the time-dependent behavior of the coefficients in the expression for the gravitational potential, $A(t)$ and $C(t)$, is drawn from $e(t)$ as detailed in Table 1, immediately below.

Table 1:   Time-Dependent Coefficients of the Gravitational Potential $\mathrm{\Phi }(\overrightarrow{x},t)=A(t){\varpi }^2+C(t)z^2,$ where, it is understood that the eccentricity of the spheroid, $e(t)$, varies with time.
Oblate Spheroid	Prolate Spheroid
$\frac{A(t)}{\pi G\rho (t)}$ $=$ $-\frac{(1-e^2)}{e^2}+\frac{(1-e^2{)}^{1/2}}{e^3}\left[{\sin }^{-1}e\right],$ $\frac{C(t)}{\pi G\rho (t)}$ $=$ $\frac{2}{e^2}-\frac{2(1-e^2{)}^{1/2}}{e^3}\left[{\sin }^{-1}e\right],$ where, $e(t)$ $=$ $\{1-(1-e_0^2)[\frac{Z(t)}{R(t)}{]}^2{\}}^{1/2},$ $\frac{\rho (t)}{{\rho }_0}$ $=$ $\frac{1}{R^2(t)Z(t)}.$	$\frac{A(t)}{\pi G\rho (t)}$ $=$ $\frac{1}{e^2}-\frac{(1-e^2)}{2e^3}\cdot \ln \left(\frac{1+e}{1-e}\right),$ $\frac{C(t)}{\pi G\rho (t)}$ $=$ $\frac{(1-e^2)}{e^3}\cdot \ln \left[\frac{1+e}{1-e}\right]-\frac{2(1-e^2)}{e^2},$ where, $e(t)$ $=$ $\{1-(1-e_0^2)[\frac{R(t)}{Z(t)}{]}^2{\}}^{1/2},$ $\frac{\rho (t)}{{\rho }_0}$ $=$ $\frac{1}{R^2(t)Z(t)}.$

Our Numerical Integration

Case of the Oblate Spheroid

Governing, Dimensionless Differential Equations

We are primarily interested in determining, for various initial values of the parameter pair $(e_0,\mathrm{\Omega })$, how rapidly an oblate spheroid collapses to an infinitesimally thin disk $(e\to 1)$, and what the radius of this disk is at the instant it forms. As has been detailed above, according to LB62 the time-dependent behavior of $R$ and $Z$ — and, hence, of $e$ — can be obtained by integrating the following pair of governing relations:

$\ddot{R}-\frac{{\mathrm{\Omega }}^2}{R^3}$	$=$	$-2AR,$
$\ddot{Z}$	$=$	$-2CZ.$

Turning both of these 2^nd-order ODEs into a pair of 1^st-order ODEs, then multiplying through by ${\tau }_{\mathrm{f}\mathrm{f}}\equiv [3\pi /(32G{\rho }_0){]}^{1/2}$ or, as appropriate, by ${\tau }_{\mathrm{f}\mathrm{f}}^2$, we have,

${\tau }_{\mathrm{f}\mathrm{f}}\cdot \frac{d({\tau }_{\mathrm{f}\mathrm{f}}\dot{R})}{dt}$	$=$	$\frac{3\pi }{32G{\rho }_0}[\frac{{\mathrm{\Omega }}^2}{R^3}-2AR]$
	$=$	$\frac{3{\pi }^2}{16}\{\frac{{\mathrm{\Omega }}^2/(2\pi G{\rho }_0)}{R^3}-[\frac{\rho (t)}{{\rho }_0}]\frac{A(t)}{\pi G\rho (t)}\cdot R\},$
${\tau }_{\mathrm{f}\mathrm{f}}\cdot \frac{dR}{dt}$	$=$	${\tau }_{\mathrm{f}\mathrm{f}}\dot{R},$
${\tau }_{\mathrm{f}\mathrm{f}}\cdot \frac{d({\tau }_{\mathrm{f}\mathrm{f}}\dot{Z})}{dt}$	$=$	$-\frac{3\pi }{32G{\rho }_0}\left[2CZ\right]$
	$=$	$-\frac{3{\pi }^2}{16}\left\{\right[\frac{\rho (t)}{{\rho }_0}]\frac{C(t)}{\pi G\rho (t)}\cdot Z\},$
${\tau }_{\mathrm{f}\mathrm{f}}\cdot \frac{dZ}{dt}$	$=$	${\tau }_{\mathrm{f}\mathrm{f}}\dot{Z}.$

Adopting the shorthand notation,

${\sigma }_0^2$	$\equiv$	$\frac{{\mathrm{\Omega }}^2}{2\pi G{\rho }_0},$
$\tau$	$\equiv$	$\frac{t}{{\tau }_{\mathrm{f}\mathrm{f}}},$

and pulling expressions for the oblate spheroid from Table 1 above, we have,

$\frac{d({\tau }_{\mathrm{f}\mathrm{f}}\dot{R})}{d\tau }$	$=$	$\frac{3{\pi }^2}{16}\{\frac{{\sigma }_0^2}{R^3}-\frac{1}{RZ}[-\frac{(1-e^2)}{e^2}+\frac{(1-e^2{)}^{1/2}}{e^3}\cdot {\sin }^{-1}e\left]\right\},$
$\frac{dR}{d\tau }$	$=$	${\tau }_{\mathrm{f}\mathrm{f}}\dot{R},$
$\frac{d({\tau }_{\mathrm{f}\mathrm{f}}\dot{Z})}{d\tau }$	$=$	$-\frac{3{\pi }^2}{16}\left\{\frac{1}{R^2}\right[\frac{2}{e^2}-\frac{2(1-e^2{)}^{1/2}}{e^3}\cdot {\sin }^{-1}e\left]\right\},$
$\frac{dZ}{d\tau }$	$=$	${\tau }_{\mathrm{f}\mathrm{f}}\dot{Z}.$

Discrete Representations

Along a uniformly segmented, discrete time grid — time step, $\mathrm{\Delta }=({\tau }_{n+1}-{\tau }_n)$ — let's define the following discretized variables:

• At integer values of $\mathrm{\Delta }$:   $R$ and $Z$
• At half-integer values of $\mathrm{\Delta }$:   $F\equiv ({\tau }_{\mathrm{f}\mathrm{f}}\dot{R});K\equiv ({\tau }_{\mathrm{f}\mathrm{f}}\dot{Z})$

The four finite-difference relations are, then:

$\frac{F_{n+1/2}-F_{n-1/2}}{\mathrm{\Delta }}$	$=$	$\frac{3{\pi }^2}{16}\{\frac{{\sigma }_0^2}{R^3}-\frac{1}{RZ}[-\frac{(1-e^2)}{e^2}+\frac{(1-e^2{)}^{1/2}}{e^3}\cdot {\sin }^{-1}e]{\}}_n,$
$\frac{R_{n+1}-R_n}{\mathrm{\Delta }}$	$=$	$F_{n+1/2},$
$\frac{K_{n+1/2}-K_{n-1/2}}{\mathrm{\Delta }}$	$=$	$-\frac{3{\pi }^2}{16}\left\{\frac{1}{R^2}\right[\frac{2}{e^2}-\frac{2(1-e^2{)}^{1/2}}{e^3}\cdot {\sin }^{-1}e]{\}}_n,$
$\frac{Z_{n+1}-Z_n}{\mathrm{\Delta }}$	$=$	$K_{n+1/2}.$

STEP #0:

• Choose time-independent values of the parameters, $(c_0/a_0)$ and ${\beta }_0$, which means that, ${\sigma }_0^2=3{\pi }^2{\beta }_0/8$ and $e_0=[1-(c_0/a_0{)}^2{]}^{1/2}$. Also set the (uniform) integration time step, $\mathrm{\Delta }$; for 201 time steps, for example, set $\mathrm{\Delta }\approx 0.005$.
  The example values shown below (inside parentheses) assume that $(c_0/a_0,{\beta }_0)=(0.90,0.0)\Rightarrow (e_0,{\sigma }_0^2)=(0.43588989,0)$, which corresponds to one of the model-evolutions presented by 📚 Lin, Mestel, & Shu (1965).

STEP #1: Initially, that is, at integration time step, $n=0$

• Set $R_0=1$, $Z_0=1$, and $e=e_0$; this means that the RHS of the first and third discrete evolution equations is known.
• Given that the configuration is collapsing from rest, we want to set ${\dot{R}}_0$ and ${\dot{Z}}_0$ to zero. This is accomplished by establishing reflection symmetry through the (time) origin, that is, by setting $F_{+1/2}=-F_{-1/2}$ and $K_{+1/2}=-K_{-1/2}$. Initially, then, the LHS of the first and third discrete equations are, respectively, $2F_{+1/2}/\mathrm{\Delta }$ and $2K_{+1/2}/\mathrm{\Delta }$.
• Use the first and third discrete equations to solve for the "R" and "Z" velocities at time step $n=0$, namely,

$F_{+1/2}$	$=$	$\frac{\mathrm{\Delta }}{2}\cdot \frac{3{\pi }^2}{16}\{{\sigma }_0^2-[-\frac{(1-e_0^2)}{e_0^2}+\frac{(1-e_0^2{)}^{1/2}}{e_0^3}\cdot {\sin }^{-1}{e}_0\left]\right\},$	$(-2.952452\times 10^{-3})$
$K_{+1/2}$	$=$	$-\frac{\mathrm{\Delta }}{2}\cdot \frac{3{\pi }^2}{16}\left\{\right[\frac{2}{e_0^2}-\frac{2(1-e_0^2{)}^{1/2}}{e_0^3}\cdot {\sin }^{-1}{e}_0\left]\right\}.$	$(-3.347849\times 10^{-3})$

• From the second and fourth discrete relations, determine the advanced coordinate positions.

$R_{+1}$	$=$	$1+\mathrm{\Delta }\cdot F_{+1/2},$	$(1.0-1.476226\times 10^{-5})$
$Z_{+1}$	$=$	$1+\mathrm{\Delta }\cdot K_{+1/2}.$	$(1.0-1.673924\times 10^{-5})$

• Determine the eccentricity at this advanced time.

$e_{+1}$

$=$

$[1-(1-e_0^2)(\frac{Z_{+1}}{R_{+1}}{)}^2{]}^{1/2}.$

$(1-0.809996797{)}^{1/2}=(0.435893568)$

STEP #2:   Repeat, in sequence for all values of $n>1$ until $Z$ passes through zero.

• Set $n\to n+1$.
• Given knowledge of the various variable values at time-step, $(n-1/2)$ and $n$, use the first and third discrete evolution relations to determine $F_{n+1/2}$ and $K_{n+1/2}$; specifically,

$F_{n+1/2}$	$=$	$F_{n-1/2}+\mathrm{\Delta }\cdot \frac{3{\pi }^2}{16}\{\frac{{\sigma }_0^2}{R^3}-\frac{1}{RZ}[-\frac{(1-e^2)}{e^2}+\frac{(1-e^2{)}^{1/2}}{e^3}\cdot {\sin }^{-1}e]{\}}_n,$	$(F_{+1/2}-5.905085\times 10^{-3})$
$K_{n+1/2}$	$=$	$K_{n-1/2}-\mathrm{\Delta }\cdot \frac{3{\pi }^2}{16}\left\{\frac{1}{R^2}\right[\frac{2}{e^2}-\frac{2(1-e^2{)}^{1/2}}{e^3}\cdot {\sin }^{-1}e]{\}}_n.$	$(K_{+1/2}-6.695907\times 10^{-3})$

• Similarly, from the second and fourth discrete relations, determine the advanced coordinate positions.

$R_{n+1}$	$=$	$R_n+\mathrm{\Delta }\cdot F_{n+1/2},$	$(0.99994095)$
$Z_{n+1}$	$=$	$Z_n+\mathrm{\Delta }\cdot K_{n+1/2}.$	$(0.999933042)$

• Determine the eccentricity at this advanced time.

$e_{n+1}$

$=$

$[1-(1-e_0^2)(\frac{Z_{n+1}}{R_{n+1}}{)}^2{]}^{1/2}.$

$(1-0.809987188{)}^{1/2}=(0.43590459)$

Results

Table 2: Collapse of a Nonrotating, Pressure-Free Oblate Spheroid
$\frac{c_0}{a_0}$	$e_0$	${\tau }_c$	$R_c$	${\dot{R}}_c$	${\dot{Z}}_c$
$0.99$	$0.141$	$$	$$	$$	$$
$0.95$	$0.312$	$$	$$	$$	$$
$0.90$	$0.436$	$$	$$	$$	$$
$0.80$	$0.600$	$$	$$	$$	$$
$0.70$	$0.714$	$$	$$	$$	$$
$0.60$	$0.800$	$$	$$	$$	$$
$0.10$	$0.995$	$$	$$	$$	$$

Key References

• 📚 D. Lynden-Bell (1962, Math. Proc. Cambridge Phil. Soc., Vol. 58, Issue 4, pp. 709 - 711): On the gravitational collapse of a cold rotating gas cloud
  NOTE … according to the new ADS listing, the authors associated with this paper should be, D. Lynden-Bell & C. T. C. Wall (Charles Terence Clegg "Terry" Wall); however, the archived article, itself, lists Lynden-Bell as the sole author while indicating that the paper was simply being communicated by Wall.

Comment by J. E. Tohline: In §II of this "1964" article, Lynden-Bell references his 1962 article with an incorrect year (Lynden-Bell 1963); within his list of REFERENCES, the year (1962) is correct, but the journal volume is incorrectly identified as 50 (it should be vol. 58).

 

• D. Lynden-Bell (1964), ApJ, 139, 1195 - 1216: On Large-Scale Instabilities during Gravitational Collapse and the Evolution of Shrinking Maclaurin Spheroids
• Classic paper by C. C. Lin, Leon Mestel, and Frank Shu (1965, ApJ, 142, 1431 - 1446) titled, "The Gravitational Collapse of a Uniform Spheroid."

See Also