UFDC Home  Search all Groups  UF Institutional Repository  UF Institutional Repository  UF Theses & Dissertations   Help 
Material Information
Subjects
Notes
Record Information

Full Text 
FLUID DYNAMICS OF CIRCUMSTELLAR MATERIAL ASSOCIATED WITH Bb STARS By Thomas Harlow Morgan A DISSERTATION PRESENTED TO THE GRADUATE COUNCIL OF THE UNIVERSITY OF FLORIDA IN PARTIAL FULFILLMENT OF THE REQUIREMENTS FOR THE DEGREE OF DOCTOR OF PHILOSOPHY UNIVERSITY OF FLORIDA 1972 To KYC and HENboth of whom pushed. ACKNOWLEDGEMENTS 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. KY Chen suggested the topic of this dissertation. Mr. W. W. Richardson drew the figures that appear here. The members of my committee, particularly Drs. KY 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. iii TABLE OF CONTENTS ACKNOWLEDGEMENTS. . . LIST OF TABLES . LIST OF FIGURES . . KEY TO SYMBOLS . . ABSTRACT . . CHAPTER I * S S S * 5* 9 5 9 S * S S . S S S . INTRODUCTION . . Discovery and Taxonomy of EmissionLine B Stars . . A Brief History. . Nomenclature . . Characteristics of Be Stars . Spectral and Luminosity Class. . Galactic Distribution, Membership in Multiple Star Systems and Incidence Among EarlyType 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. . II HYDRODYNAMICAL APPROACHES AND STEADYSTATE SOLUTIONS . . Background . . "iv Page iii ix x xii xv 1 1 1 2 3 3 3 4 5 6 6 10 10 11 12 13 15 15 * * CHAPTER The Hydrodynamical Equations and Notations. 15 Static Solutions. . 18 The More General Solutions. . 20 Form of the Circular Velocity Law 23 Final Comments ... .. ..... 27 III APPLICATION OF HYDRODYNAMICS TO THE ENVELOPES 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 NavierStokes Equations and Similarity Numbers .. .. 37 Importance of Viscous Term in the Equation Of Motion . . 41 IV MATHEMATICAL ANALYSIS OF THE HYDRODYNAMICAL EQUATIONS. . . 44 Preliminary Remarks . 44 Linearization of the Equations for Temporal and Angular Dependence. . 45 Introduction of Time and Angle Dependent Terms . 45 Selection of a SubscriptZero 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 v Page CHAPTER A Reprise. . . V COMPUTATIONS AND RESULTS . . 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 . 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 . VI EQUATION FOR ENERGY TRANSPORT . 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. ... Page 68 70 70 71 71 72 74 115 116 118 121 121 121 123 125 * , CHAPTER Page A Linearized Equation. 129 Separation of the Time and Angle Dependence .. 129 An Equation for the SubscriptOne 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 NonStatic Solutions and the Adequacy of the Equations. 138 Final Comments .. 140 VII FINAL COMMENTS . 142 Stability of the SteadyState 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 APPENDIX C . 155 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 vii Page APPENDIX F ... .. ... ...... 168 An Equation for reff . 168 BIBLIOGRAPHY ..... ........ ... 171 BIOGRAPHICAL SKETCH. ... .. 174 viii LIST OF TABLES 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 LIST OF FIGURES Figure Page 1 A Be star viewed equatoron and poleon 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 Figure Page 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 KEY TO SYMBOLS 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 steadystate 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 xii al integrated density n integer controlling the angular dependence S quantity controlling the temporal dependence 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 crosssectional flow), 22 pl, 137 Q, 50 q 44 Rs, 29 RN, 40 r, 17 Snr, 64 SN, 40 xiii 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 xiv 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 FLUID DYNAMICS OF CIRCUMSTELLAR MATERIAL ASSOCIATED WITH Be STARS By Thomas Harlow Morgan December, 1972 Chairman: Alex G. Smith CoChairman: KwanYu 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 steadystate axiallysymmetric 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 XV of linear equations for the terms with temporal and angular dependence. The known steadystate axiallysymmetric 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 radiationfluid interaction. xvi 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 steadystate flows is examined. The steadystate flows are not absolutely unstable, but the most general solutions appear to be periodic in character. xvii CHAPTER I INTRODUCTION .Discovery and Taxonomy of EmissionLine B Stars A Brief History In August 1866 Father Secchi observed HB.in emission in y Cas. Over the next halfcentury 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 earlytype emissionline 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 Btype emissionline 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 diskshaped envelope surrounding a rapidly rotating central star; this interpretation has not changed substantially since. 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 emissionabsorption features have received attention. The most recent listing of Be stars (Jaschek et al., 1971) contains 1900 objects. Nomenclature The basis of modern systems of spectral classification is the system adopted for the Henry Draper catalogue (Cannon and Pickering, 191824). 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 sequenceO, B, A, F, G, K, Mrepresented an ordering in temperature with decreasing temperatures from 0 to M. Among the HD objects typed B, those possessing hydrogenemission lines were marked with an additional lower case p for peculiarfollowing the B. Modern classification uses the notation Be for such stars, that is, stars with Btype 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 emissionline stars whose underlying spectral types are early of mid B (BOB7) 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 mainsequence (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 groupsvery 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 EarlyType 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 0B associations. Be stars are found in 0B 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 breakup velocities. Balmer emissions occur in approximately 15 percent of all Btype 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 activity. E/C and V/R variations, though strictly periodic only when they reflect orbital motions of stars in binary systems, .show much quasiperiodic 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 Origin In 1931 Struve suggested that the emissionabsorption features that are the definitive characteristic of Be stars are formed in gaseous circumequatorial rings surrounding (and ejected from) Btype central stars that are rotating at their breakup 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 r' m aa 4 tr 1r tra po+ 8 I ^ ^ e IC \ 0,'4' 0 st, SI , .C > S t a, crr r i trn  1 O1 /"^ C 0 "o 1 5 0 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 "poleon." 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 require (1) a model for the radiation field of the central star, (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 diskthat is, changes from B to Be to shell star and vice versanor does it say anything about the mechanisms that produce other spectral variations. Line Profiles The hydrogenemission 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 6level 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 steadystate 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 present. 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 are: (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 range, (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 nonuniform transfer of momentum from rapidly rotating inner regions of the central star out to the photospheric regions resulting in turn in nonuniform mass loss. Limber (1969) analyzed the shell phase of Pleione in terms of a timedependent 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 quasiperiodic 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 breakup velocity. Observations suggest that Be stars' rotational velocities lie below the breakup 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 midUV resonances in silicon, carbon, nitrogen, and sulfur, these being in high ionization states. 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) (wc/oc)1/2 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 steadystate, axiallysymmetric 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 secondorder term onto the timeindependent 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. CHAPTER II IYDRODYNAMICAL APPROACHES AND STEADYSTATE SOLUTIONS Background The steadystate axiallysymmetric 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 steadystate axiallysymmetric 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 nonzero 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 is (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 noncgs 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. coordinate 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 velocity Static Solutions Limber (1964) attempted the first hydrodynamical the material outside rapidly rotating stars. He central star to be rotating at the breakup 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 potential; 3) disk selfgravitational effects are negligible; 4) steadystate conditions prevail; 5) only axiallysymmetric 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 zero; 8) the circular velocity has the parameterized form v v O( ) a 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 usedEuler's equation. The third approximation, neglect of selfgravitational 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 solutionssolutions in which the radial velocity is zeroare,. strictly speaking, nonphysical for the same reason that static solutions are nonphysical for the solar wind (Parker, 1958). Limber extended his analysis to nonstatic 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 aJ 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 crosssectional 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 r F (r) = (c b s in the elementary solar wind theory (Parker, 1963) where b is a parameter (usually given a value near 2). Limber chose, instead, a oneparameter family of crosssectional flow terms 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 velocity 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 nearequatorial 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 crosssectional 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 crosssectional 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 breakup 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 Quasistatic motion of a small fluid element away from the star with centrifugal forces present requires that 2 v GMs W1 2 c Wc at any point in the disk; hence Vc=" "T c a !1/2c c There are four directsupport 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 outwarddirected 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 photoionization of neutral hydrogen. Maximum momentum transfer in the radial direction from photoionization occurs if the recombination following the photoionization 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 backemitted. On the basis of both this assumption and other considerations, these authors calculated an upper limit for support due to photoionization. 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 Steadystate solutions have been applied to the shell phase of Pleione (Limber, 1969) by allowing the envelope to be at any moment very near steadystate. 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 crosssectional flow parameter. CHAPTER III APPLICATION OF HYDRODYNAMICS TO THE ENVELOPES OF Be STARS 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 breakup 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 journals. 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. 30 TABLE I Characteristic Central Star and Disk Values Quantity Stellar Mass Stellar Equatorial Radius Surface Temperature Equatorial Rotational Velocity Disk Temperature Molecular Weight Average Disk Density Symbol M s R s T s V s T e Pd Value 10 M 10R 25,0000K 250 km/s 10,0000K 0.68 1012 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 form f D opf dc Tc is a conserved quantity, Dop, the Boltzmann operator, and f, the distribution function which is a solution of the integrodifferential equation D f= 0 . op 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 socalled 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 longrange 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 103 1/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 bemade 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 result. 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 singlefluid hydrodynamical description is invalid. For a gas under disk physical conditions there are two relevant relaxation timethat 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 threespecies 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 electronproton relaxation time (Tep) to the characteristic time (T ) is T 3.5 x 10+8V T3/2 ep se Tc 2irRN enAd s e d where Ad is the reduced Debye length (Spitzer, 1962). The ratio of the electronhelium relaxation time (Tehe) to the characteristic time is approximately T: V ehe Vs Tc 2~Rshe T1/2 Q S 2 Nhe e he (Zel'dovich and Raizer, 1966) where Qhe is the helium electron collision crosssection for typical disk conditions +4 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 electronproton relaxation time (Tep) to the 4 characteristic time is 4.6 x 104. The ratio which compares electronhelium relaxation time to the characteristic time is approximately 1.79 x 105. Finally one can compare he 0.4 x 101 eP 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 threeelement picture and the actual fluid is the number of neutral species. The quantity Tehe 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 Tehe.. overestimates the actual electronneutral atom relaxation time. Substitution of a smaller number for Tehe 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 form av c + + 1 + + v Vv = VPc Fc at p c c where P is the force per unit mass. The continuity C equation is SC+ V (pcc) = 0 8t c 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 NavierStokes equation, is formed. The NavierStokes 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. NavierStokes Equations and Similarity Numbers The NavierStokes 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 lefthand 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 and T = Ro/v respectively. Note that the characteristic numbers are dimensional. The characteristic numbers can be used to set up a system of dimensionless variables defined by the following relations: + + (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 NavierStokes equations in dimensionless coordinates. The result is shown below: v2o v + V GMs V + v Vv= v2 R t R 0 p R2 r O + 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 names: 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 . RN RV o 2. Froude Number. This value, referred to symbolically as FN, is 2 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 VoTo Ro The inverse of this quantity is the ratio of the time derivative of the velocity to the :inertial term. 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 1015 T5/2 XnAd d 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 108 cm2 W = 6.00 This gives n = 2.4 x 104 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 are: (1) RN = 1.2 x 10+12 (2) FH = 8 x 103 (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 NavierStokes equation to v VP GM v(s 1 P oV 00 Euler's equation. CHAPTER IV MATHEMATICAL ANALYSIS OF THE HYDRODYNAMICAL EQUATIONS 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 NavierStokes equation can be discarded. The momentum equation is GM av + + VP s + v Vv V Pt P RV 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 RV RV 00 ss will be used. The equation of continuity, P + V (p ) = o0 at 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 steadystate axiallysymmetric 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 steadystate axiallysymmetric 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 subscriptzero identifies the steadystate axially symmetric terms; the subscriptone, 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 pressure. 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 2 v +Vv = + 2 + Vpo 1 2 P V q(o(32 p p0 o p 0 o Po o and ap app 4 t + + V + Vp + 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 subscriptzero terms are themselves steady state axiallysymmetric solutions of the problem, they satisfy the equations av o + V P 1 at o o p o o r and ap 0 3^+ V (poVo) =0 The equations for the subscriptone terms become av + v V + v VO + VV = ( ) VP o o 0 p P O and 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 steadystate axiallysymmetric solutions, the subscriptzero terms, represent the gross nature of the flow. The subscriptone terms are taken to represent a refinement to this gross nature. These refinement terms would then be expected to be smallthough not negligible. The assumption will be made here that subscriptone terms are sufficiently smaller than the subscriptzero quantities that elements of the equations which include the products of subscriptone quantities can be omitted. The resulting linearized equations are i Pi 1 +v v + v Vv Vp Vp at + Vo l + V 2 1p p0 PO and api + + + t+ Po V 1 P1v Vo +o 0 P1+ V1, Vo = 0 49 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 + vW a o + W go + V zi w " 21 rz 1av y + V1 Vzo ~a 20 BS aw+ 2v v1 CO p av WO 9aw + v ao V l a zlz 00 zo z V01 1 Po aP1 _5T P1 p ? 0 w a V W1 00 W 1 ap1 p1 ap + + 2 3p o p0 + .o : Vo 1 o , 96 av + Vzo aTW v yv V + " + to 84 tVzo + Vzl  + Vzo z1 5z[ zo Vzl ap1 P,1 apo 5z p O z 2 Bz Po Po and 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 avzl + V too 9(0 v 1l 1 + po aw + 1Pl ap _9E vPl Po 1. = + + v + p 0 . + zo z zl a o z Selection of a Subscript'ze'ro Solution The subscriptzero 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 subscriptone terms. This subscriptzero solution (in dimensionless coordinates) is Vo = 0 , vso = 0 , 1 1 1 =0 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 .,pMsG while Y is given by. ..(.2..) .2 t 2e pm Vs The ratio of these two, Q/', is the Froude number (FN). Substitution of these expressions into the linearized equations for the subscriptone quantities gives avWI 27 av 1 4 = ex i 1 __ ] pi  + 3 72w 7  V7 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 and 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 subscriptzero terms adequately represent the steadystate axiallysymmetric 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 zcomponent of Euler's equation reveals that only the pressure gradient balances the zcomponent of the gravitational force, but the ratio of the former to the latter gives 2a  = 1, 8.3 x 103 3 1 qo Ro 0 r ) The dominance of the zcomponent 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 and 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 expression Vil= Vil ep 1 23 ' 2Qw the i standing for w,(, or z. Now +e 4z2 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 become yve 8ve ap 3Yz2 4" + f 2f 2e4 P at a0 1 aw 2Qw4 SI f1ve f e at W= ap and +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 subscriptone density and the subscriptone pressure to be isentropic with the isentropic constant determined from the subscriptzero quantities, then 2. _ r = 102 (vs/2 ) where c2 is the isothermal speed of sound determined from the subscriptzero 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 steadystate axially symmetric solutions, the interior boundary surface is rotating and has no radial velocity. The steadystate 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 boundaryvalue 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 subscriptone solutions. 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 and 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 subscriptone 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: +z (1) ul = vidzl z +z z 0 Under the condition requiring symmetry about the equatorial plane, the zcomponent 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 welldefined 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 ' and 18 + l a 1 .801 S u + f U+ + 2vzl(z = ) = 0 . The last term of the first equation above will be rewritten 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 61 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 2Qw 2au f *ai r= ^ (2) a B + ul + f 'u = r a U1 2 Ul + f a (3) 3 St + f a a # + + + + U(l 1 '1 = Solutions The following symbols are introduced to simplify the notation: (1) c1 Q (2) c2 = r 8a (3) L = + f p s The equations of the previous section become L Ul 2fu.= Wl p cia1 (I)l c2: 80 u=W = W 0. L 1 + 1 2 and 62 and `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 = Let 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] and 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 k3 + 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) n,k exp[ino + ik(w,n)t], ul = Cn(m) exp[in + ik(w,n)t], n,; and 1a = SnE() exp[ino + ik(w,n)t] n,k 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 n,k 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 c1 i(' + n'f) Un', 2f Cn,, 1 Sn, = 0 , in'c 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 nontrivial 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 nontrivial solutions exist, the determinant formed by the coefficients of the variables in the system above, 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 67 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 UnK X  nk and y  y =.nk nk They satisfy the system of equations c1 (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 + + ] DET r[. nc .(E +. nf) . DET where 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 steadystate axiallysymmetric 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. CHAPTER V .COMPUTATIONS AND RESULTS 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 2irT IPERI = o (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 threefold 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 2 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 and +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 rcalled Tis the value for r in the case that pl and p1 are related isentropically. The isentropic coefficient is determined from the subscriptzero quantities; it is (2w) kT r = .a 2 pm vs For an electron temperature (Te) of 10+4.K, = = 7.67 x 102 A representative value for Acalled Kcan be calculated be setting p equal to a constant; for zo equal to 0.1000, 2 A= = 3.334 x 10 . 3 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]} , and 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 nontrivial 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. 77 IPerl (h) 10  48.5 w3/2 102/ 2.0 4.0 6.0 8.0 10.0 Radius(Stellar Radii) C1  = 0.152 (2w)2 where 2Q Setting c1 = c1 one finds that 2w+ C1 1 lil= 72 ( 1 + ( 227)2 2w 3F2 The absolute value of the corresponding period is IPERI= 2 T or 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 nonzero 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 rA 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 100point 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 outsetthe 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. I w1 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. 83 4.. 