Editing SSC/Structure/BiPolytropes/51RenormaizePart2 (section)

====Treatment of the Envelope====

<table border="1" align="center" cellpadding="8">
<tr>
  <td align="center"><b>Table C2</b></td>
</tr>
<tr>
  <td align="center">[[File:ModelC Table2.png|1000px|Equilibrium Structure of Model C Envelope]]</td>
</tr>
</table>

<font color="darkgreen"><b>STEP11:</b></font>&nbsp; &nbsp; Divide the envelope into <math>(N_e+1)</math> grid lines &#8212; that is, into <math>N_e</math> radial zones &#8212; associating the first "grid line" with the location of the core/envelope interface and the last grid line with the radial location of the surface of the bipolytropic configuration; in <b>Table C2</b>, we have set <math>N_e = 20</math>.   Choosing <math>\nu \le \tilde{M}_r \le 1</math> as the principal Lagrangian coordinate, and using the available analytic expressions, assign values to the following physical quantities at each grid line:
<ul>
  <li>
Mass (see 1<sup>st</sup> column titled <font color="darkgreen">tilde M_r</font> in <b>Table C2</b>): &nbsp; &nbsp; Set <math>(\Delta m)_e = (1-\nu)/(N_e)</math>; then, for <math>n = 1 ~\mathrm{thru}~ (N_e + 1)</math>, set <math>\tilde{M}_r = \nu + (n - 1)(\Delta m)_e \, .</math>
  </li>
  <li>
First and last polytropic radial coordinates:  In the [[#VariableProfiles|above summary of various physical variable profiles]], we have provided expressions for the envelope's polytropic radial coordinate at the core/envelope interface <math>(\eta_i)</math> and at the surface <math>(\eta_s)</math>.  Evaluating these expressions for <b>Model C</b>, we find,
<table border="0" align="center" cellpadding="5">

<tr>
  <td align="right"><math>\eta_i</math></td>
  <td align="center"><math>=</math></td>
  <td align="left"><math>\biggl(\frac{\mu_e}{\mu_c}\biggr)~ \sqrt{3}  \xi_i \biggl(1+\frac{1}{3}\xi_i^2\biggr)^{-1} = 0.42286 \, ,</math></td>
</tr>

<tr>
  <td align="right"><math>\eta_s</math></td>
  <td align="center"><math>=</math></td>
  <td align="left"><math> 
\frac{\pi}{2} + \eta_i + \tan^{-1}\biggl\{
\frac{\xi_i}{\sqrt{3}} 
\biggl[ \biggl(\frac{\mu_e}{\mu_c}\biggr)^{-1}\frac{1}{\xi_i^2}\biggl(1+\frac{1}{3}\xi_i^2\biggr) - 1\biggr] 
\biggr\} = 2.67311  
\, ,</math>
</td>
</tr>
</table>
These two values appear, respectively, as the first and last numbers in the (2<sup>nd</sup>) column of <b>Table C2</b>, titled <font color="darkgreen">eta_guess</font>.
  </li>
  <li>
All other polytropic radial coordinates (see the "<math>N_e - 1</math>" numerical values that are highlighted in yellow in the column titled <font color="darkgreen">eta_guess</font> in <b>Table C2</b>): &nbsp; &nbsp; Next, we need to determine what value of <math>\eta</math> is associated with each Lagrangian fluid marker (''i.e.,'' each gridline) that lies between the interface and the surface.  Referring again to the  [[#VariableProfiles|above summary of various physical variable profiles]], we can determine analytically the value of <math>\tilde{M}_r</math> that is associated with any selected value of <math>\eta</math>, via the relation,

<table border="0" align="center" cellpadding="8">

<tr>
  <td align="right">
<math>
\tilde{M}_r(\eta)
</math>
  </td>
  <td align="center"><math>=</math></td>
  <td align="left">
<math>\mathcal{m}_\mathrm{surf}^{-1}~ \theta^{-1}_i \biggl( \frac{2}{\pi} \biggr)^{1/2} 
A\biggl[ \sin(\eta-B) - \eta\cos(\eta-B) \biggr] \, ,
</math>
  </td>
</tr>
</table>
but what we really need is to invert this relation to obtain <math>\eta(\tilde{M}_r)</math>.  Unfortunately, unlike our experience with the core, for which we were able to invert the <math>\tilde{M}_r(\xi)</math> relation to obtain an analytic prescription for <math>\xi(\tilde{M}_r)</math> &#8212; see <font color="red"><b>STEP1</b></font> above &#8212; here, we are unable to analytically invert our <math>\tilde{M}_r(\eta)</math> expression.

Instead, for each gridline we ''guessed'' the value of <math>\eta</math> which, when plugged into the <math>\tilde{M}_r(\eta)</math> expression, would give the value of <math>\tilde{M}_r</math> that was already assigned to that gridline; by trial-and-error, we revised our guess for <math>\eta</math> until the desired value of the normalized mass (see the column titled <font color="darkgreen">tilde M_r</font>) was obtained (usually, to 6 or 7 digit accuracy).  <b>Table C2</b> provides the following record:  At each gridline, our final/best iterative "guess" for <math>\eta</math> is highlighted in yellow in the (2<sup>nd</sup>) column titled <font color="darkgreen">eta_guess</font>; the value of <math>\tilde{M}_r</math> that is obtained by plugging this "best guess" value of <math>\eta</math> into the <math>\tilde{M}_r(\eta)</math> relation is recorded in the (3<sup>rd</sup>) column titled <font color="darkgreen">M_r from eta</font>; and the (4<sup>th</sup>) column titled <font color="darkgreen">error</font> shows the fractional difference between this value and the desired value &#8212; as was our goal, all measured errors are zero, to at least six significant digits.
  </li>
  <li>
Given the "best guess" value of <math>\eta</math> at each gridline, determine the associated values of <math>\tilde{r}(\eta), \tilde{\rho}(\eta), \tilde{P}(\eta)</math> &#8212; see the (5<sup>th</sup>, 6<sup>th</sup>, and 7<sup>th</sup>) columns in <b>Table C2</b> titled <font color="darkgreen">tilde r, tilde rho, tilde P</font> &#8212; using the appropriate analytic expressions for the ''Envelope'' [[#VariableProfiles|as provided above]].  For example, at the 20<sup>th</sup> gridline (associated with the first gridline just inside the surface) where our "best guess" is <math>\eta = 2.3765850</math>, we find, <math>\tilde{r} = 0.021014</math>, <math>\tilde\rho = 5.29458 \times 10^{3}</math>, and <math>\tilde{P} = 1.37706 \times 10^{4}</math>.
  </li>
</ul>


<font color="darkgreen"><b>STEP12:</b></font>&nbsp; &nbsp; As was done above, in <font color="red"><b>STEP2</b></font>, for the core, determine here for the envelope the value of <math>\tilde{M}_r/(4\pi \tilde{r}^4)</math> at each gridline in the (initially) equilibrium model; see the (8<sup>th</sup>) column of <b>Table C2</b> titled <font color="darkgreen">M/(r pi r^4)</font>.


<font color="darkgreen"><b>STEP13:</b></font>&nbsp; &nbsp; As was done in <font color="red"><b>STEP3</b></font> for the core &#8212; but, here, for the envelope &#8212; we choose to evaluate <math>\tilde{P}</math> midway between gridlines and to evaluate the pressure gradient via the (2<sup>nd</sup>-order accurate) expression,

<table border="0" align="center" cellpadding="8">
<tr>
  <td align="right"><math>\biggl[\frac{d\tilde{P}}{d\tilde{M}_r}\biggr]_n</math></td>
  <td align="center"><math>\approx</math></td>
  <td align="left">
<math>
\frac{\tilde{P}_{n + 1/2} - \tilde{P}_{n - 1/2}}{(\Delta m)_e} \, .
</math>
  </td>
</tr>
</table>
Notice that the denominator on the right-hand-side of this expression is <math>(\Delta m)_e</math> (for the envelope) rather than <math>(\Delta m)_c</math> (for the core).  And, as with our treatment of the core, the ''difference'' between the pair of discrete mid-zone values of the pressure that appears in the numerator of this right-hand-side term straddles the discrete grid in such a way that the left-hand-side pressure ''gradient'' is centered on the n<sup>th</sup> gridline.  

At the center of each ''envelope'' grid zone, we can use the exact relationship between the normalized pressure and normalized density &#8212; namely,

<table border="0" align="center" cellpadding="8">

<tr>
  <td align="right"><math>\frac{\tilde{P}_{n+1/2}}{ [\tilde\rho_{n+1/2}]^{2} }</math></td>
  <td align="center"><math>=</math></td>
  <td align="left">
<math>
\mathcal{m}_\mathrm{surf}^6 \biggl(\frac{\mu_e}{\mu_c}\biggr)^{-12} \theta^{6}_i \phi^{2}
\biggl\{
\mathcal{m}_\mathrm{surf}^5 \biggl( \frac{\mu_e}{\mu_c} \biggr)^{-9} \theta^{5}_i \phi
\biggr\}^{-2}
</math>
  </td>
</tr>

<tr>
  <td align="right">&nbsp;</td>
  <td align="center"><math>=</math></td>
  <td align="left">
<math>
\biggl(\frac{\mu_e}{\mu_c}\biggr)^{6} (\mathcal{m}_\mathrm{surf}\theta_i)^{-4} 
</math>
  </td>
</tr>

<tr>
  <td align="right"><math>\Rightarrow ~~~ \tilde{P}_{n+1/2}</math></td>
  <td align="center"><math>=</math></td>
  <td align="left">
<math>
\biggl(\frac{\mu_e}{\mu_c}\biggr)^{6} (\mathcal{m}_\mathrm{surf}\theta_i)^{-4} \biggl[\tilde\rho_{n+1/2} \biggr]^{2}
</math>
  </td>
</tr>
</table>
&#8212; to determine <math>\tilde{P}_{n+1/2}</math> from <math>\tilde{\rho}_{n+1/2}</math> or, after inversion, to determine <math>\tilde{\rho}_{n+1/2}</math> from <math>\tilde{P}_{n+1/2}</math> for all <math>1 \le n \le N_e</math>.

<font color="darkgreen"><b>STEP14:</b></font>&nbsp; &nbsp; Akin to our treatment of the core (<font color="red"><b>STEP4</b></font>), the mass contained within every spherical shell of our discrete model of the ''envelope'' is <math>(\Delta m)_e</math> and,
in an effort to satisfy the continuity equation, throughout our discrete model we will relate the gas density of each spherical shell to the bounding radii of that shell via the expression,
<table border="0" align="center" cellpadding="8">
<tr>
  <td align="right"><math>\rho_{n+1/2} = \frac{(\Delta m)_e}{(\Delta ~ \mathrm{Vol})_{n+1/2}}</math></td>
  <td align="center"><math>=</math></td>
  <td align="left">
<math>
\frac{3(\Delta m)_e}{4\pi}\biggl[ \tilde{r}^3_{n+1} - \tilde{r}^3_{n}\biggr]^{-1} \, .
</math>
  </td>
</tr>
</table>
<ul>
  <li>
Values of the normalized density computed in this manner have been recorded in the (9<sup>th</sup>) column titled <font color="darkgreen">rho_FD</font> of <b>Table C2</b>.  For example, in the first shell just inside the surface <math>(n = 20)</math>, we find <math>\tilde\rho_{n+1/2} = 2.458565 \times 10^{3}</math>.
  </li>
  <li>
We have determined the value of the normalized pressure that corresponds to each of these "Finite Difference" values of the density using the algebraic relation presented just above; their values have been recorded in the (10<sup>th</sup>) column titled <font color="darkgreen">P_FD</font> of <b>Table C2</b>.  For example, in the first shell just inside the surface <math>(n = 20)</math>, we find <math>\tilde{P}_{n+1/2} = 4.91237\times 10^{-4} [\tilde\rho_{n+1/2}]^{2} = 2.96930 \times 10^{3}</math>.
  </li>
  <li>
From our determination of <math>\tilde{P}_{n+1/2}</math> throughout the envelope, values of the normalized pressure ''gradient'' have been computed in the manner described above in <font color="darkgreen"><b>STEP13</b></font>, and have been recorded in the (11<sup>th</sup>) column titled <font color="darkgreen">(dP/dM)_FD</font> of <b>Table C2</b>.  For example, at the <math>n = 20</math> gridline, we find <math>[d\tilde{P}/d\tilde{M}_r]_{n} = -4.549 \times 10^{5}</math>.
  </li>
</ul>


<font color="darkgreen"><b>STEP15:</b></font>&nbsp; &nbsp; Throughout the envelope, compare evaluation of the finite-difference representation of the (absolute value of the) pressure gradient to an evaluation of the unperturbed and analytically prescribed profile of the quantity, <math>\tilde{M}_r/(4\pi \tilde{r}^4)</math>.  The right-hand segment of <b>Figure C1</b> provides such a comparison; actually, for reasons that will become clear later, we have multiplied both quantities by <math>4\pi \tilde{r}^2</math> before plotting.



<table border="1" align="center" cellpadding="8" width="80%">
<tr>
  <td align="center"><b>Figure C3</b></td>
</tr>
<tr>
  <td align="center">[[File:ModelC Fig3.png|750px|Finite-Difference Structure of Model C with 41 gridlines]]</td>
</tr>

<tr>
  <td align="left">
Same as <b>Figure C2</b>, except higher resolution with <math>N_e = N_c = 40</math>.
  </td>
</tr>
</table>
varies with mass shell throughout the equilibrium model, after setting <math>\mathrm{amp}~= -50</math>; this data has been recorded in the column titled "<font color="red">Alternate_FD</font> error" in <b>Table C1</b>.
  </td>
</tr>
</table>