Dimensionless Formulation

To improve numerical stability, GYRE solves the separated equations and boundary conditions by recasting them into a dimensionless form that traces its roots back to Dziembowski (1971).

Variables

The independent variable is the fractional radius \(x \equiv r/\Rstar\) (with \(\Rstar\) the stellar radius), and the dependent variables \(\{y_{1},y_{2},\ldots,y_{6}\}\) are

(12)\[\begin{split}\begin{aligned} y_{1} &= x^{2 - \ell}\, \frac{\txir}{r}, \\ y_{2} &= x^{2-\ell}\, \frac{\tP'}{\rho g r}, \\ y_{3} &= x^{2-\ell}\, \frac{\tPhi'}{gr}, \\ y_{4} &= x^{2-\ell}\, \frac{1}{g} \deriv{\tPhi'}{r}, \\ y_{5} &= x^{2-\ell}\, \frac{\delta \tS}{\cP}, \\ y_{6} &= x^{-1-\ell}\, \frac{\delta \tLrad}{\Lstar}, \end{aligned}\end{split}\]

where \(\Lstar\) is the stellar luminosity.

Oscillation Equations

The dimensionless oscillation equations are

(13)\[\begin{split}\begin{aligned} x \deriv{y_{1}}{x} &= \left( \frac{V}{\Gammi} - 1 - \ell \right) y_{1} + \left( \frac{\lambda}{c_{1} \omega^{2}} - \alphagam \frac{V}{\Gammi} \right) y_{2} + \alphagrv \frac{\lambda}{c_{1} \omega^{2}} y_{3} + \upsT \, y_{5}, \\ % x \deriv{y_{2}}{x} &= \left( c_{1} \omega^{2} - \fpigam \As \right) y_{1} + \left( 3 - U + \As - \ell \right) y_{2} - \alphagrv y_{4} + \upsT \, y_{5}, \\ % x \deriv{y_{3}}{x} &= \alphagrv \left( 3 - U - \ell \right) y_{3} + \alphagrv y_{4} \\ % x \deriv{y_{4}}{x} &= \alphagrv \As U y_{1} + \alphagrv \frac{V}{\Gammi} U y_{2} + \alphagrv \lambda y_{3} - \alphagrv (U + \ell - 2) y_{4} - \alphagrv \upsT \, U y_{5}, \\ % x \deriv{y_{5}}{x} &= \frac{V}{\frht} \left[ \nabad (U - c_{1}\omega^{2}) - 4 (\nabad - \nabla) + \ckapad V \nabla + \cdif \right] y_{1} + \mbox{} \\ & \frac{V}{\frht} \left[ \frac{\lambda}{c_{1} \omega^{2}} (\nabad - \nabla) - \ckapad V \nabla - \cdif \right] y_{2} + \mbox{} \\ & \alphagrv \frac{V}{\frht} \left[ \frac{\lambda}{c_{1} \omega^{2}} (\nabad - \nabla) \right] y_{3} + \alphagrv \frac{V \nabad}{\frht} y_{4} + \mbox{} \\ & \left[ \frac{V \nabla}{\frht} (4 \frht - \ckapS) + \dfrht + 2 - \ell \right] y_{5} - \frac{V \nabla}{\frht \crad} y_{6} \\ % x \deriv{y_{6}}{x} &= \left[ \alphahfl \lambda \left( \frac{\nabad}{\nabla} - 1 \right) \crad - V \cepsad - \alphaegv \cegv \nabad V \right] y_{1} + \mbox{} \\ & \left[ V \cepsad - \lambda \crad \left( \alphahfl \frac{\nabad}{\nabla} - \frac{3 + \dcrad}{c_{1}\omega^{2}} \right) + \alphaegv \cegv \nabad V \right] y_{2} + \mbox{} \\ & \alphagrv \left[ \lambda \crad \frac{3 + \dcrad}{c_{1}\omega^{2}} \right] y_{3} + \mbox{} \\ & \left[ \cepsS - \alphahfl \frac{\lambda\crad}{\nabla V} + \ii \alphathm \omega \cthk + \alphaegv \cegv \right] y_{5} - \left[ 1 + \ell \right] y_{6}, \end{aligned}\end{split}\]

where the dimensionless oscillation frequency is introduced as

(14)\[\omega \equiv \sqrt{\frac{\Rstar^{3}}{G\Mstar}}\sigma\]

(with \(\Mstar\) the stellar mass). These differential equations are derived from the separated equations, with the insertion of ‘switch’ terms (denoted \(\alpha\)) that allow certain pieces of physics to be altered. See the Physics Switches section for more details.

For non-radial adiabatic calculations, the last two equations above are set aside and the \(y_{5}\) terms dropped from the first four equations. For radial adiabatic calculations with reduce_order = .TRUE., the last four equations are set aside and the first two replaced by

\[\begin{split}\begin{aligned} x \deriv{y_{1}}{x} &= \left( \frac{V}{\Gammi} - 1 \right) y_{1} - \frac{V}{\Gamma_{1}} y_{2}, \\ % x \deriv{y_{2}}{x} &= \left( c_{1} \omega^{2} + U - \As \right) y_{1} + \left( 3 - U + \As \right) y_{2}. \end{aligned}\end{split}\]

Boundary Conditions

The dimensionless boundary conditions applied by GYRE are selected based on the inner_bound and outer_bound options of the &osc namelist group.

Inner Boundary

When inner_bound = 'REGULAR', GYRE applies regularity-enforcing conditions at the inner boundary:

\[\begin{split}\begin{aligned} c_{1} \omega^{2} y_{1} - \ell y_{2} - \alphagrv \ell y_{3} &= 0, \\ \alphagrv \ell y_{3} - (2\alphagrv - 1) y_{4} &= 0, \\ y_{5} &= 0. \end{aligned}\end{split}\]

(these are the dimensionless equivalents to the expressions appearing in the Boundary Conditions section).

When inner_bound = 'ZERO_R', the first and second conditions above are replaced with zero radial displacement conditions,

\[\begin{split}\begin{aligned} y_{1} &= 0, \\ y_{4} &= 0. \end{aligned}\end{split}\]

Likewise, when inner_bound = 'ZERO_H', the first and second conditions are replaced with zero horizontal displacement conditions,

\[\begin{split}\begin{aligned} y_{2} - y_{3} &= 0, \\ y_{4} &= 0. \end{aligned}\end{split}\]

Outer Boundary

When outer_bound = 'VACUUM', GYRE applies the outer boundary conditions

(15)\[\begin{split}\begin{aligned} y_{1} - y_{2} &= 0 \\ \alphagrv \alphagbc U y_{1} + (\alphagrv \ell + 1) y_{3} + \alphagrv y_{4} &= 0 \\ (2 - 4\nabad V) y_{1} + 4 \nabad V y_{2} + 4 \frht y_{5} - y_{6} &= 0 \end{aligned}\end{split}\]

(these are the dimensionless equivalents to the expressions appearing in the Boundary Conditions section).

When outer_bound = 'DZIEM', the first condition above is replaced by the Dziembowski (1971) outer boundary condition,

\[\left\{ 1 + V^{-1} \left[ \frac{\lambda}{c_{1} \omega^{2}} - 4 - c_{1} \omega^{2} \right] \right\} y_{1} - y_{2} + V^{-1} \left[ \frac{\lambda}{c_{1} \omega^{2}} - \ell - 1 \right] y_{3} = 0.\]

When outer_bound = 'UNNO' or 'JCD', the first condition is replaced by the (possibly-leaky) outer boundary conditions described by Unno et al. (1989) and Christensen-Dalsgaard (2008), respectively. When outer_bound = 'ISOTHERMAL', the first condition is replaced by a (possibly-leaky) outer boundary condition derived from a local dispersion analysis of waves in an isothermal atmosphere.

Finally, when outer_bound = 'GAMMA1', the first condition is replaced by the outer momentum boundary condition described by Ong & Basu (2020). A variant of this boundary condition is provided by the 'GAMMA2' choice.

Internal Boundaries

Across density discontinuities, GYRE applies the boundary conditions

\[\begin{split}\begin{aligned} U^{+} y_{2}^{+} - U^{-} y_{2}^{-} &= y_{1} (U^{+} - U^{-}) \\ y_{4}^{+} - y_{4}^{-} &= -y_{1} (U^{+} - U^{-}) \\ y_{5}^{+} - y_{5}^{-} &= - V^{+} \nabad^{+} (y_{2}^{+} - y_{1}) + V^{-} \nabad^{-} (y_{2}^{-} - y_{1}) \end{aligned}\end{split}\]

(these are the dimensionless equivalents to the expressions appearing in the Boundary Conditions section). Here, + (-) superscripts indicate quantities evaluated on the inner (outer) side of the discontinuity. \(y_{1}\), \(y_{3}\) and \(y_{6}\) remain continuous across discontinuities, and therefore don’t need these superscripts.

Structure Coefficients

The various stellar structure coefficients appearing in the dimensionless oscillation equations and boundary conditions are defined as follows:

\[\begin{split}\begin{gathered} V = -\deriv{\ln P}{\ln r} \qquad V_{2} = x^{-2} V \qquad \As = \frac{1}{\Gamma_{1}} \deriv{\ln P}{\ln r} - \deriv{\ln \rho}{\ln r} \qquad U = \deriv{\ln M_{r}}{\ln r} \\ % c_1 = \frac{r^{3}}{\Rstar^{3}} \frac{\Mstar}{M_{r}} \qquad \fpigam = \begin{cases} \alphapi & \As > 0, x < x_{\rm atm} \\ \alphagam & \As > 0, x > x_{\rm atm} \\ 1 & \text{otherwise} \end{cases}\\ % \nabla = \deriv{\ln T}{\ln P} \qquad \clum = x^{-3} \frac{\Lrad+\Lcon}{\Lstar} \qquad \crad = x^{-3} \frac{\Lrad}{\Lstar} \qquad \dcrad = \deriv{\ln \crad}{\ln r} \\ % \frht = 1 - \alpharht \frac{\ii \omega \cthn}{4} \qquad \dfrht = - \alpharht \frac{\ii \omega \cthn \dcthn}{4 \frht} \\ % \ckapad = \frac{\alphakar \kaprho}{\Gamma_{1}} + \nabad \alphakat \kapT \qquad \ckapS = - \upsT \alphakar \kaprho + \alphakat \kapT \\ % \ceps = x^{-3} \frac{4\pi r^{3} \rho \epsnuc}{\Lstar} \qquad \cepsad = \ceps \epsnucad \qquad \cepsS = \ceps \epsnucS \\ % \cdif = - 4 \nabad V \nabla + \nabad \left(V + \deriv{\ln \nabad}{\ln x} \right) \\ % \cthn = \frac{\cP}{a c \kappa T^{3}} \sqrt{\frac{G\Mstar}{\Rstar^{3}}} \qquad \dcthn = \deriv{\ln \cthn}{\ln r} \\ % \cthk = x^{-3} \frac{4\pi r^{3} \cP T \rho}{\Lstar} \sqrt{\frac{G\Mstar}{\Rstar^{3}}} \qquad \cegv = x^{-3} \frac{4\pi r^{3} \rho \epsgrav}{\Lstar} \end{gathered}\end{split}\]

Physics Switches

GYRE offers the capability to adjust the oscillation equations through a number of physics switches, controlled by options in the &osc namelist group. The table below summarizes the mapping between the switches appearing in the expressions above, and the corresponding namelist options.

Symbol	Option	Description
\(\alphagrv\)	alpha_grv	Scaling factor for gravitational potential perturbations. Set to 1 for normal behavior, and to 0 for the Cowling (1941) approximation
\(\alphagbc\)	alpha_gbc	Scaling factor for the \(y_1\) term in the outer gravitational potential boundary condition (the second line of eqn. 15). Set to 1 for normal behavior, and to 0 to suppress this term
\(\alphathm\)	alpha_thm	Scaling factor for local thermal timescale. Set to 1 for normal behavior, to 0 for the non-adiabatic reversible (NAR) approximation (see Gautschy et al., 1990), and to a large value to approach the adiabatic limit
\(\alphahfl\)	alpha_hfl	Scaling factor for horizontal flux perturbations. Set to 1 for normal behavior, and to 0 for the non-adiabatic radial flux (NARF) approximation (see Townsend, 2003b)
\(\alphagam\)	alpha_gam	Scaling factor for g-mode isolation. Set to 1 for normal behavior, and to 0 to isolate g modes as described by Ong & Basu (2020)
\(\alphapi\)	alpha_pi	Scaling factor for p-mode isolation. Set to 1 for normal behavior, and to 0 to isolate p modes as described by Ong & Basu (2020)
\(\alphakar\)	alpha_kar	Scaling factor for opacity density partial derivative. Set to 1 for normal behavior, and to 0 to suppress the density part of the \(\kappa\) mechanism
\(\alphakat\)	alpha_kat	Scaling factor for opacity temperature partial derivative. Set to 1 for normal behavior, and to 0 to suppress the temperature part of the \(\kappa\) mechanism
\(\alpharht\)	alpha_rht	Scaling factor for time-dependent term in the radiative heat equation (see Unno & Spiegel, 1966). Set to 1 to include this term (Unno calls this the Eddington approximation), and to 0 to ignore the term
\(\alphatrb\)	alpha_trb	Scaling factor for the turbulent mixing length. Set to the convective mixing length to include the turbulent damping term (see the Convection Effects section), and to 0 to ignore the term
\(\alphacon\)	alpha_con	Exponent in turbulent viscosity reduction factor (see the Convection Effects section)