Fluid dynamics of circumstellar material associated with Be stars


Material Information

Fluid dynamics of circumstellar material associated with Be stars
Physical Description:
xvii, 174 leaves. : ; 28 cm.
Morgan, Thomas Harlow, 1945-
Publication Date:


Subjects / Keywords:
B stars   ( lcsh )
Stars -- Spectra   ( lcsh )
Stars -- Atmospheres   ( lcsh )
bibliography   ( marcgt )
non-fiction   ( marcgt )


Bibliography: leaves 171-173.
General Note:
General Note:

Record Information

Source Institution:
University of Florida
Rights Management:
All applicable rights reserved by the source institution and holding location.
Resource Identifier:
oclc - 13989425
System ID:

Full Text



Thomas Harlow Morgan



To KYC and HEN--both of whom pushed.


I am debtor both to the Greeks, and to the barbarians;

both to the wise, and the unwise.
Romans, i14

I would like to thank my parents for their

encouragement throughout my education.

I am indebted to Dr. D. C. Swanson for first interesting

me in physics as well as to Drs. S. S. Ballard, C. F. Hooper,

Jr., and F. E. Dunnam for help and advice in the summer of

1970. Dr. F. B. Wood introduced me to variable stars, and

Dr. K-Y Chen suggested the topic of this dissertation.

Mr. W. W. Richardson drew the figures that appear

here. The members of my committee, particularly Drs. K-Y

Chen and A. G. Smith, read and criticised the rough draft.

Dr. R. S. Leacock also read sections of the rough draft,

and in addition, offered many practical suggestions. Mr.

H. E. Nutter was of great help. I am indebted to them all.










* S S S

* 5* 9 5 9 S

* S S

. S S S .


Discovery and Taxonomy of Emission-Line
B Stars . .

A Brief History. .
Nomenclature . .

Characteristics of Be Stars .

Spectral and Luminosity Class. .
Galactic Distribution, Membership
in Multiple Star Systems and
Incidence Among Early-Type Stars
Rotational Velocity. .
Spectral Variations. .

The Struve Model . .

Origin . .
Line Profiles. .. .
Hydrodynamical Models for the Flow
Spectral Variations and their
Possible Causes .
Modifications of the Struve Model.

Major Topics of this Work. .


Background . .























The Hydrodynamical Equations and Notations. 15

Static Solutions. . 18

The More General Solutions. . 20

Form of the Circular Velocity Law 23

Final Comments ... .. ..... 27

OF Be STARS .... ............. 29

Disk Dimensions .. . 29

Theoretical Foundations of Hydrodynamics. 31

Validity of the Boltzmann Transport
Equation. . . 31
Relaxation Times . 32

The Hydrodynamical Prescription . 35

The Hydrodynamical Equations. 35
Navier-Stokes Equations and Similarity
Numbers .. .. 37

Importance of Viscous Term in the Equation
Of Motion . . 41


Preliminary Remarks . 44

Linearization of the Equations for Temporal
and Angular Dependence. . 45

Introduction of Time and Angle
Dependent Terms . 45
Selection of a Subscript-Zero
Solution. . 50
Considerations Concerning Approximations
Near the Plane. . 53
Additional Relations. . 56

Integration Over z .. .... .. 58

Further Restrictions on the Flow. 58
Relations Among the Flow Variables. .. 64




A Reprise. . .


Introduction to the Computations ..

The Equation for . .

The Consequences of Imaginary Roots.
Initial Values for the Parameters.
Machine Computation of the Roots .

The Case n = 0 .

Distinctive Aspects of the
Case .
Variation with ,r and A

The Case n = 1 .

Preliminary Remarks .
The Effect of c1 ..
The Effect of A. .
Extreme Cases .

The Case n > 1 ...

General Features .
Large n and the Effects of

Sample Flow Variable Solutions

Comparison with Observations .

S. 75
* 75

S. 78

.* 78
S. 79
. 93

. 107

S. 107
c and r. 107

* 114

. 115

Actual Stellar Values for Vs .
Qualitative Aspects of the Predicted
Spectral Variations. .
Comparison with Observations Reported
in the Literature .


Form of the Equation of Energy Transport

The Distinctive Nature of the Energy
Equation . .
An Approximate Form. ..
Relative Importance of the Terms of
the Equations. ...














* ,


A Linearized Equation. 129

Separation of the Time and Angle
Dependence .. 129
An Equation for the Subscript-One
Terms. .. 131

Relationship to Previous Work. 134

Reduction to Comparable Form 134
The Question of Consistency. 137

A Reconsideration of the Approximate
Energy Equation ........... 138

Non-Static Solutions and the
Adequacy of the Equations. 138
Final Comments .. 140


Stability of the Steady-State Flows. 142

Conclusions. . 144

APPENDIX A . . 148

Exchange of a Differential and an
Integral Operator. . 148

APPENDIX B .. .. . 151

Separation of the Equations and
the k Relation . 151


Evaluation of the Determinant formed
from the Equations Relating Unk
Cnk, and k . 155
APPENDIX D . . 157

The Program CALSOL .. .. 157

APPENDIX E . . 165

Pertinent Integrals. . 165



APPENDIX F ... .. ... ...... 168

An Equation for reff . 168

BIBLIOGRAPHY ..... ........ ... 171




Figure Page
1 Characteristic Central Star and Disk
Values . . ... 30

2 The values of the Coefficients Appearing in
the Energy Equation. .. 128

3 The Program Names for Important Quantities .159

4 Main Program .. . 160

5 Subroutine ALP1 . .. 163

6 Subroutine BET 1 .. .. 164


Figure Page

1 A Be star viewed equator-on and pole-on 8

2 The coordinate system . 17

3 The radial dependence of IPERI for n = 0
and for large c .. . 77

4 The radial dependence of IPER11 and
IPER2I for small values of cl 81:

5 The radial dependence of IPERlland
IPER21for larger values of 83
6 The radial dependence of IPER3 for four
values of cl. .. . 85

7 The radial.dependence IPER11 and IPER2I
for small values of r.. . 88

8 The radial dependence of IPER11 and
IPER21for large values of r 90

9 The radial dependence of IPER31 for
small values of r . 92

10 The radial dependence of IPER31 for three
values of r . . 95

11 The radial dependence of all three periods
for large cl and r . 97

12 The radial dependence of IPER1| andIPER21
for small values of c1 and r .. 99

13 The radial dependence of the three periods
for three different values of n 102

14 The quantity IPER1I1. as a function of n 104

15 The effect of large parameter values on
the periods when n is large .. .. 106



16 The radial dependence of ul (0.) and
u~ (0) for various times. 111

17 The radial dependence of U (90) -
u,1(270) for various time . 113


In order to avoid an excessive number of symbols,

extensive use was made of subscripting. The subscripts fall

mainly into four classes. First, subscripts are used to

denote a component of a vector. For the cylindrical

coordinate system in use here (see Figure 2), these subscripts

are w,(, and z. Second, the subscript c indicates a

dimensioned quantity. Third, a capital letter subscripted o

is a characteristic dimensioning quantity. Finally, a subscript

o appearing with a lower case letter indicates a steady-state

cylindrically symmetric term, and the subscript 1 indicates

a term with temporal and angular dependence.

The abbreviations for units follow the usage of

the Manual of Style for the Astrophysical Journal. Symbols

whose usage is confined to a page or a section will not be

listed here.

Important Symbols

p density
v velocity

P,p' Pressure, dimensionless pressure

u1 integrated velocity


al integrated density

n integer controlling the angular

S quantity controlling the temporal


A More General Listing

The numbers below are the pages on which the symbols

are introduced.

A, 127

B, 127

Crn, 64 C1, 61 C2, 61

DET, 68

FN, 40 f, 56

G, 16

1, unit vector i, imaginary number

k, 63

L, 61

Ms, 30 m, proton mass

Ne, 32 Nhe' 34

IPERI, 70 IPER1 86 IPER21, 86 IPER31, 86 PN, 127

p (as a parameter in the cross-sectional flow), 22 pl, 137

Q, 50 q 44

Rs, 29 RN, 40 r, 17

Snr, 64 SN, 40


Te, 29 T,. 29 t, time

-, 64.

Vs, 29 v (as e), .58 and 133 _e, .59 and 139

x; 67

y, 67

z, 17 Z0, 59

r, 59 T, 73

n, 37

k, 123

A, 60 A, 73 A' 136 A; 136

u, 30

Pd' 29

o; 124

<, 17

T, 51

w, 17

woc' 6


Abstract of Dissertation Presented to the Graduate Council of
the University of Florida in Partial Fulfillment of the
Requirements for the Degree of Doctor of Philosophy



Thomas Harlow Morgan

December, 1972

Chairman: Alex G. Smith
Co-Chairman: Kwan-Yu Chen
Major Department: Physics and Astronomy

Temporal and angular variations in the motion and

distribution of circumstellar material about typical Be stars

are studied by means of a simple hydrodynamical approach. The

validity of any hydrodynamical study of the dynamics of the

envelope is first examined, and the problem is found to lie

within the realm of hydrodynamics for typical disk conditions.

The problem is, however, close to the limit of validity of

that discipline. Further, Euler's equation may be used as

the equation of motion throughout most of the gaseous disk

surrounding the central star.

Each flow variable- is written as the sum of a known

steady-state axially-symmetric term and an unknown term

containing the temporal and angular dependence. The squares

of these latter quantities are assumed small. Using this,

the equations of motion and continuity are reduced to a set


of linear equations for the terms with temporal and angular

dependence. The known steady-state axially-symmetric terms

are taken from the literature (Limber, D. N. 1964, Ap. J.,

140, 1391). This solution is a static isothermal one. The

geometry of the problem suggests several approximations which

may be used in the linearized equations; these approximations

introduce two parameters. Solutions can be found to this

somewhat approximated version of the linearized equations.

Since the boundary conditions have not been treated, a fully

determined set of solutions is not possible. The dependence

on both time and angle is, however, explicit.

The angular dependence enters the solutions in an

extremely simple form involving an integer n, but understanding

the temporal dependence requires computation of the periods

characteristic of the temporal variation at a given location

in the gaseous envelope and for a given choice of the two

parameters discussed above. The IBM 360/65 at the University

of Florida was used to perform these computations. Under

most conditions the characteristic periods are found to be

real. Further, their dependence on n and w is w 3/2/n where

w is the distance from the center of the star in the equatorial

plane. These calculations predict temporal variation down to

fractions of an hour, and in qualitative agreement with

relevant observations.

The third equation of hydrodynamics, the energy

equation,is examined insofar as it can be in the absence of

an effective treatment of the radiation-fluid interaction.


This approach..leads to inconsistencies with the results

of the study of the equations of motion and continuity.

It is concluded that in the absence of some accommodation

of radiation effects, no equation for energy transport can

improve the existing solutions.

Finally, the stability of steady-state flows is

examined. The steady-state flows are not absolutely

unstable, but the most general solutions appear to be

periodic in character.




.Discovery and Taxonomy of Emission-Line B Stars

A Brief History

In August 1866 Father Secchi observed HB.in emission

in y Cas. Over the next half-century Campbell, Frost, Miss

Cannon, and others observed hydrogen in emission in other

stars of what would now be called spectral type B as well

as in members of the remaining early spectral types. Miss

Cannon appears to have been the first to observe variations

in the emission features in one of these stars (Cannon, 1898).

Pickering (1912) published the first catalogue of early-type

emission-line stars; at that time there were 94.

Further work was in two directions: more emission-

line objects were discovered and individual stars were

studied in detail. The catalogue of Merrill and Burwell

(1933, 1943, 1949) contained over 200 B-type emission-line

stars. Detailed spectral studies led Struve (1931, 1942)

to interpret the spectral features characteristic of Be

stars in terms of emission and absorption in a gaseous

disk-shaped envelope surrounding a rapidly rotating central

star; this interpretation has not changed substantially


In recent years both long term variations (Limber, 1969)

and shorter term variations (Lacoarret,. 1965; Slettebak,

1969; Hutchings,: 1969 and 1970;. Hutchings et al., .1971) in

the characteristic emission-absorption features have received

attention. The most recent listing of Be stars (Jaschek

et al., 1971) contains 1900 objects.


The basis of modern systems of spectral classification

is the system adopted for the Henry Draper catalogue (Cannon

and Pickering, 1918-24). The original Harvard arrangement

of the spectra was in order of decreasing ratios of intensities

of the hydrogen Balmer lines to the intensities of a number

of other absorption lines. The Harvard observers soon

realized that a rearrangement of the original alphabetical

order into the present spectral sequence--O, B, A, F, G,

K, M--represented an ordering in temperature with decreasing

temperatures from 0 to M.

Among the HD objects typed B, those possessing

hydrogen-emission lines were marked with an additional

lower case p-- for peculiar--following the B. Modern

classification uses the notation Be for such stars, that

is, stars with B-type spectral features plus emission in

one or more of the Balmer lines. The suffix p is now

used only to indicate the presence of annomalously strong

or weak lines. If a Be star shows sharp metallic absorption

lines superposed on the Be spectrum, the star is called a

shell star.

Characteristics of Be Stars

Spectral and Luminosity Class

Inthis work the term "Be star" applied to hydrogen

emission-line stars whose underlying spectral types are

early of mid B (BO-B7) and whose luminosity classifications

are either II, IV, or V. This is somewhat at variance with

much of the literature, which considers Be stars to be

basically main-sequence (luminosity class V) objects. But

Mendoza's (1958) comprehensive luminosity class survey of

Be stars in the Merrill and Burwell catalogue simply placed

the stars in two groups--very luminous objects (I, II) and

"less" luminous objects (III, IV, V)--on the basis of a

pronounced dichotomy in his observational material. In

fact, strong doubts (Underhill, 1966) have been expressed

as to the existence of a systematic distinction between MK

luminosity classes IV and V in this early spectral region.

Absolute magnitudes of Be stars determined on the basis of

cluster or association membership appear to fall about one

magnitude above the zero age main sequence. Finally the

appellation shell will be similarly restricted although

objects with A andeven early F spectral classes are some-

times called shell objects.

Galactic Distribution, Membership in Multiple Star Systems,
and Incidence Among Early-Type Stars

Mendoza (1958) took for his sample all the stars in

the Merrill and Burwell catalogue whose spectral types fell

in the range from BOe to B2.5e and whose luminosity classes

were III, IV, or V. He found that the Be .stars populate in

a roughly uniform manner the spiral arms as delineated by

the 0-B associations. Be stars are found in 0-B associations

but are not strongly concentrated in them. They are to be

found both in galactic clusters (X Persei, Coma) and also

in binary systems (a Tau, 4 Per, V 367 Cyg).

Rotational Velocity

Be stars are the most rapidly rotating class of

stars. The most recent catalogue of observed rotational

velocities (Bernacca and Perimotto, 1971) lists 189 objects

in this spectral class whose rotational velocities have

been measured. The mean rotational velocity of the group

is 262 km/s, where the observed rotational velocity is

the true rotational velocity times sin i, i being the angle

between the axis of rotation and the observer's line of

sight. On the assumption that the axes of stellar rotation

are directed randomly in space one finds the true mean

rotational velocity (Huang and Struve, 1960) to be 344 km/s.

Traditional methods of line profile analysis used to

determine the observed rotational velocities are now believed

(Collins, 1970; Hardrop and Strittmatter, 1968) to underestimate

the actual values by as much as 40 percent. If this belief

is true the Be stars are rotating at large fractions of

their equatorial break-up velocities. Balmer emissions

occur in approximately 15 percent of all B-type spectra

(Underhill, 1966).

Spectral Variations

Spectral variations occur in most Be and shell stars.

Following the terminology of McLaughlin (1961) the spectral

variations fall into three major classes:

(1) appearance or disappearance of a shell

absorption spectrum,

(2) changes in the ratio of the intensity of an

emission line to the neighboring continuum (E/C),

(3) changes in the relative intensities of the red

and violet components of the Balmer emission

lines (V/R).

The classic example of the first variation is

Pleione (28 Tau), which has been.observed spectrographically

since 1888. Pleione was first classified as an emission-

line star, but by 1905 it showed no emission features and

possessed a normal B spectrum (with large rotational velocity).

In 1938, emission lines were present; a short time later

the star displayed a shell spectrum. By 1954 the emission

features were lost. Recent observations (Sharov and Lyutiv,

1972) indicate that Pleione is entering a new cycle of


E/C and V/R variations, though strictly periodic

only when they reflect orbital motions of stars in binary

systems, .show much quasi-periodic behavior. The time scales

for these variations may be on the order of a few years or

as short as days or fractions of a day (Slettebak, 1969;

Hutchings, 1969 and 1970; Peters, 1972).

The Struve Model


In 1931 Struve suggested that the emission-absorption

features that are the definitive characteristic of Be stars

are formed in gaseous circum-equatorial rings surrounding

(and ejected from) B-type central stars that are rotating

at their break-up velocites. Later, Struve (1942) was able

to present the model on the basis of more detailed spectral

considerations. He also proposed a circular velocity law

for the equatorial disk,

V (WoC
Vc( oc)

where Wc/Woc is the radial distance from the central axis

to any point in the disk divided by the stellar radius and

Vc (Woc) is the rotational velocity at the stellar equator.

Figure 1, which draws its inspiration from Hack (1970),

illustrates the model's ability to explain the major features

of Be spectra. The first case is a Be star whose axis of

rotation is perpendicular to the line of sight; the second,

one whose axis of rotation is parallel to the observer's line

I 0
w EO
0 4) 0 0)

00 00
.I0 0
-rl rA 4

S*ri 1O 0
>d O l
J U OU 0

m 0 ,M 0

ra 4J A
H 0 a) 0l

*0 0

%9 OR '

M0 *t4 ) HI

o Ok 0c:0
S- m H I

m aa 4

tra po+


I ^ ^
\ 0,'4' 0
SI ,

> S

t a,

r i

- 1 O1
/"^ C
1 5
C: O5

of sight.. A schematic profile of the expected emission at

one of the Balmer lines is also shown for each case. The

imposition of a small expansion velocity on the large

circular velocity gives a V/R asymmetry. That the disk is

indeed cixau\-equatorial was given credence by the work of

the Burbidges. They analyzed (Burbidge and Burbidge, 1953)

high dispersion spectra of six stars whose emission features

were like the second case above, often called "pole-on."

They found that, although a very thin layer of gas did

appear to exist above the polar regions of the photosphere,

the bulk of the gas was in or near the equatorial plane.

Another aspect of the Struve model is that it explains the

differences in Be and shell spectra as due to differences

in the density and extent of the envelope.

Although the Struve model provides a qualitative

explanation for the spectral features found in Be stars

and shell stars, detailed calculations of line profiles


(1) a model for the radiation field of the central


(2) a model which gives the physical state,.

density, and velocity of the gas at all

points in the disk,

(3) the changes in the radiation field resulting

from its passage through the gaseous disk.

The work of Mihalas (1965) and of othexsnow

provides the first requirement. The second and third

requirements. have not been so satisfactorally achieved and

are not even truly separable.

The Struve hypothesis says nothing about the

mechanisms) that cause the formation and destruction of

the disk--that is, changes from B to Be to shell star and

vice versa--nor does it say anything about the mechanisms

that produce other spectral variations.

Line Profiles

The hydrogen-emission- spectrum is a recombination

spectrum formed in a disk the various portions of which

can move relative to one another. A great deal of work

on hydrogen emission in a stationary atmosphere has been

done (Kogure, 1959a, 1959b, 1961). Sobolev (1960) and

Rothenberg (1952) have studied the formation of the hydrogen

lines in spherically symmetric envelopes expanding with

constant velocity. Marlborough (Marlborough, 1969 and 1970;

Marlborough and Roy, 1971) has calculated Ha profiles using

a 6-level hydrogen atom. He used a model for the distribu-

tion and motion of the material due to Limber (1964) which

will be discussed below.

Hydrodynamical Models for the Flow

Determination of the density and velocity field

throughout the disk is a hydrodynamical problem. In the

first such study (Limber, 1964),a steady-state axially-

symmetric disk was examined for a parameterized form of

circular velocity law and in the absence of any radial

velocity. Limber treated isothermal envelopes, polytropes,

and envelopes which have a specified temperature law. For

each case he was able to calculate the density at any

point in the envelope. Later, Limber (1967) studied the

more general case in which a radial velocity component was


Spectral Variations and Their Possible Causes

A rough division of spectral changes in Be stars

will be made here to facilitate the discussion; the divisions


(1) long term changes (time scale longer than

10 years) which characterize the appearance

and disappearance of shell or emission spectra

as in Pleione,

(2) medium term changes (1 to 10 years) particularly

V/R changes with time scales in this time


(3) short term changes ( a year or less).

There is a certain inevitability to the long term

changes; Crampin and Hoyle (1960) argued that any initial

magnetic fields present at the formation of the disk would

render it unstable in 100 years or less. Limber has

suggested that there may be a non-uniform transfer of

momentum from rapidly rotating inner regions of the central

star out to the photospheric regions resulting in turn in

non-uniform mass loss. Limber (1969) analyzed the shell

phase of Pleione in terms of a time-dependent mass flux

from the central star's equatorial regions to the disk. He

found good agreement with the observations for a flux whose

time dependence showed a slow increase to maximum followed

by a steep drop to zero.

McLaughlin (1961) reviewed what are called here

medium term variations. He observed that the most difficult

problem concerned V/R changes. McLaughlin (1962) discussed

the major attempts to explain the quasi-periodic V/R changes

of the sort seen in i Aqr,and he concluded that only a

suggestion credited to Otto Struve was consistent with the

observations. Struve's suggestion had been that the ring

or disk was elliptical in shape and that the V/R variations

result from a line of apsides rotation of the disk. There

has been no recent examination of this idea.

The most comprehensive work on short term variations

is Huang's (1972a). Huang showed that an asymmetry in the

disk will produce spectral variations observable after time

spans as short as a fraction of a day; these variations can

persist for as much as a year before the disk's differential

rotation destroys the asymmetry. Hutchings (1969, 1970)

suggested that photospheric fluctuations or "clumpiness"

in the disk may be responsible for short time changes.

Modifications of the Struve Model

It now appears that several aspects of Struve's

original model must be modified. First,. Struve's. model

requires that the equatorial regions of the photosphere

be rotating at the break-up velocity. Observations suggest

that Be stars' rotational velocities lie below the break-up

limit. Huang (1972b) recently argued that the ejection of

matter is due to both rapid rotation and a "temperature-

dependent instability." Huang feels that the "temperature-

dependent instability" may be the mechanism studied by Lucy

and Solomon (1970) in their work on stellar winds in early-

type giant atmospheres; this work argued for a mass outflow

due to radiation pressure from mid-UV resonances in silicon,

carbon, nitrogen, and sulfur, these being in high ionization


Limber and Marlborough (1968) as part of a general

discussion of the physical processes at work in Be stellar

envelopes, reanalyzed the work of Struve and Wurm (1938) in

order to show that the velocity law,

V Oc( oc)

was in better accord with the data than the circular

velocity law originally proposed by Struve.

Major Topics of this Work

A principal contention of this dissertation is that

the existing steady-state, axially-symmetric hydrodynamical

solutions explain the key elements of the Be star phenomenon.

More exact and comprehensive explanation of the body of

observational material will require solutions to the hydro-

dynamical equations which possess time and angle dependence.

The approach here will be to view the expression for each

flow variable as the superposition of a detailed but in some

sense second-order term onto the time-independent axially-

symmetric solution.

Chapter II, a literature review, discusses the

previous hydrodynamical work on Be stars. Chapter III

studies the extent to which the motion and distribution of

circumstellar matter can be treated as a hydrodynamical

problem. Chapter IV introduces the dichotomous view of the

flow discussed above into the hydrodynamical equations;

this chapter then presents the mathematical reduction of

the original set of equations for the time and angle

dependence to a set of solvable equations. The next chapter

(V) discusses calculations based on the solution to these

equations. In Chapter VI the energy equation is considered,

and Chapter VII, after a discussion of the relation of the

results of this work to turbulence and stability arguments,

summarizes the work.



The steady-state axially-symmetric solutions of the

hydrodynamical equations offer an acceptable description for

the overall phenomena presented by Be stars and shell stars.

These published solutions are the starting point for the work

which is the principal subject of this dissertation. For

this reason, and because the discussion will outline the

hydrodynamical approach to the motion of circumstellar material

in Be stars, the literature on steady-state axially-symmetric

solutions will be reviewed in some detail. After a brief

discussion of notation, the derivation and properties of

static envelopes will be discussed; then the permissible

solutions with non-zero radial velocity will be considered;

finally the evidence available for choosing among the

various possible circular velocity laws will reviewed. There

will be a short summary at the last. The order of presenta-

tion of the three main sections represents the order in which

these ideas were developed by Limber and Marlborough (Limber,

1964, 1967; Limber and Marlborough, 1968)..

The Hydrodynamical Equations and Notation

The hydrodynamical equations used here are the simplest

ones that could be chosen; they are Euler's. equation and the

continuity equation. In dimensioned coordinates the former


(v c + + so
-PC + vC v) = -Pc !-- ) -cVP
a c Gt

and the latter takes the form

+ V (pc v) = 0 .

The justification of the use of Euler's equation will

be a principal concern of the next chapter.

The coordinate system used here is shown in Figure

2. The three coordinates (w, *, z) may appear with or

without the subscript c. The presence of the subscript

signifies that the coordinate in question is dimensioned in

cgs units; without the subscript, it is a dimensionless

variable in the system of dimensionless coordinates which is

developed in the next chapter. A similar procedure is

followed for the flow variables and time: velocity v,

density p, pressure P, and time t. The equations already

presented in the first chapter are consistent with this

practice. In cases where non-cgs units are used, these units

will be explicitly stated.
In the body of this work symbols will on occasion

undergo some change in definition; for instance, a variable

symbol may be required to absorb an integrating factor. This

Figure 2

The cylindrical
use here.


-system in

system in

is done only when necessary to prevent an undue proliferation

of closely related variable names. When questions arise

about the meaning of a symbol, the best recourse is a check

of the Key to Symbols. Each change in usage is noted there

as is the first page on which the new usage occurs.

study of

took the


Static Solutions

Limber (1964) attempted the first hydrodynamical

the material outside rapidly rotating stars. He

central star to be rotating at the break-up

and made the following assumptions:

1) viscous and magnetic terms in the equation of

motion can be ignored;

2) radiation forces are either negligible or

includable through the use of an appropriately

reduced stellar mass in the gravitational


3) disk self-gravitational effects are negligible;

4) steady-state conditions prevail;

5) only axially-symmetric flows are included;

6) the Z component of the fluid velocity (perpendicular

to the equatorial plane) is zero;

7) the radial component of the fluid velocity is


8) the circular velocity has the parameterized form

v v O( )

where the new symbol a is a parameter whose

value lies between 1/2 and 1;

9) the density is an explicit function of the

pressure alone.

The first three conditions determine the nature

of the equation of motion which can be used--Euler's equation.

The third approximation, neglect of self-gravitational terms,

is quite reasonable. The validity of the first and second

conditions will be discussed in later sections. The use of

Euler's equation as the equation of motion is common to all

hydrodynamical studies of Be star envelopes. The justification

for the use of Euler's equation or, indeed, any hydrodynamical

equation under disk physical conditions will be the subject

of the next chapter.

The next group of approximations, (4) through (7),

reduces the equation of motion to the form

GM v
VP = p { V ( --s 1 1W
c c r W a
c c

where v is the circular velocity (< component of the veloc-

ity) and i is the unit vector for the radial direction in

cylindrical coordinates. Conditions (4) through (7) also

result in solutions which will identically satisfy the

continuity equation, the second of the three equations

needed to specify completely the problem. Since violation

of (5) will show up observationally as variations with time,

both (4) and (5) are good approximations to the extent that

they describe phenomena with either no time dependence or

long term time dependence. Near the equatorial plane

condition (6) is a good approximation, while condition (7)

is clearly an idealization which will be partially removed

in subsequent work.

The radial dependence of the circular velocity

law can be at least roughly determined from dilution factors

for shell absorption lines, from the widths of absorption

lines, and from the widths of emission features. It is,

then, possible to choose a parameterized form on the basis

of the observations. Calculations were done for three values

of a : 1/2, 3/4, and 1. In all cases the circular velocities

were required to match the stellar rotational velocity at the

stellar surface.

The last condition, that the density is a function

of the pressure alone, circumvents the failure to treat

the last of the three hydrodynamical equations, the equation

for energy transport, in the disk. The condition imposed

may be alternately expressed as requiring that surfaces of

constant pressure and of constant density be identical.

Limber is then able to calculate the density

throughout the disk for isothermal envelopes, for polytropes,

and for envelopes in which the temperature though constant

along any one surface of constant pressure may vary from

surface to surface.

The More General Solutions

The static solutions--solutions in which the radial

velocity is zero--are,. strictly speaking, nonphysical for the

same reason that static solutions are nonphysical for the

solar wind (Parker, 1958). Limber extended his analysis to

non-static solutions (Limber, 1967). Three assumptions of

the approach presented in the preceding section are modified;

they are (6), (7), and (9). These three are replaced by the

following less restrictive conditions:
1) only flow near the plane is considered;

2) the radial cross section of flow must be specified
(in parameterized form);

3) the temperature distribution throughout the disk
is known.

As in the static case the ideal gas law is used. Euler's

equation takes the form

av v 2P
e c = G [ 1 1 c
v vi r = GM a-J-
v)C ai jC s 8 C P d
c c c c C

where the new symbol v is the radial component of the

velocity. The equation of continuity appears in the inte-

grated form

v,(w)p(w) A(m) = vw(0o)p(wo) A(wo)

where A(w) is the radial component of the cross-sectional

flow. This form of the continuity equation is similar to

that seen in solar wind theory. The quantity A(()/A(w ) is

analogous to the term
F (r) = (c b

in the elementary solar wind theory (Parker, 1963) where b is

a parameter (usually given a value near 2). Limber chose,

instead, a one-parameter family of cross-sectional flow


A(w ) = M ( p)
oc oc

where p is a parameter. The meridional projections of the

stream lines are straight lines running from a point

W = pwoc

in the equatorial plane and sloping slowly away from the

equatorial plane.

Limber went on to establish the types of solutions

that are now permissible. Valid solutions must be subsonic

below the star's surface and approach interstellar values

at great distances from the star. The only solutions which

satisfy these two boundary conditions are those whose radial


1) first decreases with increasing w,

2) reaches a minimum as W continues to grow,

3) then increases as w continues its increase,

4) becomes supersonic at sufficient distance from

the central star.

Limber examined the range of validity of his earlier

work, the calculations for static envelopes, in the context

of his new results. He found that the density calculations

from the earlier approach agreed quite well with the new

results over most of the near-equatorial regions of the disk.

At large distances from the star the old solutions no- longer

agreed with the new calculations; the distance at which the

breakdown of the old approach occurreddepended on the exact

physical conditions, and geometric dimensions chosen for

the disk.

Form of the Circular Velocity Law

The solutions developed by Limber in his 1967 study

contain two parameterized functions, the cross-sectional flow

expression and the circular velocity. Of these two the circular

velocity is the more crucial; the solutions are fairly insen-

sitive to changes in p, the parameter which appears in the

cross-sectional term. Limber and Marlborough (1968) examined

the physical conditions and observational evidence for the

circular velocity's dependence on a throughout the disk.

The equatorial break-up velocity is only 1//2 times

the escape velocity. This means that if a small element of

matter were perturbed with a small radial velocity so that

it began to move outward, and if the element were not further

disturbed, then the element would eventually fall back onto

the star. Theavailable evidence indicates that shell material,

rather than falling back onto the star, is continually escaping

the system; there is mass loss. Limber and Marlborough call

this the problem of "support." The physical processes which

act to provide the "support" profoundly influence the radial

dependence of the circular velocity law.: Limber and Marlborough

divide the mechanisms for support into two classes: "direct"

and "centrifugal." The first type, "direct" mechanisms, are

radially directed forces which act on any fluid element to

overcome the difference between the gravitational force on that

element and the centrifugal force. These direct forces would

appear as additional terms in the equation of motion. The

"centrifugal" classification refers to mechanisms which act

to transfer angular momentum out into the disk so that the

centrifugal force always balances the gravitational force due

to the central star.

Direct mechanisms do not act to affect the angular

momentum of a small fluid element which is moving from the

stellar equator out into the disk; the angular momentum of

the fluid element is conserved. In this case

Vc`c -1

Quasi-static motion of a small fluid element away from the

star with centrifugal forces present requires that

v GMs
W1 2
c Wc

at any point in the disk; hence

Vc=" "T
c a !1/2c

There are four direct-support mechanisms available;

they are:

1) thermal support,

2) turbulent support,

3) radiative support,

4) magnetic support.

Limber was able to show that the first possibility, thermal

support, requires temperatures and densities throughout the

disk that are totally at variance with present knowledge

of the values of these quantities. Turbulent support requires

a highly turbulent flow whose effects on spectral features

should be observable but have not been reported.

Limber and Marlborough considered two possible

radiative support mechanisms, electron scattering and photo-

ionization of neutral hydrogen. The net outward-directed

radiation force on a small fluid element due to Thomson

scattering was calculated for typical disk electron densities

and was found to be negligible in comparison to the gravitational

force acting on that element. They noted that for a given

neutral hydrogen density there was an upper limit to the radial

force that could be generated by photo-ionization of neutral

hydrogen. Maximum momentum transfer in the radial direction

from photo-ionization occurs if the recombination following

the photo-ionization emits a photon backwards toward the

source of the original photon. To establish an upper limit,

Limber and Marlborough assumed that all the flux from the

central star shortward of Lyman a and incident on the disk

was absorbed and back-emitted. On the basis of both this

assumption and other considerations, these authors calculated

an upper limit for support due to photo-ionization. This value

was far less than the value needed to provide radiative

support. While both these calculations are dependent on

estimates of typical disk values for the physical parameters,

they are so much smaller than the value needed for direct

support that they may be omitted from further discussion

until such time as the best estimates of these disk quantities

change dramatically.

Finally radial support from magnetic forces was

considered. The same investigators used elementary considera-

tions to show that the field strength, H, required was

H > 75 gauss.

They then showed that fields of this sort in the disk

would disrupt it on a time scale of days. The observational.

results of.the disruption of the disk, loss of emission lines

or shell absorption lines, show much longer time scales.

Thus they concluded that magnetic fields of the required order

are not present in the disk.

Analogously there are four mechanisms which can

act to transfer angular momentum:

1) thermal viscosity,

2) turbulent viscosity,

3) radiative viscosity,

4) magnetic viscosity.

Limber and Marlborough used the term "viscous" to denote

any phenomenon which acts to transfer angular momentum. They

found that two of these possibilities, the first and third,

were insignificant under disk conditions. However, either

small scale turbulence or small (< 5 gauss) magnetic fields

could provide the needed angular momentum transfer. They

concluded that the observational evidence was consistent

with either interpretation. It is important here to note

that the amount of turbulence or the magnetic field strength

required for angular momentum transfer is far too small to

provide any direct support.

Limber and Marlborough noted that the work of Hynek

and Struve (1938) represented the only attempt to draw

quantitative inferences from observations about the form of

the circular velocity law and that this attempt had considered

only the 1/1 behavior. The data were reanalyzed to see if

the 1/w form gave a more consistent interpretation of the

data. The conclusion of this reanalysis of the data was

that the form of the circular velocity law associated with

centrifugal support was in better accord with the data than

the old Struve form.

Final Comments

Steady-state solutions have been applied to the

shell phase of Pleione (Limber, 1969) by allowing the

envelope to be at any moment very near steady-state. This

will, of course, work only for long term phenomena. Limber

suggested that the end of the shell phase represented the

result of increased matter outflow not matched by sufficient

energy flow into the disk; hence,a collapse.

Two results of Limber's studies are of sufficient

importance to the remainder of this work to bear restating.

First, over most of the disk, static isothermal envelope

calculations agree well with the more exact theory including

a radial velocity. Second, the circular velocity law is that

for centrifugal effects. The static isothermal envelopes

with a circular velocity law can be presented in closed form

and do not require the specification of a cross-sectional

flow parameter.



Disk Dimesinsns

Table 1 gives the values used here for the physical

and geometrical quantities which characterize the star and

its envelope. The mass (M ) and surface temperature (Ts)

of the star are those of BO dwarf. The radius is larger

than one might expect for such a star and represents acknowl-

edging the evidence that the star is somewhat evolved and

rotationally distorted. These quantities are those

commonly found in the literature. The equatorial rotational

velocity (Vs) represents a slight departure in that the

rotational velocity was set at just slightly more than half

the equatorial break-up velocity. The disk temperature,

which is taken to be the electron temperature (Te), the

mean molecular weight (p) of the disk gas, and the average

disk density (Pd) are typical of values found in the


From the emission and absorption features at

hydrogen and helium resonance lines one can conclude that

in a typical Be star envelope a large percentage of the

hydrogen is ionized while most of the helium is not.



Characteristic Central Star and Disk Values


Stellar Mass

Stellar Equatorial Radius

Surface Temperature

Equatorial Rotational

Disk Temperature

Molecular Weight

Average Disk Density


M s





10 M



250 km/s



Theoretical Foundations of Hydrodynamics

Validity of the Boltzmann Transport Equation

The most satisfying foundation for the hydrodynamical
equations is the Boltmann transport equation. The hydro-

dynamical equations are the result of integrals of the


f D opf dc

Tc is a conserved quantity, Dop, the Boltzmann operator,

and f, the distribution function which is a solution of the

integro-differential equation

D f= 0 .

Two assumptions are essential to the derivation of
the Boltzmann transport equation. The first, the assumption

of molecular chaos, is mentioned only for completeness;

the second, the restriction to binary forces, is open to

question in the presence of Coulomb forces. The so-called

collision integral of the Boltzmann equation includes only

the effects of binary collisions. More fundamentally, the

Boltamn, Weltansicht is that the individual particles

comprising the gas spend most of their time in free flight,

this condition being interrupted only occasionally and

briefly by collisions.

The long-range Coulomb .forces present,. when charged

particles are constituents of the gas, .generally result in

an infinite contribution from the collision term. The

underlying view of the gas motion is suspect as several

distant collisions may simultaneously interact with one

another. There is, however, one set of circumstances under

which the Boltzmann approach is still valid (Zel'dovich

and Raizer, 1966). The conditions are:

(1) the average Coulomb potential energy at the

mean separation distance is much less than

the average thermal energy;

(2) the Debye length is much greater than the

mean separation distance.

Assume that the disk gas is made of equal numbers of

electrons and protons. The ratio of the Coulomb potential

energy at mean separation to average thermal energy for

such a gas under disk physical conditions gives

(ze -= 1.2 x 10-3
N kT
e e

where Ne, e, k, and Z are the electron number density,

the electron charge, the Boltzmann constant, and the average

ionic charge. The ratio of the Debye length to average

separation is
6.90(T /N )1/2
e = 8.0
N 1/3

The first condition is nicely satisfied; the .second, less

well so but within the limit of toleration. The conclusion

to be-made is that the Boltzmann approach and, hence, the

hydrodynamical description are justified, but are near the

limit of validity. The additional species present in a

more realistic representative disk gas will not alter this


Relaxation Times

A fundamental consideration for a gaseous mixture

is the extent to which that mixture can be treated as a

single fluid, particularly to what extent a single temperature

can be used to describe the flow. The relaxation time for

collisions between two species is a measure of the time

required after an initial disturbance in one of the species

to establish equilibrium among the translational degrees

of freedom in both. For a system of the type under study

here, there is some time which is characteristic of the

time scale on which hydrodynamic variations occur. If

the ratio of the relaxation time between any two major

species of the gas and the characteristic time of the

system is equal to or greater than one, the single-fluid

hydrodynamical description is invalid.

For a gas under disk physical conditions there

are two relevant relaxation time--that between ions and

electrons and that between electrons and neutral atoms.

For purposes of studying the nature of the disk gas, that

gas will be approximated..by a three-species gas containing

electron, ionized hydrogen, and neutral helium; further,

the characteristic time (Tc) for the system will be taken

to be 27Rs/Vs

The ratio of the electron-proton relaxation time

(Te-p) to the characteristic time (T ) is

T 3.5 x 10+8V T3/2
e-p se
Tc 2irRN enAd
s e d

where Ad is the reduced Debye length (Spitzer, 1962). The

ratio of the electron-helium relaxation time (Tehe) to

the characteristic time is approximately

T: V
e-he Vs
Tc 2~Rshe T1/2 Q
S 2 Nhe e he

(Zel'dovich and Raizer, 1966) where Qhe is the helium-

electron collision cross-section for typical disk conditions
and Nhe is the helium number density. If Te is 10 K
-16 2
Qhe is approximately 5.7 x 10 cm .
For the typical physical state of the disk, the

ratio of the electron-proton relaxation time (Te-p) to the

characteristic time is 4.6 x 10-4. The ratio which compares

electron-helium relaxation time to the characteristic time

is approximately 1.79 x 10-5. Finally one can compare

he 0.4 x 101


The fist two of these three numbers indicate that the

temperature should be the same for all three components;

electrons, protons, and neutral helium.

These two conclusions should hold in the more complex

disk gas. The principal disparity between the three-element

picture and the actual fluid is the number of neutral

species. The quantity Te-he is being used as a measure

of the relaxation time for electron neutral collisions. The

relaxation time for collisions of electrons with all neutral

atoms, in the presence of additional neutral species besides

helium, can only be smaller than that for helium alone.

Therefore, the quantity Te-he.. overestimates the actual

electron-neutral atom relaxation time. Substitution of

a smaller number for Te-he in this section will only

strengthen the conclusions.

The Hydrodynamical Prescription

The Hydrodynamical Equations

The three conserved quantities used to generate the
hydrodynamical equations are mass, momentum, and energy; these

give an equation of continuity, an equation of motion,

and an energy equation, in that order. However, quantities

such as the heat flux and the deformation tensor which

appear in the equations are defined by integrals that

contain the distribution function. Indeed, the quantities

that are identified as flow variables (the pressure, density

and fluid velocity) are also so defined.

If the transport equation is solved by a successive

approximation technique in which the distribution function

is written as a series of terms generated by increasingly

higher degrees of approximation, a set of hydrodynamical

equations may be formed at each stage of the approximation.

The "zeroth" approximation, in which only the initial term

for the distribution function appears,is a locally Maxwell-

Boltmann distribution; the resultant hydrodynamical

equations are those for an ideal inviscid fluid. The

equation of motion in this case, Euler's equation, has the


c + + 1 +
+ v Vv =- VPc Fc
at p
c c

where P is the force per unit mass. The continuity
equation is

SC+ V (pcc) = 0

The next approximation leads to the hydrodynamical

equations for a viscous gas. In this case the equation of

continuity is unchanged but a new equation of motion, called

the Navier-Stokes equation, is formed. The Navier-Stokes

equation will be examined below to ascertain the extent to

which Euler's equation is a good approximation for a fluid
under disk physical conditions.

Navier-Stokes Equations and Similarity Numbers

The Navier-Stokes equation in the presence of a

gravitational potential is

avC + + 1 2+
+ v V = p + v

+ /3 V(V v) GMsV (4- .) -
C r5

The quantities n and g are the first and second vicosity

coefficients, respectively. The quantities n and g are both

always positive. The second viscosity coefficient

represents effects that occur at high density or when
species with slowly excited degrees of freedom are present..

It is included only for completeness. Even if these effects

were present, as long as n > S, all the arguments presented
below are valid.

The second term on the left-hand side in the Navier-

Stokes equation is called the inertia term; the second

and third terms on the right are the viscous terms. The

last term contains gravitational effects.

Let R be a dimension characteristic of the

boundary surface, v a typical value of the fluid velocity,

o a representative density. These characteristic values
are chosen such that the associated physical quantities

(r v pc) vary from a large fraction to several times

the characteristic values, that is, each of the values

represents the order of magnitude of the quantity of which

it is a characteristic. Characteristic values for the

pressure and.time which result from these choices are

P = I v2
o 00


T = Ro/v

respectively. Note that the characteristic numbers are


The characteristic numbers can be used to set up

a system of dimensionless variables defined by the following


+ +
(1) rc = Rr
+ +
(2) vc = VoV

(3) Pc = IIo '

(4) P = PP
C o

(5) to = Tot

Note that the subscript zero quantities are dimensioned

numbers. These relations can be used to rewrite the

Navier-Stokes equations in dimensionless coordinates. The

result is shown below:

v2o v + V GMs V
+ v Vv= -v2
R t R 0 p R2 r

+ v +7+ /(3) + /3) V(V v)
II 2 V 0
Sp oRo2 o p

Division by V2/Ro yields

at V+ v =-- V -
at p R 2 r
o o

+ v + (+ t/3) V(v v)
RVo o RVIo p
000 .000

Since the characteristic parameters represent

the order of magnitude of their respective variables, the

order of magnitude of the ratio of any two terms is the

ratio of their coefficients. Thisstatement is, however,

open to question near the flow boundary. The flow variables

may change rapidly over short distances near the boundary

surface thereby causing terms involving their gradients

to be quite large. These coefficient: ratios taken together

qualitatively describe the flow. Historically they have

been called similarity numbers and given the following


1. Reynolds Number. This dimensionless number,

labeled RN here, is the inverse of the ratio

of the viscous terms to the inertial terms,

1 = n .

2. Froude Number. This value, referred to

symbolically as FN, is

FN = .
R (GM/R2)

Its inverse is the ratio of the gravitational

term to the inertial term.

3. Strouhal Number. Called SN, it is given by

the expression

The inverse of this quantity is the ratio of the

time derivative of the velocity to the :inertial


Importance of Viscous Term in the
Equation of Motion

The relative importance of the term in the Navier-

Stokes equation can be determined by calculating the

similarity numbers. To calculate these numbers, the

viscosity must be estimated and the characteristic values

(R V I ) chosen.

Suppose the gas to be a ternary mixture composed

of hydrogen ions, electrons, and neutral helium. Because

the ratio of the electron mass to that of either of the

other two species is so small, the viscosity of the mixture

is essentially determined by the viscosities of the

hydrogen ions and neutral helium (Chapman and Cowling, 1970).

For hydrogen ions under physical disk condition the

expression given by Spitzer (1962) may be used. It is

2.2 x 10-15 T5/2

which for disk temperature (Te) gives

n = 2.2 x 106.

Chapman and Cowling (1970) give an expression that can be
used to calculate helium viscosity.

5(kmheTe W1/2
82 W

where a is the molecular diameter and W is a tabulated

function (' is the mathematical number). Under disk

conditions the two quantities a and W are

a = 2.7 x 10-8 cm2

W = 6.00

This gives

n = 2.4 x 10-4

The following characteristic values are used:

(1) R0 = R

(2) Vo = Vs/2 ,

(3) Ho = pd

For viscous effects due to hydrogen the similarity numbers


(1) RN = 1.2 x 10+12

(2) FH = 8 x 10-3

(3) SN = 1 .

The helium viscosity results in a Reynolds number which

is two orders of magnitude smaller than that for hydrogen,

but this is of little import. Nor would a more accurate

calculation change the basic result that the Reynolds

number for the disk fluid is very large. The ratio of the

viscous terms to the inertial terms is so small that the

latter completely dominate the former. Removal of the

viscous terms reduces the Navier-Stokes equation to

v -VP GM v(s 1
P oV

Euler's equation.



Prelimiinary Remarks

In the work which follows magnetic effects will be

ignored and any radiation effects will be assumed either

entirely negligible or such that they may be included through

a small adjustment in the gravitational term. On the basis

of the discussion in the last chapter, the viscous terms in

the Navier-Stokes equation can be discarded. The momentum

equation is

av + + VP s
+ v Vv -V

The coordinates, time, and the flow variables are dimensionless

and related to the normal dimensioned quantities through the

relation of the preceding chapter. In the rest of this work

the symbol

GMs GMs(2i)2
o 2 2
00 ss

will be used. The equation of continuity,

-P + V (p ) = o0

provides an additional relation between pand v.

The full hydrodynamical prescription requires the

inclusion of a third equation, that for energy transport.

No effort will be made in this chapter to introduce such

an equation. Failure to study the energy equation will

necessitate, at some point, the specification of an additional

condition linking two of the flow variables. Introduction of

this new relation reduced the number of unknowns from five to

four. Euler's equation, which contains no viscous terms, is

not valid in the boundary layer between the bulk of the

circumstellar matter and the stellar photosphere. No analysis

of the flow in this region has been made, nor will any be

attempted here.

Linearization of the Equations for Temporal
and Angular Dependence

Introduction of Time and Angle Dependent Terms

As the steady-state axially-symmetric solutions to

the hydrodynamical equations appear to explain much of the

Be star phenomenon, each flow variable will be each divided

into two parts: a steady-state axially-symmetric part and

a term containing both angular and temporal dependence.

The notation is shown below:

(1) v = v (W,z) + V1(,(,zt)

(2) p = po(wz) + pl(m,#,z,t)

3) P = po(w,z). + pl(W, ,z,t)

The subscript-zero identifies the steady-state axially-

symmetric terms; the subscript-one, the: terms containing

the temporal and angular dependence. In order to avoid

confusion with the dimensioning parameter Po used earlier,

a lower case p is used with subscript in the case of


These pairs of functions are substituted into the

equation of motion to give

aav 44 4
o + + + + Vo + Vo +V
at + + v + v Vv0 + v .v
St at o o 1 o o 1

+- + 1 + 1 11o
v +Vv = + 2 + Vpo 1 2 P V- q(o(32
p p0 o p
0 o Po o


ap app 4
-t + + V + V--p + Vo Vpo + 1 O + Po p V PV+
at at vo 0 1 0 0 o o

+ pV v1 + pl1V vo 1+ p v1 = 0

The denominator in the term

1 V(po + p)

has been approximated by first expanding in the series

-1 2. 3*
( + x) = 1 x x x + .... ,

where x is pl/p o and then excluding terms in the third

and in higher powers of x. Subsequent conditions placed on

the flow will limit the problem to cases where

P1+ Po >0 .

With the exception of this approximation for the denominator,

these equations are exact.

Since the subscript-zero terms are themselves steady-

state axially-symmetric solutions of the problem, they satisfy

the equations

o + V P 1
at o o p o o r


ap 0
3^+ V (poVo) =0

The equations for the subscript-one terms become

av- + v V + v VO + VV = ( ) VP

o o

0 p


P"VviP + v 4Vp + Vp
-_ o v + lV P V 1 + o vpl + 1 V

+ v- Vp1 = 0

The underlying physical picture presented thus far

has been that the steady-state axially-symmetric solutions,

the subscript-zero terms, represent the gross nature of the

flow. The subscript-one terms are taken to represent a

refinement to this gross nature. These refinement terms

would then be expected to be small--though not negligible.

The assumption will be made here that subscript-one terms are

sufficiently smaller than the subscript-zero quantities that

elements of the equations which include the products of

subscript-one quantities can be omitted. The resulting

linearized equations are

i Pi 1
+v v + v Vv Vp Vp
at + Vo l + V 2 1p


api + + +
-t+ Po V 1 P1v Vo +o 0 P1+ V1, Vo = 0


These must be rewritten as scalar equations

+ a V
3 + Vwl

.avvo. a- lv. I V 3 V aV
o- + o o
-T- o T- + v-W a o + W go

+ V zi w "
21 rz

-y- + V1

Vzo ~a
20 BS


2v v1
CO p

WO 9aw

+ v ao V l a
zlz 00 zo z V01



p ?

w a

V W1 00

1 ap1 p1 ap
+ + 2 3p
o p0

+ .-o :

Vo 1

o ,

av +


v yv V
+ -" +
to 84

+ Vzl -- + Vzo
z1 5z[ zo

Vzl ap1 P,1 apo
5z p O z 2 Bz


P1V PvJl P91 p aVo 0Po
+ P l e+ o + + 3W a0
S W WO 5W 1 -r-- + W1

apl 3 o 1_ P lazo
(~o & + +1 o ) + a
(v~o g-- VPI -** + V z

+ V
too 9(0

v 1l 1
+ po -aw +


vPl Po 1. =
+ + v + p 0 .
+ zo z zl a o z

Selection of a Subscript-'ze'ro Solution

The subscript-zero solution which will be used is

that for a static isothermal envelope (Limber, 1964). This

solution is chosen for the following reasons:

(1) over most of the disk, the isothermal envelopes

are in agreement with calculations based on

more advanced treatments (See Chapter II);

(2) the solution is in closed form;

(3) the simplicity of this solution results, in

turn, in a comparatively simple set of expressions

for the subscript-one terms.

This subscript-zero solution (in dimensionless coordinates)


Vo = 0 ,

vso = 0 ,
1 1 1

P =exp { 2 +) 1
Q. W (W + z)1)

1 1 1
Po = exp {- G (- 2 2 1 ) 2

The quantity Q is defined

kR T
s e

while Y is given by.

..(.2..) .2 t
pm Vs

The ratio of these two, Q/', is the Froude number (FN).
Substitution of these expressions into the linearized
equations for the subscript-one quantities gives

avWI 27 av 1 4 = ex i 1 __ ] pi
- -+ 3- 72--w 7 -- V-7 v l C r I 30)

1 p. + 22 1 e [ 1
+w r r- 77 Q W P ) 8
. 1 + {exp[- )]-} { exp[ -( ,)

at -3 W1 + w7 3/2 W Io pW r Do
2I 2( 1 1 a exp[ 1 1 1 1

+ pl exp[ r {exp[ )] I

z 2- .Zl {exp 1 1 1 P

and w r a37 0 W r J


Iexp[- I I +v- exp[ I. I
atL Ww r' J W wlaw LQMIw JW

exp 1 1 v 2r ap 1 exp[ 1 17 I

+{ep ]} = 0 .

While the formal division of the flow into two

elements is quite general, the division has physical significance

only if what are called here the subscript-zero terms adequately

represent the steady-state axially-symmetric part of the flow.

The static isothermal solutions do not adequately represent

such flow at large distance from the equatorial plane. The

solutions to the above equations are useful only for regions

near the plane. It will be the practice in the rest of the

chapter to limit the discussion to the region Izi /i < 0.1

(quite a conservative value). Such a limitation may seem a

drastic approximation; in fact, it is not. This "thin"

region is a solar diameter in thickness near the stellar

surface and larger farther out. Additionally, examination

of the z-component of Euler's equation reveals that only the

pressure gradient balances the z-component of the gravitational

force, but the ratio of the former to the latter gives

2a -- = 1, 8.3 x 10-3
3 1 qo
Ro 0 r- )

The dominance of the z-component of the gravitational term

suggests that the disk material is highly concentrated
toward the plane.

Considerations Concerning Approximations Near the Plane

For z < w the term 1/r3 becomes

1 1 a 2 -3/2 1 3 z 3
1 -( ) -- 1 -r + .

similarly, the 1/r term is

S)21 ]
r w 2

All derivatives of known functions should be
performed before approximations such as those described
above are included. This has been the practice here. On
substitution of the appropriate expressions for 1/r3 and
1/r, the set of four equations becomes

av + 2 4av_ 12
l 2w Vl 4 z2 ap
3 pexp 2Q2 2

+3P1 z2 2
+- c~4) [ exp, )]
2Q W 2Qm-

+. 2 [exp (Z)l l -- [exp (* )
at w31 w2 32 22
PVz 2 2Vzl .1 1 .Iz 2


S+ [exp (-a)]+ 4 exp (-
I 2Qw 2Qw 2Q33

+ exp a- 3 1 exp

Z1 exp (- ) exp -v = 0.
Q 22Qw 2Q
VV=v z2 z 2-

The form: of the last and most complicated of the set
of four equations suggests a change of variable which will
simplify the equation. Define the variable v1 by the

Vil= Vil ep 1 23 '

the i standing for w,(, or z. Now

+e 4-z2 vl 3z2 -z2
1- vI = exp 2Q Vl ) -- 2m4 v1 exp( 2 Q

+ [exp ( ) + [ exp -z2~ )a
2Qw 2Qw

zv 2
1 exp( 2 )
Qw 2Q

Comparison of this result with the terms in the last of the
set of four equations above shows that the last equation can
be written

api 2 p e0
7+ 227 r + _---+ v =-1e

Study of the first three members of the set, which were
derived from Euler's equation, indicates:
(1) that no component of vI appears in a partial
derivative with respect to either w or z,
(2) that the exponentials, when they appear in any
term, come before all operators.
2 3
On multiplication by exp (-z2/2Qw ) and with the introduction
of the variable v the first three equations of the set

yve 8ve ap 3Yz2
-4"- + f 2f 2-e4 P
at a0 1 aw 2Qw4

SI f1ve f e
at W= ap


+e. ap e ZP
at .l az as z

Note that the function f used above is

f 27r

The exponential in the definition of ve makes the function

strongly peaked toward the equatorial plane.

Additional Relations

There are four equations, but five unknowns. As

mentioned in the beginning of this chapter the missing

equation is that for energy transport. This equation

requires more physical information about the nature of the

flow than is needed for either the equation of motion or

the continuity equation and it is mathematically more

complex than either. The traditional approach has been to

circumvent the need for such an equation by specifying

in some manner the relationship between any two of the five

unknowns; this will also be done here. The density and

pressure will be assumed to be related by a condition of

the form

Pl = rp1

where r is a parameter. As an illustration, suppose the

relationship between the subscript-one density and the

subscript-one pressure to be isentropic with the isentropic

constant determined from the subscript-zero quantities,

2. _
r = 10-2
(vs/2 )

where c2 is the isothermal speed of sound determined from

the subscript-zero terms.

The boundary conditions for viscous hydrodynamics

require that both the tangential and normal components of

the velocity on either side of the boundary between flow

and bounding surface be equal. In the steady-state axially-

symmetric solutions, the interior boundary surface is

rotating and has no radial velocity. The steady-state

solutions do not represent a true use of boundary conditions,

for the velocity dependence is in effect assumed. Further,

the flow is taken to be Eulerian up to the boundary. An

adequate treatment of the boundary-value problem requires:

(1) treating the flow in the boundary layer,

(2) describing the behavior of both the tangential

and normal components of the velocity

over the boundary surface, particularly

dependence on time and angle.

There is no adequate treatment of the boundary layer flow,

and the observations provide little information about the

temporal or angular dependence of the stellar surface

velocity. For these two reasons boundary conditions are

of little use in the study of the nature of the subscript-one


To this point the four equations which together

comprise Euler's equation and continuity are

av av apj 3.,z2
+ f Vr + f P
at 2f vl r 2 PQw

av fv1 favl r
01z + f)P1 P1
at + D2O a 1

Bvzl av ap -. YzP1
-+ f Zi: = -r


P l p1 Val aVl 1 av
_+f + 1 + + i + + = 0

Note that the e superscripts have been dropped. The symbol

vl now contains the exponential term, and will for the rest

of this chapter.

Integration Over z
Further Restrictions on the Flow
The subscript-one flow will be assumed to possess
symmetry about the equatorial plane; that is, the flow below

the plane is the mirror image of the flow above the plane.

This condition requires that the flow .variables, with the

exception of Vzl, must be odd functions of z.

The set of four equations which describe the flow

will be integrated with respect to z over the interval

[-z + zo] where zo is a function of w alone. New variables

must be introduced; they are:
(1) ul = vidzl



Under the condition requiring symmetry about the equatorial

plane, the z-component of the integrated velocity (u l) is

identically zero. The integration with respect to z and the

introduction of new variables require that operator reversals

of the type

+z (W) +z (W)
f a g(w,z) dz = g(,) dz,
-zo(i) a d -Z (w)

where g(w,z) is a well-defined function, be performed. The

justification for such changes in the order of operators

appears in Appendix A.

These two steps result in the three equations

SU 1 3 o 2
u f u1 2ful= r D, + z 2 p dz
t Wl Wl + qw -z

.t +'" l + f l U '


18 + l a 1 .801
S u + f U+ + 2vzl(z = ) = 0 .

The last term of the first equation above will be


f + z2pl dz=Af + p1 dz
-Z -z
0 0

where the quantity A may be considered defined by this rela-

tion. A is, strictly speaking, a function, but it will be

treated as a constant parameter. The validity of this

approximation can be tested; the solutions should change

only slowly with changes in A.

The last term of the fourth equation, 2vzl(o0),

represents mass flow perpendicular to the equatorial plane;

a quantity generally believed small. It is a reasonable

approximation to set this term to zero. Such a step, while

not crucial, does simplify the calculation.

The third term in the continuity equation, which

contains a partial derivative with respect to w, should

be dominated by the neighboring terms which contain 0 and t

partial derivatives. This term will therefore' be dropped


from further consideration,as is a similar.term in the first

equation of the set.

The result of all these steps is shown below:

'8(1 )1
at W1

+ f u
8 0

S 2f Ul 4- 1 .
1 2f 1 2Q4 1

2au f *ai r= ^
(2) a B + -ul + f 'u = --r a
U1 2 Ul + f a

(3) -3

+ f a
a #

+ +
+ +

1 -'1- =


The following symbols are introduced to simplify the


(1) c1 -Q

(2) c2 = -r

(3) L = + f p s

The equations of the previous section become

L Ul 2fu.=
Wl p


c2: 80
u=W = W -0.

L 1 +
1 2




`i 1 .1.a1
La+ + -- + = ..

This set of simultaneous differential equations must

now be reduced to a set in which there is a separate equation

for each independent variable. This involves much differentia-

tion and manipulation. These reductions are shown in Appendix
B. The results of these manipulation are that all three

equations share the common form

L c 2 c 2c f cf
3 2 2L + + ( + wf2 ) L +( -2 1 } x
84 2w 2w

(variable) = 0

where "variable" may be either u l/1 uw l or *

Each of these three equations has the form

F 1 ) y =

y = exp [a(w)) + O(w)t]

where a(w) and 5(w) are as yet unspecified functions. Since

exp Ia(Mw) + 0B()t] = a(w) exp [a((w) + B(w)t]


t exp [a(m)w + (w)t] = 0 () exp [ta(W)) + 0(o))t],

one has that

F (I, ,

y = F(a(W), p(W),w)y

A solution (y) exists only for a(w) and (w) such that

Fla(), B(W),) = 0.

Unfortunately there are infinitely many such (a(w), 3(w))pairs.
The temporal and angular dependence of the flow

variables takes the form exp[ikt + inf] where the substitutions

a = iE

S= in

have been made. The quantity E is, in effect, a frequency.

For comparison with observations the associated quantity
2wt/1 which is the period associated with such a frequency
is more useful. Since and (# + 2w) are physically the

same point, n must be an integer. The equation which k and
n must satisfy is derived in Appendix B. This equation is

3 n2
k-3 + 3wnf .2 2+ -f2 + f2 cl
fi2 i4)

c2n3f 2c2fn 3 3 1cnf
+( --- wf3n + n3f3 ) = 0
w W 2w

For a given n and w there will be in general, three values

for 1.

Relations Among the Flow Variables

The temporal and angular dependence of each flow

variable takes the form

exp [ino + iI(w,n)t]

The quantity n is an integer, and the E(w,n) are roots of a

polynomial equation whose coefficients are functions of w

and n; The general solution for the problem may be written

Ul = Uk(W)

exp[ino + ik(w,n)t],

ul = Cn(m) exp[in + ik(w,n)t],


1a = SnE() exp[ino + ik(w,n)t]

The relationships among the three quantities, Unp, Cn~ and

Sn', determine the relationship among the flow variables.

Substitution of the general solution into the

original set of three equations yields

[i( + nf) U 2f C- nS { exp[inn + ik(w,n)t]} = 0,
"._ [i (k + nf) Un. -nk 2f Cn4- S4
nk 4

[%rU inC
SU- + + i(i + nf) SnV { exp[ino + ik(w,n)t]} = 0,
n,k 2

rfUnk inc n
-- -f + i(I + nf) nk S ij{exp[ino + iF(w,n)t]} = 0

These three equations all require summations over n and j.
What is desired is a set of relations among the three
quantities (Un S r, C n) for a given n and k.
Each equation is multiplied by the quantity

exp[in'O + ik'((,n')t] .

Next, double integration of the form

+1? +7
Lim -- f dt f dO exp[-inO ikt]
T o 4wT
-T -7

is performed on each equation (and the order of summation and
integration reversed). One finds

i(' + n'f) Un', 2f Cn,, -1 Sn, = 0 ,

Un' + i(k + n'f) Cn'k' S -n = '

Un,' n Cn' + i(E' + n'f) Sn'J = 0

In the following pages the primes will be surpressed.

These expressions form a system of homogeneous linear

equations. There is no a priori reason to assume that any

non-trivial solutions exist. If solutions do exist, they will

constitute at best a single infinity of solutions. Such a

case can only provide unique values for the ratios of two

of the unknowns with respect to the third.

If non-trivial solutions exist, the determinant

formed by the coefficients of the variables in the system


i(k + nf) -'2f -cl/4

f/2 i(k + nf) -in c2/m

1/W in/a i (k + nf)

should be zero. The details of the calculation of the

determinant are left to Appendix C but the results are

pleasantly familiar. The condition that the determinant

be zero is a polynomial in k which must be zero. This

polynomial is the same as that used to calculate k. Therefore


the determinant is identically zero.

There will be a single infinity of solutions if

the rank of the matrix from which the determinant above was

formed is two; then the equations can be solved for the

ratios of two of the unknowns with respect to the third.

The rank of this matrix is, indeed, two (See Appendix C).

The ratios which will be determined are

X -


y -
y =.nk

They satisfy the system of equations

-(k + nf) x 2fy = -A

fx inc2
f- + i(k + nf)y = --

x iny = i(E + nf)

These two ratios are

ri.(kO + nf) 2fnc2
x= -- 4 + + ]


r[. nc .(E +. nf) .



DET = [-( + nf) 2 + f2]

A Reprise

The analysis of this chapter is valid only to the

extent that the underlying physical picture, discussed in

the second section of this chapter, is a reasonable representa-

tion of Be stars. This requires the dichotomous view of

the flow presented there to be a real one. The steady-state

axially-symmetric effects must dominate the flow. The

elements of the flow which depend on time and angle must

be secondary.

Further,solutions were obtained only when:

(1) known functions appearing in the differential

equations were approximated by formulae valid

only for z < w;

(2) the flow variables were assumed to possess

symmetry about the equatorial plane;

(3) two parameters, r and A, were introduced;

(4) the equations were integrated with respect

to z (and new dependent variables introduced);

(5.). two terms in the .continuity, equation and one

in the equation of motion which appear to be

of secondary importance were dropped.

The solutions to this set of equations are proportional

to exp[inj + ik(w,n)t], where n is any integer and the values

of k for any n and w are the roots of a polynomial equation.

The general solution for the flow can be determined to

within the infinite set of multiplicative factors Sn- which

appear in the integrated density. Removal of the S,

indeterminance awaits a boundary value analysis.



Introduction to the Computations

Both the angular and the temporal dependence of the

solutions in the: previous chapter are controlled by the

integers n and the quantities 1E, respectively. The calculation

for I, which depends on n, w, A, and r, emerges as the

principal difficulty. Accordingly, the three principal

concerns of this chapter will be:

1) the nature of the roots of the equation for k,

2) the means used to calculate them,

3) the results of the computations for various

choices on which the ? equation depends.

As the solutions cannot be fully evaluated in the

absence of a proper boundary treatment, these results will

constitute much of what can be learned from the analysis.

In most of the discussion the period


(generally in hours) will be used.

The first major section of this chapter is devoted

to the first two topics above, the nature of the roots and

the procedure used to compute them. The results of these

computations, the third major concern of this chapter,

naturally group themselves into three divisions:

(1) calculation of the pertinent quantities

in the absence of any angular dependence

(n = 0),

(2) temporal behavior in the presence of the

simplest angular dependence (n = 1),

(3) the trends of the results for more complex

angular behavior (n < 1).

The three central sections of the chapter reflect this

three-fold division. In each the dependence of the

results on w, A, and r is discussed. The penultimate

section presents sample calculations of the flow

variables which illustrate the qualitative aspects of the

solutions. The last section of the chapter compares

the results with relevant observations.

The Equation for k

The Consequences of Imaginary Roots

The equation for k, derived in Chapter IV is

R3 +c3nn 22 2 2 c1l
W3. + nOfk2 +f2 + 3wn2f2 4 )
.W W

c2n3f 2c2nf 3 c1nf
+ ( f3n + n~ 4-) = 0
w ( 2w

Being an equation of the third degree, this equation has

either three real roots or one real and two complex roots.

If there are complex roots, one is the complex conjugate

of the other. The appearance of complex roots had immediate

and serious consequences for the methodology of Chapter IV.

One of these roots.must of necessity have a negative imaginary

part which will appear in the solutions as an increasing

exponential term, exp (+IklIt), where kIki is the absolute

value of the imaginary part. Flow variables which increase

exponentially in time will, at some time, violate the condition

under which linearization of the flow variables is valid.

It is possible to choose the three multiplicative factors

(Snj, Cnj, Un ) in such manner that these increasing terms

do not contribute to the final solution, but such choices

are suspect. Only for conditions in which all the roots

are real can the linearization procedure, and hence all

thereafter, be considered valid.

Initial Values for the Parameters
The quantities r and A, which appear in the

coefficients of the equation for IE (in c, and c2' respectively),

enter the problem through the relations

P = rP1

+z +z
z2 P dz = A o 01 dz
-Z -z
0 0

respectively. The calculations, whose results are described

in later sections,will be performed over a whole range of

values for each of these quantities; however, an initial

or representative value for each quantity must be determined.

The representative value for r--called T--is the

value for r in the case that pl and p1 are related

isentropically. The isentropic coefficient is determined

from the subscript-zero quantities; it is

(2w) kT
r = .-a
2 pm

For an electron temperature (Te) of 10+4.K,

= = 7.67 x 10-2

A representative value for A--called K--can be

calculated be setting- p equal to a constant; for zo equal

to 0.1000,

A= = 3.334 x 10 .

Machine Computation of the Roots

Although the mathematics required to determine the
values of k at any point in the disk is not difficult, it

is time consuming. This calculation must be performed at

many points and for a number of different. values for the

quantities n, A, and r; machine computation is necessary

to perform as many calculations as are needed. A program

was devised which will perform the following:

1. Input. The program reads in the values which
are to be assumed by the physical and geometrical quantities

appearing in the coefficients of the k equation. The program

must also be given initial values for A and r.

2. Roots of the W equation. On the basis of internal
control statements, the program establishes a sequence of

A, r, n, and w combinations. For each of these combinations,
the program calculates first the coefficient and then the

roots of the equation for k. Each root is searched for an

imaginary part; if the root is real, the quantity 27/k (in

hours) is calculated.

3. Solutions for the flow variables. When given SnR
as a function of w, the program uses the formalism developed

in Chapter IV to calculate both Cnj and Unk. The quantities

Re{Un-(w) exp[in( + ik(w,n) t]}

Re{Cnr (w) exp[in + ik(w,n)t]} ,


Re{S -(w) exp[in< + ik(w,n)t]}

are calculated by the program throughout the disk for

various choices of angles and times.

4. Readout. The program prints the results of the

calculations in tabular form.

The actual program steps are contained in Appendix D.

The Case n = 0

Distinctive Aspects of the n = 0 Case

The k values here are for solutions with no angular

dependence. In the absence of dependence, the I equation

takes the form

w3 (2f2 + ) 0

the three roots are
2 c1 1/2
k = 0, + (f + -V)

The first value corresponds to a trivial solution. The

last expression is real for all reasonable values of cI.

Variation with w, r and A

The expression for the non-trivial roots depends on

w and c1 but not on c2. Now

Figure 3

The dependence of IPERI on
w is shown for the case n = 0
and c1 = 36.122.




-- 48.5 w3/2


2.0 4.0 6.0 8.0 10.0
Radius(Stellar Radii)

---- = 0.152



Setting c1 = c1 one finds that

2w+ C1 1
lil= 72 ( 1 + ( 227)2


The absolute value of the corresponding period is



IPERI= 48.5 w3/2

For large c1 the behavior of IPERI near the stellar surface

(o = 1) departs slightly from 3/2 behavior, as would be

expected from the full expression for IEI but the outer

regions of the disk still show the 3/2 behavior. These

features are shown in Figure 3.

The Case n = 1

Preliminary Remarks

The n = 0 case cannot be used for a general discussion.

Exploratory calculations: indicated that the n = 1 case

illustrates most of the characteristic features of the

general non-zero n case, yet possesses certain unique

features themselves worthy of investigation. Accordingly,

the most extensive series of computations were performed

for n = 1.

On the basis of both these preliminary calculations
and physical considerations the range of A was chosen

to be [0.000, 0.200] and the range of F was taken to be

[0.001, 1.000]. Calculations were made for each of the 150

r-A pairs; these pairs uniformly spanned the combined ranges.

In a later series of calculations r was allowed to take values

as high 10.00. In subsequent discussion, the quantity c1,

through which Aenters the equation for k, will be used

rather than A; however, r continues in use since c2 = -F

Initial studies also indicated that an evenly spaced 100-point

grid spanning the range [1.0, 11.0] in the radial coordinate

was sufficient to study the variation of the calculated

quantities throughout the disk.

The dependence of the solutions to the k equation
on c1, F, and w will be the subject of the next three

subsections. One important result should be mentioned at

the outset--the roots of the k equation under almost all

circumstances studied are real.

The Effect of c

Rather than E, the absolute value of the period

Figure 4

The quantities IPER11 and jPER21
are plotted against radial
coordinate for two different values of
c,. Curves a and b are IPER21 and
PER I respectively, when c1 = 0.361;
c and d, when c = 2.709. In all
cases r = 0.075.



Figure 5

The dependence of IPER I and
IPERil on w is shown for two
different values of c,. Curves
a and b are IPER21and IPER I,
respectively, when c = 18. 061;
c and d, when cl = 36.122. In
all cases r = 0.075.