SSC/Stability/BiPolytropes/PlannedApproach

From jetwiki
Jump to navigation Jump to search

Marginally Unstable (nc,ne) = (5,1) Bipolytropes

Overview

file = Dropbox/WorkFolder/Wiki edits/EmbeddedPolytropes/CombinedSequences.xlsx --- worksheet = EqSeqCombined2
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 1n6. (Another version of this figure includes the isothermal sequence.) On each sequence for which n3, 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, γg=(n+1)/n, 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, (nc,ne)=(5,1). For a given choice of μe/μc — 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, ξi, from zero to infinity. Figure 2, whose content is essentially the same as Figure 1 of this separate chapter, shows how the fractional core mass, νMcore/Mtot, varies with the fractional core radius, qrcore/R, along sequences having seven different values of μe/μc, as labeled: 1 (black), ½ (dark blue), 0.345 (brown), ⅓ (dark green), 0.316943 (purple), 0.309 (orange), and ¼ (light blue).

When modeling bipolytropes, the default expectation is that an increase in ξi along a given sequence will correspond to an increase in the relative size — both the radius and the mass — of the core. This expectation is realized along the Figure 2 sequences that have the largest mean-molecular weight ratios: μe/μc = 1 and ½. But the behavior is different along the other five illustrated sequences. For sufficiently large ξi, the relative radius of the core begins to decrease; along each sequence, a solid purple circular marker identifies the location of this turning point in radius. Furthermore, along sequences for which μe/μc<13, eventually the fractional mass of the core reaches a maximum and, thereafter, decreases even as the value of ξi continues to increase; a solid green circular marker identifies the location of this maximum mass turning point along each of these sequences. (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

  1. Earlier Planned Approach
  2. Headscratching
  3. Succinct Discussions
  4. 51Models
Figure 2: Equilibrium Sequences of Bipolytropes
with (nc,ne)=(5,1) and Various μe/μc

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 (nc,ne)=(0,0). 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 (nc,ne)=(5,1) 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. Astron. Soc. Australia, Vol. 6, no. 2, pp. 222 - 226). 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 (nc,ne)=(1,5). In an accompanying chapter, we describe in detail how 📚 Murphy & Fiedler (1985b) obtained these equilibrium bipolytropic structures and detail some of their equilibrium properties.

Here are the steps we initially plan to take:

  • Governing LAWEs:
  • 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, μe/μc=1/4. 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 (γc,γe)=(6/5,2) and initially for the fundamental mode of oscillation.

Virial Analysis

In an accompanying chapter we have examined the viral stability of (nc,ne)=(5,1) bipolytropes — see also an associated discussion of the free-energy of these configurations. We present, here, an overview of this separate discussion in an effort to provide a broader perspective to the focus of this current chapter.

Various Energy Expressions

Drawing from our accompanying detailed discussion, we have,

𝔊

=

Wgrav|core+𝔖A|core+Wgrav|env+𝔖A|env

 

=

Wgrav|core+[23(γc1)]Score+Wgrav|env+[23(γe1)]Senv.

Recognizing that,

q

riR=ηiηs,

ν

McoreMtot,

ξ3,

χeq

=

(π2336)1/2(νq3)2(1+i2)3(qi)5,

we have,

χ[WgravEnorm]core

=

[(324)(qi)5(νq3)2(1+i2)3]eq[i(i483i21)(i2+1)3+tan1i]

[WgravEnorm]core

=

χeqχ(3825π)1/2[i(i483i21)(i2+1)3+tan1i].

Next, given that,

bη=ηs

       and       

bηq=ηi=3(μeμc)i(1+i2)1,

we have,

[WgravEnorm]env

=

χeqχ(μeμc)3(123π)1/2A2[6bηx3sin[2(bηxB)]4bηxsin2(bηxB)]q1

 

=

χeqχ(μeμc)3(123π)1/2A2{[6ηs3sin[2(ηsB)]4ηssin2(ηsB)][6ηi3sin[2(ηiB)]4ηisin2(ηiB)]}.

Also, given that,

qaξ1/2

=

ξi3,

we have,

(𝔖AEnorm)core

=

23(γc1)(χχeq)33γc{(3827π)1/2[tan1[aξ1/2q]aξ1/2q(1aξq2)(1+aξq2)2]}

 

=

23(γc1)(χχeq)33γc{(3827π)1/2[tan1(ξi3)(ξi3)(1ξi2/3)(1+ξi2/3)2]}

Finally, then, we can write,

(𝔖AEnorm)env

=

23(γe1)(χχeq)33γe{(μeμc)3A2(3225π)1/2[2bηxsin[2(bηxB)]]q1}

 

=

23(γe1)(χχeq)33γe(μeμc)3A2(3225π)1/2{[2ηssin[2(ηsB)]][2ηisin[2(ηiB)]]} .

Free Energy and Its Derivatives

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

𝔤𝔊Enorm

=

(χeqχ)𝔴+[23(γc1)](χχeq)33γc𝔰core+[23(γe1)](χχeq)33γe𝔰env,

where,

𝔰core

(3827π)1/2[tan1(ξi3)(ξi3)(1ξi2/3)(1+ξi2/3)2],

𝔰env

(μeμc)3A2(3225π)1/2{[2ηssin[2(ηsB)]][2ηisin[2(ηiB)]]},

𝔴

(3825π)1/2[i(i483i21)(i2+1)3+tan1i]

 

 

+(μeμc)3(123π)1/2A2{[6ηs3sin[2(ηsB)]4ηssin2(ηsB)][6ηi3sin[2(ηiB)]4ηisin2(ηiB)]}.

The first derivative is,

d𝔤dχ

=

+(χeqχ2)𝔴2χeq[(χχeq)23γc𝔰core+(χχeq)23γe𝔰env].

And the second derivative is,

d2𝔤dχ2

=

2(χeqχ3)𝔴2χeq2[(23γc)(χχeq)13γc𝔰core+(23γe)(χχeq)13γe𝔰env] .

What to Expect for Equilibrium Configurations

In equilibrium we should set d𝔤/dχ=0 and χ=χeq. In which case we should expect to find that,

0

=

1χeq[𝔴2(𝔰core+𝔰env)].

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

χeq22d2𝔤dχ2

=

𝔴[(23γc)𝔰core+(23γe)𝔰env].

Then, the transition from stable to unstable configurations occurs when d2𝔤/dχ2=0, that is, when,

𝔴

=

[(23γc)𝔰core+(23γe)𝔰env]

0

=

(43γc)𝔰core+(43γe)𝔰env

𝔰core𝔰env

=

(43γe)(43γc).

For example, if we set γc=6/5 and γe=2, this implies,

𝔰core𝔰env

=

2(418/5)=5.


Free-Energy Stability Evaluation

Here we pull together excerpts from several different H_Book Chapters in which we have presented, from several different perspectives, an analysis of the free-energy of bipolytropes.

  • In one chapter, using purely analytic techniques, we have derived expressions that detail the structural properties of bipolytropes having (nc,ne)=(5,1). Among these are analytic expressions for various terms that make up the free-energy expression: 𝔰core, 𝔰env, 𝔴core, 𝔴core, and PiVcore. Equilibrium model sequences are defined by fixing the ratio, μe/μc, then varying the radial location, ri, of the core-envelope interface; note that the volume of the core is, then, Vcore4πri3/3.
  • [Virial Equilibrium] In a subsection of this same chapter titled, Equilibrium Condition: Global, we have shown that a statement of virial equilibrium — obtained by setting the first derivative of the free-energy expression to zero — is,

    (2S+W)tot=2(Score+Senv)+(Wcore+Wenv)

    =

    0.

    In another subsection of this same chapter titled, Equilibrium Condition: In Parts, we showed that for each bipolytropic equilibrium structure, the statements

    2Score+Wcore=3PiVcore

        and    

    2Senv+Wenv=3PiVcore,

    also hold separately. Therefore, for every equilibrium configuration we should expect the CASE1 expression (see Table XXX) to precisely sum to unity.
  • [Marginally Unstable Model] Near the bottom of this same chapter, in a subsection titled, Stability Condition, we point out that the model along each sequence that is marginally (dynamically) unstable — obtained setting the second derivative of the free-energy expression to zero — is identified by the configuration for which,

    2(γeγc)𝔰core+(𝔴core+𝔴env)(γe43)

    =

    0.

    Therefore, along each equilibrium sequence, the marginally unstable model can be identified by the configuration for which the CASE2 expression (see Table XXX) precisely sums to zero. Immediately above, in a subsection titled, What to Expect for Equilibrium Configurations, we have shown that this same marginally unstable model can be identified by the configuration for which the CASE3 energy ratio, 𝔰core/𝔰env, has a value that is precisely 5. And, thirdly, as highlighted in our accompanying Tabular Overview, this same marginally unstable model can be identified by the configuration for which the CASE4 expression precisely sums to zero.


Table XXX: Properties of Marginally Unstable Bipolytropes Having

(nc,ne)=(5,1) and (γc,γe)=(65,2)

Determined from Free-Energy Arguments

μeμc ξi Rsurf* qrcoreRsurf νMcoreMtot ρcρ¯ Pi 𝔰core 𝔴core 𝔰env 𝔴env 3PiVcore CASE1 CASE2 CASE3 CASE4
1 2.41610822 2.8049 0.59520261 0.68306067 16.3788 0.039116848 4.446748782 - 6.606135366 0.889349762 - 4.066061722 2.287362198 1 0 5.0000000 2.1 × 10-8
12 4.1853093 8.8058 0.328419479 0.70131896 354.089 0.003126324 5.76978580 - 10.58931853 1.153956968 - 3.258165567 0.95025163 1 6.3 × 10-8 5.0000002 0
0.345 7.64325 44.116 0.119714454 0.52700045 2.85 × 104 0.000116533 6.230343527 - 12.24495934 1.1246068658 - 2.707865028 0.215727713 1 6.4 × 10-8 5.0000002 0
13 8.548103 59.643 0.099032423 0.47901529 6.30 × 104 6.1337 × 10-5 6.261548334 - 12.36425897 1.252309682 - 2.663457063 0.158837699 1 0 4.9999999 6.0 × 10-8
0.316943 10.7441565 108.14 0.068655205 0.38238387 2.93 × 105 1.6252 × 10-5 6.301810768 - 12.52005323 1.260362204 - 2.604292714 0.083568307 1 0 4.9999998 2.0 × 10-7
0.309 12.77156 166.06 0.053145011 0.31696879 8.70 × 105 5.8905 × 10-6 6.318902171 - 12.58692884 1.26378042 - 2.57843634 0.050875500 1 1.9 × 10-8 5.0000001 0

CASE1

16PiVcore[(2𝔰core+𝔴core)(2𝔰core+𝔴core)]

CASE2

2(γeγc)𝔰core+(𝔴core+𝔴env)(γe43)

CASE3

𝔰core𝔰env

CASE4

(43γc)𝔴core+(43γe)𝔴env+32(γcγe)PiVcore


The left-hand panel of Figure 3 is identical to Figure 2, above. It displays in the qν parameter space, the behavior of (nc,ne)=(5,1) bipolytropic equilibrium sequences that have, as labeled, seven different values of the ratio of mean-molecular-weights, μe/μc. Using a numerical root-finding technique, we have determined where the virial stability condition, 𝔰core/𝔰env=5, is satisfied along each of these sequences — as well as along a number of additional equilibrium sequences. Key properties of each of these identified models have been recorded in Table 1 of an accompanying discussion; see also an associated discussion of the free-energy of these configurations. Pulling from this tabulated data, the solid-red circular markers that appear in the right-hand panel of Figure 3 identify where this virial stability condition is satisfied along the separate equilibrium sequences while the accompanying dashed red curve identifies more broadly how the qν parameter space is divided into stable (below and to the right) versus unstable (above and to the left) regions.

Figure 3

dynamical stability in qNu plane

In what follows we use a complementary — and more quantitatively rigorous — approach to evaluating the stability of equilibrium models, and contrast the results of that analysis with the virial-analysis results presented graphically here in Figure 3.

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

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

dxdξ|ξ=0=0;

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

dlnxdlnξ

=

α+ω2γg(14πGρc)ξ(θ')         at         ξ=ξs.

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

dlnxdlnξ|s

=

α+ω2R3γgGMtot.

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:

Ω2

ω2[R3GMtot]

 

=

ω2[34πGρ¯]=ω2[34πGρc]ρcρ¯=3ω2(nc+1)[(nc+1)4πGρc]ρcρ¯

 

=

3ω2(nc+1)[an2ρcPcθc]ρcρ¯=3γ(nc+1)ρcρ¯[an2ρcPcω2θcγ].

For a given radial quantum number, k, the factor inside the square brackets in this last expression is what 📚 Murphy & Fiedler (1985b) refer to as ωk2θc. Keep in mind, as well, that, in the notation we are using,

σc2

3ω22πGρc

σc2

=

(2ρ¯ρc)Ω2=6γ(nc+1)[an2ρcPcω2θcγ]=6γ(nc+1)[ωk2θc].

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

dlnxdlnξ|s

=

Ω2γgα.

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 (nc,ne)=(0,0) bipolytropes; specifically, we will draw from STEP 4: in the Piecing Together subsection. Following the discussion in §§57 & 58 of 📚 Ledoux & Walraven (1958), the proper treatment is to ensure that fractional perturbation in the gas pressure (see their equation 57.31),

δPP

=

γx(3+dlnxdlnξ),

is continuous across the interface.

Reaffirmation:   In our introductory discussion of the eigenvalue problem, we adopted the following expression for the time-dependent pressure,

P(m,t)

=

P0(m)+P1(m,t)=P0(m)[1+p(m)eiωt].

In this expression, P0(m) is the function that details how the unperturbed pressure varies with Lagrangian mass shell (m), and P1=δP(m)eiωt traces the variation of the pressure away from its equilibrium value at each mass shell. The time-dependent and spatially dependent behavior of P1 has been separated, with δP carrying information about the function's spatial dependence. Furthermore, we have adopted the shorthand notation,

p(m)δP(m)P0(m).

We can just as well use r0(m) to tag the (initial, unperturbed location of the) Lagrangian mass shells, in which case we write,

p(r0)=δP(r0)P0(r0),     and, similarly,    d(r0)=δρ(r0)ρ0(r0),     and,     x(r0)=δr(r0)r0.

These three spatially dependent quantities — p,d, and x — are related to one another via the set of linearized governing relations, namely,

Linearized
Equation of Continuity
r0dxdr0=3xd,

Linearized
Euler + Poisson Equations
P0ρ0dpdr0=(4x+p)g0+ω2r0x,

Linearized
Adiabatic Form of the
First Law of Thermodynamics

p=γgd.

Using the third of these expressions to replace d in favor of p in the first expression, we find that,

r0dxdr0

=

3xpγg

p

=

γgx[3+dlnxdlnr0],

which is identical to the pressure-perturbation expression used by 📚 Ledoux & Walraven (1958) and referenced above. As they state, the function, p=δP/P0, should be continuous across the core-envelope interface.

That is to say, at the interface (ξ=ξi), we need to enforce the relation,

0

=

[γcxcore(3+dlnxcoredlnξ)γexenv(3+dlnxenvdlnξ)]ξ=ξi

 

=

γe[γcγe(3+dlnxcoredlnξ)(3+dlnxenvdlnξ)]ξ=ξi

dlnxenvdlnξ|ξ=ξi

=

3(γcγe1)+γcγe(dlnxcoredlnξ)ξ=ξi.

In the context of this interface-matching constraint (see their equation 62.1), 📚 Ledoux & Walraven (1958) state the following:   In the static (i.e., unperturbed equilibrium) modeldiscontinuities in ρ or in γ might occur at some [radius]. In the first case — that is, a discontinuity only in density, while γe=γc — the interface conditions imply the continuity of 1xdxdξ 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 γ.

The algorithm that 📚 Murphy & Fiedler (1985b) used to "… [integrate] through each zone …" was designed "… with continuity in x and dx/dξ being imposed at the interface …" Given that they set γc=γe=5/3, their interface matching condition is consistent with the one prescribed by 📚 Ledoux & 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

d2xdr02+[4r0(g0ρ0P0)]dxdr0+(ρ0γgP0)[ω2+(43γg)g0r0]x=0

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 (nc,ne)=(5,1), it makes sense to adopt the normalizations used when defining the equilibrium structure, namely,

ρ*

ρ0ρc

;    

r*

r0[Kc1/2/(G1/2ρc2/5)]

P*

P0Kcρc6/5

;    

Mr*

Mr[Kc3/2/(G3/2ρc1/5)]

We note as well that,

g0

=

GM(r0)r02

 

=

G[Mr*ρc1/5(KcG)3/2][r*ρc2/5(KcG)1/2]2

 

=

GMr*(r*)2[ρc3/5(KcG)1/2].

Hence, multiplying the LAWE through by (Kc/G)ρc4/5 gives,

0

=

d2xdr02+[4r0(g0ρ0P0)]dxdr0+(ρ0γgP0)[ω2+(43γg)g0r0]x


0

=

d2xdr*2+[4r*ρc2/5(KcG)1/2(g0ρ0P0)]dxdr*+ρc4/5(KcG)(ρ0γgP0)[ω2+(43γg)g0r0]x

 

=

d2xdr*2+{4r*ρc2/5(KcG)1/2GMr*(r*)2[ρc3/5(KcG)1/2][ρcρ*P*Kcρc6/5]}dxdr*+ρc4/5(KcG)[ρcρ*γgP*Kcρc6/5]{ω2+(43γg)GMr*(r*)2[ρc3/5(KcG)1/2]ρc2/5r*(GKc)1/2}x

 

=

d2xdr*2+{4r*Mr*(r*)2[ρ*P*]}dxdr*+(1Gρc)[ρ*γgP*]{ω2+(43γg)GρcMr*(r*)3}x

 

=

d2xdr*2+{4r*(ρ*P*)Mr*(r*)2}dxdr*+(ρ*P*){ω2γgGρc+(4γg3)Mr*(r*)3}x

 

=

d2xdr*2+{4(ρ*P*)Mr*(r*)}1r*dxdr*+(ρ*P*){2πσc23γgαgMr*(r*)3}x.

Profile

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


Table 1: Radial Profile of Various Physical Variables

Variable

Throughout the Core
0ξξi

Throughout the Envelope
ηiηηs

Plotted Profiles

ξi=0.5

ξi=1.0

ξi=3.0

r*

(32π)1/2ξ

(μeμc)1θi2(2π)1/2η

 

ρ*

(1+13ξ2)5/2

(μeμc)θi5ϕ(η)

P*

(1+13ξ2)3

θi6[ϕ(η)]2

Mr*

(23π)1/2[ξ3(1+13ξ2)3/2]

(μeμc)2θi1(2π)1/2(η2dϕdη)

In order to obtain the various envelope profiles, it is necessary to evaluate ϕ(η) and its first derivative using the information presented in Step 6, above.


Therefore, throughout the core we have,

ρ*P*

=

(1+13ξ2)1/2;

Mr*r*

=

(23π)1/2[ξ3(1+13ξ2)3/2](2π3)1/21ξ=2ξ2(1+13ξ2)3/2;

and, throughout the envelope we have,

ρ*P*

=

(μeμc)θi1ϕ(η)1;

Mr*r*

=

(μeμc)2θi1(2π)1/2(η2dϕdη)[(μeμc)1θi2(2π)1/2η]1=2(μeμc)1θiη(η2dϕdη).

NOTE on 15 May 2019:  Prior to this date the last RHS expression had an incorrect exponent on η. It previously (incorrectly) read,

Mr*r*

=

2(μeμc)1θiη(η2dϕdη).

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,

0

=

x+r*x+𝒦x,

where,

x

=

dxdr*

      and      

x

=

d2xd(r*)2.

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

x

x+x2δr*

      and      

x

x+2xj+x(δr*)2.

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

x+2xj+x(δr*)2

=

r*[x+x2δr*]𝒦xj

x+2xj+x

=

δr*2r*[x+x](δr*)2𝒦xj

xj+1[1+(δr*2r*)]

=

[2(δr*)2𝒦]xj[1(δr*2r*)]xj1.

In what follows we will also find it useful to rewrite 𝒦 in the form,

𝒦(σc2γg)𝒦1αg𝒦2.

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

{4(ρ*P*)Mr*(r*)}

      ,      

𝒦1

2π3(ρ*P*)

      and      

𝒦2

(ρ*P*)Mr*(r*)3.

Special Handling at the Center

In order to kick-start the integration, we set the displacement function value to x1=1 at the center of the configuration (ξ1=0), then draw on the derived power-series expression to determine the value of the displacement function at the first radial grid line, ξ2=δξ, away from the center. Specifically, we set,

x2

=

x1[1(n+1)𝔉(δξ)260]n=5=x1[1𝔉(δξ)210].

Special Handling at the Interface

Integrating outward from the center, the general approach will work up through the determination of xj+1 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 x+ 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,

x'i|core

12δr*{x+xi1}

 

=

xi12δr*+12δr*{[2(δr*)2𝒦]xi[1(δr*2r*)]xi1}[1+(δr*2r*)]1

 

=

12δr*{[2(δr*)2𝒦]xi[1(δr*2r*)]xi1[1+(δr*2r*)]xi1}[1+(δr*2r*)]1

 

=

12δr*{[2(δr*)2𝒦]xi2xi1}[1+(δr*2r*)]1

Conversely, as viewed from the envelope, if we assume that we know xi and x'i, we can determine the amplitude, xi+1, at the first zone beyond the interface as follows:

x

xi+12δr*x'i|env

xi+1[1+(δr*2r*)]

=

[2(δr*)2𝒦]xi[1(δr*2r*)][xi+12δr*x'i|env]

xi+1[1+(δr*2r*)]+[1(δr*2r*)]xi+1

=

[2(δr*)2𝒦]xi+[1(δr*2r*)]2δr*x'i|env

xi+1

=

[112(δr*)2𝒦]xi+[1(δr*2r*)]δr*x'i|env

Splitting Analysis Into Separate Core and Envelope Components

Core:

Given that, 2π/3r*=ξ, lets multiply the LAWE through by 3/(2π). This gives,

0

=

d2xdξ2+{4(ρ*P*)Mr*(r*)}1ξdxdξ+32π(ρ*P*){2πσc23γgαgMr*(r*)3}x.

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

x+2xj+x(δξ)2

=

ξ[x+x2δξ][3𝒦2π]xj

x+2xj+x

=

δξ2ξ[x+x](δξ)2[3𝒦2π]xj

xj+1[1+(δξ2ξ)]

=

[2(δξ)2(3𝒦2π)]xj[1(δξ2ξ)]xj1.

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

[dxdξ]interface

=

12δξ{[2(δξ)2(3𝒦2π)]xi2xi1}[1+(δξ2ξ)]1.

Envelope:

Given that,

(μeμc)θi2(2π)1/2r*=η,

let's multiply the LAWE through by (2π)1θi4(μe/μc)2. This gives,

0

=

d2xdη2+{4(ρ*P*)Mr*(r*)}1ηdxdη+12πθi4(μeμc)2(ρ*P*){2πσc23γgαgMr*(r*)3}x.

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

x+2xj+x(δη)2

=

η[x+x2δη](μeμc)2[𝒦2πθi4]xj

x+2xj+x

=

δη2η[x+x](δη)2(μeμc)2[𝒦2πθi4]xj

xj+1[1+(δη2η)]

=

[2(δη)2(μeμc)2(𝒦2πθi4)]xj[1(δη2η)]xj1.

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,

xi+1

=

[112(δη)2(μeμc)2(𝒦2πθi4)]xi+[1(δη2η)]δη[dxdη]interface.

Interface

If we consider only cases where γe=γc, then at the interface we expect,

dlnxdlnr*

=

dlnxdlnξ=dlnxdlnη

r*dxdr*

=

ξdxdξ=ηdxdη

dxdr*

=

(2π3)1/2dxdξ=(μeμc)θi2(2π)1/2dxdη.

Switching at the interface from ξ to η therefore means that,

[dxdη]interface

=

13(μeμc)1θi2[dxdξ]interface.

If, however, we want to consider values for the adiabatic index that are different in the two regions, we have to follow the above-outlined guidelines, that is,

dlnxdlnη|i

=

3(γcγe1)+γcγe(dlnxdlnξ)i

[dxdη]i

=

xiηi{3(γcγe1)+γcγeξixi[dxdξ]i}

 

=

3xiηi(γcγe1)+γcγe13(μeμc)1θi2[dxdξ]i.

Eigenvectors for Marginally Unstable Models with (γc, γe) = (6/5, 2)

We now have the tools in hand to identify the eigenvectors — that is, various radial eigenfunctions and the corresponding eigenfrequency for each — associated with various modes of oscillation in (nc,ne)=(5,1) bipolytropes. Which models should we examine?


In our accompanying review of the bipolytrope stability analysis presented by Murphy & Fiedler (1983b), our primary objective was to show that we were able to match their results quantitatively. We therefore set μe/μc = 1 — the only μ-ratio that they considered — and picked values of the core-envelope interface radius, ξi, that were listed among their set of chosen models. For a fixed value of ξi, we integrated the relevant LAWE from the center toward the surface for many different eigenfrequency (σc2) guesses until an eigenfunction was found whose behavior at the surface matched with high precision the physically justified surface boundary condition.


With the above virial stability analysis in mind (see especially Figure 3), here we have chosen to focus on models that reside along six of the equilibrium sequences that have already been analytically identified, above — specifically, the sequences for which μe/μc = 1, ½, 0.345, ⅓, 0.309, and ¼ — and to examine oscillation modes under the assumption that,

γc=(nc+1)nc=65,       and       γe=(ne+1)ne=2.

Fundamental Modes

We decided to examine, first, whether any model along each sequence marks a transition from dynamically stable to dynamically unstable configurations. We accomplished this by setting σc2 = 0, then integrating the relevant LAWE from the center toward the surface for many different guesses of the core-envelope interface radius until an eigenfunction with no radial nodes — i.e., an eigenfunction associated with the fundamental mode of radial oscillation — was found whose behavior at the surface matched with high precision the physically desired surface boundary condition. We were successful in this endeavor. A marginally unstable model was identified on each of the six separate equilibrium sequences.

Equilibrium Properties of Marginally Unstable Models

Table 2 summarizes some of the equilibrium properties of these six models. For example, the second column of the table gives the value of the core-envelope interface radius, ξi, associated with each marginally unstable model. The table also lists: the value of the model's dimensionless radius, Rsurf*, the key structural parameters, q & ν, and the central-to-mean density associated with each model; and in each case the dimensionless thermal energy (𝔰) and dimensionless gravitational potential energy (𝔴) associated, separately, with the core and the envelope. Note that, once the pair of parameters, (μe/μc,ξi), has been specified, we can legitimately assign high-precision values to all of the other model parameters because they are analytically prescribed.

Table 2: Properties of Marginally Unstable Bipolytropes Having

(nc,ne)=(5,1) and (γc,γe)=(65,2)

Determined from Integration of the LAWE

μeμc ξi Rsurf* qrcoreRsurf νMcoreMtot ρcρ¯ 𝔰core 𝔴core 𝔰env 𝔴env
1 1.6686460157 2.139737 0.53885819 0.497747626 8.51704656 3.021916335 -3.356583022 1.47780476 -5.642859167
12 2.27925811317 5.146499 0.306021732 0.401776274 63.29514949 4.241287117 -6.074241035 4.284931508 -10.97819621
0.345 2.560146865247 9.554041 0.185160563 0.234302525 209.7739052 4.639705843 -7.125754184 11.72861751 -25.61089252
13 2.582007485476 10.120558 0.176288391 0.218241608 230.4125398 4.667042505 -7.200966267 13.15887139 -28.45086152
0.309 2.6274239687695 11.464303 0.158362807 0.184796947 279.0788798 4.722277318 -7.354156963 17.1374434 -36.36528446
14 2.7357711469398 15.895632 0.118924863 0.11071211 430.0444648 4.84592201 -7.70305421 37.84289623 -77.67458196


As was expected from our above discussion of virial equilibrium conditions, we found that to high precision for each of these equilibrium models,

(𝔴core+𝔴env)+2(𝔰core+𝔰env)

=

0.

However, contrary to expectations, in no case did we find that 𝔰core/𝔰env=5. That is to say, we found that none of the models lies on the (red-dashed) curve in the qν parameter space that separates stable from unstable models as defined by our above free-energy-based stability analysis. The left-hand panel of Figure 4 shows this (red-dashed) demarcation curve; for all intents and purposes, it is a reproduction of the right-hand panel of Figure 3, above — turning-point markers have been removed to minimize clutter, the equilibrium sequences have been labeled, and the horizontal axis has been extended to unity in order to include a longer portion of the μe/μc=1 sequence. The orange triangular markers that appear in the right-hand panel of Figure 4 pinpoint where each of the Table 2 "marginally unstable" models resides in this qν plane. Clearly, all six of the orange triangles lie well off of — and to the stable side of — the red-dashed demarcation curve. This discrepancy, which has resulted from our use of two separate approaches to stability analysis, will be discussed further and gratifyingly resolved, below.

Figure 4
Marginally unstable models

Eigenfunction Details

Here we examine some of the properties of the fundamental-mode eigenfunctions that we have found are associated with marginally unstable, (nc,ne)=(5,1) bipolytropes.

Figure 5

Example eigenvector

Consider the model on the μe/μc=1 sequence for which σc2=0; key properties of this specific equilibrium model are enumerated in the first row of numbers provided in Table 2, above. Figure 5 shows how our numerically derived, fundamental-mode eigenfunction, x=δr/r0, varies with the fractional radius over the entire range, 0r/R1. By prescription, the eigenfunction has a value of unity and a slope of zero at the center (r/R=0). Integrating the LAWE outward from the center, through the model's core (blue curve segment), x drops smoothly to the value xi=0.81437 at the interface (ξi=1.6686460157q=rcore/Rsurf=0.53885819). Our numerical integration of the LAWE showed that, at the interface, the logarithmic slope of the core (blue) segment of the eigenfunction is,

{dlnxdlnr|i}core={dlnxdlnξ|i}core

=

0.455872.

Next, following the above discussion of matching conditions at the interface, we determined that, from the perspective of the envelope, the slope of the eigenfunction at the interface must therefore be,

{dlnxdlnr|i}env={dlnxdlnη|i}env

=

3(γcγe1)+γcγe{dlnxdlnξ|i}core=1.47352.

Adopting this "env" slope along with the amplitude, xi=0.81437, as the appropriate interface boundary conditions, we integrated the LAWE from the interface to the surface, obtaining the green-colored segment of the eigenfunction that is shown in Figure 5. The amplitude continued to steadily decrease, reaching a value of xs=0.38203, at the model's surface (r/R=1). At the surface, this envelope (green) segment of the eigenfunction exhibits a logarithmic slope that matches to eight significant digits the value that is expected from astrophysical arguments for this marginally unstable (σc2=0) model, namely,

dlnxdlnη|s=[(ρcρ¯)σc202γe(34γe)]=1.


Key Reminder: We were able to find an eigenfunction whose surface boundary condition matched the desired value — in this particular case, a logarithmic slope of negative one — to this high level of precision only by iterating many times and, at each step, fine-tuning our choice of the equilibrium model's radial interface location, ξi before performing a numerical integration of the LAWE.


The discontinuous jump that occurs in the slope of the eigenfunction at the interface results from our assumption that the effective adiabatic index of material in the core (γc=6/5) is different from the effective adiabatic index of the envelope material (γe=2). In an effort to emphasize and more clearly illustrate the behavior of this fundamental-mode eigenfunction as it crosses the core/envelope interface, we have added a pair of dashed line segments to the Figure 5 plot. The red-dashed line segment touches, and is tangent to, the blue segment of the eigenfunction at the location of the core/envelope interface; it has a slope,

dxd(r/R)|i=xi(ri/R){dlnxdlnr|i}core

=

0.455872(0.814370.53885819)=0.68895.

On the other hand, the purple-dashed line segment touches, and is tangent to, the green segment of the eigenfunction at the location of the core/envelope interface; it has a slope,

dxd(r/R)|i=xi(ri/R){dlnxdlnr|i}env

=

1.47352(0.814370.53885819)=2.22691.

For comparison purposes, the eigenfunction shown in Figure 5 has been presented again in Figure 6, along with several other of our numerically derived eigenfunctions, but in Figure 6 the plotted amplitude has been renormalized to give a surface value — rather than a central value — of unity.


In Figure 6 we show the behavior of the fundamental-mode eigenfunction for each of the marginally unstable models identified in Table 2. In the top figure panel, each curve shows — on a linear-linear plot — how the amplitude varies with radius; in the bottom figure panel, the amplitude is plotted on a logarithmic scale. On each curve, the black plus sign marks the radial location of the core-envelope interface; in the bottom panel, these markers are accompanied by the values of ξi that are associated with each corresponding model (see also the second column of Table 2). Each eigenfunction has been normalized such that the surface amplitude is unity. In the top panel, the value of the central amplitude of the eigenfunction that results from this normalization is recorded near the point where each eigenfunction touches the vertical axis. (In each case, the value provided on the plot is simply the inverse of the value of xs given in Table 3, below.)

file = Dropbox/WorkFolder/Wiki edits/BiPolytrope/LinearPerturbation/FaulknerBipolytrope1.xlsx --- worksheet = Mode0Ensemble
file = Dropbox/WorkFolder/Wiki edits/BiPolytrope/LinearPerturbation/FaulknerBipolytrope1.xlsx --- worksheet = Mode0Ensemble
Figure 6: Eigenfunctions Associated with the Fundamental-Mode of Radial Oscillation

in Marginally Unstable Models having Various μe/μc

Eigenfunctions for Marginally Unstable Models

Notice that, especially as they approach the surface, the "envelope" segments of these six marginally unstable eigenfunction appear to merge into the same curve, irrespective of their value of the ratio of mean molecular weights. Note as well that the discontinuous jump that occurs in the slope of each eigenfunction at the radial location of the core/envelope interface — resulting from our choice to adopt a different adiabatic index, γg, in the core from the one in the envelope — becomes less and less noticeable for smaller and smaller values of the ratio of mean molecular weights.

Is There an Analytic Expression for the Eigenfunction?

After noticing that, in Figure 6, the envelope segments of all of the marginally unstable eigenfunctions merge into the same curve, we began to wonder whether a single expression — and, even better, an analytically defined expression — would perfectly describe the eigenfunction. We had reason to believe that this might actually be possible because, in pressure-truncated polytropic configurations, we have derived analytic expressions for the marginally unstable, fundamental-mode eigenfunctions of both n=5 and n=1 systems.

Very quickly, we convinced ourselves that a parabolic function does indeed perfectly match the "core" segment of each displayed eigenfunction. Specifically, throughout the core (0ξξi),

xP|core

1ξ215

dxPdξ|core

2ξ15

dlnxPdlnξ|core

2ξ215[(15ξ2)15]1=2ξ2(15ξ2).


The envelope segment posed a much greater challenge. In the context of our discussion of Radial Oscillations of n = 1 Polytropic Spheres, and in an accompanying Ramblings Appendix chapter we have detailed some trial derivations that are mostly blind alleyways. Twice — once in January, 2019 and again (independently) in April 2019 — we have analytically demonstrated that the following appears to work for the envelope:

Given that the Structural Properties of the envelope are described by the Lane-Emden function,

ϕ

=

a0[sin(ηb0)η]

Qdlnϕdlnη

=

[1ηcot(ηb0)],

the relevant LAWE is satisfied by the fractional displacement function,

xP

=

3c0Qη2,

where, c0 is an arbitrary scale factor.

Note that,

dxPdη

=

3c0η2[ηcot(ηb0)+ηcot2(ηb0)]6c0Qη3

dlnxPdlnη=ηxPdxPdη

=

η33c0Q{3c0η2[ηcot(ηb0)+ηcot2(ηb0)]6c0Qη3}

 

=

1Q{[η2ηcot(ηb0)+η2cot2(ηb0)]2Q}

 

=

1Q[η2(1Q)+(1Q)2]2

 

=

1Q[η2Q+Q2]2

 

=

η2Q+Q3.

But, as far as we have been able to determine (as of 16 April 2019), this analytic displacement function does not match the displacement function that has been generated through numerical integration of the LAWE (see the light-green segment of the eigenfunction displayed above in Figure 5). It remains unclear whether (a) the numerical integration is at fault, (b) we are imposing an incorrect slope at the core-envelope interface, or ( c) we are misinterpreting how to compare the two separately derived (one, numerical, and the other, analytic) envelope eigenfunctions.

Other Modes

file = Dropbox/WorkFolder/Wiki edits/BiPolytrope/LinearPerturbation/FaulknerBipolytrope1.xlsx --- worksheet = mu100plot
file = Dropbox/WorkFolder/Wiki edits/BiPolytrope/LinearPerturbation/FaulknerBipolytrope1.xlsx --- worksheet = mu100plot
Our Determinations for Marginally Unstable Model Having μe/μc=1
 

[ξi=1.6686460157]

Mode σc2 Ω2σc22(ρcρ¯) xi dlnxdlnr*|i xsurf dlnxdlnr*|surf rR|1 1MrMtot|1 rR|2 1MrMtot|2 rR|3 1MrMtot|3
core env expected measured
1
(Fundamental)
0.00 0.00 +0.81437470 -0.455872 -1.473523 +0.3820 -1 -0.999999992 n/a n/a n/a n/a n/a n/a
2 2.51513333 10.7107538 0.20482050 -7.09124 -5.4547441 - 0.9962 4.355376917 4.35537692 0.64133 0.3502 n/a n/a n/a n/a
3 5.72371888 24.3745901 -0.14269277 +8.046019 +3.627611 +0.9308 11.18729505 11.18729506 0.4837 0.5864 0.842 0.0854 n/a n/a
4 10.3458476 44.0622916 -0.20845197 -0.6949966 -1.61699793 -1.1443 21.03114578 21.03114577 0.3939 0.7154 0.6902 0.2777 0.9115 0.0284

Our determination of eigenvector for mu_ratio = 1  Our determination of multiple eigenvectors for mu_ratio = 1

Reconciliation

Figure 7: Conflicting Instability Regions
Marginally unstable models

Figure 7, shown here on the right, is identical to the right-hand panel of Figure 4, as displayed above. From a standard global, free-energy analysis — such as the one summarized above — we have determined that the red-dashed curve shown in the right panel of Figure 6 divides the qν plane into dynamically stable (below and to the right) and unstable (above and to the left) regions.


Variational Principle

Setup

Let's follow the guidelines of the variational principle. Instead of starting with the form of the LAWE given above, namely,

Adiabatic Wave (or Radial Pulsation) Equation

d2xdr02+[4r0(g0ρ0P0)]dxdr0+(ρ0γgP0)[ω2+(43γg)g0r0]x=0

we will start with a form that is more amenable to the variational principle, namely,

0

=

ddr0[r04γP0dxdr0]+[ω2ρ0r04+(3γ4)r03dP0dr0]x.

ASIDE: Let's show that these two expressions are equivalent. Remembering that,

dP0dr0

=

g0ρ0=(GMrr02)ρ0,

the second expression becomes,

0

=

r04γP0d2xdr02+γdxdr0[4r03P0r04g0ρ0]+[ω2ρ0r04(3γ4)r03g0ρ0]x

 

=

r04γP0{d2xdr02+[4r0g0ρ0P0]dxdr0+ρ0γP0[ω2(3γ4)g0r0]x}.

Hence, we must multiply the first expression through by r04γP0 in order to obtain the second expression.


From above, we realize that multiplying the second expression through by (Kc/G)ρc4/5 gives,

0

=

r04γP0{d2xdr*2+[4(ρ*P*)Mr*(r*)]1r*dxdr*+(ρ*P*)[2πσc23γgαgMr*(r*)3]x}

 

=

(r*)4[Kc1/2/(G1/2ρc2/5)]4γP*[Kcρc6/5]{d2xdr*2+[4(ρ*P*)Mr*(r*)]1r*dxdr*+(ρ*P*)[2πσc23γgαgMr*(r*)3]x}

 

=

[Kc3G2ρc2/5](r*)4γP*{d2xdr*2+[4(ρ*P*)Mr*(r*)]1r*dxdr*+(ρ*P*)[2πσc23γgαgMr*(r*)3]x}.

That is, multiplying the second expression through by, (Kc/G)ρc4/5G2ρc2/5/Kc3=G/(Kc2ρc2/5) , should give a desirable, totally dimensionless version of the LAWE. Remembering that,

ρ*

ρ0ρc

;    

r*

r0[Kc1/2/(G1/2ρc2/5)]

P*

P0Kcρc6/5

;    

Mr*

Mr[Kc3/2/(G3/2ρc1/5)]

dP*dr*

=

Mr*ρ*(r*)2

;    

Enorm

=

[Kc5G3]1/2

let's try it.

0

=

GKc2ρc2/5{ddr0[r04γP0dxdr0]+[ω2ρ0r04+(3γ4)r03dP0dr0]x}

 

=

GKc2ρc2/5{(Kc2ρc2/5G)ddr*[(r*)4γP*dxdr*]+[ω2(Kc2G2ρc8/5)ρcρ*(r*)4+(3γ4)(KcGρc4/5)Kcρc6/5(r*)3dP*dr*]x}

 

=

ddr*[(r*)4γP*dxdr*]+[(ω2Gρc)ρ*(r*)4+(3γ4)(r*)3dP*dr*]x.

Now, guided by the accompanying summary, if we multiply through by 4πxdr* and integrate over the entire volume, we obtain the governing variational relation, namely,

0

=

0rcore*4πxd[(r*)4γcP*dxdr*]+0rcore*[(ω2Gρc)ρ*(r*)4+(3γc4)(r*)3dP*dr*]4πx2dr*

 

 

+rcore*R*4πxd[(r*)4γeP*dxdr*]+rcore*R*[(ω2Gρc)ρ*(r*)4+(3γe4)(r*)3dP*dr*]4πx2dr*

 

=

[4πx(r*)4γcP*dxdr*]0rcore*0rcore*4π[(r*)4γcP*(dxdr*)2]dr*0rcore*[(3γc4)Mr*ρ*]4πx2r*dr*

 

 

+[4πx(r*)4γeP*dxdr*]rcore*R*rcore*R*4π[(r*)4γeP*(dxdr*)2]dr*rcore*R*[(3γe4)Mr*ρ*]4πx2r*dr*+0R*[(ω2Gρc)ρ*(r*)4]4πx2dr*

 

=

0rcore*x2(dlnxdlnr*)2γc4π(r*)2P*dr*+0rcore*(3γc4)x2(Mr*r*)4πρ*(r*)2dr*[4πx2(r*)3γcP*(dlnxdlnr*)]0rcore*

 

 

rcore*R*x2(dlnxdlnr*)2γe4π(r*)2P*dr*+rcore*R*(3γe4)x2(Mr*r*)4πρ*(r*)2dr*[4πx2(r*)3γeP*(dlnxdlnr*)]rcore*R*+0R*4π(ω2Gρc)ρ*(r*)4x2dr*.


Energy Normalization:

(γ1)dUint

=

4πr2Pdr

 

=

4π[Kc1/2G1/2ρc2/5]3Kcρc6/5(r*)2P*dr*

 

=

4π[Kc5G3]1/2(r*)2P*dr*

Enorm

[Kc5G3]1/2

dWgrav

=

(GMrr)4πr2ρdr

 

=

(GMr*r*)4π(r*)2ρ*dr*[KcGρc4/5]ρc[Kc3/2G3/2ρc1/5]

 

=

(Mr*r*)4π(r*)2ρ*dr*Enorm

Hence, the dimensionless governing variational relation becomes,

0R*(2π3)σc2(r*)2x2dMr*

=

0rcore*γc(γc1)x2(dlnxdlnr*)2dUint*0rcore*(3γc4)x2dWgrav*+[4πx2(r*)3γcP*(dlnxdlnr*)]0rcore*

 

 

+rcore*R*γe(γe1)x2(dlnxdlnr*)2dUint*rcore*R*(3γe4)x2dWgrav*+[4πx2(r*)3γeP*(dlnxdlnr*)]rcore*R*

 

=

2γc30rcore*x2(dlnxdlnr*)2dStherm*(3γc4)0rcore*x2dWgrav*+4πxi2(rcore*)3γcPi*{dlnxdlnr*|i}core

 

 

+2γe3rcore*R*x2(dlnxdlnr*)2dStherm*(3γe4)rcore*R*x2dWgrav*4πxi2(rcore*)3γePi*{dlnxdlnr*|i}env

 

=

2γc30rcore*x2(dlnxdlnr*)2dStherm*(3γc4)0rcore*x2dWgrav*

 

 

+2γe3rcore*R*x2(dlnxdlnr*)2dStherm*(3γe4)rcore*R*x2dWgrav*4πxi2(rcore*)3Pi*[γc{dlnxdlnr*|i}coreγe{dlnxdlnr*|i}env]

 

=

2γc30rcore*x2(dlnxdlnr*)2dStherm*(3γc4)0rcore*x2dWgrav*

 

 

+2γe3rcore*R*x2(dlnxdlnr*)2dStherm*(3γe4)rcore*R*x2dWgrav*4πxi2(rcore*)3Pi*[3(γeγc)],

where,

dMr*

=

4π(r*)2ρ*dr*,

dUint*

=

1(γ1)[4π(r*)2P*dr*]=[23(γ1)]dStherm*,

dWgrav*

=

(Mr*r*)4π(r*)2ρ*dr*.

Or, for inclusion in our accompanying Tabular Overview,

(2π3)σc20R*(xr*)2dMr*

=

γc(γc1)0rcore*x2(dlnxdlnr*)2dUint*(3γc4)0rcore*x2dWgrav*

 

 

+γe(γe1)rcore*R*x2(dlnxdlnr*)2dUint*(3γe4)rcore*R*x2dWgrav*

 

 

+32(γcγe)xi2Pi*Vcore*.

Implementation

Table 3

μeμc ξi Core Envelope Virial
Integrals over

dStherm*
Integrals over

dWgrav*
Integrals over

dStherm*
Integrals over

dWgrav*
Numerical

[2𝔰tot|𝔴tot|1]
Analytic
𝔰core
Numerical
𝔰core
TERM1 Analytic
𝔴core
Numerical
𝔴core
TERM2 Analytic
𝔰env
Numerical
𝔰env
TERM3 Analytic
𝔴env
Numerical
𝔴env
TERM4
1 1.6686460157 3.021916335 3.021921 0.116389175 -3.356583022 -3.35666 -2.649752079 1.47780476 1.47791 1.0720821 -5.642859167 -5.642820 -1.91142893 0.000020
12 2.27925811317 4.241287117 4.241410819 0.440878529 -6.074241035 -6.074317546 -4.150731169 4.284931508 4.28547195 1.44651932 -10.97819621 -10.97847622 -0.92598634 0.000057
0.345 2.560146865247 4.639705843 4.6399114 0.6794857 -7.125754184 -7.125854025 -4.5487829 11.72861751 11.730381 1.51410084 -25.61089252 -25.6115597 -0.4496513 0.000097
13 2.582007485476 4.667042505 4.667254935 0.700414598 -7.200966267 -7.201068684 -4.57274936 13.15887139 13.1608467 1.51408246 -28.45086152 -28.45153761 -0.4170461 0.000101
0.309 2.6274239687695 4.722277318 4.722504339 0.744964507 -7.354156963 -7.3542507 -4.61961058 17.1374434 17.1399773 1.51055838 -36.36528446 -36.36591543 -0.3524855 0.000110
14 2.7357711469398 4.84592201 4.846185027 0.857001395 -7.70305421 -7.703178009 -4.7163542 37.84289623 37.8479208 1.47966673 -77.67458196 -77.67408155 -0.2194152 0.000128

𝔰core

0rcore*dStherm*

          ;          

TERM1

0rcore*x2(dlnxdlnr*)2dStherm*

𝔴core

0rcore*dWgrav*

          ;          

TERM2

0rcore*x2dWgrav*

𝔰env

rcore*R*dStherm*

          ;          

TERM3

rcore*R*x2(dlnxdlnr*)2dStherm*

𝔴env

rcore*R*dWgrav*

          ;          

TERM4

rcore*R*x2dWgrav*

NOTE: In all integrals, the fractional radial-displacement function, x, has been normalized
to unity at the center of the spherical model.

[σc2]μe/μc=1VP

=

32π𝐓𝐄𝐑𝐌𝟓{2γc3[𝐓𝐄𝐑𝐌𝟏](3γc4)[𝐓𝐄𝐑𝐌𝟐]+2γe3[𝐓𝐄𝐑𝐌𝟑](3γe4)[𝐓𝐄𝐑𝐌𝟒]4πxi2(rcore*)3Pi*[3(γeγc)]}

 

=

0.036312577×[0.0931113401.059900832+1.4294428(3.82285786)1.782200484×(2.4)]=0.036312577×[0.00819]

 

=

0.0002973.

[σc2]μe/μc=1/2VP

=

0.020704186×[0.3527028231.660292468+1.928692426(1.85197268)1.029029184×(2.4)]=0.020704186×[0.003405]

 

=

0.0000705.

[σc2]μe/μc=0.345VP

=

0.019377703×[0.543588561.81951316+2.01880112(0.8993026)0.68747441×(2.4)]=0.019377703×[0.003585]

 

=

0.0000695.

[σc2]μe/μc=1/3VP

=

0.019467503×[0.5603316781.81099744+2.018776613(0.8340922)0.658354814×(2.4)]=0.019467503×[1.6022030521.580051555]

 

=

0.00043123.

[σc2]μe/μc=0.309VP

=

0.019781992×[0.5959716061.847844232+2.01407784(0.704971)0.60905486×(2.4)]=0.019781992×[1.4671762141.461731665]

 

=

0.0001077.

[σc2]μe/μc=1/4VP

=

0.021344909×[0.6856011161.88654168+1.972888973(0.4388304)0.499262049×(2.4)]=0.021344909×[1.2107788091.198228918]

 

=

0.000267876.


Table 4

μeμc ξi Analytic
𝔪tot
Numerical
𝔪tot
TERM5 xi rcore* Pi* {dlnxdlnr*|i}core {dlnxdlnr*|i}env [σc2]VP
1 1.6686460157 4.818155928 4.818145 13.14874521 0.814374698 1.153014872 0.139506172 +0.455871977 +1.47352 0.0002973
12 2.27925811317 9.020985415 9.021084268 23.0612702 0.653665497 1.574940686 0.049058481 +1.059668912 +1.835801347 0.0000705
0.345 2.560146865247 17.41399388 17.4141672 24.63990825 0.563043202 1.769031527 0.030957085 +1.552125296 +2.131275177 0.0000695
13 2.582007485476 18.8449906 18.8451614 24.52624951 0.555549156 1.78413696 0.029889634 +1.600041467 +2.16002488 0.0004312
0.309 2.6274239687695 22.61541791 22.6155686 24.13633675 0.539776219 1.815519219 0.027798189 +1.705239188 +2.223143513 0.0001077
14 2.7357711469398 39.12088278 39.1208075 22.36902623 0.501037082 1.890385851 0.023427588 +1.991720516 +2.39503231 0.0002679

𝔪tot

0R*dMr*

          ;          

Isphere*

0R*(r*)2dMr*

          ;          

TERM5

0R*(r*)2x2dMr*

NOTE: In the TERM5 integral (as elsewhere), the fractional radial-displacement function, x,
has been normalized to unity at the center of the spherical model.

Revised Free-Energy Analysis

If we set x = constant in the variational principle relation, we have,

(2π3)σc20R*(r*)2dMr*

=

(3γc4)0rcore*dWgrav*(3γe4)rcore*R*dWgrav*4π(rcore*)3Pi*[3(γeγc)]

(2π3)σc2Isphere*

=

(43γc)𝔴core+(43γe)𝔴env+Pi*Vcore*[32(γcγe)],

Notice the similarity between this last expression and the pair of expressions — numbered and — that arise in the context of pressure-truncated polytropes.

See Also


Tiled Menu

Appendices: | VisTrailsEquations | VisTrailsVariables | References | Ramblings | VisTrailsImages | myphys.lsu | ADS |