
Citation 
 Permanent Link:
 https://ufdc.ufl.edu/UF00097423/00001
Material Information
 Title:
 Thermal effects in stellar pulsations
 Creator:
 Pesnell, William Dean, 1957 ( Dissertant )
Buchler, J. Robert ( Thesis advisor )
Bailey, Thomas L. ( Reviewer )
Hunter, James H. ( Reviewer )
Kumar, Pradeep ( Reviewer )
 Place of Publication:
 Gainesville, Fla.
 Publisher:
 University of Florida
 Publication Date:
 1983
 Copyright Date:
 1983
 Language:
 English
 Physical Description:
 iv, 129 leaves : ill. ; 28 cm.
Subjects
 Subjects / Keywords:
 Approximation ( jstor )
Astrophysics ( jstor ) Cepheids ( jstor ) Eigenvalues ( jstor ) Eigenvectors ( jstor ) Hydrogen ( jstor ) Ionization ( jstor ) Luminosity ( jstor ) Opacity ( jstor ) Stellar pulsations ( jstor ) Dissertations, Academic  Physics  UF Physics thesis Ph. D Pulsating stars ( lcsh )
 Genre:
 bibliography ( marcgt )
nonfiction ( marcgt )
Notes
 Abstract:
 A SturmLiouville operator, embedded in the linearized radial
pulsation equations, is used to analyze the thermal effects of stellar
pulsations. Asymptotic expansions are utilized to demonstrate the nature
of these eigenfunctions and to define the thermal response time of a
stellar model. The eigenvalues and eigenvectors for models of a Cepheid
and Beta Cephei variable are given showing the differences in the two
classes of variables. The time scales introduced agree with earlier ones
in their common use of predicting the location of possible driving
zones. In addition, a way of displaying the temperature variations in
the hydrogen ionization zone without the customary "spike" is
demonstrated, and a new class of pulsations, called "sudden" modes, is
introduced.
 Thesis:
 Thesis (Ph. D.)University of Florida, 1983.
 Bibliography:
 Bibliography: 125128.
 General Note:
 Typescript.
 General Note:
 Vita.
 Statement of Responsibility:
 by William Dean Pesnell.
Record Information
 Source Institution:
 University of Florida
 Holding Location:
 University of Florida
 Rights Management:
 Copyright [name of dissertation author]. Permission granted to the University of Florida to digitize, archive and distribute this item for nonprofit research and educational purposes. Any reuse of this item in excess of fair use or other copyright exemptions requires permission of the copyright holder.
 Resource Identifier:
 030301751 ( AlephBibNum )
11353743 ( OCLC ) ACK4888 ( NOTIS )

Downloads 
This item has the following downloads:

Full Text 
THERMAL EFFECTS IN STELLAR PULSATIONS
B '
WILLIAl [DEAJ PESF.ELL
.A, ISSlF.R'aTION t R'LSLITEL' TO rhE iCrADLTATE SCHOOL
OF THFE lNIVE'ii' ,F FitLORIDL'A
id P[A.TIAiL ILFILLfLtJT OF THE PFKELilR.IHE NTS
FOF THE [EFGR.EE Of( L)CTO'n Of Pti lLOSu'piu
LObfjLLKSIT'i OF FLORIDA
ACKNOWLEDGEMENTS
I would like to thank Dr. J. Robert Buchler for his assistance
during the past four years as well as his great patience. In addition, I
would like to thank A. N. Cox for insight into the physical nature of
stellar pulsations and Dr. Oded Regev for useful suggestions and
discussions.
Several other people have assisted in parts of the work. They
include Dr. Robert Coldwell, who assisted in the writing of the programs
that are inchludic and, Sharoni Eullivant, who demonstrated the use of the
urd pro.ceSsor. I oulid also iike to thank Marie Jose Goupil for reading
the r..gh draft .f t.c z.an'iscript.
.:,rue *f the co.niputations in this work were performed at the Los
ij.mo) NaLrional L borat:,r' Los Alamos, New Mexico, while I was a Summer
'7raduatE Fesearch .As.istsnt there in the summer of 1982. The bulk of the
comput~tion were d.,e at the University of Florida, primarily under the
IUSIC sIs ,'teu, for which I ad grateful.
TABLE OF CONTENTS
ACKNOWLEDGMENTS ... ............ ........ ................... i
ABSTRACT ... ........ ...... ....... .... ..... .............. iv
CHAPTER
I INTRODUCTION.....................................
II BACKGROUND .......................................3
III THEORY...........................................7
IV THE THERMAL OPERATOR ............................10
Asymptotic Expansions...........................16
Asymptotic Thermal Eigenfunctions...............22
V NUMERIC RESULTS.................................33
BW Vulpeculae...................................39
Cepheid model...................................41
VI DISCUSSION ......................................74
VII CONCLUSIONS....................................91
APPENDICES
A THE ANALYTIC EQUATION OF STATE ..................93
B THE INITIAL MODEL INTEGRATOR ...................103
BIBLIOGRAPHY ...... .......................................... 125
BIOGRAPHICAL SKETCH ...................................... 129
iii
Abstract of Dissertation Presented to the Graduate School
of the University of Florida in Partial Fulfillment of the
Requirements for the Degree of Doctor of Philosophy
THERMAL EFFECTS IN STELLAR PULSATIONS
By
WILLIAM DEAN PESNELL
AUGUST 1983
Chairman: Dr. Jean Robert Buchler
Major Department: Physics
A SriurmLi,':u. lie operator, embedded in the linearized radial
pulSation equations, is ujed to analyze the thermal effects of stellar
pul aioni. ampctoric exF.panions are utilized to demonstrate the nature
of these eicenf.,nicoioni .and, to define the thermal response time of a
scell3r moiel. The eiten.,alues and eigenvectors for models of a Cepheid
and et5 Crcphel variable are given showing the differences in the two
ilasse of varlt3rle'. TrIe tie scales introduced agree with earlier ones
in L eir ci m1i:n ue of predicting the location of possible driving
zor.e:. in additiolr, uwa., of displaying the temperature variations in
the hydrogen loniiatlion ione without the customary "spike" is
demont rated, ani a. new class of pulsations, called "sudden" modes, is
introduced.
CHAPTER I
INTRODUCTION
The use of time scales to interpret phenomena found in stars has
been a favorite tool in astrophysics, (see Cox and Ciuli 1968; Saio and
Wheeler 1982; and Shibahashi and Osaki 1981). Primarily, they are used
to simplify the governing equations by ignoring time dependence that
are much longer or averaging those that are much shorter than the
subject under study. Examples of this include the suppression of
hydrodynamic and possibly thermal phenomena when examining nuclear
evolution (Iben 1965), nuclear evolution effects in stellar pulsations
(Christy 1966), and hydrodynamic phenomena when studying thermal flashes
(Henyey and Ulrich 1972). These time scales arise by recasting the
various relationships into nondimensional forms and defining the
factors multiplying the time derivatives as position dependent time
scales. This produces a good approximation for the discussion of some
attributes of pulsations.
However, we have been using asymptotic formulations to study
stellar pulsations (Buchler 1978), and when we attempted Lo include Lhe
outer layers (Buchler and Regev 1982a), it was found necessary to
understand the thermal structure in a more quantitative way. Ocher
authors have found modes that are dominated by thermal eftects (Wood
2
1976; King 1980; and Saio and Wheeler 1982) which have been called
"strange modes." In this work, a mathematical formalism for studying
modes with these characteristics is given.
We shall clarify what the thermal response of a star is and show
effects in pulsations that result. A SturmLiouville operator, embedded
in the linear, nonadiabatic radial pulsation equations, is used to give
what will be called the thermal eigenfunctions of a star. It will be
shown that under physically realistic assumptions, the eigenvalues are
real and represent decays back to the initial state, here assumed to be
one of thermal equilibrium.
Asymptotic forms for the eigenvalues and eigenfunctions are shown,
and several realistic models have their spectra calculated. A discussion
of the physical requirements for the quasiadiabatic approximation for
the stability coefficient is given, along with a way of evaluating it
without the introduction of cutoffs in the integration variable. The
lowest order pulsation operator is shown to have a more complicated
structure than the simple adiabatic operator, when the outer layers are
included. This operator, and the accompanying stability coefficient, are
derived for a general star and compared to earlier work. A new class of
modes, dominated by thermal effects, is shown to exist.
CHAPTER II
BACKGROUND
The theoretical study of stellar pulsations can be divided into
several general areas. First is the adiabatic wave equation and the
discovery that the periods predicted by this theory agree with
observations (Ledoux and Walraven 1958; Cox, J.P. 1980). Here the entire
star is assumed to react faster via an acoustic process than thermally,
implying that stellar pulsations are basically sound waves propagating
in the envelope. In much of a star's mass, adiabaticity is a very good
assumption and the periods calculated will be correct as long as this is
true. The stability analysis of an adiabatic mode uses the quasi
adiabatic approximation (Ledoux and Walraven 1958; Ledoux 1965),
although the integrals must have a cutoff introduced to prevent the
outer layers from dominating the result. These calculations demonstrate
that the quasiadiabatic interior would damp any oscillations unless an
active source of driving was present. Several authors (Cox 1960 and
1958; Zhevakin 1963 and his earlier papers) showed that the
thermodynamic work derivable from the ionization of once ionized helium
could provide the necessary driving for the classical Cepheid variables.
This driving mechanism requires that some of the star be nonadiabatic
and a finely tuned relationship between the location in the star .:. che
4
ionization zone and the transition from adiabatic to nonadiabatic
behavior be satisfied. Predictions of this theory are narrow instability
strips in the HR diagram (Cox, J.P. 1980), with sharply defined blue
edges where the evolving star first becomes unstable.
Other classes of stars are driven by similar mechanisms. Both
radial and nonradial oscillations in the ZZ Ceti variables are driven
bY hydrogen ionization i(Sarrfield et al. 1982; Winger, Saio and
Robinion ]9:). Hot wnite dwarfs, in particular PG1159035, have been
.hown to be jntable to radial pulsations in Starrfield, Cox and Hodson
(19,0), arn this analysis has been extended in Starrfield et al. (1983),
to include nonradial motion and to show the driving is due to the
ionization of carbon and *:.ygen. Variables such as the RR Lyrae, BL
Herculae, :. Scuti and W VirgLnis classes, are driven primarily by helium
ionization although hydrogen does participate in the destabilization
(,XC J.P. 19~,!.
imnce crie des5tablilzinrg echanism was known, computer codes were
developed that followed tre nonlinear behavior of pulsating stars
.Cnhristy 19o; King et al. 1964), and reproduced many observed features
of the Cepheid and KR Lyrae variables. Details such as the asymmetry of
the light and velocity cur'.es and the phase shift of the luminosity
ma3.imurn ith respect to the radius minimum were explained using large
land expensive) codes. A number of new problems arose from this work,
two major ones being tne mass discrepancy in Cepheids (Cox, A. N. 1980a)
and the lack or drivlini in beta Cephei variables (Osaki 1982; Lesh and
Aizcrnman 19v'). in the auttept to understand these difficulties, fully
nonadiabatlc, but linear, stability analyses were developed (Baker and
5
Kippenhahn 1962; Castor 1971; and Saio, Winget, and Robinson 1983).
These analyses have the advantages of being much easier to use and less
expensive than the nonlinear codes. They are easily modified to include
more physics, such as convection (see, for example, Baker and Gough
1979).
The linear analysis shows only whether a given model is self
exciting or damping. An unstable model in the linear regime is a
candidate for further study in a hydrodynamic code. However, the use of
the linear, nonadiabatic eigenvectors in a perturbation expansion is
limited due to their nonperiodic nature. To understand how non
linearities affect an oscillation, without using a hydrodynamic computer
code, a purely real spectrum of eigenvalues is required in first order.
The adiabatic eigenfrequencies satisfy this requirement, but treat the
outer layers in an extremely poor fashion.
As the bulk of the star is quasiadiabatic, Buchler (1978)
developed a formalism using twotime asymptotic techniques modeled after
Cole (1968) and Kuzmak (1959). This method links the linear analysis
with the nonlinear calculations by providing a consistent expansion in
powers of the amplitude (Buchler and Perdang 1979). The requirement of
periodicity is inherent in an expansion of this type, but the payoff is
the capability to study a pulsation amplitude as it starts at an IniLtal
value and evolves to a final state. Additional problems that can be
studied include modal selection, mode switching, and the feedback of the
pulsation on the thermal evolution of the star, all of which are open
questions in this field (Regev and Buchler 1981; Simon, Cox and ho~d:n
1980).
In the simplest version of the theory, the quasiadiabatic results
are recovered in the zero amplitude limit, but at the cost of
introducing a cutoff in the integrals. When this was done, it became
obvious that these layers were being incorrectly treated and that a more
reliable method needed to be developed. A straightforward way of
including these layers was presented in Buchler and Regev (1982a), where
the structure of the eigenvectors was determined by the condition of
constant luminosity above the cutoff point in the model. While giving
good results for stars of the Beta Cephei class (Pesnell, Regev and
Buchler 1982), cooler stars, like Cepheids, were not treated much better
than with the simpler adiabatic approximation. The present work was
undertaken to quantify the thermal behavior of pulsations. The results
presented are encouraging and explain or clarify several phenomena found
in stellar pulsations.
CHAPTER III
THEORY
The structure and evolution of stars is governed by a set of four
equations: the conservation of momentum, energy, and mass, and a flux
equation. Other relationships are constituitive, such as the equation of
state, opacity (see Appendix A), and nuclear energy generation rates.
Written in a Lagrangian form which uses the interior mass as the
independent variable, these equations are (see Clayton 1968, ch. 6;
Appendix B below)
b r3
= 4ip (31)
Sm
d2 r Gm 4r2 g(r,S) (32)
d t2 r2 4rr2 g(rS) (32)
d t2 r2 a
TdS q(p,T) T k(r,S) (33)
d t )8 m
L(m) = (4r22 ac T(34)
3< s m3
The symbols are defined here to be
m = mass beneath the point of interest
8
r = radius from the center of the star to this point
p = density
T = temperature
P = pressure, P = P(p,T)
S = entropy, S = S(p,T)
L = luminous flux
K = opacity, K = K(p,T)
q = nuclear energy generation rate, q = q(p,T).
We have assumed radiative heat transfer in the diffusion
approximation and that all constitutive quantities are given by analytic
expressions.
To study this system, the equations are solved with all time
derivatives set to zero to find an initial model in hydrostatic and
thermal equilibrium. The equations are then linearized in r and S and
combined to l w.'e trw operator equations with two unknowns, 6r and 6S.
Written in operator form these are
i *.:r + "6fS = A.6r B*6S (35)
d k \
S (i + +r6S C.6r + D6S. (36)
it i~ well inoJn that, with the appropriate boundary conditions,
the A operator i1 of the SturmLiouville class and yields a symmetric
matrix operator (Ca(tor 1971; Buchler and Regev 1982b).
The operator A is the adiabatic pulsation operator (6S=0) which can
be written as (see Cox, J. P. 1980)
d Ir Ppr ]} + {r3d [(3r14)P] pr } 0
dr r r dr dt2 r
and in the Schrodinger form of this equation (Morse and Feshbach 1953,
ch. 6), the new independent variable has the dimensions of time and
represents the acoustic travel time from a point in the star to the
surface and is given by
R, Fr (r) P(r)/2 R
tacous (r)= r dr p(r) = r dr/c(r) (37)
where c(r) is the local adiabatic sound speed and R* is the equilibrium
radius of the star. The fundamental radial pulsation period is roughly
given by t acous(Rcore), the sound travel time from the core to the
surface. As the acoustic time corresponds to the time for a dynamic
response to occur, t acous(r) can be identified as the dynamic time scale
acous
as a function of position in a star.
The nice properties of A come from the fact that it is of second
order in the spatial derivatives (the time dependence is assumed to be
exp( iWt )). Looking at the D operator, the same second order
derivatives are present. We now show that there is a representation
where D has the same properties as A, and a natural definition of the
thermal time scale will result.
CHAPTER IV
THE THERMAL OPERATOR
For an optically thick medium, the thermal operator D can be
written in a continuum form as
d6S 6S dS 1 d 6S/
d 6d 6S d
T = q(pe ) qT C d m m)(4 T) + Lj) dI d
(41)
/I d L
Here q (ln) deq = q(p ,T ), i.e. an initial model in thermal
a la\e nT/ 'd m nq eq
balance, T = T1MT C =() and all derivatives are evaluated at
\ P v \ Tp
the equilibrium densities and temperature. This form for D has been
linearized at constant density (radius) and is strictly valid at points
deep in the star so that T > 1.
We asume a time dependence of the form exp( ot ), and note that
the operator Is of second order in the spatial derivative. Therefore, an
irLtegracirin fctor, p exists that transforms D into a selfadjoint
operator. liil factor can be written as
m0 d
,im' = e (42)
witerc r, Is .a reference mass, here chosen to be more, the mass of the
10
11
core where nuclear burning occurs. Defining a new eigenvector as
d In L
= 6S/C and X = ___ eq, equation 41 becomes
v dm
T d cT d Inr d I
T Co = [X(4 KT qT) ] + [ mdl m ]
eq
(43)
Appropriate boundary conditions (which make 43 selfadjoint) are
d 6L/L b d
d m a m(m,) + b (m,) = 0 (44a)
and
F(m ) = 0. (44b)
The surface boundary condition (44a) is equivalent to the assumption
that at the surface there is no radiation incident from the outside
(Castor 1971). In this form, 44a shows that the thermal readjustment
time of the outer surface is infinitely short. As the weight function on
the lefthand side is positive definite, and the function multiplying
the derivative on the righthand side is negative definite, 44 is of
the SturmLiouville class of operators (Boyce and DiPrima 1969, ch. 11).
Comparing 43 to the standard SturmLiouville equation, since
In 1 > 0 and pTC /L > 0 all eigenvalues are real and
d m v eq
positive. With each eigenvalue on is associated an eigenvector Fn whose
number of nodes corresponds to the numeric position of the {On} in an
increasing ordering of the eigenvalues starting with a0. In addition,
12
the eigenvectors are nondegenerate, orthogonal with respect to the
weight function pCvT/ Le and can be normalized so that ( m, is the
total mass of the star)
J p T C Si dm = 6...
0 v 1 J L 1J
eq
Due to their SturmLiouville character, these eigenvectors form a
complete set for describing the 6S readjustment in a star. In
combination with the adiabatic radial wave functions (Cox, J. P. 1980),
a complete set for describing spherically symmetric disturbances of a
sAir is found, provided that equations 31 through 34 are always valid.
To understand the structure of the eigenvectors, the approximation
of neglecting nuclear burning can be made. This reduces the equation to
a version suitable for discussing thermal effects in a stellar envelope.
Without nuclear burning the initial luminosity is constant and the
analysis is halted at a mass ( more ) outside the region where this
assumpction i not valid. is the structure of the n 's will show, this is
equi.alent co as5dmanf wee are looking at "short" thermal time scales
chat do not probe the core. Radial stellar pulsations are an envelope
phenomenon, io this aisumption is justified.
I:electcin nuclear burning and defining L, to be the static
lurninoit c entcri&, the envelope from below, equation 43 becomes
ST L T n d + (45)
Fllo and Pekers d m d 4
Fciloiare Ledou.: and Pekeris (1941), we multiply eq. 45 by n and
integrate over the interval m = more m m.. A selfadjoint system
satisfies a minimum principle and any well behaved eigenvector can be
used to estimate the lowest eigenvalue, a0. Choosing S0(m) = 1 for all
m, this estimate is
m, d IT m* d T
J LI(m) dm J dm
m0 dmm m0 d m
S= L L
0 m
a p TC dm
T dm (46)
If the thermal time is defined as the reciprocal of o0 this gives the
overall thermal time of the envelope as
t (47)
th L m d KT
pI dm
m0 dm
There is an apparent mathematical inconsistency with this
expression. When the opacity derivative, KT, is constant in the
envelope, the thermal time is undefined. However, comparing 46 to the
expression for the adiabatic radial pulsation fundamental frequency
(Ledoux 1965),
S 4xr2 dm [ ( 3 rl 4 ) P ] dm
WO m=
0 r2 dm
dP
they have the same form only if = 0 so that the formula
dm
dP d 1
for tth is incomplete. If = 0 and d = 0, then w0= 0, but the
th d m dm 0
variations of P are large in a star, varying many orders of magnitud=,
while rF can vary only from 1 to 5/3 and the neglect of the pressure
14
derivative would not be correct. The inclusion of nuclear burning in 47
remedies this situation and the net thermal time for the star is then
C =< (47')
th,nuc m d KT
L* 0f [(4 T q) d ] dm
which allows for the luminosity variations and shows that the apparent
infinite value of tth from 47 when KT is constant in the star is not
real. The term allowing for nuclear burning will always give a finite
value for tth,nuc When KT is constant in the envelope, the assumed
form of & is a poor approximation and gives an unrealistic result. The
posicI.'lt of t is difficult to guarantee, but a simple case would
th.nuc d K
d InL T
require tht  qT > 0,  m> 0, and  < 0 all of which can
te satisfed in most of a star's mass, although in some nuclear burning
regions, qT can be much greater than zero, and the opacity variations
are .i,:j3ssed below. The inclusion of nuclear burning as it relates to
scability of Ftellar models is discussed in Hansen (1978) and Gabriel
(19;'). in this work, this area will not be discussed as the
simplii.:acion that result from the neglect of the core make the
interpretation of the results easier.
the opacity is of the Kramers' form ( K T3.5 < 0 ) near the
surface but inside the ionization zones. Near the core, electron
sca terling ill dominate ( KT >0), so that KT decreases going
outward in the s r. Therefore, except in a small region near the
ionizatlion one=, tth will have only positive contributions in the
inregral. The mass of the ionization zones is small compared to the
total mass of a scar, so that the negative contributions to the integral
15
will not cause a negative overall time scale. In addition, each point in
the star is weighted by the ratio
m 4 Tm) T()
P(m) T(mo
which gives the interior regions far more importance than the outer
layers, further reducing the influence of the negative portions of the
integrand. Physically tth represents the time necessary to remove the
p weighted internal energy by the radiative luminosity present,
allowing for the position (and therefore temperature) dependence of the
opacity derivative.
Another thermal time can be found from equation 33 as the time
necessary to remove the material thermal energy at the static luminosity
present in the envelope. Writing this in a integral form (Saio, Winget
and Robinson 1983) yields
t ,p m J T C dm/L,. (48)
Comparing the two definitions, there are several differences. The value
tth,p is counting only the thermal energy in the matter that makes up
the star. In equation 47, the formula emphasizes the radiation
( p ~ T4 ) and the method of energy transport, as well as giving the
interior stellar regions more weight than the outer layers. These
differences are minor and are overshadowed by a major shortcoming
of tth. The definition in 48 can be generalized to calculate a thermal
time as a function of mass by replacing the lower limit in the Incegral
by the desired mass. Values for tth,p found are always positive and can
be interpreted as time scales. However, when the same technique is done
with tth, the value returned is not positive definite due to the
ionization zones near the surface. Mathematically, this generalization
means that ((m') = 1, m'>m and F(m') = 0, m'
function at m. The derivatives of S in 45 will yield delta functions at
m, adding contributions to tth that are not easily understood. To
alleviate these difficulties, tth shall be redefined in terms of the
independent variable that results from analyzing eq. 43 by standard
asymptotic methods.
Asymptotic Expansions
Following the technique outlined in Morse and Feshbach (1953, ch.
6) we make the substitutions:
y(x) = [ 2 (dl n
m TC /2
1 _v d InT 2 d
J m6 L. d m
mj TTCv d InT f/2dm
mO m [ d m
k = J2a ,
2 a
1 d2 d v /4 2 L d KT
'dm TCv2 ]
(d 1 1/TC] dx2 C2 TCv dm
\d m/ v
to transform eq. 43 to
2
d  + ( k2 w(x) ) y = 0. (49)
dx
The units of J are (time)/2 so a natural definition of the overall
thermal time of a stellar envelope is
tth= 2. (410)
We note that the integrating factor p does not appear in this
definition. As p is an artificial quantity introduced for mathematical
convenience, this is physically appealing.
From the equilibrium model (B4 in Appendix B below) there is a
relationship between L, and dlnT/dm. Inserting this into eq. 410 gives
m, 3 C r
th [ J ( ) 2dm ]2 (411a)
t 0 (4tr2)2 ac T3
or, using the continuity equation (31), RO = r(m0) and R* = r(m,),
R, 3 C K 1,
tth R ( v )1/2 pdr ]2. (411b)
0 a c T3
This definition shows the effect of both the internal energy that must
be dissipated and the method of dissipation. Rewriting 411(b) in a
slightly different way
R 3 K T C 1
tth =[ R CT )2 2p dr]2 (1lc)
R0 c aT4
or, defining Ug = T C ,
R U /r
t = [ I (Kpgas )dr ]2 (412)
th c R0 Urad
we have a interesting definition for tth. The photon diffusion length is
present ( 1/Kp ), so increasing the optical thickness increases the
relaxation time. The 3/c term can be interpreted as the diffusion
velocity of the photon; however, this is discussed below. From the ratio
of internal energies, it is seen that in the central regions where this
ratio is large, the contributions to tth increase. Near the surface,
where the density and opacity approach zero and the temperature becomes
constant, the radiation dominates and the integrand becomes small.
If the approximation of a monatomic ideal gas with black body
radiation is made, a simple relationship exists for the ratio
U / U Defining n to be the reciprocal of the mean molecular weight
gas rad
(or the molar density in moles/gr), the pressure of the gas can be
written as
a 4
P = P + P = nRpT + a T,
gas rad 3
or, introducing p (Chandrasekhar 1967,ch. 6),
P = pP and P = (lP)P.
rP rad
19
The internal energy of the gas can be written as
3 4
U = U + U = nRT + aT /p,
gas rad 2
so that
3 3
Ug = P /p B P/
gas 2 gas 2
and
Urad = 3Prad/p = (1P)P/p.
"rad rad
With these assumptions, the ratio becomes
U gas/Urad =1/2 P/(0 ).
gas rad 2P1 
This mixture has a thermal time given by
tth = R [ p P/(1p) 2dr }2. (413)
If P is constant in the envelope
t =3 R (/2dr ]2.
th 2c 1 p Rf
20
Here we see that the thermal time is a diffusion length
1= [ R (Kp) 2dr ]2 (415)
0
divided by a diffusion velocity
v = c (1 p)/, (416)
where c is the speed of light in vacuo.
The diffusion velocity must be less than c; therefore, the value
of 3 must be greater than 0.4. This formula shows that in stars where
0 i' of order I, the thermal times are large as it takes a long time for
Spnoton to mo.e through the material. In a star that is radiation
dominated, the thermal time is short as the photons can move freely
[lirroujh the enelope. These interpretations are qualitative, depending
on the i5sumed nature of U They are suitable for estimating whether
gas
a tcar ujrnergoe; adiabatic pulsations or moves in some way that shows
large thermal effects. A star with a very tenuous envelope, and
cliereforc mall p, would be expected to show large thermal effects.
Star. of this cpe include the R Corona Borealis class (King 1980; Saio
and wheeler 1982', as well as a Cygnus (Lucy 1976; Pesnell 1981).
RecturnLir. to the original form of tth, a position dependent time
scii e car be defined by replacing the lower limit by the desired radius.
The res;llr nt formula represents the thermal travel time from that point
tIo cti Surface and can be written as
th(r R, 3 [ j C /2 p dr ]2. (417)
r a c T3
In figure 41, we compare this thermal time with tth,p, also
generalized to a position dependent form as discussed above (Cox, J.P.
1980; Saio, Winget, and Robinson 1983). The time scales are plotted
versus the logarithm of temperature (with the surface of the star at the
left of the diagram) for a model of BW Vulpeculae and in figure 42, for
a generic Cepheid, discussed in chapter 5. Except for the constant
separation, it is readily seen that the two functions agree quite well
in the interior but diverge from one another at the surface. This effect
is a consequence of tth,p containing only the material thermal energy,
which for a hot star is quite considerable and the drop to zero at the
surface is due only to the integral formulation for tth,p. In the
equation for tth, the method of energy transport reduces to the grey
atmosphere approximation at the surface, and the thermal time must
therefore approach zero more smoothly for any class of star.
This disagreement does not change the normal use of the thermal
time. Several classes of variable stars are driven by ionization zones,
including Cepheids (Cox, J.P. 1980) and ZZ Ceti (Starrfield et al. 1982;
Winget, Saio, and Robinson 1982). For this driving mechanism to
function, the ionization zone must coincide with the point in the star
where the transition from (quasi)adiabatic (constant entropy) motion to
highly nonadiabatic (constant luminosity) motion occurs (Cox and
Whitney 1958). This transition zone is found at the radius, Rtr'
where tth,p(Rtr) ~ IO, the thermal time is of the same order as tIe
th,p tr 0
22
period of the pulsation. The expression for tth is derived from a
frequency, and the location of this region would be where
1,
o(Rtr ) = tth(Rtr) 0= 2n/H0, which is analagous to critical damping.
The constant separation in the figures is approximately 2n; so both
formulations predict the same location for the transition zone, provided
this region is below the atmosphere.
Asymptotic Thermal Eigenfunctions
Once the normalized coordinate has been defined, the asymptotic
form for the eigenvalues and eigenfunctions of eq. 49 can be
3calculated. Here, asymptotic implies that k2 = J2 o > w(x), for all
values cof .:; conSidEred. Due to the surface effects, this is not a good
assuptcion, but ai Important feature in the eigenfunctions can be seen
in trhi appro:imation; so, it is discussed.
RTe asyivpcotic form of eq. 49 is the harmonic oscillator and has
boundary, conditions
Yk'll) + b' d(l) = 0 (418a)
d
and
,0) = 0 (418b)
where a' and b' are a and b from eq. 44 transformed to the x coordinate
jatem. a tmuming ,I.;i = A sin( kx + 6 ) replaces the outer boundary
23
condition by a transcendental equation of the form
tan( k + 6 ) = b'/a' k (419)
In the limit k > 1 this equation has roots that correspond to half
integer multiples of n (Butkov 1968, ch. 9), k = ( 2 + 1 ) ,
S ( ( + ) n )2/J2, and we shall let 6 = 0 .
The physical eigenvector can be found from the asymptotic
expression as
(x) (x) [12 TC v/d inT II
n(x) = Yn(X) / [ 2 dln (420)
n Yn L:)( Ld m ) ( )
In Figure 43, we plot the x variations for the Beta Cephei model.
For log(T)<5.7, the shape of the eigenvector is determined solely by the
denominator in 420. In this star the integrating factor decreases
outward and the net effect is an eigenvector that peaks at the surface.
(Figures 44 and 45). The minimum in i at log(T)~4.6 is due to a
maximum in KT at that temperature caused by the second ionization stage
of helium. The shoulder in y at this temperature is also from this
minimum.
Moving to a Cepheid model, the situation changes considerably. In
Figure 46, the thermal coordinate is plotted, and we see that x~l
until log(T) >4.5, well interior to the hydrogen ionization zone. When
the asymptotic (n=10) eigenvector is plotted (Figure 47), a large peak
occurs at log(T)=3.9. This variation cannot be due to the cosine term in
420 as the argument is essentially 1.0; so, this feature must be due to
24
the denominator, whose major variations come from the integrating
factor. A large minimum in p can be seen in Figure 48, where the
logarithm of p is plotted. This minimum corresponds to the maximum
in and is caused by the large positive value for KT in the cool side
of the hydrogen ionization zone. The integrating factor has a
4 T
temperature dependence of roughly T and in regions where
KT > 4, p increases outward instead of decreasing as it does elsewhere.
This structure corresponds to the peak in the temperature eigenvector in
a standard linear, nonadiabatic pulsational stability analysis. As the
numeric results will demonstrate, any ionization zone will cause a local
peaking of the eigenvector, but a region where KT > 4 causes an effect
that is best removed when calculating thermal effects in these stars,
and plotting the resultant eigenvectors. There are secondary minima
in p at log(T)~4.2 and 4.7 due to the ionization stages of helium. Any
effects these have on the eigenvector y (Figure 47) are completely
hidden by the peak from the hydrogen.
However, it is not true that k2 = J2o > w(x) everywhere. At the
surface, T>TO and < >0, so that w(x) increases outward. While the
interpretation of the "bumps" in the eigenfunctions as maxima in KT
will be .alid in the actual numeric results, the frequencies and
eigenfiJnct ions are not accurately represented in these expressions. With
crts knowliege in hand we now move on to solving the equations
nIuaericall, and show the thermal modes for several stars.
m
oa
4O
0
N
t)
P,
bo >
^ I
HO
p r
,c
4J' CO
*ll 
Hr
60'6 ILi E S SCE
*H
C,
0
a)
* H
0
0
H
04
o >O
l
0
*1
a)
4 0
.e4
.C
0
aJc
n)
Fr
h   "L t4 4
4 *t' *.' '44 i,,C
~~_~ ~
(1U
w
r_ H
tc
C
g;
a
3
0
S0
r1
U
C
A
06_ __ 0 9 : a ho r
C
JJ
C
/ "''
N
+ ~. 
G I
U
r)
C
0 0
*3 '
z
a
00:1 06:L 00: L)LO 0:0 0!0 0:0 0:*
31Ui~lOH00 ~UHYFI a
28
a
0
bo
b 0
>
4
O "9 02 "L on '
o3
t0
4
3,^ 4 C friCL C CC
o
4 r
m
aC
Ul
C
SJ ,ll ] AHr 3.]. i LS
a
LL 09 0 l1 rh 0Il ? 1 9
I
ji. I 1ii
10
ci
0
c0
H
'0
a)
E9
'0
JU
;i
' 12 H
0
oJ
0
o
e
*r
.
1
.i
l o
tT s"0a Go'61 00" 1 00 6 0 n" 0"
l)
I;
I 
*n
Cr
0
',
a
0
r1
I U
ca
o
0
a)
E
A '
0 o7 0L 02 09 oh' s0
^^icI 4
I
CHAPTER V
NUMERIC RESULTS
In this section the eigenvalues and eigenvectors for two stellar
models will be discussed. In each case, the initial model was found by
using the model builder described in Appendix B, with mass, luminosity,
effective temperature, and composition appropriate to the class of star
the model represents. One model is of BW Vulpeculae, a hot supergiant of
the Beta Cephei class of variables, and the other is a generic Cepheid.
To find the matrix operator, equations 33 and 34 are
differentiated at constant radius (density) and written in a form
resembling that of Castor (1971):
d 6S
T d 6S = AK2(T6S). (51)
d t
Although the matrix AK2 is not symmetric, it is tridiagonal and a
numeric analog of the integrating factor does e
Applying this transformation to AK2 yields the D operator as a real,
symmetric, tridiagonal matrix. The diagonal transformation m*riy is
defined as
1/2 (1) = 1
/2
2 (1) P 2 (,1) [ AK2(I1,3)/AK2(I,I) 2, I=2,N,
and the matrix D is found by using /2 as a similarity transformation
D = 4/2 AK2 L '2
so that the spectrum of D is identical to that of AK2. In a convective
region of a star, one procedure is to neglect the Lagrangian
perturbation of the convective flux (Saio, Winget and Robinson 1983).
1/2
This will have no direct effect on the definition of p/2 as the two
matrix elements, AK2(I1,3) and AK2(I,1), are multiplied by the same
fraction of the luminosity being carried by the radiation. Any gross
features of 1/2 such as peaks in ionization regions, will still be
present in a model where convection is permitted, although they may be
extended over a larger region due to the smaller temperature gradient.
The transformation has one side benefit, the matrices analogous to AG2
and AK1 in Castor (1971), are much closer in magnitude to the elements
of the symmetric matrices A and D and have smaller variations as a
function of position in the model, resulting in better numeric accuracy.
Real eigenvalues are guaranteed by the symmetric, tridiagonal form
of the matrix operator (Ralston 1967). For negative definite eigenvalues
(representing decays) the diagonal elements must be negative. This can
be seen in a simple example.
A negative diagonal element means that increasing the entropy
temperaturee) in a zone increases the luminous flux out of the zone
35
more than it decreases thi flux entering it. This stabilizes the zone as
the additional energy is lost and the temperature decays to its initial
value. A simple example can be used if the luminosity gradient had a
local relationship. Equation 51 would then be written as
dS
= a(S S0
so that
S(t) SO = Sat
If a>0, the entropy decays, but if a<0 the entropy increases. The latter
situation is possible in a nuclear burning core but not in the envelope.
To preserve the SturmLiouville character of the continuum operator
in 45 in the matrix operator, the offdiagonal elements must be
positive. This condition is composed of two separate requirements. The
luminosity exiting a zone must increase if the temperature of the zone
above is decreased, and the luminosity entering a zone must also
increase if the temperature of the zone below is increased. Both of
these are easily satisfied if the opacity is independent of temperature
or decreases with temperature, ( KT < 0 ). From the form of the
continuum operator (45), only if KT > 4 can violations of thEes
conditions occur and these regions are precisely where the interpoiatlon
formulae for the luminosity are essential. In the initial model the
temperature ratio is not allowed to exceed a given nuriber in che
ionization zones. This process limits the variation of [he opacity 3a
36
well and raises the upper bound for positive matrix elements. A first
order expansion of the superdiagonal matrix element shows that for
KT < 1/ln(T(I1)/T(I))~17 ,
these matrix elements are always positive. We have never found any
negative offdiagonal matrix elements in any calculated envelope model.
This process also limits the variations of the matrix elements so that
the sum of the offdiagonal elements is less than the diagonal elements,
permitting the use of Geschgorin's Theorem (Wilkinson 1965).
The boundary conditions are incorporated into the matrix operator
by giving certain elements of the matrix special functional forms.
Assuming the outer surface to be radiative and given by the grey
atmosphere formula is equivalent to assuming that this surface has an
infinitely short thermal response time. Therefore, the analytic form of
this boundary condition is that of eq. 44, even though the diffusion
approximation is not very accurate here. For numeric computations this
expression is unsuitable, as it is valid at a point exterior to the
actual surface of a model, and we have assumed that the exiting
luminosity is given by
L(N+I) = 4tI R2(N+) T4 + 2
eff 4 3
Whin this formula is linearized and inserted into the operator equation,
the outer equation reduces to one of first order and closes the system.
The inner boundary condition is 6S/C > 0 and is included by zeroing
V
any elements in the matrix that refer to 6S/C (0). Again, this reduces
the appropriate equation (here the innermost zone) to one of first
order.
As the thermal time of the outer surface is zero, and the overall
1
thermal time ( ~ 00 ) is normally longer than a pulsation period, the
spectrum of D must contain eigenvalues which bracket the radial
fundamental adiabatic pulsation frequency. Here is the information that
we seek regarding the transition zone. All parts of a star that can
relax faster through a thermal process than acoustically will move in a
highly nonadiabatic fashion. These regions of a star will be called the
"sudden" regions and correspond, in some stars (see below), to the outer
layers where d L/ = 0 in a fully nonadiabatic stability analysis.
d m
The inner regions of a star have thermal response times longer than the
fundamental period ( o < 0 ) and are referred to as quasiadiabatic.
Eigenvalues for this matrix can be found by a variety of means. For
the low order modes iterating on the outer boundary condition is
possible (Castor 1971), but the structure of the high order modes is
such that this method fails. A general eigenvalue routine (Smith et al.
1976) is used to find all eigenvalues of D and the necessary
eigenvectors are found by standard Gaussian back substitution. The
eigenvectors of D and AK2 are related by the same similarity
transformation as the matrices and the more desirable set can be
calculated, depending on the effective temperature i T ). of the star
under study. If Teff > 15,000 K, then the hydrogen lonlzal1,,i, zone is
absent and the En's give a good representation of cne rhermal modes.
When Teff < 15,000 K, this ionization zone causes a spike to be pre.ent
eff
38
in Fn as shown in the section on asymptotic eigenvectors in Chapter 4.
To remove this feature, the eigenvectors of D are calculated and
displayed. This results in a much clearer picture of what is occurring
in the surface regions of the star.
We cannot look at the long decay rate modes as we have removed the
core and the effect it has on the eigenvectors. The lower modes can have
significant amplitude in the core, but the higher modes are effectively
isolated from the core even more so than the radial pulsation modes.
mhit .zauieo the dlftfcultv roced earlier, the number of nodes in a
caliculated eienv.'ector is no longer an indicator of it's position in the
pectrum: the i:ze of tre nuoDrs causes the computer to underflow and
the interior p.:rtion of the elaen.'ector is set to zero! Physically, this
impilee thltc cre e'.ter.r canr react to perturbations on much shorter
time icjles than the interior. Ho.weer, the outer layers are linked
quite strjonell, t tcer interior through the low order modes. When a heat
pulse is lnIerted at the baje of the envelope, the entire model is
affected. The outer la.ers react to maintain thermal balance with the
neu lui;n.: .itL from the heat pholie. The eigenfunctions of the diffusion
equJatli.n e>xhibst an infincte propagation velocity for perturbations
Iliore and fel.hbach 1953, ch. 7), causing this behavior. The asymptotic
cl.genfInctiuonl are trigir)onetric wjth an amplitude that increases
outward and we can therefore speak of the thermal time of a point in the
.tar a the reciprocal of the largest eigenvalue whose eigenvector has
Aigntilcant ampittude at that p,ition. This may seem ambiguous, but as
thermal mjadei canr be used to e pand the pulsation equations in the
therral tiCesZ, tniS identtficatrion will enable us to study the thermal
39
effects in an oscillation. In particular, the quasiadiabatic stability
coefficient and pulsations that are dominated by thermal effects can be
defined.
Concentrating on those modes with have eigenvalues comparable to
that of the radial fundamental pulsation frequency, we now examine two
stellar models.
BW Vulpeculae
The physical parameters for this star are plotted in Figure 51 and
the eigenfunctions 6r/r, 6p/p and 6L/L for the radial, adiabatic,
fundamental pulsation mode are shown in Figure 52. The results of a
linear, nonadiabatic radial pulsation analysis (LNA) are shown in Table
51. To find a realistic model for this star, it was necessary to
integrate into the nuclear burning core slightly, but the thermal modes
discussed below are not altered by this.
4 1
The fundamental frequency (w0) is 3.9 10 sec and the four
thermal eigenvectors with eigenfrequencies closest to this are plotted
in Figure 53.
Important features are the cutoff near log(T)=5.2 and the two
peaks. The peak near log(T)=4.6 is due to the second ionization of
helium causing a local maximum in KT. At log(T)~5.1, a ecoriar, peak
associated with the first due to the use of the Rosseland .eain opacity.
The weighting function in this mean has a maximum at a temperature
roughly four times that where the ionization occurs and cause a3 ma limum
in KT at this temperature. This peak has been studied by tellingwerf
40
(1978) as well as Saio and Cox (1980), as a possible driving mechanism
for these stars.
The luminosity variations ( 6L/L ) for each mode are shown in
d 6L/L
Figure 54. Near the surface dm = 0 showing the domination of the
short thermal times in this region. Moving inward (increasing the
temperature), the luminosity variation starts to decrease at log(T)=5.0,
oscillates due to the nodal structure of , and drops to zero by
log(T)=6.0. These modes are completely isolated from the nuclear burning
core (log(T)>6.5) as stated earlier.
Figures 55 and 56 show the real and imaginary parts of the
entr pc.p variations ai calculated from the LNA pulsational stability
inali t. The pcakj at log(T)=4.6 and 5.1 as well as a sharp cutoff
abofv 1lo,(Ti=5.3 car. be seen in both components. Each feature has an
iana Ioa in the Lhcrml L eigenvectors in Figure 53. The modulus of the
luLmnosit .arijatcion from this analysis is shown in Figure 57. The
d ,:L L
change from d = 0 is present, albeit at a slightly higher
cemperatire chan it, Figure 54. Some of this is due to contamination
wicn low order cnearmil modes, but the significant amplitude in 6L/L (up
to iot>T,1=1.)l ii due to the adiabatic density fluctuations in this
rec ion.
Tis i.model sihos that the shape of 6S/C in an LNA pulsational
stability analvyis 1i determined by the thermal modes of that star. As
always, [he model is stable and, at present, we do not find any new ways
of dcscablitzinr IL. The structure of the eigenvectors serves to
demonrrace and clArlaify E[e physical reasons behind a transition zone
ana the spilttins of the itar in two regions, adiabatic (log(T)>5.4) and
instantaneous thermal balance (or sudden) at log(T)<5.4. More
interesting behavior can be seen when a model of a Cepheid variable is
examined.
Cepheid model
This star is a generic Cepheid whose characteristics were obtained
from A. N. Cox (1980a). The results of the LNA pulsational analysis and
summary of the input data are found in Table 52. A graph of the
physical parameters is found in Figure 58 and a graph of the adiabatic
eigenvectors 6r/r, 6p/p and 6L/L is shown in Figure 59. The model is
unstable in the first three pulsation modes and corresponds to a model
of the star DL Cassiopelae. The results presented here are illustrated
very nicely by this star and are similar to those of other Cepheids
calculated.
In Figure 510 and 511 the four thermal modes that bracket the
pulsation frequency 0 are plotted. These are the eigenvectors of AK2,
and we see almost no detail (or the individual modes) due to the spike
in the hydrogen ionization zone. The origin of this spike was discussed
in the section on asymptotic eigenvectors, and as it eliminates any
other features present, we choose to remove it from the eigenvectors.
This is accomplished by using the eigenvectors of D, as shown in Figure
512. These vectors have the same eigenfrequencies as those in Figure 5
10, but much more detail can be seen. The bumps at log(T)=3.82 and 4.12
are due to the hydrogen ionization zone. They are the lower slopes of
the spike, and the minimum at log(T)=4.0 is a result of modulation
42
effect of the integrating factor being a maximum at this point. The
helium ionization zone is visible at log(T)=4.6 as an inverted peak, but
the second peak at log(T)=5.1 is absent. From Figure 46, where the
thermal coordinate for this star is plotted, we see that the value for x
starts to drop away from one at this point. In the interpretation of
features for log(T')<.7. the neglect of the trigonometric function
in r 1t valid. Dut interior to this, the actual nodal structure
of r will start to d>rmin3ie.
n
Examining the luminosity variations which correspond to Figure 512
(Figure 51). the iriv.irig mechanism responsible for destabilizing this
tlar is visible .se below). To accentuate this, the variations
for o < h.' are sh,;n In Figure 514 and 515 and a > w0 in Figure 516
and 5!;. In the low frequency results, the luminosity looks much like
that of BW v*l ( Figure 5); increasing the frequency to those in Figure
5lr, ani 517, the dccrea3e in 6L/L near log(T)=4.6 has become larger
than any other part of the vector. Here is the physical reason for the
drn ing zone of this jtar.
For a scar to be unstable in a pulsation mode, a region where
d UL
d < mJiit exitr (i :a, J. P. 1980). This is seen if the work
d
intteral is wricr In itr e form (Aizenman and Cox 1975; Buchler and Regev
1982 b,).
[ d m. 6L/
) f 1 U) () (d L/L dm .
3 P d m
The quan.,tity r3 I is always positive and for the fundamental mode
< (jas  iricreae. outward) (see Figure 52); thus, when
p r
d 8L/L > 0 the star is stabilized and for d L/ < 0 the star is
dm dm
d 8L/L
unstable. In the BW Vul model (Figure 57) dm is zero in the
atmosphere and 6L/L generally decreases with increasing temperature. The
thermal modes have the same basic luminosity behavior (Figure 54),
showing that this star is stable. In the Cepheid model, the pulsational
6L/L (Figure 526) shows a definite region near log(T)~4.6 where
d 6L/L < 0 This outward decrease of 6L/L cannot be due to the
dm
adiabatic density fluctuations as they will always have 6L/L increasing
outward. The thermal modes show a change in structure as the decay rate
is increased through the value a~". The change in the shape of 6L/L in
this region is the cause of the destabilization of the CGpheid.
Calculations of the BW Vul model with o/u < 1 and w/o < 1 do not show
the same change. Examining Figure 518 and 519 where the thermal modes
with /l~1lO, and Figures 520 and 521 where /o~l10, we do not see the
same change in the structure of the thermal modes. This leads to the
conclusion that the star's thermal structure prevents this type of
driving mechanism from functioning.
The location of the transition zone in each model should be
discussed as well. The BW Vul model shows a definite change from
adiabaticity to sudden behavior at log(T)=5.2. The Cepheid shows a less
dramatic change, as the thermal eigenvectors do not show an extremely
rapid cutoff at the transition. However, unlike BW Vul, the Cepheid has
d 6L/L
a region where d L < 0 that begins near the transition. Also, it
dm
would seem more difficult to establish a transition zone for this model
d 6L/L
as d is not zero until well above the helium ionization. In fact,
dm
the luminosity variations from the pulsation analysis (Figure 12bi do
44
not have a vanishing derivative until above the hydrogen ionization
zone. It is difficult to use the timescale arguments to locate the
transition zone for such a model. The use of a w0 to partition the
thermal modes is, in this case, a fairly reliable way to locate this
region. For o < wo the thermal mode has a 6L/L vector that is, in
general, increasing outward and has only a small dip in the helium
ionization zone. When a > w0, 6L/L is dominated by this region,
until on > w when the amplitude of the eigenvector (6S/C ) becomes
small enough that the mode is isolated from that region of the star.
To see the benefits of the integrating factor, the real and
imaginary parts of the LNA pulsational analysis are shown in Figures 5
22, 523, 524 and 525. In Figures 522 and 523, the eigenvectors of
the AK2 matrix are shown, with the "spike" due to the hydrogen
ionization zone dominating the graphs. When the eigenvectors of D are
examined in Figures 524 and 525, the variations have been greatly
reduced. Lb 522, (6S/C )real varies by two orders of magnitude and the
imajinar', part changes by almost three. The eigenvectors of D vary by
less than one order of magnitude and the details of the helium
iont at3ion an other structures can be seen.
The eigenvectors shown in this section have some numeric noise that
*jc3cr'c: coouient. The system that is solved here is identical to that
used in tre L1IA stability analysis, and any noise in the thermal modes
is also present in that analysis. The spikinesss" of the eigenvectors
preenc in Figure 514 can be eliminated by putting more zones in the
model here tr250) but this puts the eigenvalue routine at a loss for
aolutlons. The envelope mass can be decreased, but at the expense of the
accuracy in the pulsation modes. The basic features of the thermal modes
are well represented and accurately calculated. Interior of the
transition zone, the adiabatic density fluctatiuons will overwhelm the
small luminosity changes due to the entropy variations.
The optimal numeric solution to an eigenvalue problem is found when
the matrix elements are approximately equal in magnitude. The
transformation of the mechanical variable introduced by Castor (1971),
6r XO = /DM2 6r does this for the adiabatic wave equation. In a model
where the shell mass is a geometric progression increasing inward, this
transformation is equivalent to each zone having the same acoustic
traversal time. The sound speed is roughly proportional to the square
root of the temperature and the thickness of the shell is roughly
proportional to the mass; so, the integrand in 37 has equal
contributions from each zone.
This transformation has the unfortunate side effect of decreasing
the numeric niceness of the thermal operator, D. The thermal diffusion
velocity (416) decreases as the sound speed increases. To optimize the
matrix D, the mass of each zone should decrease inward, the opposite
direction required of A. This tradeoff between the adiabatic and thermal
operators will be important only if the mode is dominated by thermal
effects. For quasiadiabatic modes, the interior is dominated by the
adiabatic density fluctuations and the inaccuracy in the thermal modes
is unimportant. We are working on ways of circumventing this difficulty
so that both operators are represented as accurately as possible.
I
/ /
/
/
S/
I.
I
ILI
OO
o
J)
0
ua
41
., ) H
IO
60
ca.
u0
n 3 r
.H
r
I :I"IL I WIl.lnl I II '. 'I" l IL II
47
0U
( t
HU
cu
oo inu
cc
/ 2 2VT LSD '0 C6,0 9911 AD 2 61 
i C **"
f m ,0
I / M)^
I I =r *r
1 / ..CT U 
1 / Tl ^
1~, /u nl
1 / U t
\ / ,
I / a
I~ / 1
S~OL3 nN31" 1 tJI0
,n
Oh'o LEO L2iU
A3/SU
nC
., C *
r o
0C
*S C
m HI
1
ucc
Crd
o m
~3
m I
nl
Cd '
a)
BOl
uv I ~u ;
49
0
EO
M)
S lI
ul <
00
c'a
) l
SI
0
O EI

" I tt0 El
ri
". 01
H U
l 0)
1 1 1 + 2 0H
0'"h [r4 0' hy2 00'? ts!i 6 09!0 20!) L 0 0 r 0
1' ig
*ri H
>g
0
41
4
00
4 U
o
0 ,
u a
00 0
02
tHl
0
S0
041
'r
aP
Ul 4
0 3H
i *
a a)
0
CO
so
(u
0 7
0
0
a *
,8
m
o
0g,
C
,i
, ,
5.
14
'I. l I II
0a
0
*u
00
:3 4J
OC
00
0r
oa
o
41
O H
0 U
t C
0 1
a,
40
0''
.HB
r1
11 901 ONU
~
SC,
1
o
0
Ul
vO
.H
,
U 1
4
0
, M
u) H
IJ0
U)
n 
0U
'4
CL (
1111 1
01
.,
ar
m
Cl
40
w 0
C) *H
z0 0
o cC
0''0
A/ CO
C )r
C.c
EI"2eh C5ELt OCCCE Ch'LEE LZ'PE8 COCIh Ee'C)I 29'Chi OCC'96 BI'Ch *' '1
A/SO La
C
L n
0
to
S0i
C
., 3"
)E
w) C
4i
0 ,0
0 1
o CGO
.I
CC
1]v 0:0
WEC
h;E LI E OE : hhf2 1012 OL:T EE I
A3!/s
09*a ES2O E['
ij..
T ,
.';,
Elu i
C
la)
o 1
10'
i 0
co
Sca
U ca
) U
El t>
II l
n 
a) o
O; 41
<"".
C;
 I
i 1. I I I i~
, '
"'*
.,',1 *
1 I
r ~
.L .
1
~J 1
? ii i
rr
I'
L" 
L'
;2I
,; '
_~nrc
CuoE S92'E ;3 C6UI h'l
 H ~ t 
OL O 90O0 o0' 22'
A l.'q
1r
U)
, u
41
C*'1 '
o
0 >r i
El
.. I
44
C)O
o C
Ca) .
.ci
oaC
K r1
.T _
3; ~
IT~7~ii
~ ..:;. ~~. ~. .I: I
I
__ I. r
_ ~~r 
;. .
;e
T
~~?;,
~;~
~; I
=I~
~i~.r;r.~~_
c
~~~~
___~_~_~
~
60
10
!D
*I
m I
O w
*H1
, "
C C
I.I.
S
ct
+ 4
CP 4 1cI
cc
4.Icc
'C
m
s ^ .; or
I~oe I'
i'ii n
i *r
(" a g 5
/ '' 5
i a n
/ ~~ t ^
[ r~in n'in o
*d,
61
en
u
0 u
7 .0
14 C*
I_
. .,
a
0  __
: "
C ) I
6h' V 0" 29'2 912 1 Ou I 99"0 F.'o 1'0 6 414'a F.E
A31SO
LOa)
CT*H tU
is si
Z 3 t
I M p
(1; 01
U
u *
(U 0 /
in e !
6h'E~~~~~~. (0' I9Eh' EI 19L h[ Od hd LI
h3 >l^
62
ai o
P c
00
SO
'V
S 0
*H
a)
o
S .
, wI O.9 aT iqsg g r i a)Qn} )s*Lr? ra "a rven I
a
^ a 73
( 
U~
r 1
*=G__^1;'.. 1
"~"'" ^s'3s^"w ^^ ^.^ ^, __ a j >
T3> ")^' 0; 4' n
C.o U
U a
C, , Ir
WU *
4J
3
't
: Ul,
La 0
SI
a) 0 l
0,II
11 1 Ira
El" U 
Y^s~f'^^" ' '
""l/'n "
; *'
L .
G'O z[' E9'0
K
T [f;n
1/C1
CO
1:
i
d) M
w I
0_
0 0
0
*H'
0 !
cI
Sv
0 4
U
4C
Cj
uT! r
r *
* irl
~tt~
.0
o~~i
 .. 
i; I H1
Lh '* o~L
o n
>
0
S o0
14 0
> In
S 0 C
4 
I
S0 '
H )
 '^ H *
Lh'O g9cD 92 a'
Al/SO
so' 9 I '=
66
0
o
Sr
co
I ,
c .
S '
JI
C"C
rl
4 "oa
0
0
a)
0H
0)0
co
0 0
4.j
00

4
,,
CC
04(0
0 n
Cv
9C 1 bL*6 W12 9 4(h's woi0 L'21 f
C )
0) .
S4,
68
0a
w 0
00
a
.Ja
a
o
a o
I u 0
00
a um
o"3
A 1/SeJ t s u 1 2 W cc 'EL 2 a
*rl
no
69
a, o
0
mi
0 0
a
C :
Sa0)
O
o pc
C,'
ca
al,
a,>
S 0 Ifl On "0 [j 90" 0E0 '00 I0
f (
m a~
\ S m
(Z~ ~ ~ ~ QS OIY [^?1 l
02*0
?0'0 20'0
I n ll I
0 0
a)
a)
0 u
CO
CO
C
50
J C
am
O "
0 r,
C
La
m w
E o
Q
'Cq
I
H .
0
* H C
C1 O
HI
0"
ajc
a' 1
CZ
LO'O
f 4 t'c _______
0/C' )C'C $4 Ci/' fC)
71
(o
0 0
e
0 o
0
0 
I 0
h L 2 ..,[ :
_00
// "' L
4J 44V
N.' c6
<^  ft
:
Table 51
BW Vulpeculae
M=11.0 M
sun
log(L/Lsun)=4.23
Tff= 25120 K
oI= 0 193 ~n=1.9 104
nIl= 0d147 '1=2.3 103
II2= 0 118 ,2=1.5 102
Composition X=0.70, Y=0.28, Z=0.02
Envelope mass = 5.02 M
sun
Table 52
DL Cassiopelae
M=5.69 M
sun
log(L/Lsun )=3.57
T = 5848 K
110= 7 87 0 = 6.8 102
nI= 5 64 n1= 1.9 01
nl2= 4d28 2= 6.2 102
Composition X=0.70, Y=0.28, Z0.02
Envelope mass = 2.76 M
sun
CHAPTER VI
DISCUSSION
To understand the utility of the thermal eigenvectors, we shall
look at various forms that the linearized pulsation equations (35 and
36) acquire when the n's are used to expand them in the ratio of
thermal time to pulsation period. Eliminating 6S one can rewrite 35 and
36 as an integrodifferential equation:
w2 6r = ( A + B i1 C ) 6r. (61)
i) D
we have an eifenvalue problem for w2 and 6r. Eliminating 6r 35 and 3
, can b. written in a form suitable for studying the secular behavior of
i.^: = (D + C B) 6S. (6la)
W2 A
The i,. 0 cerm is a resolvent that can be expanded in the complete set
or eieenveccors of D. Adopting the Dirac (or braket, see Davydov 1965,
ch. 5) notation of quantum mechanics, this set is { 0 ,I n > }
.n.d equaijon r l becomes
[n >< Snl
02 6r = ( A + B C ) 6r .
n
Next, the mechanical eigenvector is
and eigenvectors of A, { 21 yi >
adiabatic effects on a mode that is
expansion will be written as
expanded in the set of eigenvalues
}. We would like to examine the non
primarily adiabatic, and this
6r E I y > = I i > + ji j I yj >
w = ). + Ai
I
where
Inserting this expansion and multiplying eq. 62 by < y give
(1 + Yjiajl 2 ) 02 = < yilAl Yi > +
< Y, IBI 5n >< n ICI yi>
n i) + o
n
< Yj BI Sn >
j n(j i + o
n
(62)
"i2 I yi > = A I yi > .
< yi IBI n >< n ICI Yj >
i + 0 j
n
< Yk IBI >< n CI Yj >
kj' iw + a
(63)
This expression is rather cumbersome, but following the normal quantum
mechanical perturbation scheme, we evaluate the change in the
eigenfrequency using the initial eigenvector. Choosing the adiabatic
mode as the starting point, the quasiadiabatic correction can be
written as
< IBI n >< nICI Yi >
t2_a002 = n (64)
1 n id) + o
n
which is Identical cit equation (26) in Castor (1971). The normal quasi
adiabtiL: approximation is a /w. / 1 for all possible n, reducing 64
0 1
' = < Yi BC y. > (65)
I I = I.
L
II iI'n
Tri; e:rpresio:.n fr lth quasiadiabatic correction can be further
as3lzed it .1V C I. Then w2 = w + 2w.Aw and eq. 65 can be rewritten
1 1
AW = I1 iBCl yi >
2 W2
1
a result identical to the Ti (equation 8, in Buchler 1983). Here, Aw is
purely imaginary, as all quantities on the right are real, and generally
shows that the star is stable (iAw < 0) as the diagonal elements of B
are positive and the diagonal elements of C are negative.
This recapitulation of the quasiadiabatic analysis serves to point
out that by knowing the structure of the thermal modes, a better
approximation can be made. Partitioning the thermal modes into two sets,
S = { o n n = 1, NT : O i < I }
and
S2 = { ,n n = NT+ : 0 ia < 1 },
eq. 62 can be rewritten as an eigenvalue problem assuming an adiabatic
starting point only to determine NT
2 I y > = A y > +
BT B n >< n n ICI y >
N BI u >< I CI Y >
l iu + 0
n=NT+1 n
The second term on the right is the quasiadiabatic correction. rno .iLh
an appropriate cutoff. However, the cutoff is not in the spatial
78
variable, as the integrals implied by the notation extend over the
entire mass of the star. The correct cutoff is to stop the sum at the
thermal eigenvalue whose magnitude is just below the pulsation
frequency, as the sum over the set Sl is doing. For this reason, the
second term should be more correct in the calculation of this part of
the stability coefficient.
The third term produces a correction in the lowest order, and
therefore cannot be treated as a perturbation. To see this effect, the
denominators of the sums are expanded in a binomial series and eq. 66
becomes
w2 I y > = A y > +
1 T BIn Sn IC y >
it I + a /io +
n=l n
N BI &n >< X n I CI y >
a (1 + ) (67)
n=N +1 n n+
or
)2 1 y > = Al y > +
1 nT a
i T n >(1 )< n IC y > +
n=l
N
I BI n >(1 )< X n IC! y >. (68)
n=NT+1 n
The first term is the quasiadiabatic correction already discussed
anld cr.e second is the effect of the outer layers on the pulsation. The
expansion is in /a n, and in the limit l/on >0, or very short thermal
times; this sum becomes
NI
I B Cn >o n IC y >.
n=NT+1 n
This term represents a zeroorder quantity and should not be neglected
when the initial eigenvalue is calculated. If the thermal modes had a
thetafunction behavior at the transition zone, the set S2 would form a
complete set in and of itself for describing the thermal response of the
1
outer layers of the star on time scales shorter than w The actual
eigenvectors have a fairly sharp fall at this point; so, a first try at
1
including these terms is to write the sum over states as D' where D'
is the outer [(N NT)*(N NT)] submatrix of D that includes only those
zones with I > NT Substituting D' for the second sum over states gives
2 I y > = A y > +
B n >< n ICI y > 
n=1
BI(D')IIC y > + iwBI(D')2CI y >. (69)
This equation can be expanded in the ratio of the proper time
scales; the second term on the right al= o/o0 and the last two
in a2= 0/0 Writing a=al=a2 to give the simplest expression and
arranging the equations in a hierachy of terms by order in a g i:
0 Y > = ( A B(D') C ) Y0 > (610a)
Sy1 > + 2w0w Y0 > = ( A B(D')C)I y > +
l 0Bl ><^ C yg>+
%=T B n n CI YO > +
0 n=l
iwoB(D')2 C y >; (610b)
where
I y > = YO > + a I Yl > + a2 y2 > +
and
w = 0W + aw + a2w2+ O(a3).
To eliminate secular terms in the equation for yl >, the scalar
product < yO I y, > is set to zero. For convenience we shall set
=1. This gives the initial eigenvalue and the first order
ztabililt c.,efficient as
..= < yO IA B(D')1CI Yo > (611a)
1 NT n ICI y > +
2w6 n=l
1/2< Y IB(D')2CI YO>. (611b)
At any point the matrix D' can be replaced by the sum over states that
it appro ii,.actes
N
0 I(A + X BI n > < n C) yO >, (612a)
n=NT+1 n
S xNT < Y0 IBI n >n X n 1CI y0 > +
2w20 n=l1
N
1/2 < Y IB > < n C yo>. (612b)
n=NT+1 02
Equations 611a and 611b agree with equations 25 and 39 respectively in
Buchler and Regev (1982b), who derived them from physical arguments
regarding the constant luminosity of the outer layers. It should be
d 6L/L
noted that the assumption of = 0 is not valid over the entire
dm
section of the Cepheid model above the transition zone. The sum over
states is valid, and comparing the thermal mode 6L/L to that of the LNA
analysis (Figure 516), we see that the functional dependence of the
LNA 6L/L is reproduced by the thermal modes.
Using the infinitely sharp cutoff limit, the eigenvectors of 610a
can be calculated in a straight forward way. Defining the error in yO >
as
E = y(I) r(I) 122.
where 6r is the real part of the mechanical eigenvector given b7 the LNA
stability analysis and both eigenvectors are normalized :o units in
their respective scalar products. The error in the eigenfrewuencv is
defined as
E = 1W0 LNAI'LNA
These quantities are plotted in Figure 61 against the location of the
transition zone, log( Ttrans ). A sharp minimum is seen at
log(Ttrans ) = 5.23, corresponding to the location of the actual
trans
transition zone as predicted by the luminosity eigenvector plotted in
Figure 56. In this model, the change is rather abrupt, falling a factor
of 4 when log( Ttrans ) changes by 0.13 or +0.04 The change w0 is not
as dramatic, reflecting the "minimum" principle that a quasiadiabatic
oscillation will follow.
When the same quantities are plotted for the Cepheid (Figure 62),
it is obvious that the agreement is not much better than using the
adiabatic eigenvectors. The error is reduced at a point where the
transition zone is located, but not by a factor of 10 as in BW Vul. This
results from the more complicated behavior of 6L/L in the Cepheid (see
Figure 516). The thermal modes have a large peak in the ionization
zones due to the integrating factor's temperature dependence but do not
have the rapid decrease inward of those in BW Vul. Therefore, the
approcimation of replacing the sum over states by the submatrix D' is
not a rooa one for this class of star. As the sum over states contains
the necesr ary information about the luminosity variations in the
eitericr, including the integrals as in 612a should give much closer
agree.~ienc. This problem is discussed in Pesnell and Buchler (1983).
Comparing the initial eigensystem 612a to the adiabatic wave
equation produces three differences:
83
1) There is, at present, no guarantee of either the reality or
positivity of the eigenvalues of either 611a or 612a.
2) The eigenvectors < yo I and I yo > are not necessarily
related by a transposition, much less being equal. The matrix
A B(D') C is not symmetric but is of an upper Hessenberg
form whose only advantage is being diagonally dominant. The
eigenvectors come as right and left hand vectors and must be
calculated separately.
3) There is no proven relationship between the number of nodes
in I Yo > and the numeric position of w02 in an increasing
1
ordering of the spectrum of A B(D') C .
However, for stars which are predominantly adiabatic, as most
pulsators are, the eigenvalues are close to the adiabatic frequencies,
the second part of the operator having little or no effect on A. The
transposed eigenvector has a serious numeric problem. In the region of
the transition zone, the transposed matrix has a structure such that a
discontinuity appears in < y0 The use of eq. 612a instead of 611a
will spread this feature over several zones, removing the discontinuity.
When the transposed eigenvector is calculated using the LNA analysis,
the feature is present but is much smoother. A change in the method of
calculation is indicated, not a revision of 612a.
A final observation regarding 612a concerns the limit Ni >i.. A
a normal star is almost adiabatic ( NT >N), the opposice liiEc must
describe an abnormal star. For the condition NT .>0 o be satisfied,
the nonadiabatic regions of the star must be ectensl'c, coverlri, all of
the star that is oscillating. These conditions may be found in the R
84
Corona Borealis stars or other highly evolved stars where the envelopes
are tenuous and have low gas to total pressure ratios. In these stars,
trying to recover the adiabatic limit in one obvious way will lead to
some rather strange behavior of w2.
Introducing an arbitrary scale factor, X into equations 35 and
36, which represents the ratio of the pulsation period to some overall
thermal time, they become
w26r = A'6r + B.6S (613a)
iw6S = X(C6r + D6S). (613b)
For most stars, when +O the adiabatic limit is recovered as 6S is
forced to 0. However, in the integrodifferential form, splitting the
thermal modes as before,
M2 I y > = A y > +
NT BI >
w 1 + Xo /i +
n=1 n
N B >< I Y > (614)
n=N+1 ni + Xn
Tne .alue of f:_ may depend on X for most stars, but for the radiation
domninsce. .acmospheres, the thermal response time is very short for all
mode= who=e cutott in amplitude is above the normal transition zone.
I.Ising the ideal eas as before, the atmosphere should have a low density
85
gradient and opacity gradient so that the diffusion length (415) is
large. If the outer modes have o n/ > 1 for all n > N then N is
almost constant during this discussion.
In the limit of X >0 ( NT constant), eq. 614 is expanded using
the binomial theorem as before,
I2 y > = A y > +
N
BI n > < n IC y > + O(X) (615)
n=NT+1 n
the B(D') C term does not disappear. For highly nonadiabatic
envelopes, this term can dominate and change the spectrum so that it no
longer reduces to the adiabatic case. This would, in part, explain the
modes that have been named "strange" modes (King 1980; Saio and Wheeler
1982).
This explanation of the strange modes does not require that a mode
be near (in some sense) to an adiabatic mode. There are no theorems that
prove that all of the eigenvalues of 612a are real, and starting from
an adiabatic mode, the nonadiabatic corrections could move the
frequency far away from the initial value. It is IncerestLrni chat
if NT can be considered constant, the adiabatic limlt is never recovered
in some models. This property is just what a strange mode show.
Modes that possess the property of NT being independent of shall
be called "sudden" modes. They are characterized by 3a ondiciori of
thermal balance existing over a large fraction of the mechanical
eigenvector's amplitude. The quasiadiabatic scabilic, coefficient would
86
be a poor approximation for these modes, as the second term in 612a
will dominate w1. We must stress that the condition of thermal balance
does not imply that the luminosity has a zero derivative in these modes.
0
o
60
0
0
0
04
Tl
>4
a)
0 c
o
S' 1 
on o h o L 't o 2 O 2 0 
^^. a) a
,l~^
i., .i : .. *., ... . . ,.: *4 os co
i1, .i II I 1 ,111
0
r4
m
) P.
0c
w C
1 I I
) 0
CON
Cra
0 >
CN
X C.
4 44
o)
1 0
JC
4O
,za)
v4
s e
096 iPgo L.d h ii Z o 1 ' 9;
I' Jr
14
II
to
0
a ,.
0
) *
H
a)
41
BO
c 0
001
0
C Z
H
0
01
tw
0 11
\ 4J (U
,7u
*44
0 O
0.
'., 4 1
*il o
H 4
n) n
90
cl
to
0
5 A
ro
Q) 0
4j
mC
OC
0
m0o
0, a)
4 FJl 
4c
^ 41
0il
^ (1)
CONCLUSIONS
We have shown that the thermal structure of stellar pulsations can
be represented by the eigenfunctions of the diffusion equation. As this
equation is of second order in the spatial derivative, it can be mapped
onto a SturmLiouville operator, with the concomitant theorems regarding
the reality of the eigenvalues and the form of the eigenfunctions. This
transformed operator defines a thermal time scale through a new
independent variable with the units of (time) 2 This time scale agrees
with that of earlier authors in the interior and predicts the same
location for the transition zone in a star.
The thermal modes that arise from the eigenvalue problem can be
used to analyze the thermal effects in a first order perturbation
expansion. The traditional adiabatic wave operator is replaced by an
operator that includes the knowledge that the outer layers are in
d 6L/L
thermal balance, and that this does not imply that d 6 = 0 in these
dm
layers.
Partial spectra of the thermal modes are presented for two
realistic stellar models, along with diagrams of the eigenvectors.
Mathematically consistent expansions of the pulsation equations are
presented that give the quasiadiabatic stability coefficient and the
92
analogous expression when the outer layers are correctly accounted for.
It is shown that a new type of mode, the "sudden", mode can be defined
using these expressions.
The driving mechanism in Ocpheids has features in the thermal modes
that do not show up in a star of the Beta Cephei class. This casts doubt
on the viability of ionization driving being responsible for these
stars' variability.
In an aesthetic vein, a way of displaying the entropy or
temperature variations for stars that normally have a "spike" in the
hydrogen ionization zone is given, showing much of the detail in this
region without the introduction of an arbitrary cutoff. This
transformation also increases the numeric accuracy of solutions in the
linear nonadiabatic pulsation stability analysis by limiting the
variation of the elements of the coupling matrices.
APPENDIX A
THE ANALYTIC EQUATION OF STATE
For the purposes of these calculations, an equation of state has
been developed, based in part on a routine supplied by Arthur Cox (Cox,
A. N. 1980b). In this appendix the various formulae and constants used
are discussed.
As our analysis requires the entropy of an ionizing gas, we first
calculate the free energy (F) as a function of the specific volume (v),
temperature (T), and composition (X,Y,Z). All thermodynamic functions
can be evaluated by derivatives of F (Zel'dovich and Razier 1966, ch.
III):
entropy: S = ( ,N
pressure: P = 
internal energy: E = T2 ("F)
N
In order to calculate F, the approximation of an ideal, lrnizing
gas is adopted with the characteristics of the five species listed in
Table AI. Values for the mass weighting of the elements has been
93
94
adapted from Cox and Tabor (1976), and the atomic masses and ionization
potentials are from Novotny (1973). For completeness, a term is added to
include the effects of black body radiation. This is assumed to be that
appropriate to isotropic radiation so that
F v
rad A
The 1inized fraction of ezac element is determined by minimizing
the free energy at the gi.en temperature and specific volume. This
procedure leads to a Saha equation that is solved using a Newton Method
to a relatl.e accuracy of i part in 102 for the electron density. When
the electron density iis no,n, the reciprocal of the mean molecular
weight fin moles'gram) 1i
n = 1,'. = :,I. 1.0'797 + Y.'4.0026 + Z/21.02447 + ne.
The relative abunrdance of each ionization stage, ni, is calculated from
the ratios of the appropriate partition functions, and the free energy
15 wurittei as a .su o.er these ,pecies.
F = I r n(n uI ,n + n 2.5048 + ln(vT3/2)

where the constant term is:
1 + In(2xm k) ln(NA) 3 ln(h),
2 p A
mi is the species gram molecular weight, and ui is the temperature
independent partition function that is assumed to the ground state
degeneracy.
The expressions for the entropy (ergs/gram/K), P (dynes/cm2), and E
(ergs/gram) are
S = R{ nl n(m32ui/n ) + n[ 1.0048 + ln(vT3/2) ] + vT3
P = nRT/v + T
E = nRT + avT4 + ni
Derivatives of the pressure and the internal energy are found by
analytic differentiation of the above formulae.
The subroutine NEOS performs the calculations required to find the
thermodynamic quantities as outlined above. A listing of this routine
follows.
r)
e n>
w0 P
'& E
UHU< E
1F> U
CP U 
.n W _
C E P
Lj.JI E 0 I 0
> r I 0* 0
F. lI H O1
1 x W; W Q wz pq 11 O
zI JU I I C)  Il a 
4 ) HI m :l I g* 0= H W )W
4 4ail 4 z rZCNJ
I/fllllu l/] (1 IC PTWL
3 E i0 I H nr
 S.>iiZ. H I N ulcr
L 3j H EI.4.. 3 oiN
N E  I rrc
u >.>C I . o
2 11 I a rn w inu M r w m < 0
r C u w 1 ratsci.xr :j . > <* In :F u
I H i 3 nO (COP5 (
* > aS z::z...! 1 Hi HEi II* 0 US ua3 ei nc
>I ~lLrJJJ E g E(IEDa4l IY H H. H E4I=
* V It UHUC )iU O4E
1 ~ Cc toC* C "* eumC c x T rU H
1I ,2..ccCa.. ij C) Q NrCariU 0. .3Z
i C '5l Ow fiyllui g N nr M EN (0O
= =C =i viH E a ......n.*. .... i 1I Ez EU;.
," ,,n H I I QNNNc* 0 01. a*0 oa z 0
 i*A ^' ^^. J rl j C NEUC .11 (/2
* 'C a'^CLS .2zlC ^ t I HIli NN N I E(oC
r* C 4 i' 1 ,E** ( flo Ul C(J
* v Je *T .'. *. H *C .) = I VN a 4 OW
* 0. "t czi.LuLj
3: 0 INN.(CHN* H II1w HU liii Q 11
XLLa %_ 11 i i 11 1i MI H H i iH C U NHNaC ZN 50CC
Y Cli us~ 4 iC 4 4i a P I C.H E24W4. 00
I *
>u I *
I *
I*
I*
.JUUUtJUU'jUUUtt)JUU'JUUUUUUUUUUUUC0 ) UUUUUU JUUUUU

Full Text 
PAGE 1
THERMAL EFEECTS IN STELLAR PULSATIONS By WILLIAM DEAN PESNELL A DISSERTATION PRESENTED TO THE GRADUATE SCHOOL OF THE UNIVERSITY OF FLORIDA IN PARTIAL FULFILLMEOT OF THE REQUIREMENTS FOR THE DEGREE OF DOCTOR OF PHILOSOPHY UNIVERSITY OF FLORIDA 1983
PAGE 2
ACKNOWLEDGEMENTS I would like to thank Dr. J. Robert Buchler for his assistance during the past four years as well as his great patience. In addition, I would like to thank A. N. Cox for insight into the physical nature of stellar pulsations and Dr. Oded Regev for useful suggestions and discussions. Several other people have assisted in parts of the work. They include Dr. Robert Coldwell, who assisted in the writing of the programs that are included and Sharon Bullivant, who demonstrated the use of the word processor. I would also like to thank Marie Jose Goupil for reading the rough drafts of the manuscript. Some of the computations in this work were performed at the Los Alamos National Laboratory, Los Alamos, New Mexico, while I was a Summer Graduate Research Assistant there in the summer of 1982. The bulk of the computations were done at the University of Florida, primarily under the MUSIC subsystem, for which I am grateful. ii
PAGE 3
TABLE OF CONTENTS ACKNOWLEDGMENTS ii ABSTRACT iv CHAPTER I INTRODUCTION 1 I I BACKGROUND 3 III THEORY 7 IV THE THERMAL OPERATOR 10 Asymptotic Expansions 16 Asymptotic Thermal Eigenf unctions 22 V NUMERIC RESULTS 33 BW Vulpeculae 39 Cepheid model 41 VI DISCUSSION 74 VII CONCLUS IONS 91 APPENDICES A THE ANALYTIC EQUATION OF STATE 93 B THE INITIAL MODEL INTEGRATOR 103 BIBLIOGRAPHY 125 BIOGRAPHICAL SKETCH 129 ill
PAGE 4
Abstract of Dissertation Presented to the Graduate School of the University of Florida in Partial Fulfillment of the Requirements for the Degree of Doctor of Philosophy THERMAL EFFECTS IN STELLAR PULSATIONS By WILLIAM DEAN PESNELL AUGUST 1983 Chairman: Dr. Jean Robert Buchler Major Department: Physics A SturmLiouville operator, embedded in the linearized radial pulsation equations, is used to analyze the thermal effects of stellar pulsations. Asymptotic expansions are utilized to demonstrate the nature of these eigenf unctions and to define the thermal response time of a stellar model. The eigenvalues and eigenvectors for models of a Cepheid and Beta Cephei variable are given showing the differences in the two classes of variables. The time scales introduced agree with earlier ones in their common use of predicting the location of possible driving zones. In addition, a way of displaying the temperature variations in the hydrogen ionization zone without the customary "spike" is demonstrated, and a new class of pulsations, called "sudden" modes, is introduced. iv
PAGE 5
CHAPTER I INTRODUCTION The use of time scales to Interpret phenomena found in stars has been a favorite tool in astrophysics, (see Cox and Giuli 1968; Saio and Wheeler 1982; and Shibahashi and Osaki 1981). Primarily, they are used to simplify the governing equations by ignoring time dependences that are much longer or averaging those that are much shorter than the subject under study. Examples of this include the suppression of hydrodynamic and possibly thermal phenomena when examining nuclear evolution ( Iben 1965), nuclear evolution effects in stellar pulsations (Christy 1966), and hydrodynamic phenomena when studying thermal flashes (Henyey and Ulrich 1972). These time scales arise by recasting the various relationships into nondimensional forms and defining the factors multiplying the time derivatives as position dependent time scales. This produces a good approximation for the discussion of some attributes of pulsations. However, we have been using asjnnptotic formulations to study stellar pulsations (Buchler 1978), and when we attempted to include the outer layers (Buchler and Regev 1982a) , it was found necessary to understand the thermal structure in a more quantitative way. Other authors have found modes that are dominated by thermal effects (Wood 1
PAGE 6
2 1976; King 1980; and Saio and Wheeler 1982) which have been called "strange modes." In this work, a mathematical formalism for studying modes with these characteristics is given. We shall clarify what the thermal response of a star is and show effects in pulsations that result. A SturmLiouville operator, embedded in the linear, nonadiabatic radial pulsation equations, is used to give what will be called the thermal eigenf unctions of a star. It will be shown that under physically realistic assumptions, the eigenvalues are real and represent decays back to the initial state, here assumed to be one of thermal equilibrium. Asymptotic forms for the eigenvalues and eigenf unctions are shown, and several realistic models have their spectra calculated. A discussion of the physical requirements for the quasiadiabatic approximation for the stability coefficient is given, along with a way of evaluating it without the introduction of cutoffs in the integration variable. The lowest order pulsation operator is shown to have a more complicated structure than the simple adiabatic operator, when the outer layers are included. This operator, and the accompanying stability coefficient, are derived for a general star and compared to earlier work. A new class of modes, dominated by thermal effects, is shown to exist.
PAGE 7
CHAPTER II BACKGROUND The theoretical study of stellar pulsations can be divided into several general areas, first is the adiabatic wave equation and the discovery that the periods predicted by this theory agree with observations (Ledoux and Walraven 1958; Cox, J, P. 1980). Here the entire star is assumed to react faster via an acoustic process than thermally, implying that stellar pulsations are basically sound waves propagating in the envelope. In much of a star's mass, adiabaticity is a very good assumption and the periods calculated will be correct as long as this is true. The stability analysis of an adiabatic mode uses the quasiadiabatic approximation (Ledoux and Walraven 1958; Ledoux 1965), although the integrals must have a cutoff introduced to prevent the outer layers from dominating the result. These calculations demonstrate that the quasiadiabatic interior would damp any oscillations unless an active source of driving was present. Several authors (Cox 1960 and 1958; Zhevakin 1963 and his earlier papers) showed that the thermodynamic work, derivable from the ionization of once ionized helium could provide the necessary driving for the classical Cepheid variables. This driving mechanism requires that some of the star be nonadiabatic and a finely tuned relationship between the location in the star of the 3
PAGE 8
4 ionization zone and the transition from adiabatic to nonadiabatic behavior be satisfied. Predictions of this theory are narrow instability strips in the HR diagram (Cox, J. P. 1980), with sharply defined blue edges where the evolving star first becomes unstable. Other classes of stars are driven by similar mechanisms. Both radial and nonradial oscillations in the ZZ Ceti variables are driven by hydrogen ionization (Starrfield et al. 1982; Winget, Saio and Robinson 1982). Hot white dwarfs, in particular PG1159035, have been shown to be unstable to radial pulsations in Starrfield, Cox and Hodson (1980), and this analysis has been extended in Starrfield et al. (1983), to Include nonradial motion and to show the driving is due to the ionization of carbon and oxygen. Variables such as the RR Lyrae, BL Herculae, 6 Scuti and W Virginis classes, are driven primarily by helium ionization although hydrogen does participate in the destabilization (Cox, J. P. 1980). Once the destabilizing mechanism was known, computer codes were developed that followed the nonlinear behavior of pulsating stars (Christy 1964; King et al. 1964), and reproduced many observed features of the Cepheid and RR Lyrae variables. Details such as the asymmetry of the light and velocity curves and the phase shift of the luminosity maximum with respect to the radius minimum were explained using large (and expensive) codes. A number of new problems arose from this work, two major ones being the mass discrepancy in Cepheids (Cox, A. N. 1980a) and the lack of driving in Beta Osphei variables (Osaki 1982; Lesh and Aizenman 1978). In the attempt to understand these difficulties, fully nonadiabatic, but linear, stability analyses were developed (Baker and
PAGE 9
5 Kippenhahn 1962; Castor 1971; and Saio, Winget, and Robinson 1983). These analyses have the advantages of being much easier to use and less expensive than the nonlinear codes. They are easily modified to include more physics, such as convection (see, for example. Baker and Gough 1979). The linear analysis shows only whether a given model is selfexciting or damping. An unstable model in the linear regime is a candidate for further study in a hydrodynamic code. However, the use of the linear, nonadiabatic eigenvectors in a perturbation expansion is limited due to their nonperiodic nature. To understand how nonlinearities affect an oscillation, without using a hydrodynamic computer code, a purely real spectrum of eigenvalues is required in first order. The adiabatic eigenf requencies satisfy this requirement, but treat the outer layers in an extremely poor fashion. As the bulk of the star is quasiadiabatic, Buchler (1978) developed a formalism using twotime asymptotic techniques modeled after Cole (1968) and Kuzmak (1959). This method links the linear analysis with the nonlinear calculations by providing a consistent expansion in powers of the amplitude (Buchler and Perdang 1979). The requirement of periodicity is inherent in an expansion of this type, but the payoff is the capability to study a pulsation amplitude as it starts at an initial value and evolves to a final state. Additional problems that can be studied include modal selection, mode switching, and the feedback of the pulsation on the thermal evolution of the star, all of which are open questions in this field (Regev and Buchler 1981; Simon, Cox and Hodson 1980).
PAGE 10
6 In the simplest version of the theory, the quasiadiabatic results are recovered in the zero amplitude limit, but at the cost of introducing a cutoff in the integrals. When this was done, it became obvious that these layers were being incorrectly treated and that a more reliable method needed to be developed. A straightforward way of including these layers was presented in Buchler and Regev (1982a), where the structure of the eigenvectors was determined by the condition of constant luminosity above the cutoff point in the model. While giving good results for stars of the Beta Cephei class (Pesnell, Regev and Buchler 1982), cooler stars, like Cepheids , were not treated much better than with the simpler adiabatic approximation. The present work was undertaken to quantify the thermal behavior of pulsations. The results presented are encouraging and explain or clarify several phenomena found in stellar pulsations.
PAGE 11
CHAPTER III THEORY The structure and evolution of stars is governed by a set of four equations: the conservation of momentum, energy, and mass, and a flux equation. Other relationships are constituitive , such as the equation of state, opacity (see Appendix A), and nuclear energy generation rates. Written in a Lagrangian form which uses the interior mass as the independent variable, these equations are (see Clayton 1968, ch. 6; Appendix B below) ^= 4up (31) m dr Gm /^^Pzon ft o^ = 4Tir<^ r Â— = g(r,S) (32) o o m d t^ r^ T^. The symbols are defined here to be m = mass beneath the point of interest 7
PAGE 12
8 r = radius from the center of the star to this point p = density T = temperature P = pressure, P = P(p,T) S = entropy, S = S(p,T) L = luminous flux K = opacity, k = k(p,T) q = nuclear energy generation rate, q = q(p,T). We have assumed radiative heat transfer in the diffusion approximation and that all constitutive quantities are given by analytic expressions. To study this system, the equations are solved with all time derivatives set to zero to find an initial model in hydrostatic and thermal equilibrium. The equations are then linearized in r and S and combined to give two operator equations with two unknowns, 6r and 6S. Written in operator form these are 77= (^Â•^^(^^^A5rB.6S (35) d 6S _/5 k\.^^ +(14). 6S = C.6r + D.6S. (36) d t \d r/ V5 S It is well known that, with the appropriate boundary conditions, the A operator is of the SturmLiouville class and yields a symmetric matrix operator (Castor 1971; Buchler and Regev 1982b).
PAGE 13
9 The operator A is the adiabatic pulsation operator (6S=0) which can be written as (see Cox, J. P, 1980) dl l^,fp''[^^ II * I'^d? [(3^i*>^] ^^' ^ 1 ^ dt^' and in the Schrodinger form of this equation (Morse and Feshbach 1953, ch. 6) , the new independent variable has the dimensions of time and represents the acoustic travel time from a point in the star to the surface and is given by t (r) = / di acous r*^ r^(r) P(r) p(r) = ^/ dr/c(r) (37) where c(r) is the local adiabatic sound speed and R^ is the equilibrium radius of the star. The fundamental radial pulsation period is roughly given by t (R ), the sound travel time from the core to the acous core surface. As the acoustic time corresponds to the time for a dynamic response to occur, t (r) can be identfied as the dynamic time scale ' acous ^ as a function of position in a star. The nice properties of A come from the fact that it is of second order in the spatial derivatives (the time dependence is assumed to be exp( itot )), Looking at the D operator, the same second order derivatives are present. We now show that there is a representation where D has the same properties as A, and a natural definition of the thermal time scale will result.
PAGE 14
CHAPTER IV THE THERMAL OPERATOR For an optically thick medium, the thermal operator D can be written in a continuum form as ^ d 6S , Â„ . 6S d d L AC u 1 tAI d 6S/C T /\ / / \ ^S , ^ , . /d lnT\ V L(m)(4 KÂ„) pr+ Uta) eq T C^^ 6q y d m y dm (41) Here q =( ^ "^  , 3 Â— ^ = q( P ,T ), i.e. an initial model in thermal ^T \dlnTy ' d m ^ ^eq eq balance, k = L "!1 j , C ={r^) Â» and all derivatives are evaluated at the equilibrium densities and temperature. This form for D has been linearized at constant density (radius) and is strictly valid at points deep in the star so that t; > 1. We assume a time dependence of the form exp( at ) , and note that the operator is of second order in the spatial derivative. Therefore, an integrating factor, \x , exists that transforms D into a selfadjoint operator. This factor can be written as rÂ°,, s d InT , Â„ J (4 K ) dm m^"^ T d m li(m) = e " (42) where ra^ is a reference mass, here chosen to be m , the mass of the core 10
PAGE 15
11 core where nuclear burning occurs. Defining a new eigenvector as d In L I = 6S/C and \ = Â— ; ^^ , equation 41 becomes ^ V dm d K_ , An Al p cÂ„. c [x(4 .,,,)7^ k 3^ [ .(VS'J*dm. eq ' (43) Appropriate boundary conditions (which make 43 selfadjoint) are and i^^ = a ^(m,) + b 4^(m^) = (44a) dm * dm " C(mQ) = 0. (^^b) The surface boundary condition (44a) is equivalent to the assumption that at the surface there is no radiation incident from the outside (Castor 1971). In this form, 44a shows that the thermal readjustment time of the outer surface is infinitely short. As the weight function on the lefthand side is positive definite, and the function multiplying the derivative on the righthand side is negative definite, 44 is of the SturmLi ouville class of operators (Boyce and DiPrima 1969, ch. 11). Comparing 43 to the standard SturmLiouville equation, since A1 > and uTC /L > , all eigenvalues are real and V en /d InfT v eq positive. With each eigenvalue a is associated an eigenvector ^^ whose number of nodes corresponds to the numeric position of the {a^} in an increasing ordering of the eigenvalues starting with a^. In addition.
PAGE 16
12 the eigenvectors are nondegenerate, orthogonal with respect to the weight function jiC T/ L ^ , and can be normalized so that ( m^ is the total mass of the star) / fi T C I. l. dm = 6. ., eq Due to their SturmLiouviile character, these eigenvectors form a complete set for describing the 6S readjustment in a star. In combination with the adiabatic radial wave functions (Cox, J. P. 1980), a complete set for describing spherically symmetric disturbances of a star is found, provided that equations 31 through 34 are always valid. To understand the structure of the eigenvectors, the approximation of neglecting nuclear burning can be made. This reduces the equation to a version suitable for discussing thermal effects in a stellar envelope. Without nuclear burning the initial luminosity is constant and the analysis is halted at a mass ( m ) outside the region where this ' core ^ assumption is not valid. As the structure of the ^ 's will show, this is equivalent to assuming we are looking at "short" thermal time scales that do not probe the core. Radial stellar pulsations are an envelope phenomenon, so this assumption is justified. Neglecting nuclear burning and defining L^ to be the static luminosity entering the envelope from below, equation 43 becomes L^ vn dm n dm I d m J dm (45) Following Ledoux and Pekeris (1941), we multiply eq. 45 by C and
PAGE 17
13 integrate over the interval m^= m < m < m . . A selfadioint system *' core * J ^ satisfies a minimum principle and any well behaved eigenvector can be used to estimate the lowest eigenvalue, aÂ„. Choosing Cr)(Â™) 1 for all m, this estimate is m^ d <^ m^ d K^ J i^^"") TITdm Â„ J ^i T^rdm ^ dm m^'' d m oÂ„ = Lj. : = L^ r * / H T C dm ^^ ,, ,m Â•" *^ v (46) If the thermal time is defined as the reciprocal of aÂ„ , this gives the overall thermal time of the envelope as _ f, ,. t^, = r . (47) th ^ m^ d K^ * I a Â—: dm m ^ dm There is an apparent mathematical inconsistency with this expression. When the opacity derivative, <Â„, is constant in the envelope, the thermal time is undefined. However, comparing 46 to the expression for the adiabatic radial pulsation fundamental frequency (Ledoux 1965), ^ ''^^^ dT t ( 3 r^4 ) P ] dm o/ r'^ dm d P they have the same form only if Â— Â— = , so that the formula d P "^ ^1 for t^, is incomplete. If Â— = and Â— = 0, then cjo_= 0, but the th '^ dm dm variations of P are large in a star, varying many orders of magnitude, while r can vary only from 1 to 5/3 and the neglect of the pressure
PAGE 18
14 derivative would not be correct. The inclusion of nuclear burning in 47 remedies this situation and the net thermal time for the star is then th.nuc m^ d k L* 0^ I't^^^ K^q^) ^Â— ] dm which allows for the luminosity variations and shows that the apparent infinite value of t^, from 47 when k_ is constant in the star is not th T real. The terra allowing for nuclear burning will always give a finite value for t . . When k_ is constant in the envelope, the assumed th,nuc T form of C is a poor approximation and gives an unrealistic result. The positivity of t , is difficult to guarantee, but a simple case would '^ ' th.nuc J d K require that 4 kÂ„ qÂ„ > 0, Â— ; > 0, and Â—r Â— < , all of which can T T dm dm be satisfied in most of a star's mass, although in some nuclear burning regions, q can be much greater than zero, and the opacity variations are discussed below. The inclusion of nuclear burning as it relates to stability of stellar models is discussed in Hansen (1978) and Gabriel (1972). In this work, this area will not be discussed as the simplifications that result from the neglect of the core make the interpretation of the results easier. The opacity is of the Kramers' form ( k ~3.5 < ) near the surface but inside the ionization zones. Near the core, electron scattering will dominate, ( k Â— >0), so that kÂ„ decreases going outward in the star. Therefore, except in a small region near the ionization zones, t , will have only positive contributions in the integral. The mass of the ionization zones is small compared to the total mass of a star, so that the negative contributions to the integral
PAGE 19
15 will not cause a negative overall time scale. In addition, each point in the star is weighted by the ratio / s ~ I T(,m; \ T a(m) ' ' ' Â» T( which gives the interior regions far more importance than the outer layers, further reducing the influence of the negative portions of the integrand. Physically t represents the time necessary to remove the J. weighted internal energy by the radiative luminosity present, allowing for the position (and therefore temperature) dependence of the opacity derivative. Another thermal time can be found from equation 33 as the time necessary to remove the material thermal energy at the static luminosity present in the envelope. Writing this in a integral form (Saio, Winget and Robinson 1983) yields ^h n= J *T Cdm/L.. (48) 'th,p m Comparing the two definitions, there are several differences. The value t^j^ is counting only the thermal energy in the matter that makes up the star. In equation 47, the formula emphasizes the radiation 4 ( IJ. ~ T ) and the method of energy transport, as well as giving the interior stellar regions more weight than the outer layers. These differences are minor and are overshadowed by a major shortcoming of t^^. The definition in 48 can be generalized to calculate a thermal time as a function of mass by replacing the lower limit in the integral
PAGE 20
16 by the desired mass. Values for t , found are always positive and can ^ th,p be interpreted as time scales. However, when the same technique is done with t , , the value returned is not positive definite due to the th ionization zones near the surface. Mathematically, this generalization means that CCni') = 1> m'>m and C(ni' ) 0Â» ni'
PAGE 21
17 to transform eq. 43 to ,2 I + ( k2 w(x) ) y = 0. (49) d X 1/The units of J are (time)^ , so a natural definition of the overall thermal time of a stellar envelope is t^h= J2. (410) We note that the integrating factor [i does not appear in this definition. As ^ is an artificial quantity introduced for mathematical convenience, this is physically appealing. From the equilibrium model (B4 in Appendix B below) there is a relationship between L^ and dlnT/dm. Inserting this into eq. 410 gives m^ 3 C K i= [ / ( )'2dm ]2 (4lla) ^^ ' ""O (47ir2)2 ac T or, using the continuity equation (31), R^ = r(m_) and R^^ = r(m^). K^ J '^ 1/ ^h = t J ( Â— ^Â—y^p^]'(^iib) a c T^ This definition shows the effect of both the internal energy that must be dissipated and the method of dissipation. Rewriting 4ll(b) in a slightly different way R. 3 K T C 1/ c aT^
PAGE 22
18 or, defining U = T C , ' Â° gas V R^ U "0 "rad ^h i^y<'^^>'''^^' <*'Â» we have a interesting definition for t . . The photon diffusion length is present ( l/
PAGE 23
19 The internal energy of the gas can be written as 3 4 U = U , + U = T nRT + aT /p. 'as rad 2 so that U = ^ P /p = 4 P P/P gas 2 gas "^ 1 '^ and "rad = 3^ad/P = (1P)P/PWith these assumptions, the ratio becomes U /U , = 1/2 P/ ( 1 p) . gas rad ^ ^ ^ This mixture has a thermal time given by R 1 ^h = If U^ * t
PAGE 24
20 Here we see that the thermal time is a diffusion length ^* 1/ 1 = [ u/ (
PAGE 25
21 R^ 3 K C t,u(r) = [ / \ ^/^2p dr ]2. (417) r a c T^ In figure 41, we compare this thermal time with t , , also ^ '^ th,p' generalized to a position dependent form as discussed above (Cox, J. P. 1980; Saio, Winget, and Robinson 1983). The time scales are plotted versus the logarithm of temperature (with the surface of the star at the left of the diagram) for a model of BW Vulpeculae and in figure 42, for a generic Oepheid, discussed in chapter 5. Except for the constant separation, it is readily seen that the two functions agree quite well in the interior but diverge from one another at the surface. This effect is a consequence of t containing only the material thermal energy, tn , p which for a hot star is quite considerable and the drop to zero at the surface is due only to the integral formulation for t^, . In the th,p equation for t , the method of energy transport reduces to the grey atmosphere approximation at the surface, and the thermal time must therefore approach zero more smoothly for any class of star. This disagreement does not change the normal use of the thermal time. Several classes of variable stars are driven by ionization zones, including Gepheids (Cox, J. P. 1980) and ZZ Ceti (Starrfield et al. 1982; Winget, Saio, and Robinson 1982). For this driving mechanism to function, the ionization zone must coincide with the point in the star where the transition from (quasi)adiabatic (constant entropy) motion to highly nonadiabatic (constant luminosity) motion occurs (Cox and Whitney 1958). This transition zone is found at the radius, R , where t , (R ) ~ 11^, the thermal time is of the same order as the th,p tr 0'
PAGE 26
22 period of the pulsation. The expression for t is derived from a frequency, and the location of this region would be where a{R ) = t (R ) ~ cjj = 2Ti/n , which is analagous to critical damping. The constant separation in the figures is approximately In; so both formulations predict the same location for the transition zone, provided this region is below the atmosphere. Asymptotic Thermal Eigenfunctions Once the normalized coordinate has been defined, the asymptotic forms for the eigenvalues and eigenfunctions of eq. 49 can be calculated. Here, asymptotic implies that k^ = J^ a >w(x), for all values of x considered. Due to the surface effects, this is not a good assumption, but an important feature in the eigenfunctions can be seen in this approximation; so, it is discussed. The asymptotic form of eq. 49 is the harmonic oscillator and has boundary conditions a'y(l) + b' ^(1) = (418a) d x and y(0) = (418b) where a' and b' are a and b from eq. 44 transformed to the x coordinate system. Assuming y(x) = A sin( kx + 6 ) replaces the outer boundary
PAGE 27
23 condition by a transcendental equation of the form tan( k + 6 ) = b'/a' k (419) In the limit k > 1 , this equation has roots that correspond to halfinteger multiples of n (Butkov 1968, ch. 9) , k = ( j + i ) 't. a = ( ( y + 1 ) n )2/j2, and we shall let 6 = . The physical eigenvector can be found from the asymptotic expression as 5^(,, .,_^Cx)/[.^^(V?yf'*. (420, In Figure 43, we plot the x variations for the Beta Cephei model. For log(T)<5.7, the shape of the eigenvector is determined solely by the denominator in 420. In this star the integrating factor decreases outward and the net effect is an eigenvector that peaks at the surface. (Figures 44 and 45). The minimum in \x at log(T)"'4.6 is due to a maximum in k at that temperature caused by the second ionization stage of helium. The shoulder in y at this temperature is also from this minimum. Moving to a Cfepheid model, the situation changes considerably. In Figure 46, the thermal coordinate is plotted, and we see that x~l until log(T) >4.5, well interior to the hydrogen ionization zone. When the asymptotic (n=10) eigenvector is plotted (Figure 47), a large peak occurs at log(T)=3.9. This variation cannot be due to the cosine term in 420 as the argument is essentially 1.0; so, this feature must be due to
PAGE 28
24 the denominator, whose major variations come from the integrating factor. A large minimum in \i can be seen in Figure 48, where the logarithm of j. is plotted. This minimum corresponds to the maximum in 5 and is caused by the large positive value for k in the cool side of the hydrogen ionization zone. The integrating factor has a 4K^ temperature dependence of roughly T , and in regions where < > 4, \i increases outward instead of decreasing as it does elsewhere. This structure corresponds to the peak in the temperature eigenvector in a standard linear, nonadiabatic pulsational stability analysis. As the numeric results will demonstrate, any ionization zone will cause a local peaking of the eigenvector, but a region where k > 4 causes an effect that is best removed when calculating thermal effects in these stars, and plotting the resultant eigenvectors. There are secondary minima in (1 at log(T)~4.2 and 4.7 due to the ionization stages of helium. Any effects these have on the eigenvector y (Figure 47) are completely hidden by the peak from the hydrogen. However, it is not true that k^ = J^o > w(x) everywhere. At the surface, T Â— >Tq and k Â— >0, so that w(x) increases outward. While the interpretation of the "bumps" in the eigenfunctions as maxima in ic will be valid in the actual numeric results, the frequencies and eigenfunctions are not accurately represented in these expressions. With this knowledge in hand we now move on to solving the equations numerically and show the thermal modes for several stars.
PAGE 29
25 e c p.
PAGE 30
26 a Â•H 0) & W1
PAGE 31
27 r, o d^'O Oh'O 0Â£'0 Jltl'HQHOOa 'lUWySHi
PAGE 32
28 yOlDbJ 0NIlBy03iMI
PAGE 33
29 H il
PAGE 34
30 llOO Dhlj QOtf'
PAGE 35
31 0) o
PAGE 36
32 H 00
PAGE 37
CHAPTER V NUMERIC RESULTS In this section the eigenvalues and eigenvectors for two stellar models will be discussed. In each case, the initial model was found by using the model builder described in Appendix B, with mass, luminosity, effective temperature, and composition appropriate to the class of star the model represents. One model is of BW Vulpeculae, a hot supergiant of the Beta Cephei class of variables, and the other is a generic Cepheid, To find the matrix operator, equations 33 and 34 are differentiated at constant radius (density) and written in a form resembling that of Castor (1971): T ^T^ = AK2Â«(T6S). (51) d t Although the matrix AK2 is not symmetric, it is tridiagonal and a numeric analog of the integrating factor does exist (Wilkinson 1965), Applying this transformation to AK2 yields the D operator as a real, symmetric, tridiagonal matrix. The diagonal transformation matrix is defined as V9 /2(1) = 1 33
PAGE 38
34 /2(i) = /2 (i_i) [ AK2(I1,3)/AK2(I,1) P, 1=2, N, Vo and the matrix D is found by using i ^ as a similarity transformation D = \i^ AK2 (i 2 so that the spectrum of D is identical to that of AK2. In a convective region of a star, one procedure is to neglect the Lagrangian perturbation of the convective flux (Saio, Winget and Robinson 1983). This will have no direct effect on the definition of \i "^ as the two matrix elements, AK.2( 11,3) and AK2(I,1), are multiplied by the same fraction of the luminosity being carried by the radiation. Any gross V9 features of [i '^ , such as peaks in ionization regions, will still be present in a model where convection is permitted, although they may be extended over a larger region due to the smaller temperature gradient. The transformation has one side benefit, the matrices analogous to AG2 and AKl in Castor (1971), are much closer in magnitude to the elements of the symmetric matrices A and D and have smaller variations as a function of position in the model, resulting in better numeric accuracy. Real eigenvalues are guaranteed by the symmetric, tridiagonal form of the matrix operator (Ralston 1967). For negative definite eigenvalues (representing decays) the diagonal elements must be negative. This can be seen in a simple example. A negative diagonal element means that increasing the entropy (temperature) in a zone increases the luminous flux out of the zone
PAGE 39
35 more than it decreases thÂ» flux entering it. This stabilizes the zone as the additional energy is lost and the temperature decays to its initial value. A simple example can be used if the luminosity gradient had a local relationship. Equation 51 would then be written as TT=^^^o> so that S(t) Sq = S^e^*^ If a>0, the entropy decays, but if a<0 the entropy increases. The latter situation is possible in a nuclear burning core but not in the envelope. To preserve the SturmLiouville character of the continuum operator in 45 in the matrix operator, the offdiagonal elements must be positive. This condition is composed of two separate requirements. The luminosity exiting a zone must increase if the temperature of the zone above is decreased, and the luminosity entering a zone must also increase if the temperature of the zone below is increased. Both of these are easily satisfied if the opacity is independent of temperature or decreases with temperature, ( k < ), From the form of the continuum operator (45), only if < > 4 can violations of these conditions occur and these regions are precisely where the interpolation formulae for the luminosity are essential. In the initial model the temperature ratio is not allowed to exceed a given number in the ionization zones. This process limits the variation of the opacity as
PAGE 40
36 well and raises the upper bound for positive matrix elements. A first order expansion of the superdiagonal matrix element shows that for K^ < l/ln(T(Il)/T(I))~17 , these matrix elements are always positive. We have never found any negative offdiagonal matrix elements in any calculated envelope model. This process also limits the variations of the matrix elements so that the sum of the offdiagonal elements is less than the diagonal elements, permitting the use of Geschgorin's Theorem (Wilkinson 1965). The boundary conditions are incorporated into the matrix operator by giving certain elements of the matrix special functional forms. Assuming the outer surface to be radiative and given by the grey atmosphere formula is equivalent to assuming that this surface has an infinitely short thermal response time. Therefore, the analytic form of this boundary condition is that of eq. 44, even though the diffusion approximation is not very accurate here. For numeric computations this expression is unsuitable, as it is valid at a point exterior to the actual surface of a model, and we have assumed that the exiting luminosity is given by L(Wl) = 4Tta r2(N+1) T^^^ ^ (. x +  ) When this formula is linearized and inserted into the operator equation, the outer equation reduces to one of first order and closes the system. The inner boundary condition is 6S/C Â— > and is included by zeroing
PAGE 41
37 any elements in the matrix that refer to 5S/C (0). Again, this reduces the appropriate equation (here the innermost zone) to one of first order. As the thermal time of the outer surface is zero, and the overall thermal time ( ~ a^ ) is normally longer than a pulsation period, the spectrum of D must contain eigenvalues which bracket the radial fundamental adiabatic pulsation frequency. Here is the information that we seek regarding the transition zone. All parts of a star that can relax faster through a thermal process than acoustically will move in a highly nonadiabatic fashion. These regions of a star will be called the "sudden" regions and correspond, in some stars (see below), to the outer layers where Â— rÂ—^ = in a fully nonadiabatic stability analysis. dm The inner regions of a star have thermal response times longer than the fundamental period ( a < wÂ„ ) and are referred to as quasiadiabatic. Eigenvalues for this matrix can be found by a variety of means. For the low order modes iterating on the outer boundary condition is possible (Castor 1971), but the structure of the high order modes is such that this method fails. A general eigenvalue routine (Smith et al. 1976) is used to find all eigenvalues of D and the necessary eigenvectors are found by standard Gaussian back substitution. The eigenvectors of D and AK2 are related by the same similarity transformation as the matrices and the more desirable set can be calculated, depending on the effective temperature ( T ,), of the star under study. If T ^^ > 15,000 K, then the hydrogen ionization zone is absent and the C 's give a good representation of the thermal modes. When T ,, < 15,000 K, this ionization zone causes a spike to be present
PAGE 42
38 in ^ as shown in the section on asyiiiptotic eigenvectors in Chapter 4. To remove this feature, the eigenvectors of D are calculated and displayed. This results in a much clearer picture of what is occurring in the surface regions of the star. We cannot look at the long decay rate modes as we have removed the core and the effect it has on the eigenvectors. The lower modes can have significant amplitude in the core, but the higher modes are effectively isolated from the core even more so than the radial pulsation modes. This causes the difficulty noted earlier, the number of nodes in a calculated eigenvector is no longer an indicator of it's position in the spectrum; the size of the numbers causes the computer to underflow and the interior portion of the eigenvector is set to zero! Physically, this implies that the exterior can react to perturbations on much shorter time scales than the interior. However, the outer layers are linked quite strongly to the interior through the low order modes. When a heat pulse is inserted at the base of the envelope, the entire model is affected. The outer layers react to maintain thermal balance with the new luminosity from the heat pulse. The eigenf unctions of the diffusion equation exhibit an infinite propagation velocity for perturbations (Morse and Feshbach 1953, ch. 7), causing this behavior. The asymptotic eigenfunctions are trigonometric with an amplitude that increases outward and we can therefore speak of the thermal time of a point in the star as the reciprocal of the largest eigenvalue whose eigenvector has significant amplitude at that position. This may seem ambiguous, but as thermal modes can be used to expand the pulsation equations in the thermal times, this identification will enable us to study the thermal
PAGE 43
39 effects in an oscillation. In particular, the quasiadiabatic stability coefficient and pulsations that are dominated by thermal effects can be defined. Concentrating on those modes with have eigenvalues comparable to that of the radial fundamental pulsation frequency, we now examine two stellar models. BW Vulpeculae The physical parameters for this star are plotted in Figure 51 and the eigenf unctions 6r/r, 6p/p and 6L/L for the radial, adiabatic, fundamental pulsation mode are shown in Figure 52. The results of a linear, nonadiabatic radial pulsation analysis (LNA) are shown in Table 51. To find a realistic model for this star, it was necessary to integrate into the nuclear burning core slightly, but the thermal modes discussed below are not altered by this. 4 1 The fundamental frequency ((jl)Â„) is 3.9 * 10 sec , and the four thermal eigenvectors with eigenf requencies closest to this are plotted in Figure 53. Important features are the cutoff near log(T)=5.2 and the two peaks. The peak near log(T)=4.6 is due to the second ionization of helium causing a local maximum in k . At log(T)~5.1, a secondary peak associated with the first due to the use of the Rosseland mean opacity. The weighting function in this mean has a maximum at a temperature roughly four times that where the ionization occurs and causes a maximum in K at this temperature. This peak has been studied by Stellingwerf
PAGE 44
40 (1978) as well as Saio and Cox (1980), as a possible driving mechanism for these stars. The luminosity variations ( 6L/L ) for each mode are shown in Figure 54. Near the surface Â—r. = , showing the domination of the dm short thermal times in this region. Moving inward (increasing the temperature), the luminosity variation starts to decrease at log(T)=5.0, oscillates due to the nodal structure of C Â» and drops to zero by log(T)=6.0. These modes are completely isolated from the nuclear burning core (log(T)>6.5) as stated earlier. Figures 55 and 56 show the real and imaginary parts of the entropy variations as calculated from the LNA pulsational stability analysis. The peaks at log(T)=4.6 and 5.1 as well as a sharp cutoff above log(T)=5.3 can be seen in both components. Each feature has an analogy in the thermal eigenvectors in Figure 53. The modulus of the luminosity variations from this analysis is shown in Figure 57. The change from Â— ; = is present, albeit at a slightly higher dm temperature than in Figure 54. Some of this is due to contamination with low order thermal modes, but the significant amplitude in 6L/L (up to log(T)=6.0) is due to the adiabatic density fluctuations in this region. This model shows that the shape of 6S/C in an LNA pulsational stability analysis is determined by the thermal modes of that star. As always, the model is stable and, at present, we do not find any new ways of destabilizing it. The structure of the eigenvectors serves to demonstrate and clarify the physical reasons behind a transition zone and the splitting of the star in two regions, adiabatic (log(T)>5.4) and
PAGE 45
41 instantaneous thermal balance (or sudden) at log(T)<5.4. More interesting behavior can be seen when a model of a Cepheid variable is examined. Cepheid model This star is a generic Cepheid whose characteristics were obtained from A. N, Cox (1980a). The results of the LNA pulsational analysis and summary of the input data are found in Table 52. A graph of the physical parameters is found in Figure 58 and a graph of the adiabatic eigenvectors 6r/r, 6p/p and 6L/L is shown in Figure 59. The model is unstable in the first three pulsation modes and corresponds to a model of the star DL Cassiopeiae. The results presented here are illustrated very nicely by this star and are similar to those of other Cfepheids calculated. In Figure 510 and 511 the four thermal modes that bracket the pulsation frequency CjO_ , are plotted. These are the eigenvectors of AK2 , and we see almost no detail (or the individual modes) due to the spike in the hydrogen ionization zone. The origin of this spike was discussed in the section on asymptotic eigenvectors, and as it eliminates any other features present, we choose to remove it from the eigenvectors. This is accomplished by using the eigenvectors of D, as shown in Figure 512. These vectors have the same eigenf requencies as those in Figure 510, but much more detail can be seen. The bumps at log(T)=3.82 and 4.12 are due to the hydrogen ionization zone. They are the lower slopes of the spike, and the minimum at log(T)=4.0 is a result of modulation
PAGE 46
42 effect of the integrating factor being a maximum at this point. The helium ionization zone is visible at log(T)=4.6 as an inverted peak, but the second peak at log(T)=5.1 is absent. From Figure 46, where the thermal coordinate for this star is plotted, we see that the value for x starts to drop away from one at this point. In the interpretation of features for log(T)<4.7, the neglect of the trigonometric function in ^ is valid, but interior to this, the actual nodal structure of t will start to dominate, n Examining the luminosity variations which correspond to Figure 512 (Figure 513), the driving mechanism responsible for destabilizing this star is visible (see below). To accentuate this, the variations for a < 0) are shown in Figure 514 and 515 and a > co in Figure 516 and 517. In the low frequency results, the luminosity looks much like that of BW Vul (Figure 54); increasing the frequency to those in Figure 516 and 517, the decrease in 5L/L near log(T)=4.6 has become larger than any other part of the vector. Here is the physical reason for the driving zone of this star. For a star to be unstable in a pulsation mode, a region where ^ ^^^^ < must exist (Cox, J. P. 1980). This is seen if the work d m integral is writen in the form (Aizenman and Cox 1975; Buchler and Regev 1982b), = nJ (^31) (^) (j^) dm . The quantity T^1 is always positive and for the fundamental mode Â— ^ < (as Â— increases outward) (see Figure 52); thus, when P r
PAGE 47
43 '^ > the star is stabilized and for Â—3 < the star is dm dm unstable. In the BW Vul model (Figure 57) Â— ^ Â— is zero in the atmosphere and 6L/L generally decreases with increasing temperature. The thermal modes have the same basic luminosity behavior (Figure 54), showing that this star is stable. In the Cepheid model, the pulsational 6L/L (Figure 526) shows a definite region near log(T)~4.6 where Â— < . This outward decrease of 6L/L cannot be due to the d m adiabatic density fluctuations as they will always have 6L/L increasing outward. The thermal modes show a change in structure as the decay rate is incresased through the value a~(i). The change in the shape of 6L/L in this region is the cause of the destabilization of the Cepheid. Calculations of the BW Vul model with o/w < 1 and co/o < 1 do not show the same change. Examining Figure 518 and 519 where the thermal modes with o/u)~10, and Figures 520 and 521 where Wa~10, we do not see the same change in the structure of the thermal modes. This leads to the conclusion that the star's thermal structure prevents this type of driving mechanism from functioning. The location of the transition zone in each model should be discussed as well. The BW Vul model shows a definite change from adiabaticity to sudden behavior at log(T)=5.2. The Cepheid shows a less dramatic change, as the thermal eigenvectors do not show an extremely rapid cutoff at the transition. However, unlike BW Vul, the Cepheid has a region where Â—: < that begins near the transition. Also, it dm would seem more difficult to establish a transition zone for this model as rÂ— is not zero until well above the helium ionization. In fact, d m the luminosity variations from the pulsation analysis (Figure 526) do
PAGE 48
44 not have a vanishing derivative until above the hydrogen ionization zone. It is difficult to use the timescale arguments to locate the transition zone for such a model. The use of a Â« a)_ to partition the thermal modes is, in this case, a fairly reliable way to locate this region. For a < a)_ , the thermal mode has a 6L/L vector that is, in general, increasing outward and has only a small dip in the helium ionization zone. When a > oo., 6L/L is dominated by this region, until a > a)Â„ , when the amplitude of the eigenvector (5S/C ) becomes small enough that the mode is isolated from that region of the star. To see the benefits of the integrating factor, the real and imaginary parts of the LNA pulsational analysis are shown in Figures 522, 523, 524 and 52 5. In Figures 522 and 523, the eigenvectors of the AK2 matrix are shown, with the "spike" due to the hydrogen ionization zone dominating the graphs. When the eigenvectors of D are examined in Figures 524 and 525, the variations have been greatly reduced. In 522, ( 6S/ C ) , varies by two orders of magnitude and the imaginary part changes by almost three. The eigenvectors of D vary by less than one order of magnitude and the details of the helium ionization and other structures can be seen. The eigenvectors shown in this section have some numeric noise that deserves comment. The system that is solved here is identical to that used in the LNA stability analysis, and any noise in the thermal modes is also present in that analysis. The "spikiness" of the eigenvectors present in Figure 514 can be eliminated by putting more zones in the model (here N=250) but this puts the eigenvalue routine at a loss for solutions. The envelope mass can be decreased, but at the expense of the
PAGE 49
45 accuracy in the pulsation modes. The basic features of the thermal modes are well represented and accurately calculated. Interior of the transition zone, the adiabatic density fluctatiuons will overwhelm the small luminosity changes due to the entropy variations. The optimal numeric solution to an eigenvalue problem is found when the matrix elements are approximately equal in magnitude. The transformation of the mechanical variable introduced by Castor (1971), 6rÂ»XO = /DM2 6r , does this for the adiabatic wave equation. In a model where the shell mass is a geometric progression increasing inward, this transformation is equivalent to each zone having the same acoustic traversal time. The sound speed is roughly proportional to the square root of the temperature and the thickness of the shell is roughly proportional to the mass; so, the integrand in 37 has equal contributions from each zone. This transformation has the unfortunate side effect of decreasing the numeric niceness of the thermal operator, D. The thermal diffusion velocity (416) decreases as the sound speed increases. To optimize the matrix D, the mass of each zone should decrease inward, the opposite direction required of A. This tradeoff between the adiabatic and thermal operators will be important only if the mode is dominated by thermal effects. For quasiadiabatic modes, the interior is dominated by the adiabatic density fluctuations and the inaccuracy in the thermal modes is unimportant. We are working on ways of circumventing this difficulty so that both operators are represented as accurately as possible.
PAGE 50
46 > nJ ;:^ ? ^. en
PAGE 51
47 3
PAGE 52
48
PAGE 53
49
PAGE 54
50
PAGE 55
51
PAGE 56
52
PAGE 57
53 91 fi If.] iP.'02U'Z39'ti(l)DOl ONa lOMQ (jUdO>3D1 ICiHi.!) DOT >> u
PAGE 58
5h
PAGE 59
55 Â«
PAGE 60
56
PAGE 61
57
PAGE 62
58
PAGE 63
59
PAGE 64
60 B o f, n rvj 59'06 tc'os iiTq otif^t; /.(.'ic srÂ•^^
PAGE 65
61 S r^
PAGE 66
62
PAGE 67
63 ts u I
PAGE 68
64 r.fo ^ro
PAGE 69
65
PAGE 70
66 o
PAGE 71
67 TÂ•hi;UfsjÂ£6"<;''l( lOfX) iltiJ'A) A3/':.G o
PAGE 72
68 O TD e c c
PAGE 73
69
PAGE 74
70 c
PAGE 75
71 U 0) CO
PAGE 76
72 Table 51 BW Vulpeculae M=11.0 M sun log(L/L )=4.23 sun T ..= 25120 K eff nQ= 0^193 , t1q=1.9 * 10~^ n^= 0^147 , Ti^=2.3 * lO"^ n^= 0^118 , Ti2=1.5 * 10"^ Composition X=0.70, Y=0.28, Z=0.02 Envelope mass = 5.02 M sun
PAGE 77
73 Table 52 DL Cassiopeiae M=5.69 M log(L/L )=3.57 sun T3ff= 5848 K nQ= 7^87 , tIq = 6.8 * lO"^ n^= 5^64 , r^= 1.9 * 10~^ n^= 4^28 , 11^= 6.2 * lO"^ Composition X=0.70, Y=0,28, Z=0.02 Envelope mass = 2.76 M '^ sun
PAGE 78
CHAPTER VI DISCUSSION To understand the utility of the thermal eigenvectors, we shall look at various forms that the linearized pulsation equations (35 and 36) acquire when the ^ 's are used to expand them in the ratio of thermal time to pulsation period. Eliminating 6S one can rewrite 35 and 36 as an integrodif f erential equation: 0)2 6r = ( A + B Â—^ Â— C ) 6r. (61) lU) Â— u we have an eigenvalue problem for ur and 6r. Eliminating 6r , 35 and 36 can be written in a form suitable for studying the secular behavior of a star: ia)6S = (D + C Â— B) 6S. (6la) /,l2 _ A The ioj D term is a resolvent that can be expanded in the complete set of eigenvectors of D. Mopting the Dirac (or braket, see Davydov 1965, ch. 5) notation of quantum mechanics, this set is { ~o Â» C > } Â» and equation 61 becomes 74
PAGE 79
75 0)2 6r = ( A + ^^ B .^" ^ / C ) 6r . (62) Next, the mechanical eigenvector is expanded in the set of eigenvalues and eigenvectors of A, { w? ,  y . > } . We would like to examine the nonadiabatic effects on a mode that is primarily adiabatic, and this expansion will be written as 6r H I y > = I y. > + I a. I y. > and 0) = (0. + Au where (0.^ y. > = A y. > . 1 ' 1 ' ^1 Inserting this expansion and multiplying eq. 62 by < y  give (1 + ^^jmI'^jI^ ) 0)2 = < y.A y. > + < y. B l X I C y.> y 1 ' ' n n ' ' ^ 1 ^n ici) + a n . < y. B I X E, C y. > J Tl J 10) + o
PAGE 80
76 : T a. ) + ICO + a J n ^k ^j ^ \ : ^ a. . loj + a J (63) This expression is rather cumbersome, but following the normal quantum mechanical perturbation scheme, we evaluate the change in the eigenfrequency using the initial eigenvector. Choosing the adiabatic mode as the starting point, the quasiadiabatic correction can be written as < y, B I X ^ C y > c,2 _ ^2 ^ I i . ^ ^(64) X ^n 10) + a which is identical to equation (26) in Castor (1971). The normal quasiadiabatic approximation is a /w. < 1 , for all possible n, reducing 64 to 2 _ ^2 =1 < Â„_ iBCI y. > (65) OJ (Ji)' If = 7^ < y. BC y. > 1 10). ^1 ' ' 1 as y I C X C 1=1. ^n ' n n ' This expression for the quasiadiabatic correction can be further analyzed if Ao)/oj < 1. Then o)^ Â« o)? + 2o).Ao) and eq. 65 can be rewritten
PAGE 81
77 Aw = 7'< y. BC y. >, wr a result identical to the y. (equation 8, in Buchler 1983). Here, Aw is purely imaginary, as all quantities on the right are real, and generally shows that the star is stable (lAco < 0) as the diagonal elements of B are positive and the diagonal elements of C are negative. This recapitulation of the quasiadiabatic analysis serves to point out that by knowing the structure of the thermal modes, a better approximation can be made. Partitioning the thermal modes into two sets. and h = t o^, l^, a= 1, N^ : 0^0,. < 1 S, = I a , C , n = N + 1, N : w,/a < 1 }, eq. 62 can be rewritten as an eigenvalue problem assuming an adiabatic starting point only to determine NÂ„ u)2 I y > = a y > + ^ , io) + n=l n N B I X & C y > y :2Â— 2 . (66) M ^1 1(1) + a n=N^+l n The second term on the right is the quasiadiabatic correction, now with an appropriate cutoff. However, the cutoff is not in the spatial
PAGE 82
78 variable, as the integrals implied by the notation extend over the entire mass of the star. The correct cutoff is to stop the sum at the thermal eigenvalue whose magnitude is just below the pulsation frequency, as the sum over the set Si is doing. For this reason, the second term should be more correct in the calculation of this part of the stability coefficient. The third term produces a correction in the lowest order, and therefore cannot be treated as a perturbation. To see this effect, the denominators of the sums are expanded in a binomial series and eq. 66 becomes 0)2 I y > = a y > + 1 h Â«l ^n>< ^n 1^1 y> 10) ^ 1 + a /io) n=l n N B C^ X 5^ C y > ^ ^ ^^ a (1 + io)/a ) ^^^^ n=N +1 n n 0)2 I y > = A y > + ^ f BU XI ^ X C C y > + 10) '' , ' n 10) n ' ' N, I n=l I B C^Xl iii^X C^ Id y >. (68) n=N^+l n The first term is the quasiadiabatic correction already discussed and the second is the effect of the outer layers on the pulsation. The
PAGE 83
79 expansion is in (x^l a , and in the limit w/o Â— >0, or very short thermal '^ n n times; this sum becomes N 1 y Bl 5 >^ I C y >. n=N^+l n This term represents a zeroorder quantity and should not be neglected when the Initial eigenvalue is calculated. If the thermal modes had a thetafunction behavior at the transition zone, the set S2 would form a complete set in and of itself for describing the thermal response of the outer layers of the star on time scales shorter than o) . The actual eigenvectors have a fairly sharp fall at this point; so, a first try at including these terms is to write the sum over states as D' , where D' is the outer [(N N )*(N N )] submatrix of D that includes only those zones with I > N . Substituting D' for the second sum over states gives 0)2 I y > = a y > + 1 \ JÂ—Y^\l >< I C y > i CO ^ ' n ^n ' ' ^ n=l B(D')"^C y > + ia)B(D') C y >. (69) This equation can be expanded in the ratio of the proper time scales; the second terra on the right a = o/uÂ„ and the last two in a = ojÂ„/a . Writing a=a =a to give the simplest expression and arranging the equations in a hierachy of terms by order in a give
PAGE 84
80 oo2 I y^ > = ( A B(D') ^C ) Yq > (6lOa) 0)2 I y^ > + 20)^0)^ y^ > = ( A B(D') ^C) Y^ > + N ^ I B\ I X C C y> + i(i)_ ^ ' n n ' ' '0 n=l IWqBCD') ^C Yq >; (6lOb) where I y > = I Yq > + a I y^ > + a2 I y^ > + ... and 0) = u)_ + ao) + a2co + O(a^), To eliminate secular terms in the equation for  Yi >, the scalar product < Yq I yj^ > is set to zero. For convenience we shall set =l. This gives the initial eigenvalue and the first order stability coefficient as 0)2 = < y^ A B(D')"^C Yq > (6lla) 1=;^ r + liiit. n=l V2< Yq B(D')~^C Yq>. (6llb) At any point the matrix D' can be replaced by the sum over states that it approximates
PAGE 85
ol = < Yq \iA+ I bU^ > ^< ^^ C) Yq >, (612a) n=NÂ„+l n T 2oj^ n=l n=NÂ„+l a^ T n Equations 61 la and 61 lb agree with equations 25 and 39 respectively in Buchler and Regev (1982b), who derived them from physical arguments regarding the constant luminosity of the outer layers. It should be noted that the assumption of ; = is not valid over the entire dm section of the Cepheid model above the transition zone. The sum over states is valid, and comparing the thermal mode 6L/L to that of the LNA analysis (Figure 516), we see that the functional dependence of the LNA 6L/L is reproduced by the thermal modes. Using the infinitely sharp cutoff limit, the eigenvectors of 6lOa can be calculated in a straight forward way. Defining the error in Iyq > as E = { Ii I yd) 6r(I) 2f/2 where 6r is the real part of the mechanical eigenvector given by the LNA stability analysis and both eigenvectors are normalized to unity in their respective scalar products. The error in the eigenf requency is defined as
PAGE 86
82 w ' LNA' LNA These quantities are plotted in Figure 61 against the location of the transition zone, loe( T )Â• A sharp minimum is seen at ' Â° trans ^ log(T ) = 5.23, corresponding to the location of the actual '^ trans ' transition zone as predicted by the luminosity eigenvector plotted in Figure 56. In this model, the change is rather abrupt, falling a factor of 4 when log( T ) changes by 0.13 or +0.04 . The change w_ is not * trans ' Â» Â•' Â° as dramatic, reflecting the "minimum" principle that a quasiadiabatic oscillation will follow. When the same quantities are plotted for the Cepheid (Figure 62), it is obvious that the agreement is not much better than using the adiabatic eigenvectors. The error is reduced at a point where the transition zone is located, but not by a factor of 10 as in BW Vul. This results from the more complicated behavior of 6L/L in the Cepheid (see Figure 516). The thermal modes have a large peak in the ionization zones due to the integrating factor's temperature dependence but do not have the rapid decrease inward of those in BW Vul. Therefore, the approximation of replacing the sum over states by the submatrix D' is not a good one for this class of star. As the sum over states contains the necessary information about the luminosity variations in the exterior, including the integrals as in 612a should give much closer agreement. This problem is discussed in Pesnell and Buchler (1983). Comparing the initial eigensystem 61 2a to the adiabatic wave equation produces three differences:
PAGE 87
83 1) There is, at present, no guarantee of either the reality or positivity of the eigenvalues of either 6lla or 612a. 2) The eigenvectors < y^  and  yÂ„ > are not necessarily related by a transposition, much less being equal. The matrix A B(D') C is not symmetric but is of an upper Hessenberg form whose only advantage is being diagonally dominant. The eigenvectors come as right and left hand vectors and must be calculated separately. 3) There is no proven relationship between the number of nodes in I y^ > and the numeric position of oo in an increasing ordering of the spectrum of A B(D') C . However, for stars which are predominantly adiabatic, as most pulsators are, the eigenvalues are close to the adiabatic frequencies, the second part of the operator having little or no effect on A. The transposed eigenvector has a serious numeric problem. In the region of the transition zone, the transposed matrix has a structure such that a discontinuity appears in < yQ . The use of eq. 612a instead of 6lla will spread this feature over several zones, removing the discontinuity. When the transposed eigenvector is calculated using the LNA analysis, the feature is present but is much smoother. A change in the method of calculation is indicated, not a revision of 61 2a. A final observation regarding 61 2a concerns the limit N Â— >0. As a normal star is almost adiabatic ( N Â— >N) , the opposite limit must describe an abnormal star. For the condition N Â— >0 to be satisfied, the nonadiabatic regions of the star must be extensive, covering all of the star that is oscillating. These conditions may be found in the R
PAGE 88
84 Corona Borealis stars or other highly evolved stars where the envelopes are tenuous and have low gas to total pressure ratios. In these stars, trying to recover the adiabatic limit in one obvious way will lead to some rather strange behavior of w^. Introducing an arbitrary scale factor, X , into equations 35 and 36, which represents the ratio of the pulsation period to some overall thermal time, they become u)26r = A'6r + BÂ»6S (613a) iu)6S = \(CÂ»6r + D'6S). (613b) For most stars, when \K) the adiabatic limit is recovered as 6S is forced to 0. However, in the integrodif ferential form, splitting the thermal modes as before , (i)2 I y > = A y > + X A Â«l ^n><^n \^\y> ^ 1 + Xa /ico n=l n N B IX I C y > X I Â—. 2_,j . (614) n=N^+l n The value of N may depend on X for most stars , but for the radiation dominated atmospheres, the thermal response time is very short for all modes whose cutoff in amplitude is above the normal transition zone. Using the ideal gas as before, the atmosphere should have a low density
PAGE 89
85 gradient and opacity gradient so that the diffusion length (415) is large. If the outer modes have a /co > 1 for all n > NÂ„ , then N is almost constant during this discussion. In the limit of X Â— >0 ( N constant), eq. 614 is expanded using the binomial theorem as before. 0)2 I y > = A y > + N 1 T B C > ^ < ^ C y > + 0(\) , (615) n=N^+l n the B(D') C term does not disappear. For highly nonadiabatic envelopes, this term can dominate and change the spectrum so that it no longer reduces to the adiabatic case. This would, in part, explain the modes that have been named "strange" modes (King 1980; Saio and Wheeler 1982). This explanation of the strange modes does not require that a mode be near (in some sense) to an adiabatic mode. There are no theorems that prove that all of the eigenvalues of 61 2a are real, and starting from an adiabatic mode, the nonadiabatic corrections could move the frequency far away from the initial value. It is interesting that if N can be considered constant, the adiabatic limit is never recovered in some models. This property is just what a strange mode shows. Ilodes that possess the property of N being independent of X. shall be called "sudden" modes. They are characterized by a condition of thermal balance existing over a large fraction of the mechanical eigenvector's amplitude. The quasiadiabatic stability coefficient would
PAGE 90
86 be a poor approximation for these modes, as the second term in 612a will dominate co. . We must stress that the condition of thermal balance does not imply that the luminosity has a zero derivative in these modes,
PAGE 91
87
PAGE 92
88
PAGE 93
89
PAGE 94
90 r' GOVo(r 0) 13
PAGE 95
CONCLUS IONS We have shown that the thermal structure of stellar pulsations can be represented by the eigenfunctions of the diffusion equation. As this equation is of second order in the spatial derivative, it can be mapped onto a SturmLiouville operator, with the concomitant theorems regarding the reality of the eigenvalues and the form of the eigenfunctions. This transformed operator defines a thermal time scale through a new independent variable with the units of (time)'2 . This time scale agrees with that of earlier authors in the interior and predicts the same location for the transition zone in a star. The thermal modes that arise from the eigenvalue problem can be used to analyze the thermal effects in a first order perturbation expansion. The traditional adiabatic wave operator is replaced by an operator that includes the knowledge that the outer layers are in thermal balance , and that this does not imply that ^ ^ = in these layers. Partial spectra of the thermal modes are presented for two realistic stellar models, along with diagrams of the eigenvectors. Mathematically consistent expansions of the pulsation equations are presented that give the quasiadiabatic stability coefficient and the 91
PAGE 96
92 analogous expression when the outer layers are correctly accounted for. It is shown that a new type of mode, the "sudden", mode can be defined using these expressions. The driving mechanism in Cepheids has features in the thermal modes that do not show up in a star of the Beta Cephei class. This casts doubt on the viability of ionization driving being responsible for these stars' variability. In an aesthetic vein, a way of displaying the entropy or temperature variations for stars that normally have a "spike" in the hydrogen ionization zone is given, showing much of the detail in this region without the introduction of an arbitrary cutoff. This transformation also increases the numeric accuracy of solutions in the linear nonadiabatic pulsation stability analysis by limiting the variation of the elements of the coupling matrices.
PAGE 97
APPENDIX A THE ANALYTIC EQUATION OF STATE For the purposes of these calculations, an equation of state has been developed, based in part on a routine supplied by Arthur Cox (Cox, A. N. 1980b). In this appendix the various formulae and constants used are discussed. As our analysis requires the entropy of an ionizing gas, we first calculate the free energy ( F) as a function of the specific volume (v), temperature (T) , and composition (X,Y,Z). All thermodynamic functions can be evaluated by derivatives of F (Zel'dovich and Razier 1966, ch. Ill): entropy: S = (y^ ' v,N pressure: P = \^j^^ 0/5 (F/T)\ internal energy: E = T^l Â— :Â—z Â—  . \ A.N In order to calculate F, the approximation of an ideal, ionizing gas is adopted with the characteristics of the five species listed in Table A1. Values for the mass weighting of the elements has been 93
PAGE 98
94 adapted from Cox and Tabor (1976), and the atomic masses and ionization potentials are from Novotny (1973). For completeness, a term is added to include the effects of black body radiation. This is assumed to be that appropriate to isotropic radiation so that rad 3 The ionized fraction of each element is determined by minimizing the free energy at the given temperature and specific volume. This procedure leads to a Saha equation that is solved using a Newton Itethod 12 to a relative accuracy of 1 part in 10 for the electron density. When the electron density is known, the reciprocal of the mean molecular weight (in moles/gram) is n = 1/ti = X/1. 00797 + Y/4.0026 + Z/21. 02447 + n^. The relative abundance of each ionization stage, n^ , is calculated from the ratios of the appropriate partition functions, and the free energy is written as a sum over these species. F = RT { y. n,ln(m^^^u./n.) + n[ 2.5048 + ln(vT^^^) ] } ' ''1 i 1 1 1 'vT'* 3 where the constant term is ;
PAGE 99
95 1 +1ln(2ran k) ln(N.) 3 ln(h), 2 p A m. is the species gram molecular weight, and u^ is the temperature independent partition function that is assumed to the ground state degeneracy. 2 The expressions for the entropy (ergs/gram/K) , P (dynes/cm ), and E (ergs/gram) are S = r{1^ n.ln(m^^^u./n^) + n[ 1.0048 + ln(vT^^^) ] } + ^ vT^ a 4 P = nRT/v + ^ T E = Y nRT + avT** + l^ n^ x^ Derivatives of the pressure and the internal energy are found by analytic differentiation of the above formulae. The subroutine NEDS performs the calculations required to find the thermodynamic quantities as outlined above. A listing of this routine follows.
PAGE 100
96
PAGE 101
97
PAGE 102
98 Q 2! < 1.5 to OW ^ Wfct \ ,Â»Â» mÂ«c * ^. CN CU> CNCN X wca \* Hca uw > Â«>i + "s o Ott ~^ u pmm tJ I a* as: Â•Â•2; \ (JW C/10 ^ MfcJ U4"Â— CJ to H U* '> UJ i 23 Qj ,^+CQWW \m\ CN w CO w + "Â— ^tNOJOh*; U ^fr* C5C5fNN,\HH + H \ tn o X** .Â—Â— .wcQÂ•)^Â• CJ3 c/j ^^ cncj <Â«. xo^x::srNb^r^l^>4 4W ,.U WQ3 OXWI lfrHU^ s,.cj (MtJ ^^w 3c * u w lij xa oÂ»< wt/j+ crjw <: r3 X * "Â— Â— ' I * Â•Â» s .Â«:wca wu >^5 *Â— Â• o* * "if^ro uicj^iNtuÂ— mÂ«j: Â«Â• w js cn uij i_> u ^' Â— cnH"Â— \ ^2= m a H U S >X Oi* II \MWt3 Wpq UO"Â— XC^^^'a QtSOJUM Qa X z: H II II II II II II II >J II II II II II iS II ti MO Â«Â— 3(x4nsÂ»Â«OXrS XIS13 r(J UUUUCDIil HXH>^^^^E iriCH aofoooioi CO Â»m w n
PAGE 103
99 H 134 o
PAGE 104
100 en c^ w * tn ?> Â«: Â» (N H rOrtJtH H in (N ^OlO U OlTIQ >T^. Â•Of < OQ I r\HÂ»(N^jtQj Qinaao\H+ ca 3 to .m ^ ^ Â»+0=>3\ H c3 4vO I 'vO i.n rnm + * uun iaS * IS] + X :3 ro frn^fNZCN ^^ zrm I Â«Â»;> i uij Â«fct tsi* X* o >)'Z> % 0)cn >Â— Â»>w>iwa Ln WW + t<:oi +1Â— Teta< h lUiau):^.Mrsxaunx ro^o* tt* H Â• CN OO t Â• Kl "' I W** Si* M rO i t  rÂ« > \W* (N CJ => Â•HQ''Â— ^tnO ::^Â«CWWWÂ« O W03i O EH Oi WJ 2= W I . 3" * OIZ O O* O vO Â•Â» OiO Â•< W < C PfNit: isj tsun 0^ Â— cnohJo a^HtnoJsoo J Â» nil ii li'H^too tota n* to ^oz >!Â»Â«C3 DD fNZ:^ M* CNZi^ :3 < 25 Â•Â» H =3  II >Cir>'0!SIO II II SCittf'N^. MQIinpHIIII H^ll M 05 MÂ»rMrrMH  ill' M 11 Oil UH I utHW :d w fciP3D5a5C3 II i o in ii n i i n n tH PtH"Â— II iitiÂ«ioii II II ^MEHÂ» II Hbf ^ Â»<;Â«:Â«:Â«Â«; II acorinutaui (N OCUCUOjltiH^ 00="CK H OO O H en lAj W OS>S>S*S> ^^0'ESOUftiOH>HeH^OUf^OÂ« U UNN>'!>'!=50nnn3t3n)D!=3D>PÂ»HHHE' o o u U UUU OO
PAGE 105
101 a . Â• . < o o o wj ana X4 ^^ Ul Oi Ul EÂ» T rr,Â•a! Â— Â» Â» Â» H 3 Z J5 II II II M a: ^ e" U Â«: h: Â» or* CLi _5 > tn tOÂ»o* <_) rn 1X4 c/j CN I 3,^ tH o (/I O H I*< lilKTÂ— O^D 2:3Â» 3 ^J*+ * H aq W W ftJCCHfN Â» 31"^ Â» 23PÂ»25 HÂ— '3 hJ 2 fM*  > :Â» > W^NH in CtiMrt: Q 115 CdZ* ^ tH riHiX M M 3 3:S>M\ rMZ cCr^H Â«s:fcÂ»Â«x:5HMJia(i o *H^ I >tMsÂ»w:>0 hJ:3co f Â»coz: Â«:rn_* * * i aq ,. ^^ur^rN Â•a Â»z %Z5 jssn CO p^ uj Â•Â— '> rNCMmi^Â«i:ocN<^Utoa3 Tiji Â»3:Â« urn ^la; e rÂ» HO* . 3OUQJ CflHHX>^>HtsJ>3XSHNiaLn OZ<=S cqw* I _jg ^(sj _jr\^ tnQO:3Â»Â— 03 M OCX'' CC Â•Â— 6iH *?!;C3Â«j: r/5 O Â« tÂ« rN Â» irQ>Â— Â» H> I tOV5U II II II II II II II II II fn Â»C3W C50CH U OO OO OO as Â•r') V. U4 3<_)(lJ II U txjrÂ«Â»; UltviUM M .XI Â•Â—33 rS3 wrn i;ij Oi'"! MQS pq.Â—  h1,_ Â— ..,^ Â— ^.^^^^^^^.COZ II Oi II liJ^Â— Â— ~' H mÂ— Â».Â— Â«irzTr^z f\j ii:^ 3.Q(N^L)vor*r^cno> on rso** to sOÂ— o"Â— Â• vr>Â«Â— LDcro II WH zo Â— *Â— ~*~*~'ll *Â«:z II z; o >Â£i ~H *^ O N H>0 II SZZZZZZZZ IIO^HMII l 03 ?; UJ< uJkC WÂ«3; II II II xujH 3 tT* HlHHcHt'HHctlriX 3hHWtl W:3 'O MS Oi tH SUiH S Uj fri'C C=Jt3 z; ,CN(Nt=iu;03KC5iJ5 33C3a5tirt: ^aesti (^c" < M a3 O M u; u it = U HHHUUu ens: u HHQjCuajOiCXtOjtiiCuaiCOWQ cjmdj qs q 3=cuc/5 3eiA,co3Cvii/} < to cu o o o oa t*cNicN ro ct Â»Â— O rO TO u u uuuuu u u uuu
PAGE 106
TABLE A1 PHYSICAL PROPERTIES OF THE COMPONENTS OF THE GAS f/
PAGE 107
APPENDIX B THE INITIAL MODEL INTEGRATOR For the present calculations, an initial envelope model gives an adequate description of the variation of the physical parameters as a function of position in a star. To construct an initial model we must solve the equations governing stellar structure with time derivatives set equal to zero. With the appropriate boundary conditions, a solution can be found by integrating inward from the surface until the desired envelope mass is reached. In this appendix, the equations, boundary conditions and integration techniques are discussed. A spherically symmetric model in hydrostatic and thermal balance is governed by the equations d r^ conservation of mass Â— Â— = Â— ; Â— rr (B1) d m 4ti/3 d P G M ,^ Â„, conservation of momentum Â— Â— = r K'i~^) dm ,4 4ur conservation of energy : Â— = (B3) dm 4 T/ X // 2.2 ac d T energy transport LCm) = (4nx ; Â— :; ;; Â— ^' '^ 3k d m 103 (B4)
PAGE 108
104 where they have been written in a Lagrangian form that uses the interior mass as the independent variable. The energy transport is assumed to be purely radiative and the opacity is a fit due to Stellingwerf ( 1975a, b). We divide the model into N mass shells, or zones, with each variable assigned to the zone center or zone interface as shown in Figure B1. The mass of shell I (zonecentered mass) is called DMl(I), and the mass associated with interface I is DM2(I+1) = V2 (DM1(I)+DM1(I+1)). A mass derivative of a zonecentered quantity (pressure etc.) is interface centered and viceversa. All quantities are seen to be correctly centered with the exception of the opacity in equation B8. The presence of the zonecentered opacity in the equation for the interfacecentered luminosity necessitates the use of special 4 interpolation schemes to correctly center the Â— ^ Â— term. Several *^ ^ K d m authors have published formulae for this purpose and we have chosen that of Stellingwerf (1975a) for our calculations. The specific form of the difference equations is such that they are compatible with the nonlinear hydrodynamic code developed at Los Alamos National Laboratory (King et al. 1964; Cox, Brownlee , and Eilers 1966). This work forms the basis for further work in nonlinear behavior of models and this compatibility enables us to compare our results with those of fully nonlinear models. With this convention for the masses and derivatives, equations B1 through B5 can be written as
PAGE 109
105 Rd+l)^ R(I)^ _ v(I) .Â„_.. DM1 (I) 4ti/3 Â• ^ ^ P(I+1) P(I) ^ G M(I+1) DM2(I+1) "4, R(i+i)4 (B6) and L(IM) L(I) ^ _ DM1 (I) 4 L(I+1) = (4tiR(I+1)^)^ ^ <^>. (B8) J K a m The boundary conditions at the surface are P=Pq and T=TQ=0.841Tg^Â£ , where Pq is the pressure due to an external atmosphere and T ^^ is the effective temperature of the model. As these conditions are valid at a point beyond the actual model surface, they must be incorporated into the differenced equations by giving some elements of B6 and B8 special forms. In the momentum equation (B6), the pressure derivative is rewritten as P(Wl) P(N) _ P(N) . Qs DM2(Nfl) y Â• When f .= 1 the outer pressure (P(IÂ«1)=Pq) is zero, but there is a better approximation for f (see King et al. 1964). If one assumes that an infinite number of zones extends outward from the surface with each zone mass a constant fraction of the one interior, f . can be written as A f = " ^ \ (B10) A a + 1
PAGE 110
106 where a = DM1(N)/DM1(N1) is the constant mass ratio. Values for f are normally 0.10.2 . This modification moves the zero pressure boundary condition to the outermost zone as desired. For a consistent treatment, radiation pressure is added to this "atmospheric" contribution. When the boundary condition T = Tq is inserted into the luminosity equation (B7), an analogous expression can be derived for the temperature of the outermost mass layer. Starting with equation B4 , and assuming the luminosity and radius are constant, a grey atmosphere relationship for T(N) is found: t'*(N) = T^^^ I ( T + 2/3 ) (B11) where X = Jk dm/4rtR(Nfl)2 (B12) Given T , DMl(N), R(N+1), and the opacity, T(N) can be evaluated so that the boundary condition is satisfied. Along with the assumed form for the outermost mass shell, there are two additional relationships that need to satisfied to give a physical model. One of these is that the radius evaluated at the effective temperature be the photospheric radius given by the black body formula L = 4ti0 R J ^ /.. . (B13) photo eff This is satisfied by introducing an auxiliary expression
PAGE 111
107 R(Wl) = R . ( 1 + P ) (B14) photo and varying P so that R( T(m) = T^^p = Rp^^^^ The second condition requires the optical depth of the outermost 3 2 shell be at a chosen value, xÂ„ , where 10 < x^ < 10 . This is achieved by varying the outer zone mass until x( eq. B12 ) = x^ . Both of conditions are satisfied to a relative accuracy of one part in one million. Given eq. B5 through B8 and boundary conditions B9 , B11, Br12, and B14, we have a complete set of equations to describe the envelope. To integrate this set, the mass, luminosity, effective temperature, and composition are read in, guesses for DMl(N), f^ , and P are made, then values for P(N) and T(N) are found from B9, B6 and BU. This zone is automatically in thermal balance and v(N) is changed to satisfy hydrostatic balance in a one dimension Newton method. The optical depth of the zone, x(N) , is compared to x^ , and if they do not agree, the mass of the zone is modified following another Newton method until x(N) = xÂ„ to the desired accuracy. After the outer zone is converged, eq. B5 is used to find the inner radius R(N). The mass for zone N1 is a constant multiple of DMl(N), and this constant can be varied so that T(N1)/T(N) < y, where Y is normally 1.06. In a gas whose dominant component is hydrogen, the electron density changes rapidly during the ionization of hydrogen. This
PAGE 112
108 Increase in electrons corresponds to a rapid increase in the opacity as well. A star whose surface temperature is greater than 15,000 K will not have this difficulty as the ionization of helium does not produce as rapid an increase in the electron concentration. By limiting the temperature gradient as above, the opacity and opacity derivatives are not permitted to vary so much that they introduce further inaccuracies into the calculations. Once DM1 (N1) is known, the hydrostatic pressure necessary to keep DM1 (N) in suspension is given by eq. B6. A first guess for T(Nl) and v(Nl) is made and they are corrected by a two dimension Newton method. Convergence of T and v to the desired accuracy ends the iteration in this zone. This process of guessing (if necessary) the mass ratio, evaluating the hydrostatic pressure and zeroing the functions P P. , I and IL t, o~ L ^ ' eos hydro' ' eq.Bo star is repeated for each succesive interior zone until the desired envelope mass has been integrated. A subsidary loop controls the radius relationship. Using the initial guess for p , eq. B5 to B8 are integrated to a point where the temperature exceeds the effective temperature. By linear interpolation in the mass shells, the radius at T ^^ is found and compared ' eft photo eff photo' phot true, p is changed using a secant method, R(Nfl) is recalculated, and to R . , . While the condition  R(T ..) R H^^^l / ^r^hoin> ^Â° ^Â® ohoto eff photo photo the integration restarted.
PAGE 113
109 At the end of the integration, all equations are satisfied to the desired accuracy (here one part in one million.) This integration has been extremely stable for all physically viable models and no traps or inconsistencies are known to exist. The routine ATMHAT and support routines for the model Duilder are listed below.
PAGE 114
110 (JtSJ M Â« HO bi Â• Ha xaao W m HM OM B3hJ tqat DÂ« tOH O VO EH CkSO OS Â» Cfc(25 HCQ OH C< MS a w Â»^Z Oi 003 OK Â• OiQO OOlQj PSQ o< w ZOXM whs:: Â« HOatH .kM CGCMOSfcJ sjEtcu oca %H C3Â«(N 'Â»Mcn Â» >Â— H *KH M0 OH 'scqcKiu *i Â— bdU Â»CC1E>HUC3U *Â«i:u CuiÂ«:<.,DOO Â»s:ort<oÂ»o H Â»>jotoH Â»ciqBucaoÂ— CON VOÂ— QOIEU t,o&4SUQ Â» Â»,. Â— Â«T< w pj H \0 \H < XWNrHca uwssH cu TaLnoortyiw\33(aa;oH CQCQ33CQauOHUCOHHO ooooooooooooz ss:Â£3ESs::s:;s::x^is:z:ss OOOOOOOOOOOOH UUUUUUUUUUUUQ Â•
PAGE 115
Ill 25 04 O ^ Q M 2: w cu a: H V H S 53 o a MM MOKOco C5W a SSOi HH * w \ (Â»<>Â• H OO [J< Â«5P3 r>M (N,. HVJ WO P3S:: O tn PS Â«Â« H&< +0 *QO (( USn KiK rZjZ O"Â— H W ^W CJ^^: ooa Hn acM DU naH hw 11 P4< H"* M 0> HW (N# WH <Â« CL,0 a: Â•)Â»Â• 00 w ^~ 03 osu CnH * M Qj . tti Q II H CmOj ew Cm usa cao Zftfe* Uta W HD H wo wo H Â• (r> Â«0 X ?: h:i )*Â• wcn CQU WW \H w_ V3 e: < o3> ^^ xÂ»Â«: Olio zw UM ^* H* E VJ"Â« was :aH < M H O > HÂ«< E,Â» a."Â— O HO4 Hcn > CJP3 Â«WIH H Qfc. OW kjH MM w\ Â— a:w o <Â« M*II ttf Â— \ OiH tOZ! aw * rH Sh ^11 H< CJQ =rS O OÂ— CB > OH =3 H HQ ^J HH < ^ C5 a t/IU HU Oi A O 03 E HH 3 UH Â— ,. tne)tÂ» a M n SjÂ«;hq 0 C MOM O < E1 Â« 0tninrt:30 fMK UOS 3 Â•Â«: 03 S CU Â— HOWKHH o sw ,^ o a > OH Hw WO wz ,j3csi c^DH 03 z en m OQtnHEa: M Â»i5 mu wwcnH . > rO(MH'J O" a CHO O JÂ«<< Q r^rÂ«Â— I iias t ^jatN w < . akoh < or~ kxiJwu O I I cnOQ^OI 00 t sen Q W HH"* 2: QrHWHWfci QS Â— ~+OOUD3 tSCn, Â— .O >3t?a3 ww n tt Â»3;Q3cnMÂ«co ou q Â«Â«: Â•I > II II Â»E> M u 10 =e 11 o uuPx 3333 ifi 03 siMcn H :3Â«: US o Â«s:e= ans W cuij^^rt:Â«:DK^s:*' Â«:m ioch O tt>Â«WWQ3S<<Â«:QjW WWCfit. S H OESO''C't^HO H QCHQQW3HUÂ«H HHEM H tn UMnu^^caccju UUUU UUUUU UUU UUUUUUUUU
PAGE 116
112
PAGE 117
113 Ot n S3 H H a. s Q g H a, C3 J w 2: Qj Q * w ej <=> O 3 dQ ^ W * H Â» ^ H a;Q H\KQrÂ— M 2 2! M * O 03 rOJOO Â— ' * U Â« O M Â» ca Q4 CÂ»^ Z5 s a 5Â» CO Q < O Oj ,^ O HO Q4 .rH 04"Â— H tt, Q o # ^Q II Q40 Â» ^o s scosj: ^ rZ. a o 00 "ins HCJ3 CO riN. u Â«:Â»Â«jxa x: 03 o r 3^" ^ ^ ?. "^ ^'^ ==> 3Â«:Â«Q Qin H (X. . a, H* ^ Q^CN Q4 HÂ« to OH^M VN w ^ < \Q HO ^ \_ 01 in^itNQ. a o M H > C^ M ^ OH tl. H .Mr\* < J* w . ^ S'lT' S "" '^ ^ ^"^ " Â• u3o^+cN^ :<:^ ae Â• "^ >9j J^ Â•Â• \ seu ZH cNrc3Â»ca < * n Â» r^CO^* O. rufc, TatO^* NltX 02!QÂ«^Â— iCD  + )Â».0303ro fH ^iWo +*cacan a nae* aiTrjÂ«cQ rri Q5 "Â»03,,Â«Â— Oi O WHH03O Â• "Â—ca ,Â— ."^ Cli icC ginDH iÂ» H Â— C3wnÂ«(x^ to^ .W3Â»!)Â• Q3Â«MQfc< 03 OCbHD^H cJHO.mrt O _J4U5*^H 05 UfcHairo< KZ iniC* .oaH . *. b \ Â»*H* Â» Â— ,^ \ ,^ **e* Hw Â•r.caD M ,^=r o > ^S"^^ "^tl ^HHHHkJtcaocs,. WH ZOElD.HEt00 HQi W W w+Â„^.Q^:^ .3WC. .w Q +^Â« ^"^ 'llÂ§lÂ«WÂ§Â§H ^11 H Â£ ^S Stl Â°^^ HÂ«:QHHQ^ Ts ccM WW EH => era r^ n Â— ECO^MÂ— MMrcjirlÂ— UÂ— Ci II StH.JHÂ— >M SoSaqS Â— ll^^ g^ M COWQCKCdÂ— CU XSQ+ .HCOII . t/J H WQCCWÂ— Oj W 25 OXH^JOJCO rWtU rM SS" i,rs ^=='^^S"SÂ«22S2'^5" Â„r2, ^^ ^E^ <5miiqiiSwSs Â« g MH II a,""H tJCCHQHM M II OiÂ— M Qj S Â«HWQMlOM HS3 > OM^II fet M Â—OHS S. OH K 3 E "S Ei lO Â«;^L0CO Â— HHhI '^)'~Â— ~~^Ha5tOtO ^ HS Q Â— <: JÂ— Â»'M33Â«
PAGE 118
114 h3 00 Vi z ^M \{ tn CqpJQWO HtÂ«QO O4 > CQDUCmH .SDW Q CO g;MtnQ>' D to fN O Â«< HUÂ« O MtnW \ h^ B Q Â• 33 Mtucqwca tn QU (N ifl U ^ 0* DOOKW 1 E<<00 "' HHUCaOa CXU Â»H &Â• U OfMUU Â»<Â»2 <Â«8fc< >jHa:pt, M G Mro W _.^ cua HW so MUOM H H+ H ^~ Â»Â— . O iJ H EWU WZUQ I rtJH Â— . rÂ»530 Q <,> CBOnK Z! o< O S H H * oinuUQ 10 rtOcoHWo caaawjg o ^Â» tow o < M + fn(NÂ«si^ 12! >H Q tow H3HH * MO cGtoHxsa ODatoHO to to^* xoa C3t3o e \uuc3 Â• Â•Â«: ^'<>^ MWZk^rs Z HhJ W tOpviCjJ lO tO "'SHH * t^A >i QeaCOH^SHX tOM< >J O I M II 33tt C3HIJ'p^ zz Â«: io m ^ Hn< ai tv. 11 to tnnuueitH 110 (KEco tn< oiMuasto Â« H +11Â— w 11 to 11 to fN totoÂ«:Â«: II (N II to DjOhOcs se w h < as Â»lou II Â— aitO W tl Â«etO^Â— Â— ^Â— OifcH Qj r"Â— C3ZiDW=3Â« WW MQ H H Â»JE COM 03 On Pu tJ rf S= O WÂ«<;IX( W&jUkJOil WEtiQQ^>,_J
PAGE 119
115 'N OOOtU Vi MMSCU M II MM H kJ II i1hJ II II Z Ml) Â— >^^^ U HUHMMUUlJ^fci CO M > H H Â•< EÂ» M Â« (3 Â« H P4 s w Q Z kC W >) o E> w 61 I) D u ll Â•a: u Et
PAGE 120
116
PAGE 121
117 Â•TQ en Â• <Â•Â« Â«i^ Â»ata II ^ o ^ Â»ho x: 3\ yo ^ Qll Â• u\ i* ^ Â•to II w < Â» fn M II o> Z Ct,m M Q Â— t/JMco o (ri Â» j: sÂ«: a. MOW NO W H D< Q SH <Â» U "r" Â» Â» o o OM JU X tXII Dm [i4 rsi i) iÂ» Â«. Â— H Â•* Cj HC1CO wr%ir, Q Ns Â« ia= x:^H ow Â» &iQ\iu a w aq,.j^oi Ht^*HSi^ cua, ^H o t^ '^03Â» 7ia\ CU Â• 0'~Lr)B> H U \ O^MUOH. artnoWDs Â» Â• h >j o Â»Â«Â— cj33 Â« cu,Â— , THQSZt/lTr* Z33> OOiCO <Â» Â— VÂ— . axiioErf Â»'^ cuwa ^cQQjwi a*Â— .Qjtn?>'r^i: II t e > Â» Â»h mm&i Â»" tn Â» Â«:o Â«i w uÂ— Â»:3" Â• Oi""Â— H a4CN en *o!/i Â• Â«3:Ha3 oÂ—* *ci^ m* o Â«Â»cij a: w ^j B^ II u t vonrvi * >^h Â•XCaOCSW H cnZOlW Â— C5MW cur'' livoo'~o>a2W rriM a cli %qjÂ« > + >< vj Â»v)3pw M so >^* Â•=: Â« M^_ H\ xiinM fcEt Qaou Â«o H H X^***M Â— aiÂ«j; earn ro* Oi Â« HWW Â»30U II Â»Â«Â« 3M rj; ai< Â»Wrr)^^ ^^ccj Â•Â— Â»Â«s; OH=5XHtn hJ H^ ^J wo D3t<:OcgmÂ»TÂ— . 25Â»H MLr)Â«3:Â«etou >Â— oo Â«3:xa:s <: m,^r \ mcd ujQÂ»u5EuK H* ;>caH \% n^.,>_,^ j^^^ coHcsMO^xiSH o<: tamSfj "Coxh iOiOiH .5Â»s: rH Ea3 RjW EC30< OiOQO OTO C3 05 O tj Â«(JW [IJH.B CqcH HC3 M SH*Â— MHSVÂ—^Qj WQ3 xwM3:ttNs>r" HOW a:Â»MU33 Â« k ^EmO^Evoz Â» 2:h tnU3rQQ25Qj U4 25 55 C E H O O OH Q cy O) Â»^S E ~hj " Hm " HoSSS E O"Â— o *o 20S H o 00Â«j:a(CuajcuzaE OooujÂ«ckjw<:>_jwzoe h '^ s"^ s*^ W(^H U UUQC3C3C5a3WWH U UUHU KiCBU Â«Â« Wli,Htn Q ^
PAGE 122
118 Â• in Ln t * o I QO H Qj Q ino 2s Â» rn insr o '33 LnÂ«Â— M Â• O tvO HMon r* r>o UMos 1<Â«r~ D \ u3 con Q><^OH Â• I ^ OHrcam CO qvd a:Q\oi3" Â« >x)in MCU E * I ^Â«* ot fccco a I k o Â»Â«S Â•>( >Â— vx)0 \Q S Â» W % cr r>^(^J Â• Hrn H3'\H .OOWO wia o=ij\p~ Â•orÂ» "O I HO 2 Â» Â» LpMincOWQ 02:551 O<^^'^*~CL^00D3O Â»Â«Â»:Â«: Â•ucaocN *rMrtn Â• CH 0< Â«."a:ulÂ«i;CT> .Hf^Â» I toco M *Â«< Â» .SIlLn:^* "O Â«j>;Â£;w tut* *carMUÂ»Â»ooq "oai ojosrcses Â•M^Joacsin couw ca Â•jjrotn'QHocN 33 U\>CN Â«. tCO Â» Â• Â• JhJCu \ o"*a<^'^coaÂ»o <3Â«:co H\s rf^u "Oiz Â» Â•Â• WHO tOaDO:li<3TJa3 MNeQ^OOl/i nt<] 5riH2iÂ«a ooo ^'~ ootsJOPt^o:au2:CMa:ti,cGOOW25iÂ»E ftj ooo< ^mcN"* (jSlZQM HM2;HUHM3H3fc*Pt,QSWPaH H UUUO Q Q O t/1 O
PAGE 123
119 tn
PAGE 124
120
PAGE 125
121 o a ^ <: o >i VO W Â« ^ Q H > >< <Â» 03 Q 03 .^ Ol CQ Oj O Â« Q tl vo PH Â» S5 Â— 33 Oi Tu I) << cq ,. Â— Â» tÂ« * Â« H o^ a z pq Qj '^ Â» H O H Â«~O % Â» H Â« OiTko ca w Pu sÂ» C4S "Â— O > CU Q Kro O) fri Q \ CQ ^oÂ«a:oÂ«;a: Â« w Â« h^q h q Â—^DH^ Â»H Q4 z # r)* t^ CQ EÂ— )l*Q Â•05 H O M MM a Q UÂ«Â— Â•* H U W ,^ fc>n::>.HO Â«Â— aa z a, , . cattoM oq v^o %Â»fcW^ o Q!iiSna\uMMv;*a en Â»rÂ» Â— ^O Dj ^tiS 2: 6MH Â•Â» Â» ra aM ^Â«Â— CMO X* cuvO CN*Hocao2 ino en hh e" * h W'rU'SQ om oarEfi^QvO Â•Â» Â• Â» rsCS O f>H= W * .TOQ ttf W O^W 03 M CQ QinXCOWtq M M Â•Â•Â» Â«WHQ3<Â• OiOjOj <; N,t5,Â«rNCQ25 kiJS CUCD*. ,.4a&4a3 *WM Â— . U M^ ^JHOTH M M* M Â•O OMOOU >*M OOjCVUOSISIH n M OiQ W* Q OiOiÂ— T^W>^ 330* Â•* VJ fcMOi O rH > Â• TMMrn < ^tCfNca* OMCSHCTiZlM O MOiH "" W Â» X W^M UMl5J II <; Oi CU Â«: ciCiKjmwz^ciMM t Â«s; >^Q II .OQZH Hiieiomu m* hh ia c3h+ 11 vhmuieum'Â— *u tn Â» cutaua q m Â— Â«s;tÂ«H h iÂ»KHÂ«:fcj,Â» :3 w* m Â— n3cQajFiucacuajtnH< w or^^HOQ 25 .^3 11 11 Mfon m w ca HCufa,^ OWWWEt> CUh1 33 WQMCU II II M Â—.lIHOill ,Jr: Oi OI II *Â— 3 <; ZOaOiM iCSKKZ^TSz II OfÂ»H H 3 Â•^C II II fe H "W 0Â»Mu ptj SIZ2: H M Â»^OOOOOOII II MZ HO^M M^ ^.Q M3:mÂ—  Â— caCwM II \[ Uu TsaESCIS^E ttj klQtn >6HÂ«tM rM<; HCa^ttll MOi^JM M W KijsawKiajMMtn*^ Â— Hicawj'Â— UQfcHÂ— ' eÂ—^ hmosm ;5 ^hhz; _jmmw QOoocooQjOiODMfe rt:QD3Cijaaa3W> rnrua xkehmo puzwo ^b^c^x UCJUUUUH>!^MM U aiMQQmu UCJQ UUCJUOU UWQ5CJ UCupgO^ '0 OS uuu uuu u u uuu
PAGE 126
122
PAGE 127
123 II Â•"Â— ^ ^ H Â«Â•a; fcti o a H CQ H a SB ^ W Q O ^^ O Zi en U 2,^ (4 *^ I ~ ^Â», yj Â»^ W ^ TW ^ an < a Â«Â»; ez an Q H O H U QO ^ ^03 \ CD Oi _JS H H ttlÂ»^ Z orQ . ^^ ^ (3 on Â•0=tT&t oÂ»< i^ HQfjww <;tt, z: o ta^^ OS n Â«Â« *~WN, Â•HcQoq El Q u oo k'I ca ttJ o* Â«3:..Â£u^., Â— I H + II c> Â» o t/^a CQ nzEiHQHo tow M H Â«a ,^ OCO oo. Â«i:*. I t^ww I wo C'Z I ,_ HO 2Â«CKNOO EW Til ,11 cn03 Z HcnCQQHQ'>OfN >Â— ca C3 O S:3toow zaa o Â»" C3 1*^ fcfci ^w Â— Â« zcÂ» ii 11^3 :y uiui 'O'NOWWen 'Â— + Â•>**cncQ'y)tOlOWHM SÂ« H W II Â» (NTrÂ» Ota*Â— CJ ca CdTOS"Â— TQD W Q iO *'*Â— li^rniirns; II O II HU II II a! II <*rr)z  o Z P3 ZtnZ'HH ""Mil HO OW M II (N H MÂ«eQ3W << Â«t:oaÂ«a:oH < 11 ii tan < hoh Â•scn h MEOtttnEE HHÂ«!Mc<^Â«i;tt Â— a oj Â— EH "^M tt z Â«: H z zttHMocowQ WOWW0 0HWUCDCJCi
PAGE 128
FIGURE B1 VARIABLE ASSIGNMENT IN THE MASS MESH ZONE Il
PAGE 129
BIBIOGRAPHY Alzenman, M. L. and Cox, J. P. 1975, Astrophysical Journal, 195, 175. Baker, N. H. and Kippenhahn, R. 1962, Zeitschrift fur Astrophysik, 54, 114. Baker, N. H. and Gough, D. 0. 1979, Astrophysical Journal, 234, 232. Boyce, W. E. and DiPrima, R. C. 1969, Elementary Differential Equations and Boundary Value Problems (New York: John Wiley and Sons). Buchler, J. R. 1978, in Nonlinear Equations in Physics and Mathematics ed. A. 0. Barut and F. Calgero (Dordrecht: Reidel), p. 85. Buchler, J. R. 1983, Astronomy and Astrophysics, 118, 163. Buchler, J.R. and Perdang, J. 1979, International Journal of Quantum Chemistry, Symposium 13, 181. Buchler, J.R. and Regev, 0. 1982a, Astronomy and Astrophysics, 114, 188. Buchler, J. R. and Regev, 0. 1982b, Astrophysical Journal, 261, 301. Butkov, E. 1968, Mathematical Physics (Reading, Massachusetts: Addison Wesley) . Castor, J. I. 1971, Astrophysical Journal, 166, 109. Chandrasekhar, S. 1967, An Introduction to the Study of Stellar Structure (New York: Dover Publications). Christy, R. P. 1964, Reviews of Modern Physics, 36, 555. Christy, R. P. 1966, Astrophysical Journal, 144, 108. Clayton, D. D. 1968, Principles of Stellar Evolution and Nucleosynthesis (New York: McGrawHill). Cole, J. D. 1968, Perturbation Methods in Applied Mathematics (Waltham, Massachusetts: Blaisdell). Cox, A. N. 1980a, Annual Reviews of Astronomy and Astrophysics, 18, 15. Cox, A. N. 1980b, private communication. Cox, A. N. and Tabor J. E. 1976, Astrophysical Journal Supplement, 31, 271. 125
PAGE 130
126 Cox, A. N. , Brownlee, R. R. and Eilers, D. D. 1966, Astrophysical Journal, 144, 1024. Cox, J. P. 1958, Astrophysical Journal, 127, 194. Cox, J. P. 1960, Astrophysical Journal, 132, 594. Cox, J. P. 1980, Theory of Stellar Pulsations (Princeton, New Jersey: Princeton University Press). Cox, J. P. and Giuli , R. T. 1968, Principles of Stellar Structure (New York: Gordon and Breach). Cox, J. P., and Whitney, C. A. 1958, Astrophysical Journal, 127, 561. Davydov, A. S. 1965, Quantum Mechanics , 2nd edition, trans. D. ter Haar (Oxford, England: Pergamon Press Ltd.). Gabriel, M. 1972, Astronomy and Astrophysics, 18, 242. Hansen, C.J. 1978, Annual Reviews of Astronomy and Astrophysics, 16, 15. Henyey, L. G. and Ulrich, R. K. 1972, Astrophysical Journal, 173, 109. Iben, I., Jr. 1965, Astrophysical Journal, 141, 993. King, D. S. 1980, Space Science Reviews, 27, 519. King, D. S., Cox, J. P., Eilers, D. D. and Davey, W. P. 1964, Astrophysical Journal, 182, 859. Kuzmak, G. E. 1959, Journal of ^plied Mathematics and Mechanics, 23, 730. Ledoux, P. 1965, in Stellar Structure , ed. L, H. Aller and D. B. McLaughlin (Chicago: University of Chicago Press), p. 499. Ledoux, P. and Pekeris, C. L. 1941, Astrophysical Journal, 94, 124. Ledoux, P. and Walraven, Th. 1958, Handbuch der Physik, 51, 353. Lesh, J. R. and Aizenman, M. L. 1978, Annual Reviews of Astronomy and Astrophysics, 34, 203. Lucy, L. B. 1976, Astrophysical Journal, 206, 499. Morse, P. M. and Feshbach, H. 1953, Methods of Theoretical Physics (New York: McGrawHill). Novotny, E. 1973, Introduction to Stellar Atmospheres and Interiors (New York: Oxford University Press).
PAGE 131
127 Osaki , Y, 1982, Proceedings of Pulsations in Classical and Cataclysmic Variable Stars , ed, C. J. Hansen and J. P. Cox (Boulder: University of Colorado and National Bureau of Standards), p. 303, Pesnell, W. D. 1981, talk presented at Los Alamos National Laboratory, August 1981. Pesnell, W. D. and Buchler, J. R, 1983, in preparation. Pesnell, W. D., Regev, 0. and Buchler, J. R. 1982, Proceedings of Pulsations in Qassical and Cataclysmic Variable Stars , ed. C. J. Hansen and J. P. Cox (Boulder: University of Colorado and National Bureau of Standards), p. 216. Ralston, A. 1967, Mathematical Methods for Digital Computers (New York: Wiley). Regev, 0. and Buchler, J. R. 1981, Astrophysical Journal, 250, 769. Saio, H. and Cox, J. P. 1980, Astrophysical Journal, 236, 549. Saio, H. and Wheeler, C. J. 1982, Proceedings of Pulsations in Classical and Cataclysmic Variable Stars , ed. C. J. Hansen and J. P. Cox (Boulder: University of Colorado and National Bureau of Standards), p. 327. Saio, H., Winget, D. E. and Robinson, E. L. 1983, Astrophysical Journal, 265, 982. Shibahashi, H. and Osaki, Y. 1981, Publications of the Astronomical Society of Japan, 33, 427. Simon, N. R. , Odx, A. N. and Hodson, S. W. 1980, Astrophysical Journal, 237, 550. Smith, B. T. , Boyle, J. M. , Dongarra, J. J., Garbow, B. S. , Dcebe, Y. , Klema, V. C. and Moler, C. B. 1976, Eispack Guide , 2nd edition (Berlin: SpringerVerlag) . Starrfield, S. G. , Cox, A. N. and Hodson, S. W. 1980, Space Science Reviews, 27, 621. Starrfield, S. G. ,Cox, A. N. , Hodson, S. W. and Pesnell, W. D. 1982, Proceedings of Pulsations in Classical and Cataclysmic Variable Stars , ed. C. J. Hansen and J. P. Cox (Boulder: University of Colorado and National Bureau of Standards), p. 78. Starrfield, S. G., Cox, A. N. , Hodson, S. W. and Pesnell, W, D. 1983, Astrophysical Journal, 268, L27. Stellingwerf , R. F. 1975a, Astrophysical Journal, 195, 441, Stellingwerf , R. F. 1975b, Astrophysical Journal, 195, 705.
PAGE 132
128 Stellingwerf , R, F. 1978, Astronomical Journal, 83, 1184. Stobie, R. S. 1969, Itonthly Notices of the Royal Astronomical Society, 144, 485. Wilkinson, J. H. 1965, The Algebraic Eigenvalue Problem (Oxford, England: Clarendon Press). Winget, D. E. , Saio, H. and Robinson, E. L. 1982, Proceedings of Pulsations in Classical and Cataclysmic Variable Stars , ed. C. J. Hansen and J. P. Cox (Boulder: University of Colorado and National Bureau of Standards), p. 68. Wood, P. R. 1976, Monthly Notices of the Royal Astronomical Society, 174, 531. Zel'dovich, Ya. B. and Razier, Yu. P. 1966, Physics of Shock Waves and HighTemperature Hydrodynamic Phenomena , ed. W. D. Hayes and R, F. Probstein (New York: Academic Press). 21ievakin, S. A. 1963, Annual Reviews of Astronomy and Astrophysics, 1, 367.
PAGE 133
BIOGRAPHICAL SKETCH The author was born on January 19, 1957, in Wilmington, Delaware. He attended Claymont High School and entered the University of Delaware in September, 1974, receiving a Bachelor of Science in physics in June, 1978. During this time, he held an Undergraduate Summer Research Grant from the National Science Foundation under the supervision of Dr. Harry Shipman. In September, 1978, he entered graduate school at the University of Florida in the Department of Physics and began working under the direction of Dr. J. Robert Buchler. For three summers, 1980, 1981 and 1982, he was employed as a Summer Graduate Research Assistant at Los Alamos National Laboratory, where he worked with Arthur N. Cox. EKiring his tenure at Florida, he has held a oneyear research assistantship and taught undergraduate physics for four years. 129
PAGE 134
I certify that I have read this study and that in ray opinion it conforms to acceptable standards of scholarly presentation and is fully adequate, in scope and quality, as a dissertation for the degree of Doctor of Philosophy. J. Robert Buchler, Chairman Professor of Physics I certify that 1 have read this study and that in my opinion it conforms to acceptable standards of scholarly presentation and is fully adequate, in scope and quality, as a dissertation for the degree of Doctor of Philosophy. Thomas L. Bailey Professor of Physics I certify that I have read this study and that in my opinion it conforms to acceptable standards of scholarly presentation and is fully adequate, in scope and quality, as a dissertation for the degree of Doctor of Philosophy. Professor of I certify that I have read this study and that in my opinion it conforms to acceptable standards of scholarly presentation and is fully adequate, in scope and quality, as a dissertation for the degree of Doctor of Philosophy. Jg^s H. Hunter Professor of Astronomy I certify that I have read this study and that in ray opinion it conforms to acceptable standards of scholarly presentation and is fully adequate, in scope and quality, as a dissertation for the degree of Doctor of Philosophy. fy^^J^ l^lvv,^/" Pradeep Kumar Assistant Professor of Physics
PAGE 135
This dissertation was submitted to the Graduate Faculty of the Department of Physics in the College of Liberal Arts and Sciences and to the Graduate School, and was accepted as partial fulfillment of the requirements for the Doctor of Philosophy. August 1983 Dean for Graduate Studies and Research

