Difference between revisions of "User:Tohline/SSC/Stability/BiPolytropes"

From VistrailsWiki
Jump to navigation Jump to search
Line 1,772: Line 1,772:




<table border="0" align="center" cellpadding="8"><tr><td align="center">
<table border="0" align="center" cellpadding="8">
[[File:Mode0EigenfunctionsCombinedLabel.png|750px|Eigenfunctions for Marginally Unstable Models]]
<tr>
</td></tr></table>
  <td align="center">
[[File:Mode0EigenfunctionsCombinedLabel.png|600px|Eigenfunctions for Marginally Unstable Models]]
  </td>
  <td align="center">
[[File:CompositeEigenfunctions.png|400px|Eigenfunctioins for multiple modes in three separate models]]
  </td>
</tr>
</table>


=Virial Analysis=
=Virial Analysis=

Revision as of 02:08, 17 December 2018

Marginally Unstable Bipolytropes

Our aim is to determine whether or not there is a relationship between (1) equilibrium models at turning points along bipolytrope sequences and (2) bipolytropic models that are marginally (dynamically) unstable toward collapse (or dynamical expansion).


Whitworth's (1981) Isothermal Free-Energy Surface
|   Tiled Menu   |   Tables of Content   |  Banner Video   |  Tohline Home Page   |

Overview

file = Dropbox/WorkFolder/Wiki edits/EmbeddedPolytropes/CombinedSequences.xlsx --- worksheet = EqSeqCombined2
Figure 1:  Equilibrium Sequences
of Pressure-Truncated Polytropes

Equilibrium sequences of Pressure-Truncated Polytropes

We expect the content of this chapter — which examines the relative stability of bipolytropes — to parallel in many ways the content of an accompanying chapter in which we have successfully analyzed the relative stability of pressure-truncated polytopes. Figure 1, shown here on the right, has been copied from a closely related discussion. The curves show the mass-radius relationship for pressure-truncated model sequences having a variety of polytropic indexes, as labeled, over the range <math>1 \le n \le 6</math>. (Another version of this figure includes the isothermal sequence.) On each sequence for which <math>~n \ge 3</math>, the green filled circle identifies the model with the largest mass. We have shown analytically that the oscillation frequency of the fundamental-mode of radial oscillation is precisely zero for each one of these maximum-mass models. As a consequence, we know that each green circular marker identifies the point along its associated sequence that separates dynamically stable (larger radii) from dynamically unstable (smaller radii) models.


In each case, the fundamental-mode oscillation frequency is precisely zero if, and only if, the adiabatic index governing expansions/contractions is related to the underlying structural polytropic index via the relation, <math>~\gamma_g = (n + 1)/n</math>, and if a constant surface-pressure boundary condition is imposed.



In another accompanying chapter, we have used purely analytic techniques to construct equilibrium sequences of spherically symmetric bipolytropes that have, <math>~(n_c,n_e) = (5,1)</math>. For a given choice of <math>~\mu_e/\mu_c</math> — the ratio of the mean-molecular weight of envelope material to the mean-molecular weight of material in the core — a physically relevant sequence of models can be constructed by steadily increasing the value of the dimensionless radius at the core/envelope interface, <math>~\xi_i</math>, from zero to infinity. Figure 2, which has been copied from this separate chapter, shows how the fractional core mass, <math>\nu \equiv M_\mathrm{core}/M_\mathrm{tot}</math>, varies with the fractional core radius, <math>q \equiv r_\mathrm{core}/R</math>, along sequences having six different values of <math>~\mu_e/\mu_c</math>: 1 (blue diamonds), ½ (red squares), 0.345 (dark purple crosses), ⅓ (pink triangles), 0.309 (light green dashes), and ¼ (purple asterisks). Along each of the model sequences, points marked by solid-colored circles correspond to models whose interface parameter, <math>~\xi_i</math>, has one of three values: 0.5 (green circles), 1 (dark blue circles), or 3 (orange circles).

When modeling bipolytropes, the default expectation is that an increase in <math>\xi_i</math> along a given sequence will correspond to an increase in the relative size — both the radius and the mass — of the core. As Figure 2 illustrates, this expectation is realized along the sequences marked by blue diamonds (<math>~\mu_e/\mu_c = 1</math>) and by red squares (<math>~\mu_e/\mu_c = </math>½). But the behavior is different along the other four illustrated sequences. For sufficiently large <math>~\xi_i</math>, the relative radius of the core begins to decrease. Furthermore, along sequences for which <math>~\mu_e/\mu_c < \tfrac{1}{3}</math>, eventually the fractional mass of the core reaches a maximum and, thereafter, decreases even as the value of <math>~\xi_i</math> continues to increase. (Additional properties of these equilibrium sequences are discussed in yet another accompanying chapter.)

The principal question is: Along bipolytropic sequences, are maximum-mass models associated with the onset of dynamical instabilities?

Planned Approach

Figure 2: Equilibrium Sequences of Bipolytropes with <math>~(n_c,n_e) = (5,1)</math>
PlotSequencesBest02.png

Ideally we would like to answer the just-stated "principal question" using purely analytic techniques. But, to date, we have been unable to fully address the relevant issues analytically, even in what would be expected to be the simplest case:   bipolytropic models that have <math>~(n_c,n_e) = (0, 0)</math>. Instead, we will streamline the investigation a bit and proceed — at least initially — using a blend of techniques. We will investigate the relative stability of bipolytropic models having <math>~(n_c,n_e) = (5,1) </math> whose equilibrium structures are completely defined analytically; then the eigenvectors describing radial modes of oscillation will be determined, one at a time, by solving the relevant LAWE(s) numerically. We are optimistic that this can be successfully accomplished because we have had experience numerically integrating the LAWE that governs the oscillation of:

A key reference throughout this investigation will be the paper by J. O. Murphy & R. Fiedler (1985b, Proc. Astr. Soc. of Australia, 6, 222). They studied Radial Pulsations and Vibrational Stability of a Sequence of Two Zone Polytropic Stellar Models. Specifically, their underlying equilibrium models were bipolytropes that have <math>~(n_c,n_e) = (1, 5)</math>. In an accompanying chapter, we describe in detail how Murphy & Fiedler obtained these equilibrium bipolytropic structures and detail some of their equilibrium properties.

Here are the steps we initially plan to take:

  • Governing LAWEs:
    • Identify the relevant LAWEs that govern the behavior of radial oscillations in the <math>~n_c = 5</math> core and, separately, in the <math>~n_e = 1</math> envelope. Check these LAWE specifications against the published work of Murphy & Fiedler (1985b).
    • Determine the matching conditions that must be satisfied across the core/envelope interface. Be sure to take into account the critical interface jump conditions spelled out by P. Ledoux & Th. Walraven (1958), as we have already discussed in the context of an analysis of radial oscillations in zero-zero bipolytropes.
  • Determine what surface boundary condition should be imposed on physically relevant LAWE solutions, i.e., on the physically relevant radial-oscillation eigenvectors.
  • Initial Analysis:
    • Choose a maximum-mass model along the bipolytropic sequence that has, for example, <math>~\mu_e/\mu_c = 1/4</math>. Hopefully, we will be able to identify precisely (analytically) where this maximum-mass model lies along the sequence. Yes! Our earlier analysis does provide an analytic prescription of the model that sits at the maximum-mass location along the chosen sequence.
    • Solve the relevant eigenvalue problem for this specific model, initially for <math>~(\gamma_c, \gamma_e) = (6/5, 2)</math> and initially for the fundamental mode of oscillation.

Review of the Analysis by Murphy & Fiedler (1985b)

As we have detailed separately, the boundary condition at the center of a polytropic configuration is,

<math>~\frac{dx}{d\xi} \biggr|_{\xi=0} = 0 \, ;</math>

and the boundary condition at the surface of an isolated polytropic configuration is,

<math>~\frac{d\ln x}{d\ln\xi}</math>

<math>~=</math>

<math>~- \alpha + \frac{\omega^2}{\gamma_g } \biggl( \frac{1}{4\pi G \rho_c } \biggr) \frac{\xi}{(-\theta^')} </math>         at         <math>~\xi = \xi_s \, .</math>

But this surface condition is not applicable to bipolytropes. Instead, let's return to the original, more general expression of the surface boundary condition:

<math>~ \frac{d\ln x}{d\ln\xi}\biggr|_s</math>

<math>~=</math>

<math>~- \alpha + \frac{\omega^2 R^3}{\gamma_g GM_\mathrm{tot}} \, .</math>


Utilizing an accompanying discussion, let's examine the frequency normalization used by Murphy & Fiedler (1985b) (see the top of the left-hand column on p. 223):

<math>~\Omega^2</math>

<math>~\equiv</math>

<math>~ \omega^2 \biggl[ \frac{R^3}{GM_\mathrm{tot}} \biggr] </math>

 

<math>~=</math>

<math>~ \omega^2 \biggl[ \frac{3}{4\pi G \bar\rho} \biggr] = \omega^2 \biggl[ \frac{3}{4\pi G \rho_c} \biggr] \frac{\rho_c}{\bar\rho} = \frac{3\omega^2}{(n_c+1)} \biggl[ \frac{(n_c+1)}{4\pi G \rho_c} \biggr] \frac{\rho_c}{\bar\rho} </math>

 

<math>~=</math>

<math>~ \frac{3\omega^2}{(n_c+1)} \biggl[ \frac{a_n^2\rho_c}{P_c} \cdot \theta_c \biggr] \frac{\rho_c}{\bar\rho} = \frac{3\gamma}{(n_c+1)} \frac{\rho_c}{\bar\rho} \biggl[ \frac{a_n^2\rho_c}{P_c} \cdot \frac{\omega^2 \theta_c}{\gamma} \biggr] \, . </math>

For a given radial quantum number, <math>~k</math>, the factor inside the square brackets in this last expression is what Murphy & Fiedler (1985b) refer to as <math>~\omega^2_k \theta_c</math>. Keep in mind, as well, that, in the notation we are using,

<math>~\sigma_c^2</math>

<math>~\equiv</math>

<math>~\frac{3\omega^2}{2\pi G \rho_c}</math>

<math>~\Rightarrow ~~~ \sigma_c^2</math>

<math>~=</math>

<math>~ \biggl( \frac{2\bar\rho}{\rho_c}\biggr) \Omega^2 = \frac{6\gamma}{(n_c+1)} \biggl[ \frac{a_n^2\rho_c}{P_c} \cdot \frac{\omega^2 \theta_c}{\gamma} \biggr] = \frac{6\gamma}{(n_c+1)} \biggl[ \omega_k^2 \theta_c \biggr] \, . </math>

This also means that the surface boundary condition may be rewritten as,

<math>~ \frac{d\ln x}{d\ln\xi}\biggr|_s</math>

<math>~=</math>

<math>~\frac{\Omega^2}{\gamma_g } - \alpha \, .</math>


Let's apply these relations to the core and envelope, separately.

Interface Conditions

Here, we will simply copy the discussion already provided in the context of our attempt to analyze the stability of <math>~(n_c, n_e) = (0, 0)</math> bipolytropes; specifically, we will draw from STEP 4: in the Piecing Together subsection. Following the discussion in §§57 & 58 of P. Ledoux & Th. Walraven (1958), the proper treatment is to ensure that fractional perturbation in the gas pressure (see their equation 57.31),

<math>~\frac{\delta P}{P}</math>

<math>~=</math>

<math>~- \gamma x \biggl( 3 + \frac{d\ln x}{d\ln \xi} \biggr) \, ,</math>

is continuous across the interface. That is to say, at the interface <math>~(\xi = \xi_i)</math>, we need to enforce the relation,

<math>~0</math>

<math>~=</math>

<math>~\biggl[ \gamma_c x_\mathrm{core} \biggl( 3 + \frac{d\ln x_\mathrm{core}}{d\ln \xi} \biggr) - \gamma_e x_\mathrm{env} \biggl( 3 + \frac{d\ln x_\mathrm{env}}{d\ln \xi} \biggr)\biggr]_{\xi=\xi_i}</math>

 

<math>~=</math>

<math>~\gamma_e \biggl[ \frac{\gamma_c}{\gamma_e} \biggl( 3 + \frac{d\ln x_\mathrm{core}}{d\ln \xi} \biggr) - \biggl( 3 + \frac{d\ln x_\mathrm{env}}{d\ln \xi} \biggr)\biggr]_{\xi=\xi_i}</math>

<math>~\Rightarrow~~~ \frac{d\ln x_\mathrm{env}}{d\ln \xi} \biggr|_{\xi=\xi_i}</math>

<math>~=</math>

<math>~3\biggl(\frac{\gamma_c}{\gamma_e} -1\biggr) + \frac{\gamma_c}{\gamma_e} \biggl( \frac{d\ln x_\mathrm{core}}{d\ln \xi} \biggr)_{\xi=\xi_i} \, .</math>

In the context of this interface-matching constraint (see their equation 62.1), P. Ledoux & Th. Walraven (1958) state the following:   In the static (i.e., unperturbed equilibrium) modeldiscontinuities in <math>~\rho</math> or in <math>~\gamma</math> might occur at some [radius]. In the first case — that is, a discontinuity only in density, while <math>~\gamma_e = \gamma_c</math> — the interface conditions imply the continuity of <math>~\tfrac{1}{x} \cdot \tfrac{dx}{d\xi}</math> at that [radius]. In the second case — that is, a discontinuity in the adiabatic exponent — the dynamical condition may be written as above. This implies a discontinuity of the first derivative at any discontinuity of <math>~\gamma</math>.

The algorithm that Murphy & Fiedler (1985b) used to "… [integrate] through each zone …" was designed "… with continuity in <math>~x</math> and <math>~dx/d\xi</math> being imposed at the interface …" Given that they set <math>~\gamma_c = \gamma_e = 5/3</math>, their interface matching condition is consistent with the one prescribed by P. Ledoux & Th. Walraven (1958).

Radial Oscillations of (nc, ne) = (5, 1) Models

Foundation

In an accompanying discussion, we derived the so-called,

Adiabatic Wave (or Radial Pulsation) Equation

LSU Key.png

<math>~ \frac{d^2x}{dr_0^2} + \biggl[\frac{4}{r_0} - \biggl(\frac{g_0 \rho_0}{P_0}\biggr) \biggr] \frac{dx}{dr_0} + \biggl(\frac{\rho_0}{\gamma_\mathrm{g} P_0} \biggr)\biggl[\omega^2 + (4 - 3\gamma_\mathrm{g})\frac{g_0}{r_0} \biggr] x = 0 </math>

whose solution gives eigenfunctions that describe various radial modes of oscillation in spherically symmetric, self-gravitating fluid configurations. Assuming that the underlying equilibrium structure is that of a bipolytrope having <math>~(n_c, n_e) = (5, 1)</math>, it makes sense to adopt the normalizations used when defining the equilibrium structure, namely,

<math>~\rho^*</math>

<math>~\equiv</math>

<math>~\frac{\rho_0}{\rho_c}</math>

;    

<math>~r^*</math>

<math>~\equiv</math>

<math>~\frac{r_0}{[K_c^{1/2}/(G^{1/2}\rho_c^{2/5})]}</math>

<math>~P^*</math>

<math>~\equiv</math>

<math>~\frac{P_0}{K_c\rho_c^{6/5}}</math>

;    

<math>~M_r^*</math>

<math>~\equiv</math>

<math>~\frac{M_r}{[K_c^{3/2}/(G^{3/2}\rho_c^{1/5})]}</math>

We note as well that,

<math>~g_0</math>

<math>~=</math>

<math>~\frac{GM(r_0)}{r_0^2}</math>

 

<math>~=</math>

<math>~ G \biggl[ M_r^* \rho_c^{-1 / 5} \biggl( \frac{K_c}{G}\biggr)^{3 / 2} \biggr] \biggl[ r^* \rho_c^{-2 / 5}\biggl( \frac{K_c}{G}\biggr)^{1 / 2} \biggr]^{-2} </math>

 

<math>~=</math>

<math>~ \frac{G M_r^*}{(r^*)^2} \biggl[ \rho_c^{3 / 5} \biggl( \frac{K_c}{G}\biggr)^{1 / 2} \biggr] \, . </math>

Hence, multiplying the LAWE through by <math>~(K_c/G)\rho_c^{-4 / 5}</math> gives,

<math>~0</math>

<math>~=</math>

<math>~ \frac{d^2x}{dr_0^2} + \biggl[\frac{4}{r_0} - \biggl(\frac{g_0 \rho_0}{P_0}\biggr) \biggr] \frac{dx}{dr_0} + \biggl(\frac{\rho_0}{\gamma_\mathrm{g} P_0} \biggr)\biggl[\omega^2 + (4 - 3\gamma_\mathrm{g})\frac{g_0}{r_0} \biggr] x </math>


<math>~0</math>

<math>~=</math>

<math>~ \frac{d^2x}{dr*^2} + \biggl[\frac{4}{r^*} ~-~ \rho_c^{-2 / 5} \biggl( \frac{K_c}{G} \biggr)^{1 / 2}\biggl(\frac{g_0 \rho_0}{P_0}\biggr) \biggr] \frac{dx}{dr*} ~+~ \rho_c^{-4 / 5}\biggl( \frac{K_c}{G} \biggr)\biggl(\frac{\rho_0}{\gamma_\mathrm{g} P_0} \biggr)\biggl[\omega^2 ~+~ (4 - 3\gamma_\mathrm{g})\frac{g_0}{r_0} \biggr] x </math>

 

<math>~=</math>

<math>~ \frac{d^2x}{dr*^2} ~+~ \biggl\{\frac{4}{r^*} ~-~ \rho_c^{-2 / 5} \biggl( \frac{K_c}{G} \biggr)^{1 / 2} \frac{G M_r^*}{(r^*)^2} \biggl[ \rho_c^{3 / 5} \biggl( \frac{K_c}{G}\biggr)^{1 / 2} \biggr] \biggl[ \frac{\rho_c \rho^*}{P^* K_c \rho_c^{6/5}}\biggr] \biggr\} \frac{dx}{dr*} ~+~ \rho_c^{-4 / 5}\biggl( \frac{K_c}{G} \biggr)\biggl[ \frac{\rho_c \rho^*}{\gamma_\mathrm{g}P^* K_c \rho_c^{6/5}}\biggr]

\biggl\{\omega^2 ~+~ (4 - 3\gamma_\mathrm{g}) \frac{G M_r^*}{(r^*)^2} \biggl[ \rho_c^{3 / 5} \biggl( \frac{K_c}{G}\biggr)^{1 / 2} \biggr] \frac{\rho_c^{2 / 5}}{r^*}\biggl( \frac{G}{K_c}\biggr)^{1 / 2} \biggr\}  x

</math>

 

<math>~=</math>

<math>~ \frac{d^2x}{dr*^2} ~+~ \biggl\{\frac{4}{r^*} ~-~ \frac{M_r^*}{(r^*)^2} \biggl[ \frac{\rho^*}{P^* }\biggr] \biggr\} \frac{dx}{dr*} ~+~ \biggl( \frac{1}{G\rho_c} \biggr)\biggl[ \frac{ \rho^*}{\gamma_\mathrm{g}P^* }\biggr]

\biggl\{\omega^2 ~+~ (4 - 3\gamma_\mathrm{g}) \frac{G\rho_c M_r^*}{(r^*)^3} \biggr\}  x

</math>

 

<math>~=</math>

<math>~ \frac{d^2x}{dr*^2} + \biggl\{ \frac{4}{r^*} -\biggl(\frac{\rho^*}{P^*}\biggr)\frac{ M_r^*}{(r^*)^2}\biggr\} \frac{dx}{dr*} + \biggl(\frac{\rho^*}{ P^* } \biggr)\biggl\{ \frac{\omega^2}{\gamma_\mathrm{g} G\rho_c} + \biggl(\frac{4}{\gamma_\mathrm{g}} - 3\biggr)\frac{ M_r^*}{(r^*)^3}\biggr\} x </math>

 

<math>~=</math>

<math>~ \frac{d^2x}{dr*^2} + \biggl\{ 4 -\biggl(\frac{\rho^*}{P^*}\biggr)\frac{ M_r^*}{(r^*)}\biggr\}\frac{1}{r^*} \frac{dx}{dr*} + \biggl(\frac{\rho^*}{ P^* } \biggr)\biggl\{ \frac{2\pi \sigma_c^2}{3\gamma_\mathrm{g}} ~-~\frac{\alpha_\mathrm{g} M_r^*}{(r^*)^3}\biggr\} x \, . </math>

Profile

Now, referencing the derived bipolytropic model profile, we should incorporate the following relations:


Table 2: Radial Profile of Various Physical Variables

Variable

Throughout the Core
<math>~0 \le \xi \le \xi_i</math>

Throughout the Envelope
<math>~\eta_i \le \eta \le \eta_s</math>

Plotted Profiles

<math>~\xi_i = 0.5</math>

<math>~\xi_i = 1.0</math>

<math>~\xi_i = 3.0</math>

<math>~r^*</math>

<math>\biggl( \frac{3}{2\pi} \biggr)^{1/2} \xi</math>

<math>\biggl( \frac{\mu_e}{\mu_c} \biggr)^{-1} \theta^{-2}_i (2\pi)^{-1/2}\eta</math>

 

<math>~\rho^*</math>

<math>\biggl( 1 + \frac{1}{3}\xi^2 \biggr)^{-5/2}</math>

<math>\biggl( \frac{\mu_e}{\mu_c} \biggr) \theta^{5}_i \phi(\eta)</math>

DenXi05.jpg
DenXi10.jpg
DenXi30.jpg

<math>~P^*</math>

<math>\biggl( 1 + \frac{1}{3}\xi^2 \biggr)^{-3}</math>

<math>\theta^{6}_i [\phi(\eta)]^{2}</math>

PresXi05.jpg
PresXi10.jpg
PresXi30.jpg

<math>~M_r^*</math>

<math>\biggl( \frac{2\cdot 3}{\pi } \biggr)^{1/2} \biggl[ \xi^3 \biggl( 1 + \frac{1}{3}\xi^2 \biggr)^{-3/2} \biggr]</math>

<math>\biggl( \frac{\mu_e}{\mu_c} \biggr)^{-2} \theta^{-1}_i \biggl( \frac{2}{\pi} \biggr)^{1/2} \biggl(-\eta^2 \frac{d\phi}{d\eta} \biggr)</math>

MassXi05.jpg
MassXi10.jpg
MassXi30.jpg

In order to obtain the various envelope profiles, it is necessary to evaluate <math>\phi(\eta)</math> and its first derivative using the information presented in Step 6, above.


Therefore, throughout the core we have,

<math>~\frac{\rho^*}{P^*}</math>

<math>~=</math>

<math>~\biggl( 1 + \frac{1}{3}\xi^2 \biggr)^{1 / 2} \, ;</math>

<math>~\frac{M_r^*}{r^*}</math>

<math>~=</math>

<math>~ \biggl( \frac{2\cdot 3}{\pi } \biggr)^{1/2} \biggl[ \xi^3 \biggl( 1 + \frac{1}{3}\xi^2 \biggr)^{-3/2} \biggr]\biggl( \frac{2\pi}{3}\biggr)^{1 / 2} \frac{1}{\xi} = 2 \xi^2 \biggl( 1 + \frac{1}{3}\xi^2 \biggr)^{-3/2} \, ; </math>

and, throughout the envelope we have,

<math>~\frac{\rho^*}{P^*}</math>

<math>~=</math>

<math>~ \biggl( \frac{\mu_e}{\mu_c} \biggr) \theta^{-1}_i \phi(\eta)^{-1} \, ; </math>

<math>~\frac{M_r^*}{r^*}</math>

<math>~=</math>

<math>~ \biggl( \frac{\mu_e}{\mu_c} \biggr)^{-2} \theta^{-1}_i \biggl( \frac{2}{\pi} \biggr)^{1/2} \biggl(-\eta^2 \frac{d\phi}{d\eta} \biggr) \biggl[ \biggl( \frac{\mu_e}{\mu_c} \biggr)^{-1} \theta^{-2}_i (2\pi)^{-1/2}\eta \biggr]^{-1} = 2 \biggl( \frac{\mu_e}{\mu_c} \biggr)^{-1} \theta_i \eta \biggl(-\eta^2 \frac{d\phi}{d\eta} \biggr) \, . </math>

Numerical Integration

General Approach

Here, we begin by recognizing that the 2nd-order ODE that must be integrated to obtain the desired eigenvectors has the generic form,

<math>~0</math>

<math>~=</math>

<math>~ x + \frac{\mathcal{H}}{r^*} x' + \mathcal{K}x \, , </math>

where,

<math>~x'</math>

<math>~=</math>

<math>~\frac{dx}{dr^*}</math>

      and      

<math>~x</math>

<math>~=</math>

<math>~\frac{d^2x}{d(r^*)^2} \, .</math>

Adopting the same approach as before when we integrated the LAWE for pressure-truncated polytropes, we will enlist the finite-difference approximations,

<math>~x'</math>

<math>~\approx</math>

<math>~ \frac{x_+ - x_-}{2\delta r^*} </math>

      and      

<math>~x</math>

<math>~\approx</math>

<math>~ \frac{x_+ -2x_j + x_-}{(\delta r^*)^2} \, . </math>

The finite-difference representation of the LAWE is, therefore,

<math>~\frac{x_+ -2x_j + x_-}{(\delta r^*)^2}</math>

<math>~=</math>

<math>~ -~ \frac{\mathcal{H}}{r^*} \biggl[ \frac{x_+ - x_-}{2\delta r^*} \biggr] ~-~ \mathcal{K}x_j </math>

<math>~\Rightarrow ~~~ x_+ -2x_j + x_-</math>

<math>~=</math>

<math>~ -~ \frac{\delta r^*}{2r^*} \biggl[ x_+ - x_- \biggr]\mathcal{H} ~-~ (\delta r^*)^2\mathcal{K}x_j </math>

<math>~\Rightarrow ~~~ x_{j+1} \biggl[1 + \biggl( \frac{\delta r^*}{2r^*}\biggr) \mathcal{H} \biggr]</math>

<math>~=</math>

<math>~ \biggl[ 2 - (\delta r^*)^2\mathcal{K}\biggr] x_j ~-~\biggl[ 1 - \biggl( \frac{\delta r^*}{2r^*} \biggr) \mathcal{H} \biggr]x_{j-1} \, . </math>

In what follows we will also find it useful to rewrite <math>~\mathcal{K}</math> in the form,

<math>~\mathcal{K} ~\rightarrow ~\biggl(\frac{\sigma_c^2}{\gamma_\mathrm{g}}\biggr) \mathcal{K}_1 - \alpha_\mathrm{g} \mathcal{K}_2 \, .</math>

The relevant coefficient expressions for all regions of the configuration are,

<math>~\mathcal{H}</math>

<math>~\equiv</math>

<math>~ \biggl\{ 4 -\biggl(\frac{\rho^*}{P^*}\biggr)\frac{ M_r^*}{(r^*)}\biggr\} </math>

      ,      

<math>~\mathcal{K}_1</math>

<math>~\equiv</math>

<math>~ \frac{2\pi }{3}\biggl(\frac{\rho^*}{ P^* } \biggr) </math>

      and      

<math>~\mathcal{K}_2</math>

<math>~\equiv</math>

<math>~ \biggl(\frac{\rho^*}{ P^* } \biggr)\frac{M_r^*}{(r^*)^3} \, . </math>

Special Handling at the Center

In order to kick-start the integration, we set the displacement function value to <math>~x_1 = 1</math> at the center of the configuration <math>~(\xi_1 = 0)</math>, then draw on the derived power-series expression to determine the value of the displacement function at the first radial grid line, <math>~\xi_2 = \delta\xi</math>, away from the center. Specifically, we set,

<math>~ x_2 </math>

<math>~=</math>

<math>~ x_1 \biggl[ 1 - \frac{(n+1) \mathfrak{F} (\delta\xi)^2}{60} \biggr]_{n=5} = x_1 \biggl[ 1 - \frac{\mathfrak{F} (\delta\xi)^2}{10} \biggr] \, .</math>

Special Handling at the Interface

Integrating outward from the center, the general approach will work up through the determination of <math>~x_{j+1}</math> when "j+1" refers to the interface location. In order to properly transition from the core to the envelope, we need to determine the value of the slope at this interface location. Let's do this by setting j = i, then projecting forward to what <math>~x_+</math> would be — that is, to what the amplitude just beyond the interface would be — if the core were to be extended one more zone. Then, the slope at the interface (as viewed from the perspective of the core) will be,

<math>~x'_i\biggr|_\mathrm{core}</math>

<math>~\approx</math>

<math>~ \frac{1}{2\delta r^*} \biggl\{ x_+ - x_{i-1} \biggr\} </math>

 

<math>~=</math>

<math>~ -\frac{x_{i-1}}{2\delta r^*} + \frac{1}{2\delta r^*} \biggl\{ \biggl[ 2 - (\delta r^*)^2\mathcal{K}\biggr] x_i ~-~\biggl[ 1 - \biggl( \frac{\delta r^*}{2r^*} \biggr) \mathcal{H} \biggr]x_{i-1} \biggr\}\biggl[1 + \biggl( \frac{\delta r^*}{2r^*}\biggr) \mathcal{H} \biggr]^{-1} </math>

 

<math>~=</math>

<math>~ \frac{1}{2\delta r^*} \biggl\{ \biggl[ 2 - (\delta r^*)^2\mathcal{K}\biggr] x_i ~-~\biggl[ 1 - \biggl( \frac{\delta r^*}{2r^*} \biggr) \mathcal{H} \biggr]x_{i-1} ~-~\biggl[1 + \biggl( \frac{\delta r^*}{2r^*}\biggr) \mathcal{H} \biggr]x_{i-1} \biggr\}\biggl[1 + \biggl( \frac{\delta r^*}{2r^*}\biggr) \mathcal{H} \biggr]^{-1} </math>

 

<math>~=</math>

<math>~ \frac{1}{2\delta r^*} \biggl\{ \biggl[ 2 - (\delta r^*)^2\mathcal{K}\biggr] x_i ~-~2x_{i-1} \biggr\}\biggl[1 + \biggl( \frac{\delta r^*}{2r^*}\biggr) \mathcal{H} \biggr]^{-1} </math>

Conversely, as viewed from the envelope, if we assume that we know <math>~x_i</math> and <math>~x'_i</math>, we can determine the amplitude, <math>~x_{i+1}</math>, at the first zone beyond the interface as follows:

<math>~x_-</math>

<math>~\approx</math>

<math>~ x_{i+1} - 2\delta r^*\cdot x'_i\biggr|_\mathrm{env} </math>

<math>~\Rightarrow ~~~ x_{i+1} \biggl[1 + \biggl( \frac{\delta r^*}{2r^*}\biggr) \mathcal{H} \biggr]</math>

<math>~=</math>

<math>~ \biggl[ 2 - (\delta r^*)^2\mathcal{K}\biggr] x_i ~-~\biggl[ 1 - \biggl( \frac{\delta r^*}{2r^*} \biggr) \mathcal{H} \biggr] \biggl[ x_{i+1} - 2\delta r^*\cdot x'_i\biggr|_\mathrm{env} \biggr] </math>

<math>~\Rightarrow ~~~ x_{i+1} \biggl[1 + \biggl( \frac{\delta r^*}{2r^*}\biggr) \mathcal{H} \biggr] ~+~ \biggl[ 1 - \biggl( \frac{\delta r^*}{2r^*} \biggr) \mathcal{H} \biggr] x_{i+1} </math>

<math>~=</math>

<math>~ \biggl[ 2 - (\delta r^*)^2\mathcal{K}\biggr] x_i ~+~ \biggl[ 1 - \biggl( \frac{\delta r^*}{2r^*} \biggr) \mathcal{H} \biggr] 2\delta r^*\cdot x'_i\biggr|_\mathrm{env} </math>

<math>~\Rightarrow ~~~ x_{i+1} </math>

<math>~=</math>

<math>~ \biggl[ 1 - \tfrac{1}{2}(\delta r^*)^2\mathcal{K}\biggr] x_i ~+~ \biggl[ 1 - \biggl( \frac{\delta r^*}{2r^*} \biggr) \mathcal{H} \biggr] \delta r^*\cdot x'_i\biggr|_\mathrm{env} </math>

Splitting Analysis Into Separate Core and Envelope Components

Core:

Given that, <math>~\sqrt{2\pi/3}~r^* = \xi</math>, lets multiply the LAWE through by <math>~3/(2\pi)</math>. This gives,

<math>~0</math>

<math>~=</math>

<math>~ \frac{d^2x}{d\xi^2} + \biggl\{ 4 -\biggl(\frac{\rho^*}{P^*}\biggr)\frac{ M_r^*}{(r^*)}\biggr\}\frac{1}{\xi} \cdot \frac{dx}{d\xi} + \frac{3}{2\pi}\biggl(\frac{\rho^*}{ P^* } \biggr)\biggl\{ \frac{2\pi \sigma_c^2}{3\gamma_\mathrm{g}} ~-~\frac{\alpha_\mathrm{g} M_r^*}{(r^*)^3}\biggr\} x \, . </math>

Specifically for the core, therefore, the finite-difference representation of the LAWE is,

<math>~\frac{x_+ -2x_j + x_-}{(\delta \xi)^2}</math>

<math>~=</math>

<math>~ -~ \frac{\mathcal{H}}{\xi} \biggl[ \frac{x_+ - x_-}{2\delta \xi} \biggr] ~-~ \biggl[ \frac{3\mathcal{K}}{2\pi} \biggr]x_j </math>

<math>~\Rightarrow ~~~ x_+ -2x_j + x_-</math>

<math>~=</math>

<math>~ -~ \frac{\delta \xi}{2\xi} \biggl[ x_+ - x_- \biggr]\mathcal{H} ~-~ (\delta \xi)^2 \biggl[ \frac{3\mathcal{K}}{2\pi} \biggr] x_j </math>

<math>~\Rightarrow ~~~ x_{j+1} \biggl[1 + \biggl( \frac{\delta \xi}{2\xi}\biggr) \mathcal{H} \biggr]</math>

<math>~=</math>

<math>~ \biggl[ 2 - (\delta \xi)^2\biggl( \frac{3\mathcal{K}}{2\pi} \biggr) \biggr] x_j ~-~\biggl[ 1 - \biggl( \frac{\delta \xi}{2\xi} \biggr) \mathcal{H} \biggr]x_{j-1} \, . </math>

This also means that, as viewed from the perspective of the core, the slope at the interface is

<math>~\biggl[ \frac{dx}{d\xi}\biggr]_\mathrm{interface}</math>

<math>~=</math>

<math>~ \frac{1}{2\delta \xi} \biggl\{ \biggl[ 2 - (\delta \xi)^2 \biggl( \frac{3\mathcal{K}}{2\pi} \biggr)\biggr] x_i ~-~2x_{i-1} \biggr\}\biggl[1 + \biggl( \frac{\delta \xi}{2\xi}\biggr) \mathcal{H} \biggr]^{-1} \, . </math>

Envelope:

Given that,

<math>~\biggl( \frac{\mu_e}{\mu_c} \biggr) \theta_i^2 (2\pi)^{1 / 2}~r^* = \eta \, ,</math>

let's multiply the LAWE through by <math>~(2\pi)^{-1} \theta_i^{-4}( \mu_e/\mu_c)^{-2} </math>. This gives,

<math>~0</math>

<math>~=</math>

<math>~ \frac{d^2x}{d\eta^2} + \biggl\{ 4 -\biggl(\frac{\rho^*}{P^*}\biggr)\frac{ M_r^*}{(r^*)}\biggr\}\frac{1}{\eta} \cdot \frac{dx}{d\eta} + \frac{1}{2\pi \theta_i^4} \biggl( \frac{\mu_e}{\mu_c} \biggr)^{-2} \biggl(\frac{\rho^*}{ P^* } \biggr)\biggl\{ \frac{2\pi \sigma_c^2}{3\gamma_\mathrm{g}} ~-~\frac{\alpha_\mathrm{g} M_r^*}{(r^*)^3}\biggr\} x \, . </math>

Specifically for the envelope, therefore, the finite-difference representation of the LAWE is,

<math>~\frac{x_+ -2x_j + x_-}{(\delta \eta)^2}</math>

<math>~=</math>

<math>~ -~ \frac{\mathcal{H}}{\eta} \biggl[ \frac{x_+ - x_-}{2\delta \eta} \biggr] ~-~ \biggl( \frac{\mu_e}{\mu_c} \biggr)^{-2}\biggl[ \frac{\mathcal{K}}{2\pi \theta_i^4} \biggr]x_j </math>

<math>~\Rightarrow ~~~ x_+ -2x_j + x_-</math>

<math>~=</math>

<math>~ -~ \frac{\delta \eta}{2\eta} \biggl[ x_+ - x_- \biggr]\mathcal{H} ~-~ (\delta \eta)^2 \biggl( \frac{\mu_e}{\mu_c} \biggr)^{-2}\biggl[ \frac{\mathcal{K}}{2\pi \theta_i^4} \biggr] x_j </math>

<math>~\Rightarrow ~~~ x_{j+1} \biggl[1 + \biggl( \frac{\delta \eta}{2\eta}\biggr) \mathcal{H} \biggr]</math>

<math>~=</math>

<math>~ \biggl[ 2 - (\delta \eta)^2 \biggl( \frac{\mu_e}{\mu_c} \biggr)^{-2} \biggl( \frac{\mathcal{K}}{2\pi \theta_i^4} \biggr) \biggr] x_j ~-~\biggl[ 1 - \biggl( \frac{\delta \eta}{2\eta} \biggr) \mathcal{H} \biggr]x_{j-1} \, . </math>

This also means that, once we know the slope at the interface (see immediately below), the amplitude at the first zone outside of the interface will be given by the expression,

<math>~x_{i+1} </math>

<math>~=</math>

<math>~ \biggl[ 1 - \tfrac{1}{2}(\delta \eta)^2 \biggl( \frac{\mu_e}{\mu_c} \biggr)^{-2} \biggl( \frac{\mathcal{K}}{2\pi \theta_i^4} \biggr)\biggr] x_i ~+~ \biggl[ 1 - \biggl( \frac{\delta \eta}{2\eta} \biggr) \mathcal{H} \biggr] \delta \eta \cdot \biggl[ \frac{dx}{d\eta} \biggr]_\mathrm{interface} \, . </math>

Interface

If we consider only cases where <math>~\gamma_e = \gamma_c</math>, then at the interface we expect,

<math>~\frac{d\ln x}{d\ln r^*}</math>

<math>~=</math>

<math>~\frac{d\ln x}{d\ln \xi} = \frac{d\ln x}{d\ln \eta}</math>

<math>~\Rightarrow ~~~ r^*\frac{dx}{d r^*}</math>

<math>~=</math>

<math>~\xi \frac{dx}{d \xi} = \eta \frac{d x}{d \eta}</math>

<math>~\Rightarrow ~~~ \frac{dx}{dr^*}</math>

<math>~=</math>

<math>~\biggl(\frac{2\pi}{3} \biggr)^{1 / 2}\frac{dx}{d\xi} = \biggl(\frac{\mu_e}{\mu_c}\biggr) \theta_i^2 (2\pi)^{1 / 2} \frac{dx}{d\eta} \, .</math>

Switching at the interface from <math>~\xi</math> to <math>~\eta</math> therefore means that,

<math>~ \biggl[ \frac{dx}{d\eta}\biggr]_\mathrm{interface}</math>

<math>~=</math>

<math>~\frac{1}{\sqrt{3}} \biggl(\frac{\mu_e}{\mu_c}\biggr)^{-1} \theta_i^{-2} \biggl[ \frac{dx}{d\xi}\biggr]_\mathrm{interface} \, .</math>

Model 10

As we have reviewed in an accompanying discussion, equilibrium Model 10 from Murphy & Fiedler (1985, Proc. Astr. Soc. of Australia, 6, 219) is defined by setting <math>~(\xi_i, m) = (2.5646, 1)</math>. Drawing directly from our reproduction of their Table 1, we see that a few relevant structural parameters of Model 10 are,

<math>~\xi_s</math>

<math>~=</math>

<math>~6.5252876</math>

<math>~\frac{r_i}{R} = \frac{\xi_i}{\xi_s}</math>

<math>~=</math>

<math>~0.39302482</math>

<math>~\frac{\rho_c}{\bar\rho} </math>

<math>~=</math>

<math>~34.346</math>

<math>~\frac{M_\mathrm{env}}{M_\mathrm{tot}} </math>

<math>~=</math>

<math>~5.89 \times 10^{-4}</math>

Here we list a few other model parameter values that will aid in our attempt to correctly integrate the LAWE to find various radial oscillation eigenvectors.

A Sampling of Model 10's Equilibrium Parameter Values

Grid
Line
<math>~\frac{r}{R}</math> <math>~\xi</math> <math>~\eta</math> <math>~\Delta</math> <math>~\phi</math> <math>~- \frac{d\phi}{d\eta}</math> <math>~r^*</math> <math>~\rho^*</math> <math>~P^*</math> <math>~M_r^*</math> <math>~g_0^*\equiv \frac{M_r^*}{(r^*)^2}</math>
25 0.12093071 0.789108         0.31480842 0.89940188 0.80892374 0.122726799 1.23835945
40 0.19651241 1.2823         0.51156369 0.74761972 0.55893525 0.473819194 1.81056130
79 0.393025 2.5646         1.02312737 0.21270605 0.04524386 2.150231108 2.05411964
79 0.393025   1.4806725 2.6746514 1.000000 1.112155 1.02312737 0.21270605 0.04524386 2.15023111 2.0541196
100 0.49883919   1.8793151 2.7938569 0.6505914 0.69070815 1.2985847 0.0247926 0.0034309 2.15127319 1.2757189
150 0.7507782   2.8284641 2.9982701 0.2149684 0.30495637 1.95443562 9.7646E-05 4.4649E-06 2.15149752 0.563246
199 0.9976784   3.7586302 3.1404305 0.00150695 0.17269514 2.59716948 1.653E-15 5.2984E-19 2.15149876 0.31896316

Our chosen (uniform) grid spacing is,

<math>~\frac{\delta r}{R} = \frac{1}{78}\biggl( \frac{r_i}{R} \biggr) \approx 0.00503878 \, ;</math>

as a result, the center is at zone 1, the interface is at grid line 79, and the surface is just beyond grid line 199.

Eigenvectors

Keep in mind that, for all models, we expect that, at the surface, the logarithmic derivative of each proper eigenfunction will be,

<math>~\frac{d\ln x}{d\ln r^*}\biggr|_\mathrm{surf}</math>

<math>~=</math>

<math>~\frac{\Omega^2}{\gamma} - \alpha \, .</math>

Also, keep in mind that, for Model 10 <math>~(\xi_i = 2.5646)</math>:

<math>~\frac{r_i}{R}</math>

<math>~=</math>

<math>~0.39302482</math>

      ,    

<math>~\frac{\rho_c}{\bar\rho}</math>

<math>~=</math>

<math>~34.3460405</math>


Our Determinations for Model 10

Mode <math>~\sigma_c^2</math> <math>~\Omega^2 \equiv \frac{\sigma_c^2}{2} \biggl( \frac{\rho_c}{\bar\rho}\biggr)</math> <math>~x_\mathrm{surf}</math> <math>~\frac{d\ln x}{d\ln r^*}\biggr|_\mathrm{surf}</math> <math>~\frac{r}{R}\biggr|_1</math> <math>~1 - \frac{M_r}{M_\mathrm{tot}}\biggr|_1</math> <math>~\frac{r}{R}\biggr|_2</math> <math>~1 - \frac{M_r}{M_\mathrm{tot}}\biggr|_2</math> <math>~\frac{r}{R}\biggr|_3</math> <math>~1 - \frac{M_r}{M_\mathrm{tot}}\biggr|_3</math>
expected measured
1
(Fundamental)
0.92813095170326 15.93881161 +85.17 8.963286966 8.963085 n/a n/a n/a n/a n/a n/a
2 1.237156768978 21.24571822 - 610 12.14743093 12.147337 0.5724 3.05E-05 n/a n/a n/a n/a
3 1.8656033984 32.0380449 +3225 18.62282676 18.6228 0.4845 1.35E-04 0.787 2.05E-07 n/a n/a
4 2.65901504799 45.66331921 -9410 26.79799153 26.797977 0.4459 2.620E-04 0.7096 1.834E-06 0.8632 1.189E-08
Match Figure 2 from MF85


For Model 17 <math>~(\xi_i = 3.0713)</math>:

<math>~\frac{r_i}{R}</math>

<math>~=</math>

<math>~0.93276717</math>

      ,    

<math>~\frac{\rho_c}{\bar\rho}</math>

<math>~=</math>

<math>~3.79693903</math>


Our Determinations for Model 17

Mode <math>~\sigma_c^2</math> <math>~\Omega^2 \equiv \frac{\sigma_c^2}{2} \biggl( \frac{\rho_c}{\bar\rho}\biggr)</math> <math>~x_\mathrm{surf}</math> <math>~\frac{d\ln x}{d\ln r^*}\biggr|_\mathrm{surf}</math> <math>~\frac{r}{R}\biggr|_1</math> <math>~1 - \frac{M_r}{M_\mathrm{tot}}\biggr|_1</math> <math>~\frac{r}{R}\biggr|_2</math> <math>~1 - \frac{M_r}{M_\mathrm{tot}}\biggr|_2</math> <math>~\frac{r}{R}\biggr|_3</math> <math>~1 - \frac{M_r}{M_\mathrm{tot}}\biggr|_3</math>
expected measured
1
(Fundamental)
1.149837904 2.182932207 +1.275 0.7097593 0.7097550 n/a n/a n/a n/a n/a n/a
2 7.34212930615 13.93880866 - 2.491 7.763285 7.763244 0.7215 0.24006 n/a n/a n/a n/a
3 16.345072567 31.03062198 +4.33 18.01837 18.01826 0.5806 0.5027 0.848 0.0541 n/a n/a
4 27.746934203 52.6767087 -9.1 31.0060 31.0058 0.4859 0.6737 0.7429 0.1974 0.8957 0.0171
Match Figure 3 from MF85


Numerical Values for Some Selected <math>~(n_c, n_e) = (1, 5)</math> Bipolytropes
[to be compared with Table 1 of Murphy & Fiedler (1985)]

MODEL Source <math>~\frac{r_i}{R}</math> <math>~\Omega_0^2</math> <math>~\Omega_1^2</math> <math>~\frac{r}{R}\biggr|_1</math> <math>~1-\frac{M_r}{M_\mathrm{tot}}\biggr|_1</math>
10 MF85 0.393 15.9298 21.2310 0.573 1.00E-03
Here 0.39302 15.93881161 21.24571822 0.5724 3.05E-05
17 MF85 0.933 2.1827 13.9351 0.722 0.232
Here 0.93277 2.182932207 13.93880866 0.7215 0.24006

Summary

dynamical stability in qNu plane

dynamical stability in qNu plane


Eigenfunctions for Marginally Unstable Models

Eigenfunctioins for multiple modes in three separate models

Virial Analysis

Various Energy Expressions

Drawing from an accompanying detailed discussion, we have,

<math>~\mathfrak{G}</math>

<math>~=</math>

<math>~W_\mathrm{grav}\biggr|_\mathrm{core} + \mathfrak{S}_A\biggr|_\mathrm{core} + W_\mathrm{grav}\biggr|_\mathrm{env} + \mathfrak{S}_A\biggr|_\mathrm{env} </math>

 

<math>~=</math>

<math> ~W_\mathrm{grav}\biggr|_\mathrm{core} + \biggl[ \frac{2}{3(\gamma_c - 1)} \biggr] S_\mathrm{core} + W_\mathrm{grav}\biggr|_\mathrm{env} + \biggl[ \frac{2}{3(\gamma_e - 1)} \biggr] S_\mathrm{env} \, . </math>

Recognizing that,

<math>~q</math>

<math>~\equiv</math>

<math>~\frac{r_i}{R} = \frac{\eta_i}{\eta_s} \, ,</math>

<math>~\nu</math>

<math>~\equiv</math>

<math>~\frac{M_\mathrm{core}}{M_\mathrm{tot}} \, ,</math>

<math>~\ell</math>

<math>~\equiv</math>

<math>~\frac{\xi}{\sqrt{3}} \, ,</math>

<math>~\chi_\mathrm{eq} </math>

<math>~=</math>

<math>~\biggl( \frac{\pi}{2^3 \cdot 3^6}\biggr)^{1 / 2} \biggl( \frac{\nu}{q^3} \biggr)^{2} \biggl( 1 + \ell_i^2 \biggr)^{3} \biggl( \frac{q}{\ell_i}\biggr)^{5 } \, , </math>

we have,

<math>~- \chi \biggl[ \frac{W_\mathrm{grav}}{E_\mathrm{norm}} \biggr]_\mathrm{core}</math>

<math>~=</math>

<math> \biggl[\biggl(\frac{3}{2^4}\biggr) \biggl( \frac{q}{\ell_i}\biggr)^{5}\biggl(\frac{\nu}{q^3} \biggr)^2 \biggl( 1 + \ell_i^2 \biggr)^{3} \biggr]_\mathrm{eq} \biggl[ \ell_i (\ell_i^4 - \frac{8}{3}\ell_i^2 - 1) (\ell_i^2 +1)^{-3} + \tan^{-1}\ell_i \biggr] </math>

<math>~\Rightarrow ~~~ \biggl[ \frac{W_\mathrm{grav}}{E_\mathrm{norm}} \biggr]_\mathrm{core}</math>

<math>~=</math>

<math>- \frac{\chi_\mathrm{eq}}{\chi} \biggl( \frac{3^8}{2^5 \pi}\biggr)^{1 / 2} \biggl[ \ell_i (\ell_i^4 - \frac{8}{3}\ell_i^2 - 1) (\ell_i^2 +1)^{-3} + \tan^{-1}\ell_i \biggr] \, . </math>

Next, given that,

<math>~b_\eta = \eta_s</math>

       and       

<math>~b_\eta q = \eta_i = 3\biggl( \frac{\mu_e}{\mu_c} \biggr) \ell_i (1 + \ell_i^2)^{-1} \, ,</math>

we have,

<math>~\biggl[ \frac{W_\mathrm{grav}}{E_\mathrm{norm}} \biggr]_\mathrm{env}</math>

<math>~=</math>

<math> - \frac{\chi_\mathrm{eq}}{\chi} \biggl( \frac{\mu_e}{\mu_c} \biggr)^{-3} \biggl( \frac{1}{2^3\pi} \biggr)^{1/2} A^2 \biggl[6 b_\eta x - 3\sin[2(b_\eta x-B)] - 4b_\eta x \sin^2(b_\eta x - B) \biggr]_q^1 </math>

 

<math>~=</math>

<math> - \frac{\chi_\mathrm{eq}}{\chi} \biggl( \frac{\mu_e}{\mu_c} \biggr)^{-3} \biggl( \frac{1}{2^3\pi} \biggr)^{1/2} A^2 \biggl\{ \biggl[ 6 \eta_s - 3\sin[2(\eta_s -B)] - 4\eta_s \sin^2(\eta_s - B) \biggr] - \biggl[ 6 \eta_i - 3\sin[2(\eta_i -B)] - 4\eta_i \sin^2(\eta_i - B) \biggr] \biggr\} \, . </math>

Also, given that,

<math>~q a_\xi^{1 / 2}</math>

<math>~=</math>

<math>~\frac{\xi_i}{\sqrt{3}} \, ,</math>

we have,

<math>~\biggl( \frac{\mathfrak{S}_A}{E_\mathrm{norm}} \biggr)_\mathrm{core}</math>

<math>~=</math>

<math> \frac{2}{3({\gamma_c}-1)} \biggl( \frac{\chi}{\chi_\mathrm{eq}} \biggr)^{3-3\gamma_c} \biggl\{\biggl( \frac{3^8}{2^7\pi} \biggr)^{1/2} \biggl[ \tan^{-1}[a_\xi^{1/2}q] - a_\xi^{1/2}q ~\frac{(1 - a_\xi q^2)}{(1 + a_\xi q^2)^2} \biggr] \biggr\} </math>

 

<math>~=</math>

<math> \frac{2}{3({\gamma_c}-1)} \biggl( \frac{\chi}{\chi_\mathrm{eq}} \biggr)^{3-3\gamma_c} \biggl\{\biggl( \frac{3^8}{2^7\pi} \biggr)^{1/2} \biggl[ \tan^{-1}\biggl( \frac{\xi_i}{\sqrt{3}}\biggr) - \biggl( \frac{\xi_i}{\sqrt{3}}\biggr)\frac{(1 - \xi_i^2/3)}{(1 + \xi_i^2/3)^2} \biggr] \biggr\} </math>

Finally, then, we can write,

<math>~\biggl( \frac{\mathfrak{S}_A}{E_\mathrm{norm}} \biggr)_\mathrm{env}</math>

<math>~=</math>

<math> \frac{2}{3({\gamma_e}-1)} \biggl( \frac{\chi}{\chi_\mathrm{eq}} \biggr)^{3-3\gamma_e} \biggl\{ \biggl( \frac{\mu_e}{\mu_c} \biggr)^{-3} A^2 \biggl( \frac{3^2}{2^5\pi} \biggr)^{1/2} \biggl[ 2b_\eta x -\sin[2(b_\eta x - B)] \biggr]_q^1 \biggr\} </math>

 

<math>~=</math>

<math> \frac{2}{3({\gamma_e}-1)} \biggl( \frac{\chi}{\chi_\mathrm{eq}} \biggr)^{3-3\gamma_e} \bigg( \frac{\mu_e}{\mu_c} \biggr)^{-3} A^2 \biggl( \frac{3^2}{2^5\pi} \biggr)^{1/2} \biggl\{ \biggl[ 2\eta_s -\sin[2(\eta_s - B)] \biggr] - \biggl[ 2\eta_i -\sin[2(\eta_i - B)] \biggr] \biggr\} \ . </math>

Free Energy and Its Derivatives

We can now rewrite the free-energy expression in the form,

<math>~\mathfrak{g} \equiv \frac{\mathfrak{G}}{E_\mathrm{norm}}</math>

<math>~=</math>

<math> ~-\biggl(\frac{\chi_\mathrm{eq}}{\chi}\biggr)\mathfrak{w} + \biggl[ \frac{2}{3(\gamma_c - 1)} \biggr] \biggl(\frac{\chi}{\chi_\mathrm{eq}}\biggr)^{3 - 3\gamma_c}\mathfrak{s}_\mathrm{core} + \biggl[ \frac{2}{3(\gamma_e - 1)} \biggr] \biggl(\frac{\chi}{\chi_\mathrm{eq}}\biggr)^{3 - 3\gamma_e}\mathfrak{s}_\mathrm{env} \, , </math>

where,

<math>~\mathfrak{s}_\mathrm{core}</math>

<math>~\equiv</math>

<math>~ \biggl( \frac{3^8}{2^7\pi} \biggr)^{1/2} \biggl[ \tan^{-1}\biggl( \frac{\xi_i}{\sqrt{3}}\biggr) - \biggl( \frac{\xi_i}{\sqrt{3}}\biggr)\frac{(1 - \xi_i^2/3)}{(1 + \xi_i^2/3)^2} \biggr]\, , </math>

<math>~\mathfrak{s}_\mathrm{env}</math>

<math>~\equiv</math>

<math>~ \bigg( \frac{\mu_e}{\mu_c} \biggr)^{-3} A^2 \biggl( \frac{3^2}{2^5\pi} \biggr)^{1/2} \biggl\{ \biggl[ 2\eta_s -\sin[2(\eta_s - B)] \biggr] - \biggl[ 2\eta_i -\sin[2(\eta_i - B)] \biggr] \biggr\} \, , </math>

<math>~\mathfrak{w}</math>

<math>~\equiv</math>

<math>~ \biggl( \frac{3^8}{2^5 \pi}\biggr)^{1 / 2} \biggl[ \ell_i (\ell_i^4 - \frac{8}{3}\ell_i^2 - 1) (\ell_i^2 +1)^{-3} + \tan^{-1}\ell_i \biggr] </math>

 

 

<math>~ +~ \biggl( \frac{\mu_e}{\mu_c} \biggr)^{-3} \biggl( \frac{1}{2^3\pi} \biggr)^{1/2} A^2 \biggl\{ \biggl[ 6 \eta_s - 3\sin[2(\eta_s -B)] - 4\eta_s \sin^2(\eta_s - B) \biggr] - \biggl[ 6 \eta_i - 3\sin[2(\eta_i -B)] - 4\eta_i \sin^2(\eta_i - B) \biggr] \biggr\} \, . </math>

The first derivative is,

<math>~\frac{d\mathfrak{g}}{d\chi} </math>

<math>~=</math>

<math> +~\biggl(\frac{\chi_\mathrm{eq}}{\chi^2}\biggr)\mathfrak{w} ~-~\frac{2}{\chi_\mathrm{eq}} \biggl[\biggl(\frac{\chi}{\chi_\mathrm{eq}}\biggr)^{2 - 3\gamma_c}\mathfrak{s}_\mathrm{core} ~+~\biggl(\frac{\chi}{\chi_\mathrm{eq}}\biggr)^{2 - 3\gamma_e}\mathfrak{s}_\mathrm{env} \biggr] \, . </math>

And the second derivative is,

<math>~\frac{d^2\mathfrak{g}}{d\chi^2} </math>

<math>~=</math>

<math> -~2\biggl(\frac{\chi_\mathrm{eq}}{\chi^3}\biggr)\mathfrak{w} ~-~\frac{2}{\chi_\mathrm{eq}^2} \biggl[(2 - 3\gamma_c) \biggl(\frac{\chi}{\chi_\mathrm{eq}}\biggr)^{1 - 3\gamma_c}\mathfrak{s}_\mathrm{core} ~+~(2 - 3\gamma_e)\biggl(\frac{\chi}{\chi_\mathrm{eq}}\biggr)^{1 - 3\gamma_e}\mathfrak{s}_\mathrm{env} \biggr] \ . </math>

What to Expect for Equilibrium Configurations

In equilibrium we should set <math>~d\mathfrak{g}/d\chi = 0</math> and <math>~\chi = \chi_\mathrm{eq}</math>. In which case we should expect to find that,

<math>~0 </math>

<math>~=</math>

<math> \frac{1}{\chi_\mathrm{eq}}\biggl[\mathfrak{w} ~-~2(\mathfrak{s}_\mathrm{core} ~+~\mathfrak{s}_\mathrm{env}) \biggr] \, . </math>

In addition, an evaluation of the second derivative should give,

<math>~\frac{\chi_\mathrm{eq}^2}{2} \cdot \frac{d^2\mathfrak{g}}{d\chi^2} </math>

<math>~=</math>

<math> -~\mathfrak{w} ~-~ \biggl[(2 - 3\gamma_c) \mathfrak{s}_\mathrm{core} ~+~(2 - 3\gamma_e)\mathfrak{s}_\mathrm{env} \biggr] \, . </math>

Then, the transition from stable to unstable configurations occurs when <math>~d^2\mathfrak{g}/d\chi^2 = 0</math>, that is, when,

<math>~\mathfrak{w} </math>

<math>~=</math>

<math> ~-~ \biggl[(2 - 3\gamma_c) \mathfrak{s}_\mathrm{core} ~+~(2 - 3\gamma_e)\mathfrak{s}_\mathrm{env} \biggr] </math>

<math>~\Rightarrow ~~~0</math>

<math>~=</math>

<math> ~(4 - 3\gamma_c) \mathfrak{s}_\mathrm{core} ~+~(4 - 3\gamma_e)\mathfrak{s}_\mathrm{env} </math>

<math>~\Rightarrow ~~~\frac{\mathfrak{s}_\mathrm{core} }{\mathfrak{s}_\mathrm{env}} </math>

<math>~=</math>

<math> ~-~\frac{(4 - 3\gamma_e) }{(4 - 3\gamma_c) } \, . </math>

For example, if we set <math>~\gamma_c = 6/5</math> and <math>~\gamma_e = 2</math>, this implies,

<math>~\frac{\mathfrak{s}_\mathrm{core} }{\mathfrak{s}_\mathrm{env}} </math>

<math>~=</math>

<math> ~\frac{2 }{(4 - 18/5) } = 5 \, . </math>

Reconciliation

In an accompanying discussion of the viral stability of <math>~(n_c, n_e) = (5, 1)</math > bipolytropes — see also an associated discussion of the free-energy of these configurations — we find the following table of models. Let's start by finding a set of eigenvectors associated with each the subset of models highlighted in yellow.

Table 1: Points Defining Stability Curve

<math>~\mu_e/\mu_c</math>

<math>~\xi_i</math>

<math>~q </math>

<math>~\nu </math>

 

<math>~\mu_e/\mu_c</math>

<math>~\xi_i</math>

<math>~q </math>

<math>~\nu </math>

1

2.416

0.5952

0.6830

0.375

6.259

0.1695

0.6054

0.95

2.500

0.5805

0.6884

0.350

7.341

0.1284

0.5439

0.90

2.594

0.5642

0.6937

0.340

8.0

0.1109

0.5083

0.80

2.816

0.5255

0.7031

8.6

0.0986

0.4799

0.70

3.109

0.4775

0.7104

0.32

10.2

0.0744

0.4038

0.65

3.296

0.4481

0.7124

0.31

12.4

0.05536

0.3264

0.60

3.523

0.4142

0.7125

0.305

14.4

0.04494

0.2772

0.55

3.809

0.3748

0.7096

0.3

17.733

0.03412

0.2186

½

4.186

0.3284

0.7014

0.295

25.737

0.02165

0.14347

0.45

4.719

0.2733

0.6830

0.291

75.510

0.00666

0.0450

0.40

5.574

0.2073

0.6429

 

 

 

 

PlotStabilityBest02.png

Figure 3: Largely the same as Figure 1, above, but a red-dashed curve has been added that separates the <math>~q - \nu</math> domain into regions that contain stable models (lying below the curve) from dynamically unstable models (lying above the curve), as determined by the virial stability analysis presented here.


In the context of our present discussion, we find that when <math>~\mu_e/\mu_c = 1/3</math>:

  • Our virial analysis shows that the model having marginal stability has the properties: <math>~\xi = 8.548</math>, <math>~q = 0.0990</math>, and <math>~\nu = 0.4790</math>.
  • In contrast, our solution to the LAWE identifies the marginally stable configuration has the properties: <math>~\xi = 2.78434248</math>, <math>~q = 0.1741</math>, and <math>~\nu = 0.2377</math>.

See Also


Whitworth's (1981) Isothermal Free-Energy Surface

© 2014 - 2021 by Joel E. Tohline
|   H_Book Home   |   YouTube   |
Appendices: | Equations | Variables | References | Ramblings | Images | myphys.lsu | ADS |
Recommended citation:   Tohline, Joel E. (2021), The Structure, Stability, & Dynamics of Self-Gravitating Fluids, a (MediaWiki-based) Vistrails.org publication, https://www.vistrails.org/index.php/User:Tohline/citation