THE UNIVERSITY OF MICHIGAN INDUSTRY PROGRAM OF THE COLLEGE OF ENGINEERING STABILITY OF NATURAL CONVECTION IN A VERTICAL SLOT Charles Mo Vest A dissertation submitted in partial fulfillment of the requirements for the degree of Doctor of Philosophy in the University of Michigan Department of Mechanical Engineering 1967 December, 1967 IP-798

To my Father and Mother ii

ACKNOWLEDGMENTS I wish to express my deepest gratitude to Professor V. S. Arpaci, chairman of the doctoral committee, for his guidance and encouragement during this research, and particularly for his generous correspondence and personal interest during the year that he was on sabbatical leave. I wish also to acknowledge my debt to Professor Chia-Shun Yih, who through teaching, writing, and personal discussion has exposed me to his exacting scholarship and insight. Professor J. A. Clark was most generous and helpful in discussing the experimental investigation, as were Professors P. S. Larsen and W. J. Yang at various stages of the research. The research support of the National Science Foundation is gratefully acknowledged, as is the assistance of the Industry Program of the College of Engineering. Finally, I wish to thank my wife for her patience, sacrifices, and encouragement. iii

TABLE OF CONTENTS ACKNOWLEDGMENTS...................................... LIST OF FIGURES...................................... NOMENCLATURE......................................... CHAPTER I INTRODUCTION AND LITERATURE SURVEY......... A. Introduction........................... B. Literature Survey..................... II THEORETICAL ANALYSIS....................... A. Introduction........................... B. The Base Flow.......................... C. Formulation of the Stability Problem... D. Solution of the Eigenvalue Problem..... E. Evaluation of the Stream Patterns...... III EXPERIMENTAL INVESTIGATION................. A. Apparatus.............................. B. Procedure.............................. C. Results and Comparison with Analysis... Page iii v..................... vii ~~~~~~~~ ~~~~~~~~ ~~~~~~~~ ~~~~~~~~ ~~~~~~~~ ~~~~~~~~ ~~~~~~~ ~~~~~~~~ ~~~~~~~~ ~~~0~~~~ ~~~~~~~~ ~~~~~~~~ ~~~~~~~* *.... *.... *.... *.... eeeoo oo~oo coco,,ooo,,*~oe eoeoo e~oeo oeeeo eoee. oeec~ eoeee eoeoe * * * 1 1 3 15 15 15 31 42 66 75 75 81 84 APPENDICES I INVISCID CONSIDERATIONS.................... II EVALUATION OF INNER PRODUCTS....................... 94 98 iv

LIST OF FIGURES Figure Page 1 Vertical Slot............................................ 16 2 Velocity and Temperature Profiles for H/d = 00............ 22 3 Vertical Temperature Gradient........................... 27 4 Velocity and Temperature Profiles........................ 28 5 Velocity and Temperature Profiles........................ 29 6 Velocity and Temperature Profiles........................ 30 7 Convergence Study of Neutral Curves for H/d = 00.......... 53 8 Convergence Study for the Taylor Problem................. 54 9 Neutral Stability Curve for H/d =.....................55 10 Effect of Prandtl Number on Critical Gr, H/d =oo.......... 56 11 Convergence Study for Finite H/d........................ 62 12 Neutral Stability Curves for Finite H/d................. 63 13 Tentative Neutral Stability Curves for Finite H/d, Pr = 1............................................... 65 14 Disturbance Stream Pattern, H/d = 00..................... 68 15 Disturbance Stream Pattern, H/d = 20..................... 69 16 Total Stream Pattern, H/d = 00, ~ = 0.1................. 71 17 Total Stream Pattern, H/d = 20, C = 0.1.................. 72 18 Total Stream Pattern, H/d = 20, E = 0.5.............. 73 19 Diagram of Experimental Apparatus........................ 77 20 Exploded View of Test Section.......................... 78 21 Test Section............................................ 80 22 Photograph of Experimental Apparatus..................... 82 23 Summary of Flow Regime Data.............................. 85 v

LIST OF FIGURES (Continued) Figure Page 24 Visualization of the Onset of Cellular Instability, H/d = 20.................................................. 86 25 26 27 28 I- 1 Streak Photograph of a Cell, H/d = 20, Ra = 4.5 x 105..... Visualization of the Onset of Cellular Instability, H/d = 33, Pr = 0.71....................................... Visualization of Turbulent Motion, H/d = 50, Ra = 5.2 x 105.......................................... Comparison of Analytical and Experimental Critical States................................................... Antisymmetical Velocity Profile........................... LIST OF TABLES.I Summary of Results.................................... 88 89 90 93 96 Table'age I 92 vi

NOMENCLATURE la thermal diffusivity C function defined by Equation (111) C wave speed of disturbance C' imaginary part of t Up constant-pressure specific heat C real part of C D operator, d/dy j width of slot F Froude number,' gravitational acceleration G Grashof number, ZAT c 3/2 H height of slot h aspect ratio, H/d ij thermal conductivity L linear operator /n parameter, Equation (24) P pressure P, Prandtl number Ra Rayleigh number, G Pt function defined by Equation (112) time T temperature U velocity component in X-direction a dimensionless velocity, U/(i/J) 0o / (/ vii

U characteristic velocity, ~T d 2/> V velocity component in Y-direction Jr dimensionless velocity, /(/Jl) Or V/U W velocity component in Z-direction?iwr dimensionless velocity, W/(v/d) or W/V / wave number of disturbance, 2TT/A / dimensionless temperature gradient, di/cX or wave number coefficient of thermal expansion perturbation parameter dimensionless temperature, (T —Tm ) /T )2 eigenvalue, Equation (118) A wave length of disturbance, or eigenvalue, Equation (113),U eigenvalue, Equation (114) kinematic viscosity p mass density, or eigenvalue, Equation (117) Wronskian, Equation I-5 amplitude function, Equation (78)'X/ stream function for disturbance viii

CHAPTER I INTRODUCTION AND LITERATURE SURVEY A. Introduction In this dissertation the stability of the natural-convective motion of a viscous fluid contained in a vertical slot having isothermal vertical walls of different temperatures is analyzed. At small values of wall temperature difference this motion is unicellular, with fluid rising adjacent to the hot wall and falling adjacent to the cold wall. It is shown that if the wall temperature difference exceeds a critical value this motion becomes unstable and a secondary, multicellular, laminar motion appears in the slot. The critical value of the wall temperature difference and the geometric form of the secondary motion are analytically predicted. These results are corroborated by a flow-visualization experiment. The physical concepts of stability are based on the experimental observation that in a given physical configuration simple laminar motions exist provided the characteristic parameter of the motion does not exceed a certain critical value. If this critical value is exceeded the motion becomes unstable with respect to small disturbances and tends to shift to a new, usually more stable, state of motion. This new state in some instances is in the form of a more complicated stationary laminar motion, which is referred to as a "stationary instability" or "secondary flow." In other instances the new state is in the form of a traveling wave. Often the viscous effects associated with such traveling waves are concentrated in a narrow "shear layer" within the flow. In such cases the waves are -1

-2 referred to as Tollmien-Schlichting waves and the instability leads to turbulence. In the analysis of stability it is recognized that, for a given configuration, the transport equations may allow a simple solution which represents a laminar flow. Such a flow is referred to here as a "base flow." The mathematical disturbance, in two or three dimensions, of this solution yields an eigenvalue problem in terms of the characteristic parameter of the flow. The eigenvalues of this problem characterize the states of neutral stability, in which the disturbance neither amplifies nor decays with increasing time. The base flow of the present problem is represented by two simple models. The first is an exact solution of the governing equations as the height of the slot becomes infinite. This is referred to as the "infiniteslot case." The second is an approximate solution of the governing equations when the slot is of finite height. This model is partially dependent upon experimental observations reported in the literature. It is referred to as the "finite case." The parameters which characterize the base flow are the Grashof number, Gr, the Prandtl number, Pr, and the aspect ratio (the height to width ratio of the slot), ho In some instances the Rayleigh number, Ra = Gr ~ Pr, is used. The disturbance is characterized by its wave speed, c, and its wave number O o Since the wave speed in this case is zero for states of neutral stability —the instability is of the stationary type —the eigenvalues of the disturbance equation are represented as curves in the Gr-Q plane which are parametrically dependent upon the Prandtl number and the aspect ratio. These curves, which separate the

-3 stable and unstable domains, are referred to as curves of neutral stability or simply "neutral curves." The point in the Gr-Q plane at which the neutral curve has a minimum with respect to Gr is referred to as the "critical point" and represents the state at which the flow first becomes unstable according to linearized theory. The present analysis is primarily concerned with the determination of the neutral curves. Experimental verification of the analytical results is obtained by simultaneous flow visualization and wall-temperature measurements. This is sufficient for the determination of the Grashof number and wave number in relatively early stages of instability. These values are compared with those predicted by the analysis at the critical point. B. Literature Survey 1o The Linearized Theory of Hydrodynamic Stability The monograph by C. C. Lin(3 is mainly devoted to the mathematical theory of the stability of parallel and nearly-parallel incompressible flows. A bibliography of 220 references representative of the analytical and experimental knowledge of hydrodynamic stability as of 1954 is included. A large section of Reference 42 is devoted to the same subject and numerous references through 1963 are given. Many other text books and reference works include similar sections. C. C. Lin(3 presents a detailed account of the theory of stability of two-dimensional parallel flows of homogeneous fluids. These papers contain detailed mathematical and physical discussions of the (44) stability of Poiseuille and Blasius flows. Tollmien discusses the stability of flow in symmetrical channels and boundary layers. His paper

-4 is essentially devoted to the study of inflectional instability of inviscid fluids. His treatment of the asymptotic solutions near the singularities which occur in these problems is of major importance in (41) the development of stability theory. Squire presents a fundamental theorem which indicates the sufficiency of considering the stability of parallel flows with respect to two-dimensional disturbances. Another large portion of the literature is devoted to stationary, or cellular, instabilities. The prototype of these is the cellular convection pattern discovered by Benard(4) in a layer of fluid heated from below. The book by Chandrasekhar( ) is largely devoted to a comprehensive study of this type of instability. A very large annotated bibliography of analytical and experimental work is included. Taylor's(43) experimental and analytical investigation of the stability of cylindrical Couette flow is of great importance in the historical context. 2. Related Stability Problems Gershuni(20) considers the stability of antisymmetric natural convection between parallel vertical plates. The stability problem is formulated in terms of both velocity and temperature perturbations. The instability is assumed to be of the stationary type and Galerkin's method is used to solve the resulting eigenvalue problem. Polynomial profiles of the fourth and second order are used to represent the disturbance stream function and temperature respectively. The results indicate a strong dependence of the critical Grashof number on the Prandtl number. In Reference 21 Gershuni extends this analysis to the case in which the parallel plates are oriented at an arbitrary angle to the vertical. In

-5 both of these papers the disturbance temperature profile contains two approximating functions while the disturbance stream function profile includes only one. This weights the solution in favor of the thermal equation and probably leads to an incorrect Prandtl number dependence. Ostrach and Maslin(34) and Ostrach(35) consider the stability of fully-developed natural convection in a vertical channel. The major protion of this work is an analysis of the stability of symmetrical naturalconvective flows. The analysis is essentially an extension of the Heisen(30) berg-Lin method of solution() to include thermal disturbances. It is found that the thermal body force has no effect on the stability of the flow for large values of the wave number —Reynolds number product. Therefore the neutral stability curve in the Grashoff number —wave number plane is independent of the Prandtl number and can in fact be inferred from the known solution for the stability of Poiseuille flow. It is further shown that this conclusion regarding the role of the thermal body force is also valid for antisymmetric natural-convective flows. The reasoning in this case, however, assumes a traveling-wave mode of instability. Gershuni and Zhukhovitskii(22) discuss the stability of antisymmetric natural convection of an electrically-conducting fluid in a vertical channel in the presence of a magnetic field. The Galerkin method is applied in a manner identical to References 20 and 21 discussed above. The stability of the flow with respect to both traveling —wave and stationary modes of instability is analyzed for the case of a transverse magnetic field. In the case of a longitudinal magnetic field it is found that only stationary instability is possible. This

-6 conclusion is based on the fact that the Galerkin method, as applied in this paper, yields only complex solutions for the Grashof number unless the wave speed is zero. Damaskos and Young(l2) investigated the stability of a model of thermally-induced flow in an induction furnace. This is essentially a symmetrical natural convection of an electrically conducting liquid in a vertical channelo The Galerkin method is applied using polynomials to represent the disturbance stream function and temperature profiles. In obtaining the secular determinant for this problem of eigenfunctions are assumed to be real. The discussion of the stability of three-dimensional boundary layers by Gregory, Stuart, and Watson(24) includes both analytical and experimental results. A large portion of this paper is devoted to the stability of the flow of a viscous fluid near a rotating disk. The governing equations of this problem are nearly analogous to the equations governing the stability of a fully developed natural convection. The equation. corresponding to the disturbance energy equation is, however, uncoupled from the equation corresponding to the disturbance momentum equation. The analysis is restricted to the inviscid limit. The analytical results are compared with the results of an experimental investigationo After the analysis of the stability of natural convection in an infinite slot, which is reported in this dissertation, had been completed, two Russian papers on the subject appeared in English translation6'4) In the first of these Birikh discusses the stability of a parallel flow having a cubic velocity profile. This is considered to

-7 represent the natural convection problem for very small Prandtl numberso The Galerkin method is applied using a set of orthogonal functions derived in Reference 5 as approximating functionso The emphasis is on determining the spectrum of eigenvalues for the problem. The first eleven eigenvalues are determined using 18 approximating functions. The solution of the Orr-Sommerfeld equation is considered to be complex in the same manner as in this dissertation. The disturbance stream pattern is evaluated and displayedo The lowest eigenvalues are found to be real (ioeo the instability is of the stationary type). The critical value of the tGrashof number i.s in agreement with that reported here within 0 7 percent. The paper by Rudakov(40) is essentially a correction to Birikh's paper for non-zero Prandtl numbers. The stream-function and temperature perturbations are expanded in powersof the wave number —Grashof number product. Although the higher-order eigenvalues exhibit a very strong Prandtl number dependence, the lowest eigenvalue, which represents the neutrally-stable mode, shows only a small (but not negligible) dependence upon the Prandtl numbero 3. Methods of Solution The technique of solution applied to the eigenvalue problem generated in this analysis is the Galerkin method. The Galerkin method is one of many "weighted-residual techniques" for the approximate solution of differential equations. Because of their formal similarity and mathematical relationship, it is useful to study variational techniques in conjunction with weighted-residual techniques, particularly in the context of hydrodynamic stability.

-8 (28) Kantorovich and Krylov devote a chapter to the development of the variational and Galerkin methods of solving, or reducing the order of, differential equations. The identicality of the Galerkin method and the Ritz variational procedure for many linear problems is discussed. Several examples, including linear eigenvalue problems, are solved in detailo Collatz(11 discusses variational techniques for the solution of self-adjoint eigenvalue problems. A number of mathematical results such as "bracketing theorems" are included. Many texts and reference works include detailed discussions of variational techniques for solving linear eigenvalue problems, but most are limited to the well-developed theory of second-order self-adjoint systems. Ames(1) devotes a chapter to weighted residual methods for differential equations. Becker's monograph(3) is a discussion of the theory and application of certain variational and weighted-residual methods for the approximate solution of boundary and eigenvalue problems arising in engineering analysis. Based upon a set of criteria which its author considers desirable in a variational method, this monograph generally recommends the use of least-square variational methods. A thorough review of the method of weighted residuals is given (19) by Finlayson and Scriven.9) Many methods which belong to this class are discussed and compared as tools for the approximate solution of engineering problems. The conclusion is reached that in general the selection of approximating functions is of greater importance than the selection of a particular weighted-residual method. For most problems none of the particular methods has a distinct theoretical advantage over the others. Since the Galerkin method is the computationally simplest

-9 of these methods, its use is implicitly recommended. A list of 187 references, many involving technical applications is included. A second paper by Finlayson and Scriven ) discusses the method of weighted residuals for the analysis of transport processes. It is shown that several "variational principles" proposed for non-selfadjoint differential systems are not true variational principles, and further that they are actually applications of the Galerkin method in a more complicated computational form. The extension of the method of weighted residuals to systems of differential equations and to vector differential equations is discussed. The application of the Galerkin method to problems in hydrodynamic stability is considered by Di Prima.(13) As examples he treats the stability of the flow between rotating cylinders and the problem of flow over a concave surface. These problems have solutions which decay exponentially rather than satisfy fixed-point boundary conditions. Pellew and Southwell(3) utilized a variational method to solve the eigenvalue problem arising in the determination of the stability of a horizontal layer of viscous liquid heated from below. The differential system describing this problem is self-adjoint and their method involves a classical extremal principle. Chandrasekhar(9'10) proposes and applies a technique for solving many non-self-adjoint systems which arise in stability analysis. This technique, although formally similar to that of Pellew and Southwell, does not have its foundation in a simple (39) extremal principle. Roberts,(39) however, shows that the technique does have a variational basis in terms of an adjoint differential system.

-10 In the problem under discussion here, a certain set of orthogonal functions which satisfy four homogeneous boundary conditions is utilized as the set of approximating functions for the Galerkin method. These functions have been used by Rayleigh(37) and others in the theory of beam (25) vibration. They are tabulated and discussed by Harris and Reid.5) In a second paper Reid and Harris(3 evaluate a number of integrals which arise when these functions are utilized for Fourier-type expansions. 4, Natural Convection in Vertical Slots and Channels A review of the literature regarding natural convection in vertical slots and channels was conducted in order to (1) aid in devising a simple but meaningful model of the base flow, and (2) search for experimental evidence of instabilities of such flows. The early investigations in this area consisted of experimental measurements of heat transfer across vertical air layers. The results of these investigations were generally presented as correlations of Nusselt number or equivalent thermal conductivity with Grashof number. Jakob(27) presents a compilation of such results from several sources. In the past six years, three experimental studies of laminar natural convection in vertical slots have been reported. These studies differ from those discussed above in that they present detailed information regarding the thermal field and/or the mechanics of the flow. The (14) first of these, a paper by Eckert and Carlson, concerns, natural convection of air in a vertical slot with isothermal walls of different temperatures. The temperature field was studied with a Mach- Zehnder interferometer. Three regimes of laminar convective motion were observed and defined. These were denoted as the conduction, transition, and

-11 boundary layer regimes. In the conduction regime, which exists for low Grashof numbers, the temperature varies linearly across the slot in a horizontal planeO Thus heat is transferred from one vertical wall to the other purely by conduction, except near the top and bottom boundaries. In the boundary layer regime, which occurs at higher Grashof numbers, the variation of temperature across horizontal planes is not linear. In this regime rather large thermal gradients occur in thin layers near the vertical walls, while the temperature gradient is essentially zero in the core between these two layers. Thus in this regime the transport of heat from one vertical wall to the other is primarily by convection. The transition from the conduction regime to the boundary layer regime is intimately related to the formation of a vertical temperature gradient in central portion of the slot. At low Grashof numbers no vertical temperature gradient exists, except near the upper and lower boundaries. Thus no vertical conduction or convection of heat occurs, and the temperature varies linearly from one vertical wall to the other. As the Grashof number is increased, however, the temperature is found to increase linearly as the centerline is traversed in the vertical direction. At large Grashof numbers the vertical temperature gradient assumes an asymptotic value of 06 times the difference in wall temperature divided by the height of the sloto Some low-frequency turbulent fluctuations were observed in the core of the slot. At high Grashof numbers wave motions were observed in the vertical boundary layers. Mordchelles- Regnier and Kaplan(32) utilized a differential interferometer to study the transition to turbulence of natural convection of carbon dioxide in a vertical slot with isothermal walls. The

-12 three laminar flow regimes and the vertical temperature gradient discussed above were also observed in this study. Over large ranges of the flow parameters the development of turbulence was found to be qualitatively similar to that which occurs in the free-convective boundary layer on a single vertical plate. Of particular relevance to the present investigation is the apparent observation of a multicellular secondary flow pattern preceding transition to turbulence when the height to width ratio was large. Caution must be observed, however, in making such inferences about the flow field from observations of the thermal field. Elder(15) presents a detailed experimental study of the thermal and flow fields in a vertical slot with isothermal walls. This study differs from all previous studies in two respects: (1) the test fluids have large Prandtl numbers (1000-2500), and (2) the flow field is visualized. The laminar flow regimes and vertical temperature gradient are qualitatively similar to those discussed above. A large amount of data regarding temperature and velocity profiles, and the vertical temperature gradient is presented. The vertical temperature gradient was observed to approach an asymptotic value of o5 to.55 times the temperature drop across the slot divided by the height of the slot. Above a critical value of the Rayleigh number a secondary flow in the form of several cells, all having the same sense of rotation, was observed. At still larger Rayleigh numbers small cells of the opposite rotation formed in the shear layers between the original cells. The related problem of natural convection between vertical isothermal plates which are open at both the top and bottom has been discussed by a number of authors. Ostrach(33) discussed the fully-developed

-13 natural convection of a viscous fluid between semi-infinite vertical isothermal plates. His analysis includes a regular perturbation scheme for evaluating the effect of viscous dissipation. His first order solution (no viscous dissipation) is identical to the conduction regime discussed above when the ambient temperature is equal to the average of the plate temperatures. This work is extended by Lietzke(29) to the case of constant wall heat flux. Heat transfer measurements are also reported. Bodoia and Osterle(7) present a finite-difference solution for natural convection between finite isothermal vertical plates having equal temperatures greater than that of the ambient. The development length was found to vary with the fourth power of the distance between the plates. Three analytical investigations deal directly with natural convection in a vertical slot with isothermal walls. Batchelor(2) presents an analysis motivated by consideration of double-paned windows and insulation gaps in buildings. His analysis utilizes an expansion in powers of the Rayleigh number. It is concluded that a conduction temperature profile and the corresponding cubic velocity profile give the correct asymptotic description of the flow as the height to width ratio approaches infinity while the Rayleigh number remains finite. Transition to a boundary layer regime is predicted for large Rayleigh numbers and finite height to width ratios. It is hypothesized that in this regime an isothermal core region of stagnant fluid exists. As indicated above, experiments later showed that there is a linear temperature increase in the vertical direction within this core.

-14 In a recent paper Gill(23) presents an analysis of the boundary layer regime for large values of the Prandtl number. An approximate solution is obtained by Carrier's modification of the Oseen technique of matched asymptotic expansions. Comparisons with Elder's experimental results are presented and the mechanics of the flow is discussed in some detail Elder ) presents numerical solutions for natural convection in a vertical slot. Of particular interest are cellular motions which were obtained just prior to divergence of the numerical solution. The cellular patterns are similar in appearance to those found experimentally(l5) although they occur at a much lower Rayleigh number. Even though the solutions are nearly unstable in the numerical sense, this result is most intriguing. 5. Experimental Methods A literature survey was not required for the design of the experiment as its general nature was dictated by the analysis. Much information was, however, gained from the directly-related investigation by Elder. (5) The dissertation of Caddell(8) presents useful information regarding the aluminum pigment suspension method of flow visualization.

CHAPTER II THEORETICAL ANALYSIS A. Introduction In this chapter two analytical models of the base flow are constructed. The first applies to the case of a vertical slot of infinite height. The second applies to the flow in the central portion of a slot of finite height. The velocity and temperature profiles given by the second model are compared with experimental data available in the literature. A detailed derivation of the linearized disturbance equations is presented. This derivation includes an extension of Squire's Theorem to the case of parallel natural convection. A brief review of available solution techniques is given, as is a description of Galerkin's method, which was used in this investigation. The solution of the governing eigenvalue problem is described and the resulting neutral curves are presented and their physical significance discussed. The stream patterns which result from the instability of the base flow are evaluated and displayed. B. The Base Flow 1o Formulation The two-dimensional natural-convective motion of a viscous liquid in the vertical slot of Figure 1 is now considered. In writing the governing equations of this flow the following assumptions are made: (1) The fluid density does not vary appreciably with pressure. The density is therefore considered to be a function of temperature only. -15

-16 f T1 I ////I/ AX TI +AT H I — d9g Y'I//// I / Figure 1. Vertical Slot.

-17 (2) The specific heat, thermal conductivity, and kinematic viscosity are constants. (3) The pressure, temperature, and density differ from their values under hydrostatic conditions by only small amounts. (4) Viscous dissipation and work of compression are neglected. Under these assumptions the Navier-Strokes (momentum), energy, and continuity equations for steady two-dimensional motion can be written as: U^ - VJU -R- X+ 72U (1) UV aV. -- ay 2 V (2) g~Uax + 3vY) O (4) where a= k/ / Cp is the thermal diffusivity. On the basis of assumption (3) the pressure, temperature, and density can be expressed as = to+ ^/ (5) T = To T' (6) p =po -t+ / (7) respectively. The subscript "o" denotes a value under hydrostatic conditions at an appropriate reference temperature, and the prime superscript denotes a small perturbation about that value. Since density is

considered to be a function of temperature only, and since T' is a small quantity, p'= ( )p The coefficient of volume expansion, f, is defined by y - aY J\aT)p. Hence the density perturbation, J, can be expressed as: PO = - X - T' (8) Now if the expressions (5) and (7) are utilized, the vertical pressure gradient term in the momentum Equation (1) can be expressed as: 1 a r 1 ar J a 1 _to' _ _ / ax X + X - aA x + O(A/2), Since to is by definition the hydrostatic pressure, Po = -op X + CONSTANT) thus P ax = - + -' o( ). Ignoring terms of quadratic or higher order in the perturbation quantities and dropping the subscript "o" since / - o, the following expression is obtained:

-19 p - 7 --- a- +? AT' (9) If the variation of temperature throughout the field is not large, /O is approximately constant. Then according to Equation (8) all derivatives of P will occur with the coefficient (. Since ( is small (of the order of 10-3 - 10-4 for most fluids) the derivatives of P can be neglected in the continuity equation. This approximation coupled with Equation (9) is equivalent to ignoring density variations except in the buoyancy-force term. This is known as the Boussinesq approximation, and the resulting equations are referred to as the Boussinesq equations. Utilizing the above results, U aV jv / +4T +t aU (10) U_ ~x V ay -,y + Pv72 V (11) 9U +. av a -+ ~ a (13) are the governing (Boussinesq) equations of the flow. In the analysis which follows, interest is concentrated on the case of slots for which H/d >>. It is the aim of this analysis to describe, in an approximate manner, the velocity and temperature fields far from the top and bottom of the slot. In this central region the flow is nearly parallel, so for the purpose of constructing simple analytical models of the flow in this region it is assumed that

-20 U- U(Y), V= O. The following dimensionless quantities are defined: V V x Y 7'-G I TJ0 (14) Pr-, h =.T-T^ 2 2 2 A- TT F = 3 where all quantities to the right of each equality sign are dimensional, and where Tm is the mean of the wall temperatures, T = T1 Tf T In terms of these dimensionless variables and parameters the equations describing the approximately-parallel flow in the central portion of the slot are: Y = - Gr + a x (15) a 2 a 2 y) (16) - = 0. (17) 2. The Infinite Slot Approximation The simplest model of natural convection in a vertical slot is obtained by considering a slot of infinite height, i.e. -00~. In this case the Equations (15), (16), and (17) can be written in the

-21 following form: 2 U -?T &T H -2 - j 2 a dY2- O x'2= 0r'Yg- -' Here dimensional velocity has been reverted to in order to determine a more natural characteristic velocity. The appropriate boundary conditions are: (~- ) = 0 (18) T(+ - )=-, T(-) --. (19) The solution of this system of equations is -6- T Y u = > 6(4 -? ) or, (20) - ='6 - 3) (21) where = -y.-j / Thus as k-0, the temperature field becomes linear, indicating that heat is transferred from one vertical wall to the other solely by conduction. The velocity field is described by a cubic polynomial. These temperature and velocity profiles are shown in Figure 2.

8 u x 103 6 4. Y -8 e.5.25.25 Y -.5 Figure 2. Velocity and Temperature Profiles for H/d = oo

-23 This simple model can also be thought of as describing two other physical configurations. The first of these is the "conduction regime" of natural convection in a vertical slot of finite height. It has been observed experimentally be Eckert and Carlson(1) and Elder(15) that if G Pt /h is not too large the natural convection in a vertical slot is virtually identical to that described by Equations (20) and (21) except near the top and bottom of the slot. Batchelor(2) estimates that the conduction regime exists for G'P-/4 < 500. When compared with experimental results this estimate appears to depend too strongly on the aspect ratio. A conservative guideline, based on Elder's(15) data for high Prandtl number fluids (Pr = 1000 - 2500) and Eckert and Carlson's(14) data for moderately small Prandtl number fluids (Pr =.73), is that the conduction regime exists for Gr.P ( 2X104 when the slot is of moderately large aspect ratio (h = 20 - 100). The second physical configuration which is described by Equations (20) and (21) is the fully-developed natural convection between two very long, parallel, isothermal, plates one of which is maintained at temperature in excess of the ambient by an amount -LAT and the other of which is cooler than the ambient by the same amount. Ostrach(33) has considered this problem using a perturbation scheme to include the effect of small viscous dissipation. His results indicate that the effect of viscous dissipation is negligible in the range of parameters of interest in the present investigation. 3. The Finite Slot Approximation An analytical model which approximately describes the natural convection in a slot of finite aspect ratio was also constructed. In

-24 addition to being a logical extension of the stability investigation of the "infinite case," analysis of this model permitted comparison of the predicted stability criterion with experimental results obtained by Elder.15) The construct this model, the flow is again considered to be parallel in the central portion of the slot. Thus Equations (15), (16) and (17) govern the problem. The aspect ratio, h, is large but finite. Differentiating Equation (15) with respect to Y and utilizing Equation (17) the relation d3U __ di y3 - t y (22) is obtained. Since it is assumed that aU=U(Y), this indicates that = (Y). The term o in Equation (16) can be neglected since h is large and Pr is not small. Then since is independent ay of X, Equation (16) indicates that a~ must also be independent of; that is the temperature must be of the form 0 = it(Y) + X 2 (y). Equation (22), however, indicates that i2(Y) must be a constant, hence the temperature can be expressed as'-= T'() + j3:, (23) where y - x'

-25 A solution of the form (23) could be strictly valid only at X=O because the imposed thermal boundary conditions do not admit vertical temperature gradients at the walls. Neglecting the term whose coefficient is, and introducing the form (23) of the temperature field, Equations (15), (16) and (17) can be combined to give g4 -d4t + 4- 4' = 0 O n 4 = Ra _44 ) (A1) where Ra = G P, is the Rayleigh number. Continuity considerations and the no-slip boundary condition require: (1) Antisymmetry of U about' = O (25) and (2) U (~ +) = 0. (26) The solution of Equations (24) through (26) subject to the thermal boundary conditions, T(tC )=- yields the following expressions: 7= K(- " ((/&))I( (y) +(xO)& C; (ry) ) (27) T=. -2/anK1(42) A (n,(y)coA)('My) + Cot ( )J4 (, )) (28) ST = ~- ~2n^ K \ ((~^/ i4 (. -Y) + C." + Y (29) + C( Y-c (IVy) - ( ) ( ) ( a29)) d a-2= -T; (30) where K = 4/ ) +cM(/2) 4 /2), (31)

-26 The derivatives f and e are included here for later reference. The above analysis does not prescribe the value of the dimensionless temperature gradient 25. To estimate J recourse is made to experimental observations reported in the literature. As discussed in Section I.B. of this dissertation, Eckert and Carlson,(14) Elder,(15) and Mordchelles-Regnier and Kaplan(32) all report the existence of a "convection regime" of laminar motion in vertical slots such as that considered here. In this regime a remarkably constant vertical temperature gradient forms along the centerline of the slot. This is a stable temperature gradient, that is the temperature increases upward. Figure 3, which is based on Elder's(5) Figure 5 shows the formation of this gradient with increasing Rayleigh number. The asymptotic values of f reported in these three references varied from 0.5 to 0.6. Based on these observations, the value of the dimensionless vertical temperature gradient is take to be /3 - 05. Figures 4, 5, and 6 show typical velocity and temperature profiles obtained by evaluating Equations (27) and (28) with Ah=0.5. Also indicated on these figures are experimental points derived from Elder's(15) measurements at X=O. Of interest are the small flow reversals indicated at large Ra/h. The value of Ra/i at which these occur, however, are beyond those of interest in the stability analysis. Gill(23) has recently considered the two-dimensional motion in such slots for values of Rs/Ah sufficiently large for boundary layer theory to be applicable. Within that framework he shows that the simple model considered above is actually the correct one at X=O.

-27 106 107 Ra Figure 3. Vertical Temperature Gradient.

x.4 0.8 I 0.1 IF-.2 -.3.4 5.3.4.3.4.5 x ro 0 -1.2 1.6 (15) Ra= 2.95 x 105 H/d= 13.7 E EXPERIMENTAL VALUES ANALYSIS Figure 4. Velocity and Temperature Profiles.

0 -.4 -.8 0 x -1.6 -I.6 0 0 I _ 0 -.2 -.3 -.4 - 5.4.5 r l S EXPERIMENTAL VALUES(15) 0 EXPERIMENTAL VALUES Ra= 6.56 x 105 H/d = 13.85 - ANALYSIS Figure 5. Velocity and Temperature Profiles.

-30-.4.5 O -I 0 x -2 1. 0 Experimental Values(15) Analysis 0 Ra = 3.61 x 106 H/d = 7.19 -3.I 0 -.1.1.5 -.4 Figure 6. Velocity and Temperature Profiles.

-31 For large Prandtl number fluids, Figure 3 indicates that the solution obtained above is valid for RF > 8 O. Gill(3) estimates a lower limit of Ra/h = 2x104for the validity of his boundary layer analysis. For very small values of Ra/h the velocity and temperature profiles (27) and (28) become asymptotically identical to the profiles of the "infinite-slot case." C. Formulation of the Stability Problem 1. Development of the Disturbance Equations In this section a differential system is developed which governs the spatial and temporal behavior of a disturbance of infinitesimal amplitude but arbitrary form which is superimposed on a parallel, two dimensional, natural-convective flow. The Boussinesq equations, which were discussed in Section B.l, can be written in the following three-dimensional form: Dt ~ -- + T7- + v t72U (32) _. =_ _ p V (33) Dt A + V DT= a VT (35) Dt aU, aV + aw (36) X a Y aZ u U, V, and W are respectively the velocity components in the X, Y, Z directions of Figure 1 (Z is positive into the plane of that figure). T and p are understood to be the deviations of temperature and pressure

-32 from their values under hydrostatic conditions at an appropriate reference temperature. In order to non-dimensionalize the above equations the following dimensionless variables are defined X = e =; au= U V uJ,~ 2 U/a Uldd w= W 4TA T bT u > (37) where all quantities to the right of each equality sign are dimensional. Tm and AT have the same meaning as in the discussion of the base flow. The velocity U is defined by - Y d2(38) U = ". (38) This characteristic velocity of fully-developed natural convection was evolved in the discussion of the base flow. In terms of these dimensionless variables the momentum, energy and continuity equations become DU Dt Dv Dt Dt _ ^+ / t - ~ _ + ~ -ap +- 7 - a 7 2 /r 2 Gr (39) (40) (41)

-33DT 2T Dt- = R Y V' (42) aW + 7 + = 0 (43) where _ Gr - ^>3= U> >(44) and = -caTds Uo.p Ra - - = P.. (45) Hence the Grashof number can be interpreted as a Reynolds number in terms of the characteristic velocity U. The Rayleigh number then represents the corresponding Peclet number. The disturbance of the base flow is accomplished by a simple regular perturbation scheme. The dependent variables of Equations (39) through (43) are expressed as power series in E, which represents an unspecified constant parameter which is small compared to unity: U + U U/ + G2 U o ~~ = -VL/ +4 _ //+ o+ ElA^ 6W' 4 2A T'- 4 h a T = T + eT' + TT &/ + *' ) P = + ~' + E +,. ~ (46) Since 2, T, and P are functions of y only, this perturbation scheme yields the following sets of equations:

-34 r~ ch - T d ~ d__ = o c12T dy = at + ax =ay.'- + aT + -0. 9?^ ax a^ I9 ff~L r J f-^v ^ -r-i-O. (47) (48) (49) 00 e e + T', G. 72,/ / 7V2. + 72 wh + o Vw' (50) (51) (52) (53) (54) 72T Since it is intended to study the stability of the base flows discussed in Section II.B., the boundary conditions (18) and (19) are imposed on the total flow. Hence;~~(~-L = 0 (55) T(~ ) = + -L (56) U/(~I) = v /(~)= U/ (~ ) = T'(~') = 0. (57) (58) Terms of quadratic or higher order in E are not considered, The zero-th order equations and boundary conditions (Equations (47) through (49), (55) and (56)) represent the base flow problem which was solved in

-35 Section II.B. The first order, or "disturbance" equations and boundary conditions (Equations (50) through (54), (57) and (58)) are the linearized governing equations of the disturbance. Because the coefficients in the set of linear differential Equations (50) through (54) are functions of / only, the set accepts solutions of the form: = (,Y) ei (2 X +R a - YC t) (-Y i(ax +iR -~td t) it' =-ttr(i) e (59) -/- - t( ~)e 2(~x *sa — Rcet) These solutions can be thought of as representing single Fourier components of disturbances of more general structure. Since the governing system is linear, a flow which is found to be stable with respect to all disturbances of this simple form is stable with respect to all infinitesimal disturbances, regardless of their particular form. The real parts of the solutions are considered to have physical significance. The parameters af and /3 are the wave numbers of the disturbance in the X and a directions respectively. They are restricted to real values in order that the disturbance will be bounded as X and approach infinity. The parameter C is the wave speed of the disturbance and is in general complex.

-36 For functions of the form t(y) e + - the following equivalences of operators are apparent: a3% - Z at =oi b = D V2=(D2- 2_-^) /(60) — J = f~ where D =. Thus introducing solutions of forms (59) into the disturbance equations and boundary conditions yields *1 (o2 AdA32) 4 = iG(2 -C)h% +DU + iP-?-jT (61) Gr(D-Yz _-2)- = i u- + UD (62) (2Ra j — 2)- = 2-(3 - -)t + / P (63) ul (~)- (-'t)T-= (_-c)T - E v (64) ZiLY7 +P +i#z3 = o (65) (I) = (2) =W()= T( ) = 0. (66) This differential system governs the spatial and temporal behavior of an infinitesimal, three-dimensional disturbance superimposed upon a parallel natural-convective flow. 2. Extension of Squire's Theorem Squire(41) proved a theorem of great importance in the study of the stability of parallel incompressible flows. This theorem states

-37 that the problem of a three-dimensional disturbance in a parallel flow is equivalent to a two-dimensional disturbance in a similar flow at a lower Reynolds number. It will now be shown that this result is also valid for the present case in which a thermal field and the resulting buoyancy force are considered. For convenience, Equation (61) is multiplied by 0' and the result is added to 3 times Equation (63). The equations to be considered are then: Gr- (D2-a,-/ )1 = z (u-c)v + D (67) + Q ( D i- t- ) ( iv) = / ( ru-C) (a +w ) + ~'D- + / (t+3 2) -_ _'T (68) G-(D —a2-3 )T = b(U-d)T + _ w (69) iZ (oU + ~) = -wD. (70) Now the following transformations are introduced: a, =Z - aU + Ju ce = c. l,," p,2 +,2 ~2 2 P' G~~~~~ ~ G,.,^~ P'T~ ~ (71) G, = A Gt T = 7 NT71 Ra = a Ra tG = r,^.

-38 Under these transformations the governing differential system becomes +'1 -+ ^(Da _ z) = ZGr G +( -^ ) r + - + iz' C- r + T (72) (P2 _ LZ)r ='Gr ( -~) + 5G Di (73) (D- z) = a( - ) + dz (74) Z ^ +D = O (75) ~(~t ) = v (+-) =T(~ 3) = o. (76) The system (72) through (76) represents the problem of a two-dimensional disturbance in a base flow with Grashof number given by G&- r G+. Since -_ +,/3 it follows that > a; hence the Grashof number for the two-dimensional problem is smaller than that for the three-dimensional problem. From this it is concluded that for the determination of the stability criterion it is sufficient to consider two-dimensional disturbances. On the basis of this extension of Squire's theorem, the twodimensional form of Equations (61) through (66) is now considered to govern the disturbance, that is r=0-O,3=O. In order to combine the x- and y- disturbance momentum equations, it is first noted that according to Equation (54) the disturbance velocity field is solenoidal. A stream function can therefore be introduced as:

-39 7/ -_?/1 - ay / 1/ - aZX (77) The stream function / can have the form io (x -d t)' = e) e2(zt) (78) Thus. (x -d t) u' = (D a) e ( - ) and (79) v;~ = -z'~eia(x-!t). The two-dimensional form of Equations (61) through (66) can then be written as A - (DP-Z)D = z'i,(-c)Dp -i'Du a (80) + ~~Pp- - T - L (v~-a 2)i' = (8i-)) D (81) R ( a2- 2)T = _z4(-d)T7 -4'z (82) = D (~ ) T(~ i) 0. (83) Upon eliminating the pressure, p, by combining Equations (80) and (81), the following form of the governing differential system is obtained:

-40 (2- r ) Q + i Grr Dz; I - -( )(D2- ^) +DT =0 (84) 2( -a')T ~ + {a -(-c)T3 = 0 (85) (P~) r = D^(i) = T(+1~) =. (86) It is observed that the first Equation (84) is the classical OrrSommerfeld equation with the addition of a buoyancy-derived term which couples it to the thermal disturbance equation. Equations (84) through (86) are a system of homogeneous ordinary differential equations subject to homogeneous boundary conditions. This constitutes an eigenvalue problem, i.e. solutions of the system will exist only for certain sets of values of a', Gr, Ra, and c. The solution of this eigenvalue problem is discussed in Section II.D. 3. The Assumption of Zero Wave Speed The disturbance velocity and temperature each have a time dependence given by exp (-! C t ), where c is in general a complex number. The imaginary part of dci jdescribes the amplification or damping of the disturbance. Since the present analysis pertains to states of neutral stability, dC = 0 If the real part of d, C1, has a non-zero value the disturbance is in the form of a traveling wave. In this problem, however, there appears to be no preferred direction for wave travel. This suggests the possibility of the occurance of a stationary instability, that is one which is characterized by cr-=O. This argument implicitly assumes that the disturbance is not confined to narrow regions of the flow, for in that case there would be preferred directions except near Y=0.

-41 Striking evidence of the existence of stationary instabilities can be found in Elder's experimental investigation of the laminar natural convection of high Prandtl number fluids in a vertical slot. At small values of the Rayleigh number the motion in the slot was found to be unicellular —rising adjacent to the hot wall and falling adjacent to the cold wall. When the critical value of the Rayleigh number ( 3 x 105 ) was exceeded a new laminar flow regime in the form of a multicellular pattern was observed. The cells all had the same sense of rotation as the base flow and were slanted with respect to the axis of the slot. At higher Rayleigh numbers smaller cells of the opposite rotation were observed to form in the shear regions between these cells. Further evidence of stationary modes of instability is found in the interferometric study of natural convection of pressurized CO2 in a vertical slot reported by Mordchelles-Regnier and Kaplan.(32) They cite evidence of cellular motions preceeding transition to turbulence when the ratio of slot height to width is large. Although such interpretation of interferograms is subject to caution, it appears likely, in view of Elder's flow visualizations, that such instabilities were in fact observed. On the basis of the antisymmetry of the flow and the observations by Elder, the wave speed, C, is taken to be zero for neutral disturbances. Further comment regarding the assumption of zero wave speed is to be found in Appendix I. This analysis may therefore be classified as a study of the stability of antisymmetric natural convection with respect to the stationary disturbances. The final form of the differential system governing the dis turbance is

-42 Mlv- 2 ^( + ^4m + i+4G/[ {2 2 U -U m + T/ = 0 (87) T Z "- t2T + iR IT' - 7 T- = 0 (88) Am~++ = 1(~1t) =T( I) = 0. (89) Since Ra = Gr * Pr, the eigenvalues of this system form some functional relationship -(G,Pr, ) = 0. (90) The determination of this relationship is the subject of the next section. D. Solution of the Eigenvalue Problem lo Available Solution Techniques Most of the classical problems of hydrodynamic stability involve either Tollmien-Schlichting waves in the plane of a parallel base flow or stationary instabilities in a plane perpendicular to, or in the absence of, a base flow. The Tollmien-Schlichting waves are governed by the OrrSommerfeld equation and are represented by complex eigensolutions of rather complicated structure. The stationary instabilities are governed by less complicated equations having real eigensolutions of fairly simple structure. The present case, which involves a stationary instability in the plane of a parallel base flow, is somewhat of a hybrid of these two types. Its eigensolutions are expected to be complex but of simple structure.

-43 Because the governing differential system of this problem is not self-adjoint, the eigenvalue problem cannot be formulated as a simple variational principle as it can in some problems of stationary instability. Chandrasekhar,(90) has devised a method for the approximate solution of non-self-adjoint systems which govern certain stationary instabilities. This method could have been extended to apply to the present problem; however, it seemed to offer no advantage, because its application would involve solving the Orr-Sommerfeld equation with an additional, nonhomogeneous, term. Various numerical schemes are available for the solution of eigenvalue problems. Since the solution of this problem is expected to be of a simple structure such schemes could no doubt have been applied. Any of the weighted residual techniques could have been used to obtain an approximate solution to the problem. One of these, the Galerkin method, was in fact utilized. 2. The Galerkin Method The Galerkin method for the approximate solution of an ordinary differential equation is as follows: A differential equation, L(-L) = - 0 (91) where L is a differential operator in the independent variable x, is posed subject to homogeneous boundary conditions at the points x = a and x = b. The approximate solution is expressed as

-44 In u (X) - a ji /i(x); (92) i-!1 where each'Ui(X) satisfies the boundary conditions. This series is substituted into Equation (91) which is then multiplied by Uji and integrated over [a? b: \ L f Qi ~(x),(x) = o. (93) Equation (93) represents a system of n algebraic equations in n unknowns whose solution gives numerical values for the coefficients in the series.(9) In the case that L() and the boundary conditions are all homogeneous, the Equations (93) lead to a secular determinant whose zeros approximate the eigenvalues of the problem. The following definition of a complete system of functions given by Kantorovich and Krylov(28) is quoted here as background for the discussion of the Galerkin method: A system of functions {(^(X)} is said to be complete if no function q (X) exists with al(X)l2J- > 0 and which is orthogonal to all the functions Q(X) simultaneously, i.e., such that aS BX) (X) = o (M=1,a,2 *), If the set of approximating functions {Ui} is a complete set in [, b] Galerkin's method has the following mathematical basis: Equation (93) states that each function Ui is rendered orthogonal to L{ (Z)} over Lb]. As n ->+O the condition Li}- O0 must be fulfilled, according to the definition of completeness, because L{U3 is orthogonal to every member of the complete set ({i}.

-45 The Galerkin method can also be applied in a similar manner to systems of differential equations. This is discussed briefly in Reference 18 and is described in the following section. No proof of convergence sufficiently general to apply to the present problem is known to the author. 3. Solution of the Infinite Slot Problem The neutrally-stable states of the natural-convective flow of a fluid in a slot of infinite height are determined by solving the eigenvalue problem (87) through (89) using the base flow given by Equations (20) and (21). Let Equations (87) and (88) be denoted respectively by L1(,VT) = 0 (94) and L2 ({, T) 0. (95) Then Galerkin's method, as used here, consists of representing the solution as series, N =C q= c t (jY) (96) i-/ N where the approximating functions, (i and T7, individually satisfy the boundary conditions (89). Each approximating function is then rendered orthogonal to the appropriate residual over L-2 J + S LI Z bTi) sdy= ~ (98) ^- - 1= z=/

-46 [^ L2 (Z] ( a,, b>T,. T,. d = 0. - 2 f St,z'=~' iC 7 (99) Since the series (96) and (97) each consist of N terms, Equations (98) and (99) represent 2N linear, homogeneous, algebraic equations in 2N unknowns. Those combinations of f, Gr, and Pr for which the determinant of coefficients of this set of algebraic equations vanishes are the approximate eigenvalue sets of the problem. For the particular base flow under consideration u = ( Y - ) (100) "// = -e (101) and T' =1 (102) Because of the antisymmetry of this base flow and the symmetric disturbance boundary conditions, the solutions of the disturbance equations are expected to have simple symmetry properties. Inspection of Equations (87) and (88), however, indicates that because U and T are odd functions of " the solutions cannot be solely even or solely odd functions of'. However, solutions can be sought which display simple symmetry properties and yet are consistent with the form of the disturbance equations. They are those which can be written as,P- E + 0 T7^ ~= o: + 2ie~ X (103) U- ~~ + ie

or $ - 0 +z E > (104) T = e t z c' where E and e denote real even functions of y, and 0 and o denote real odd functions of y. This can be seen by substituting (103) into Equations (87) through (88) and separately equating the real and imaginary parts of each equation to zero. This yields E'" - 2 /E" 4E - aGa 0+.'"0-ouO" + C3/ =0 (105) 0" - _ 2C - Ra /T'e - O} =0 (106) OV - 2 2 0 a2 0 + aG {a2[UE + "/E - E"} + e =o (107) e"- e.+' Ra' E - C>} 0. (108) Each of these equations is consistent regarding symmetry. It is inconsequential whether the form (103) or (104) is considered. The resulting eigenvalues would be identical and the solutions Tr would differ only by a phase shift of E radians in the X direction. Solutions of the form (103) are assumed in the analysis described below. Because no prior knowledge of the functional structure of the solution existed, and because they can be very systematically applied, simple sets of orthogonal functions which satisfy all the boundary

conditions were used as approximating functions. In particular, c, which must vanish along with its first derivative at the boundaries, was expanded in terms of the eigenfunctions of the system. ~"V = a- 4 (109) 2 )= &(4) = ~.0 (110) These eigenfunctions fall into non-combining even and odd groups. They can be written in normalized form as C= (Y) -,a RC /,) (111)-( and, _ 4^ A(F Y) _ ______ ^( ) ver> (1 8>) XL> (i ^); (112) /n = 1, 2) ") where a,, and R, are the roots of t(A iJ4) + t ( )} = 0 (113) and V (a M) ( 1(M ) = 0o (114) These functions have been used extensively by Chandrasekhar and others in Fourier-type expansions in hydrodynamic stability analysis. They have been tabulated and discussed by Harris and Ried.(25'38) C A and S have been referred to in the literature as a complete set of functions; however, their completeness has not been proved in any published work so far as the author can discern. The property of completeness, although central to the mathematical foundation of Galerkin's

method is often overlooked in applications, where the orthogonalization is limited to a finite number of terms anyway. T, which must satisfy only the two boundary conditions T (f(~~} was expanded in terms of C% (ye) = CB (J, ) (115) and AQl_ (e/) = (l () (116) where j~ = (2 n - 1) T (117) and ){H = 2 TT; = =, 2,. (118) Because of the symmetry considerations discussed above the solution was expanded as P(y) = 3 (a^C(Mc ) + 8 b^ ^(y)), (119) T(r ) = (d b, ) +: e cad (a)), (120) 4= 1 where bt e, b, cl, and em are real constants. Upon substituting the series (119) and (120) into the Equations (87) and (88) and applying the orthogonalization criteria (98), (99) a set of homogeneous linear algebraic equations is obtained for the constants, r, b-,, d t, and em. A necessary and sufficient condition for non-trivial solutions of this set of equations to exist is that the determinant of coefficients vanish. This criterion yields a secular equation which can be represented as

-50 (X11) (X12) (X14I) (X22) (X2) - - - (XNN),._..-. _ _ _ r =0. (121) The elemental matrix, (X,/,), from (121) is constructed is defined by which the determinant of Equation 0 Gr JrYY (X /m) = oiRa AAr^n \ O - GrE^, 0 - O/Ra jirm H/i 0 INT NMm 0 (122) C Cc, ] U Ra LL^,, / where AM = (c'X c) - 2 u2 (C /IC,) + aY4 (C, ) HW/m, = c (-^t^^ Cm ) EU - =(CL|2, +(ic, ~)-(C U)(CI) LC CM =(o1 ), ) - 2 (,, ) + 4 (S IS ) LmL = ( C0, I SJ ) A A/ =( C IT'I <C ) GGnnm = (le~n, UIiCO-T) CC = (Cvz C ) - rM2 (Cey ~CF x) NN/n, = (/ 1 ) - 2" (^11) LLF =(C^, || ) > (123) J

-51Here ( 11- _, () B() ). All of the inner products which were required for the solution of Equation (121) were evaluated by exact integration. The inner products were therefore expressible in terms of the values of the C, S, cosine, and sine functions and their derivatives at the boundary Y= - The inner products were evaluated by repeated integration by parts and use of the generating equations of the functions. The method is described in Appendix II. That appendix also includes a tabulation of a large number of inner products involving the C, S, cosine, and sine functions and powers of y ranging from zero to five. The values of the C and S functions and their first three derivatives at "- = can be obtained from the paper by Harris and Ried.(25) The secular Equation (121) was solved numerically on the IBM 7090 computer at the University of Michigan Computing Center. For specified values of the Prandtl number and the wave number, the values of the elements of the determinant in Equation (121) were calculated. The determinant was then iteratively evaluated for various values of the Grashof number. The value of the Grashof number for neutral stability at each specified wave number was thus bracketed by a sign change of the determinanto This value could be bracketed to any desired degree of (algebraic) accuracy because the computer was programmed to search for the zero of the deteminant by an "interval-halving" technique. In order to investigate the convergence of the solution, this procedure was followed for secular determinants of order 4, 8, and 12. The order was increased in steps of four so that equal numbers of

-52approximating functions for the real and imaginary parts of both the disturbance stream function and temperature were always considered. This avoided artificial weighting of the disturbance momentum or energy equation. The results of the convergence study are presented in Figure 7 (Computations for this figure were carried out for Pr = 1000). This figure indicates that a secular determinant of order 12 gives sufficient convergence in the vicinity of the critical Grashof number. This of course does not constitute a proof of convergence to the correct eigenvalue. As a further check on the method it was applied in an analogous manner to the B6nard problem and to the "narrow-gap approximation" of the Taylor problem, both of which possess solutions known to a high degree of accuracy. The solution of the Benard problem was a simple task, and even a secular determinant of order four gave a reasonably accurate answero The Taylor problem, which involves an asymmetric base velocity profile, presented a slightly more difficult task. Experience gained with the Taylor problem proved useful in the extension of this work to the "finite case." Figure 8 shows a dramatic change in the accuracy of the solution of the Taylor problem with an increase in the order of the secular determinant. The neutral curve for the stability of the flow under consideration is given in Figure 9. The critical Grashof number was found to be 7880 corresponding to a wave number of 2.65. Points on the neutral curve near the critical point were computed for a wide range of values of the Prandtl number. No significant variation of the critical Grashof number with Prandtl number was found (See Figure 10). The jump in the curve of Figure 10 between Pr = 102

ro D UI^ e \ l 1 ila I a 11 LI o 8 X 7I 6 I I I 214 O 1 2 3 4 5 WAVE NUMBER, a Figure 7. Convergence Study of Neutral Curves for H/d = X.

100 k 0 a) a) C) O r-q 0 E-i H rl O a) k PL, -2 -I Q2/Q, 0 Figure 8. Convergence Study for the Taylor Problem.

-55 6 5 It 4 I0 x b.w O( 3 OF - 0 2 0 I 2 3 4 WAVE NUMBER, a Figure 9. Neutral Stability Curve for H/d - o 4 5

8.C I _____ _ _ _______ _ _ _____ JT r 7.96 i tr) 0 X 7.92 I II I II I I iI' I I I I Pr= O. I I I I I I I I I I I I I I I _ _~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ _ 3~~~~~~~~~~~~~~~~~~~~~~~~,,,,,,, l l l~~~~~~~~~~~~~~~~~~~~~ 7.88 U 7.84 780 I!\ I 6. io-'o 10-8 o106 10-4 o102 100 102 Pr Figure 10. Effect of Prandtl Number on Critical Gr. H/d = oo.

-57 and Pr = 10-3 is due to a shift in numerical dominance of terms in the secular determinant which is inherent in this range of Prandtl numbers. This jump represents only about half of one percent and should be ignored in comparison with the general accuracy of the solution. This virtual lack of Prandtl number dependence is interpreted as indicating that the problem is of a purely hydrodynamic nature. That is, the thermal buoyancy forces play no significant role in rendering the flow unstable. The validity of this interpretation can be verified by considering a fluid having a very large thermal diffusivity. As Pr- O with Gr remaining finite, the disturbance energy equation is simplified to T' -a' 2 - = (124) subject to T (i ) = o. (125) This system has only the trivial solution T = 0. Thus the problem is governed by Orr-Sommerfeld equation alone. This was solved by the Galerkin procedure using 16 approximating functions in a series of the form (119). The critical Grashof number 7932 at / = 2.6, agreed well with the foregoing results and is indicated as an asymptote for Pr = 0 in Figure 10. (20) Gershuni() found a large variation of the critical Grashof number with the Prandtl number for this problem. His work was rather incomplete and only a small number of approximating functions (polynomials) was considered. This Prandtl number dependence was probably due to the fact that more terms were used for the disturbance temperature than for the disturbance stream function thus giving additional weighting to the energy equation when the Galerkin method was applied.

-58The lack of Prandtl number dependence is in agreement with the (34) findings of Ostrach and Maslin,(3 who considered the different, yet related, problem of Tollmien-Schlichting waves in fully-developed natural convection between vertical parallel plates. They found that the stability of the flow with respect to such disturbances was not effected by the buoyancy force. For instance, if both plates have the same temperature the base velocity is parabolic and its neutral stability curve is identical to that of plane Poiseuille flow. Birikh( ) published a paper on the "infinite case" after the present research had been completed. This work corresponds to the case of zero Prandtl number, and was solved by Galerkin's method. The critical Grashof number predicted by him agrees with that reported in this chapter to within less than one percent even though different approximating functions were used in the two studies. The spectrum of the next ten eigenvalues is also studied. Rudakov(4) extended Birikh's problem of eigenvalue spectra to the case of non-zero Prandtl numbers. He utilized an expansion for small values of the wave number —Grashof number product. He extended this solution to obtain an approximation to the neutral curve for small values of the Prandtl number (0.01 - 0.2). The critical Grashof number obtained in this manner is in reasonable agreement with that reported here, but shows a small, yet significant, increase with increasing Prandtl number. 4, Solution of the Finite Slot Problem As discussed in Section II.B, recent experimental investigations provide information regarding the formation of a vertical temperature gradient in a fluid contained in a vertical slot of finite height. This

-59 gradient rapidly approaches a constant value as the Rayleigh number is increased. This information made possible the construction of a simple analytical model of the flow in the central portion of a slot of finite height. This model is approximately valid over a large range of Rayleigh numbers and aspect ratios. The velocity and temperature fields are described by Equations (27) through (31). The study of the stability of this flow appeared to be a logical extension of the analysis described in the preceding section. An additional benefit of such a study is that the results can be compared with Elder's5) observations of the formation of multicellular secondary flows. If the height to width ratio of the slot is reasonably large, the flow in the central portion of the slot is nearly parallel. Thus the stability of a parallel flow having velocity and temperature profiles given by Equations (27) and (28) was studied. A limitation of this analysis is that the model of the flow utilizes the asymptotic value of the vertical temperature gradient,,1 = 0.5. This assumption is not valid for Rayleigh numbers below about 8 x 10 The neutrally-stable states of this motion are determined by solving the eigenvalue problem (87) through (89) using the model of natural convection in a finite slot as the base flow. The symmetry conditions discussed for the infinite slot problem are applicable to this case also. Therefore the eigenvalue problem was again solved by Galerkin's method with the solutions expanded according to Equations (119) and (120). The secular equation is given by Equation (121) and the required elements are given by Equations (122) and (123). Although this method of solution is formally identical to that of the "infinite case", the numerical procedure required for the solution

-6o of the secular equation is somewhat different. This stems from the fact that the base velocity and temperature are now functions of the Rayleigh number and the aspect ratio. Because of the complexity of U and T the inner products in Equation (123) were evaluated by numerical integration, rather than by exact integration. The procedure which was followed in searching the (6, Gr) plane for zeros of the secular determinant is: (1) A value of the Grashof number was selected (the Prandtl number and aspect ratio having been specified). (2) The velocity and temperature profiles and their derivatives at the corresponding Rayleigh number were computed at 40 equal intervals of y between 0 and 0.5. (3) The values of the C, S, cosine, and sine functions and their derivatives were computed at the same values of y. (4) The appropriate combinations of these functions, velocities, and temperatures were multiplied together and integrated by Simpson's rule to provide numerical values for the inner products in Equation (123). (5) The secular determinant was iteratively evaluated for various values of Y. In this manner the value of the wave number for neutral stability at each specified Grashof number was bracketed by a sign change of the secular determinant. The neutral curve can be traced in either the (RaY) or (GraY) planes with Pr and H/d as parameters. The neutral curves were generated for aspect ratios of 10, 20, and 50 for various values of the Prandtl number. In order to investigate the convergence of the solution of the eigenvalue problem the neutral curves were evaluated in the vicinity of the critical point using secular determinants of order 4, 8, 12, and 16. No real eigenvalues were found using a determinant of order four. This is probably because the first

-61 two C and S functions do not oscillate sufficiently to allow construction of the stream functions which, in comparison with those for the "infinite case", are relatively steep near the boundaries. The results of this convergence study, for Pr = 1000, are presented in Figure 11 a, b, c. The maximum difference between the critical Rayleigh numbers corresponding to secular determinants of order 12 and 16 is three percent, which occurs for H/d = 20. Considering the computer time required for the generation of these curves it was decided that this accuracy would suffice. Curves of neutral stability for several values of the Prandtl number and aspect ratio are shown in Figure 12. These curves were generated with a secular determinant of order 12, which was deemed sufficient to illustrate the effect of varying these parameters. Generation of neutral curves for Prandtl numbers of 12, 10, and 8 was attempted but satisfactory results were not obtained. It is possible that the use of larger secular determinants would lead to solutions for these cases, but it was decided not to invest in the amount of computer time which would be required. It should also be noted that the model of the base flow is not valid for Rayleigh numbers below about "8 x 10. It is possible that the true neutral curves for these moderate Prandtl numbers extend to Rayleigh numbers below this value, in which case the present formulation would not correctly describe the physics of the problem. As Ra/H/d becomes small the base velocity and temperature profiles of the "finite case" asymptotically approach those of the "infinite case." It is therefore anticipated that the neutral curve of the "finite case" will approach that of the "infinite case" if the aspect ratio is sufficiently large, or if the Prandtl number is sufficiently

-62w uN i 0 rl CL 1.8 1.7 1.6 3.2 3.1 3.0 0 I 2 3 4 5 (a) x L\n I 0 x UZ LT\ I' a0'ti I 0 (c) c Figure 11. Convergence Study for Finite H/d.

9 8 7 Hd = 50 6 x 0 Q: w5 3 0::,0 2 10 I I I. I. I I 0 2 3 4 5 WAVE NUMBER, a Figure 12. Neutral Stability Curves for Finite H/d

-64 small. It was indeed found that the neutral curve for Pr = 0.1, H/d = 20 essentially coincides with that generated for a slot of infinite height. Figure 13 depicts neutral curves which were generated for Pr = 1 at various aspect ratios. These curves lie in a range of the Rayleigh number in which the base flow model is not valid. They are therefore only tentative results which are not expected to be quantitatively valid. If this figure is qualitatively valid, however, the neutral curve for a Prandtl number of unity is seen to be in the form of a closed ellipselike a curve which becomes smaller as the aspect ratio is decreased. Eventually the neutral curve degenerates to a point and disappears. No roots were found for an aspect ratio of 10 for a Prandtl number of unity. It is interesting to note that Eckert and Carlson(4) observed low-frequency "turbulent fluctuations" in the core region of the natural-convective flow of air in a vertical slot. These fluctuations began occurring at a Grashof number of approximately 10 when the aspect ratio was 20 but were not observed until the Grashof number was an order of magnitude greater when the aspect ratio was 10o The stability of natural convection in a slot of infinite height was found to be independent of the Prandtl number of the fluid. The physical interpretation of this result was that the instability was hydrodynamic in origin, i.e. not effected by buoyancy forces. In the present case, however, the base velocity profile is a function of the Rayleigh number. Hence the disturbance momentum equation contains the Rayleigh number as an implicit parameter as well as the Grashof number as an explicit one. For this reason the role of the buoyancy forces cannot

-65 5 H\ i H/d =20 50 oD _T cr for Finite:D LLW 2 I (Z 0 r' 0 2 3 4 5 WAVE NUMBER, a Figure 13. Tentative Neutral Stability CurveC for Finite H/d, Pr 1. I.

-66 be determined by a simple parametric study of the states of neutral stability. In order to determine the importance of the buoyancy force, the neutral curve was generated for the "finite case" base flow using the Orr-Sommerfeld equation alone. For very small Prandtl numbers this formulation yielded a neutral curve essentially coincident with that of the "infinite case". On the other hand no states of neutral stability were indicated in the vicinity of those shown in Figure 12. It is therefore concluded that the buoyancy forces do effect the stability of natural convection in a slot of finite vertical height unless the Prandtl number is sufficiently small that the base flow approaches that of the "infinite case," iOe. when the flow is in the conduction regime. E. Evaluation of the Stream Patterns Having found the eigenvalues of the differential system (87) through (89) in order to determine the states of neutral stability, it is desirable to find the eigenfunctions in order to determine the form of the cellular disturbance motion. The disturbance pattern is completely specified by the eigenfunctions since the motion is two-dimensional. The disturbance can then be superimposed upon the base flow and the resulting total flow can be displayed graphically and compared with experimental flow visualizations. Once the eigenvalues have been determined, the algebraic equations generated by Galerkin's method can be solved, within a constant multiple, for the coefficients in the series representations of C and T o Consider the determinant in Equation (121).. If the sign of each element in the first column is changed and if any single row is then deleted, the matrix of the resulting determinant is the augmented matrix

-67 of a set of simultaneous algebraic equations. The solution of this set of equations yields the coefficients in series (119) and (120) relative to a1 1, thus determining an approximate expression for the disturbance. Such solutions were obtained by a standard numerical procedure available as a system subroutine at the University of Michigan Computing Centero The constancy of these solutions when various rows were deleted provided an indication of the (algebraic) accuracy of the eigenvalues. Since only the real part of the expression - P(- ) e (126) has physical significance, the disturbance stream function is represented by I/ = 2 (a C y() - b, $ -(y) A x) (127) 4=/ The disturbance stream functions were evaluated in this manner for neutrally-stable states in the neighborhood of the critical state for both the infinite-and finite slot problems. The disturbance stream pattern for the "infinite case" is displayed in Figure 14. The analgous stream pattern for the finite case (Ra = 3.12 x 105, H/d = 20, Pr = 1000) is displayed in Figure 15o A feature which distinguishes these cellular patterns from those such as B'enard or Taylor cells is the fact that they are tilted with respect to the planes of the boundaries. This is a result of the fact that the present problem has a base flow in the plane of the disturbances whereas the classical cellular instabilities do not. This behavior could not have been described by real eigenfunctions.

-68 Figure 14. Disturbance Stream Pattern, H'd = o

-69 Figure 15. Disturbance Stream Pattern, H/d = 20.

-7 (0 Another consequence of the base flow is apparent when Figure 14 and 15 are compared. The base velocity profile of the "finite case" is concentrated near the boundaries in comparison with that of the "infinite case." This is reflected in the distortion of the cell pattern in Figure 15 relative to the pattern in Figure 14. The stream patterns of Figure 14 and 15 represent just the disturbance. What is observed in an experimental apparatus is, of course, the superposition of the base and disturbance flows. The prediction of mode selection and amplification to finite magnitudes is outside the scope of linearized theory; however, it is known that the experimentallyobserved cellular motions of the Benard and Taylor problems are essentially identical to the critical modes predicted by linearized theory if the parameter governing the flow does not exceed its critical value by a large amount. Hence it is probable that a superposition of the cell patterns of Figure 14 and 15 and their respective base flows will produce a strewaiu pattern which qualitatively represents the actual total flow. Figure 16 represents such a total flow pattern for the "infinite case." T'his figure was obtained by normalizing both the base and the disturbance stream functions to a maximum value of unity and superimposing them as = ) + ( T (128) where E = 0.1. Figure 17 represents the analogous total flow for the "finite case" (RA = 3.12 x 105, H/d = 20, Pr = 1000). Figure 18 represents the total flow with 6 = 0.5. Although this value of E is too large for the small perturbation theory to be valid, this figure is of interest in its display of small cells of the

-71 Figure 16. Total Stream Pattern, H/d = C, & — 0.1.

- 72 - Figure 17. Total Stream Pattern, H/d = 20, f = 0.1.

-73 Figure 18. Total Stream Pattern, H/d = 20, 6= 0.5.

-74opposite rotation in the shear layers between the large cells. Such cells were observed by Elder at large Rayleigh numbers. (6) Birikh presents a stream pattern which corresponds to that of Figure 14. His pattern exhibits two small cells located inside the large cell, thus indicating more oscillation in the eigenfunctions that that found in this research.

CHAPTER III EXPERIMENTAL INVESTIGATION A. Apparatus An apparatus was constructed in which the natural onset of instability was detected by flow visualization. In order to verify the analytical results presented in Chapter II, the system was designed such that flows in the following parameter ranges could be examined with facility: (1) H/d = 20, 105 < Ra < 10, Pr > 1 (2) H/d = 50, 105 Ra < 106, Pr ~ 1 (3) Ra/(H/d) < 500 and Ra < 2xl04, 103 < Gr < 104 The first two parameter ranges correspond to two of the neutralstability curves presented for the"finite case" in Figure 12. The third parameter range is intended to produce a flow in the conduction regime, thus closely approximating the "infinite case", for which the neutral-stability curve is presented in Figure 9. The following basic criteria were specified for the system: (1) The width and aspect ratio of the slot are to be variable in accordance with these parameter ranges. (2) The depth of the test section is to be significantly larger than its width to ensure two-dimensionality of the flow. (3) One wall temperature of the slot is to be constant, but the other is to be variable. -75

-76 (4) Temperature is to be uniform over each side wall. (5) The test fluid is to be such that the LT required is sufficiently small to avoid drastic property changes, yet sufficiently large that extreme accuracy is not required in temperature measurement. Figure 19 is a schematic diagram of the apparatus. One wall of the test section is cooled by a flow of water from a municipal line. The other wall is heated with a closed distilled-water loop. Included in this loop is a 15 gallon stainless steel tank into which is inserted an immersion heater (Chromalox MTS 225A) whose energy dissipation is controlled by a Variac. The tank also contains a copper cooling coil through which cold water can be circulated. An electric stirrer is provided for agitation of the water in the tank. The distilled water is circulated through the loop by a Chempump model CFH seal-less centrifugal pump. The water leaving the cold side of the test section is circulated through a cooling jacket surrounding the pump in order to reduce heat transfer from the pump motor to the distilled water. The test section is schematically depicted in Figure 20. The slot in which natural convection is observed is formed by a rectangular plexiglass frame bounded by flat aluminum side walls of dimensions 25 x 7 x 0.25 inches. Each aluminum wall also forms one side of a chamber through which water, used to heat or cool the plate, flows. Staggered flow baffles promote uniformity of the wall temperature

TO DRAIN TANK TEST SECTION IMMERSION HEATER COOLING JACKET TO DRAIN Figure 19. Diagram of Experimental Apparatus.

-78 COOLED SIDE WALL HEATED SIDE WA HOT WATER COLD WATE R OUT PLEXIGLASS WATER IN FRAME Figure 20. Exploded View of Test Section.

-79 Plexiglass frames of nominal widths 1.25 and 0.5 inches were made in order to study the "finite case" for aspect ratios 20 and 50 respectively. One of width 0.75 inch was made for use in the "infinite case" experiment. A rectangular gasket cut from thin rubber was placed between the plexiglass and the aluminum plates in order to prevent leakage of the test fluid. The assembly was drawn together by means of threaded rods. Provision was made for the introduction and removal of the test fluid by simple gravity feed. A vertical plane in the fluid was illuminated by passing light through a slot at the top of the test section, as shown in Figure 21. This light was reflected from tracer particles suspended in the fluid. The motion of these particles could then be observed or photographed from the front of the test section. A 500-watt slide projector was used as the light source. The slit was approximately 1/8 inch wide. The beam of light diverged to a thickness of about 3/4 inch at the bottom of the slot. Temperature measurements were made with sheathed copper-constantan thermocouples imbedded in the aluminum plates. The sheaths were passed through the back of the flow chambers and were sealed there by Conax glands. The thermocouples on the heated side of the test section were located as shown in Figure 21. Only three thermocouples were used on the cooled side. They were placed along the vertical centerline of the aluminum plate. The signals from the thermocouples could be read on either a precision hand potentiometer (Leeds and Northrup Model 8662) or a printing recorder (Leeds and Northrup Speedomax Model S).

-80 LIGHT SLIT TEST SECTION PLEXIG - PLEXIGLASS THROUGH WHICH FLOW IS OBSERVED THERMOCOUPLES Figure 21. Test Section.

-81 Figure 22 is a photograph of the apparatus. Silicone oils were used as test fluids in the experiments which correspond to the "finite case." This selection was based on the following desirable properties of these fluids: (1) They are available with kinematic viscosities ranging from 0.65 cs to 10,000 cs with other properties essentially constant. (2) They are Newtonian at low shear rates. (3) Their variation of properties with temperature is small compared to many fluids. (4) Their light transmissivity is excellent. Dow Corning silicone fluid DC 200/100 was utilized in the tests at an aspect ratio of 20, and DC 200/5 was utilized in the tests at an aspect ratio of 50. Air was used for the test corresponding to the "infinite case." B, Procedure The kinematic viscosity of each silicone fluid was measured at 100~F and found to agree well with the value reported by the manufacturer. A small quantity of aluminum particles (Alcoa Aluminum pigment no. 10005) was suspended in the fluid. All cold water lines in the system were opened for at least one hour prior to each test, in order to attain a steady line temperature. The silicone fluid, with suspended aluminum particles, was introduced into the test section, and the apparatus was brought to a steady operating condition with cold water circulating through the cooling coil in the distilled water loop. This did not result in a quiescent test fluid because the distilled water was heated a few degrees by the energy dissipated

-82 Figure 22. Photograph of Experimental Apparatus.

-83 by the pump motor. Next the temperature of the hot wall of the test section was very slowly increased in order to attain quasi-steady states at various wall temperature differences. This was achieved by reducing the flow in the cooling coil and by activating the immersion heater. A typical rate of wall-temperature rise was about 0.15~F/min. At each quasi-steady state for which data was obtained, the outputs of the six thermocouples located on the vertical centerlines of the plates were measured with the precision potentiometer. Then all twelve thermocouple outputs were recorded continuously by the printing recorder while one or two time-lapse photographs were takeno The six centerline thermocouple outputs were then remeasured with the precision potentiometer. The variation of wall temperature with both time and thermocouple location was less than about 0.5~F during the time required to obtain this data. All photographs were taken on Polaroid Land type 47 film (ASA speed 3000) using a Graphic camera. Bending of the light rays due to the change with temperature of the index of refraction of the silicone fluids made it difficult to obtain a uniformly-illuminated plane. This difficulty was alleviated by placing a lens system between the slit and the test section which caused the light rays to diverge, thereby compensating for the refraction. Because of the. nature of the reflection of light from the Aluminum particles, a polarizing filter in front of the camera lens was found to increase the constrast between the streaklines and the background.

-84 The procedure followed to obtain data with air as a test fluid was identical, with the exception of the visualization technique. In this case a small amount of cigarette smoke was introduced into the top of the test section. Portions of this smoke were entrained by the turning boundary layer and stretched into long filaments circuiting around the slot along the streamlines. Photographable stream patterns could be obtained in this manner for four or five minutes before diffusion of the smoke particles caused too much clouding. C. Results and Comparison with Analysis The Rayleigh (or Grashof) number was computed for each state at which a photograph had been taken. These computations were made in terms of properties evaluated at the mean of the two wall temperatures. Each photograph was then examined and the flow was catagorized as being either parallel (in the central portion of the slot), cellular, or turbulent. The results of the examination of the data are represented in Figure 23. In this figure each of the flow classifications is indicated by a different symbol. Each horizontal row of data points represents one run, The critical Rayleigh (or Grashof) numbers can be estimated from this figure. The critical wave numbers were estimated by directly measuring the size of the cells in the photographs. The tests which were carried out at an aspect ratio of 20 are to a large extent a verification of Elder's(l5) experiment with high Prandtl number fluids. The velocities were quite low in these tests and it was virtually impossible to discern the pattern of the flow by direct visual observation. The streak photographs, however, gave a clear representation of the pattern. Figure 24 shows three photographs taken at Rayleigh numbers near the critical value. In Figure 24a

H/d = 20 Pr = 900 00 0 00OO 00 0 00 A AA 1 1 1 -L — I I I A A AA A I I I I I I 0 I 2 3 4 5 6 7 8 9 10 Ra x 10'5 H/d = 50 Pr= 77 0 0 O l I III I 0 I 2 3 VV V V V V V V v v I I 6 7 4 5 6 7 8I I 8 9 10 Ro x 10-5 I cx \n I H/d=33 0 0 A 0 0 AA Pr= 0.71 0 0 A A 1 I I i.I I I I I 0 1 2 3 4 5 6 7 8 9 10 0 PARALLEL FLOW A CELLULAR FLOW V TURBULENT FLOW Gr x I03 Figure 23. Summary of Flow Regirre Data.

Ra = 3.5 x 105 Ra = 3.78. 105 Ra = 3.96 x 105 (a) (b) (c) Figure 24. Visualization of the Onset of Cellular Instability, H/d = 20.

-87 the flow is essentially parallel and hence stable. In Figure 24b the flow has become unstable as evidenced by the appearance of a cellular secondary flow in the neighborhood of the interface between the rising and falling columns of fluid. In Figure 24c the cellular flow is somewhat stronger. Figure 25 is a view of a well-developed cell. This pattern can be compared with that of Figure 17 which was obtained analytically. Figure 26 is a similar sequence of photographs for the tests which were carried out with air in order to approximate the "infinite case."t Figure 26c can be compared with Figure 16. The bright spots in these photographs are caused by small dust particles which adhered to the walls, In both of the above cases the laminar cells were quite steady. The cells often remained stationary for five minutes or more. There did appear to be a very small upward drift of the cells. The cells appeared to be very nearly two-dimensional, as could be discerned by moving the light slit forward and backwards A stationary instability was not observed in the test at an aspect ratio of 50. Instead, the flow became turbulent at a Rayleigh number of about 3.3 x 105. It was not discerned whether the turbulence originated in the growing boundary layers near the end regions, or if it actually originated in the central portion of the slot. Figure 27a is a streak photograph of the turbulent flow at Ra = 52 x 105 taken with an exposure interval of 10 sec. Although this motion appears to be chaotic, the 3-sec exposure shown in Figure 27b bears a striking resemblance to the laminar cellular flows discussed above. Similar patterns were observed in several short time

-88 Figure 25. Streak Photograph of a Cell, H/d = 20, Ra' 4.5 x 105.

Gr = 7150 (a) Gr = 8750 (b) Gr = 9500 (c) Figure 26. Visualization of the Onset of Cellular Instability, H/d = 33, Pr = 0.71.

-90 10-Sec exposure (a) Figure 27. Visi Mot: 3-Sec exposure (b) ualization of Turbulent ion, H/d = 50, Ra = 5.2 x 105

-91duration photographs. Hence it might be conjectured that the structure of the turbulence in this flow is closely related to the cellular instabilities which have been considered in this dissertation. The experimentally-determined critical states are indicated and compared with the analytical results in Table I and in Figure 28.

-92 TABLE I SUMMARY OF RESULTS /d- Pr Racr cr H/ad Pr Experiment Analysis Experiment Analysis 20 900 3.7 x 105 3.12 x 105 35 1.85 Turbulent 50 77 at 7.05 x 105 --- 1.6 3.4 x 105 00.71 8.7 x 103 7.88 x 103 2.74 2.65 (Conduction Regime)

-93 3 2 is b O I 0 0 5 0 I 2 3 4 5 WAVE NUMBER, a I l I - 4 - 1') I 0 X 3 T 2 CRITICAL STATE \S ~ (EXPERIMENT) o - _~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ / CRITICAL STATE (ANALYSIS) i I H/d= 20 Pr 1000 I I 0 I 2 3 4 5 WAVE NUMBER, a Comparison of Analytical and Experimental Critical States. Figure 28.

APPENDIX I INVISCID CONSIDERATIONS On the basis of simple physical arguments and existing experimental evidence, it was assumed in Chapter II that stationary modes of instability exist in antisymmetric flows of the type considered here. The following remarks concerning the inviscid stability analysis of these flows provide further insight concerning such modes. It was found above that the stability of fully-developed antisymmetric natural convection was not effected by the buoyancy forces. In fact the buoyancy-force term in the disturbance momentum equation becomes identically zero as the Prandtl number approaches zero. For this reason; the essence of the stability analysis of these flows can be gleaned from a study of the Orr-Sommerfeld equation alone. The Orr-Sommerfeld equation and appropriate boundary conditions for this problem are (D2_-G2)2 -i G2 {()-c)pD}- ) - (D ) } =0 (I-1) (-2 ) = D(i(2 )= 0. (I-2) Since 2Grt is large compared with unity, in the vicinity of the critical Grashof number, a first approximation to this equation is that obtained by letting GI ~-i0, i.e. the inviscid form of Equation (I-l). (3 - C) ( a^ -/ 2M) _= 0 (1-3) -94

-95 (p( -E) = 0~ (I-4) This equation has been studied in detail by Tollmien (44 for parallel flows of the symmetrical and boundary-layer types. The major difficulty in solving this equation is that introduced by its singularity at the point or points where the wave speed and base velocity are equal. By imploying a series solution of Equation (1-3) in conjunction with considerations of the complete Orr-Sommerfeld equation, Tollmien was able to circumvent this difficulty and prove a number of important theorems regarding the existence of certain types of solutions of Equation (1-3). The following arguments concerning stationary modes of instability are based on Tollmien's analysis. For states of neutral stability the wave speed, c, is a real quantity. Equation (I-3) is then singular at each point where U = C In the case of antisymmetric velocity profiles such as those considered in this dissertation there are three such singular points if c = o (Figure I-l). These singular points are denoted as yl, y2, and y3. The Wronskian of Equation (1-3) is given by r = cTTmi Gd - r d y~ n (I-5) which is proportional to the disturbance Reynold's stress. Here Qi-^ and p9i denote the real and imaginary parts of the solution respectively. By applying the method of Frobenius to Equation (1-3) Tollmien shows that the Wronskian increases suddenly by an amount I9,| -r T upon transition through a critical point, As, in the positive direction of

-96 Y2 Y3 Figure I.1. Antisymmetric Velocity Profile.

-97 y (defined here as the direction away from the wall in each half of the channel). Since the singular point Y2 is an inflection point of the velocity profile, this quantity is zero at y2. According to the theory of ordinary differential equations the Wronskian must be constant in a region which contains no singularities, so it must be zero in (Y1, Y2) and (Y2, Y3). Hence no jump in 7 can occur at the boundaries yl and y3. This condition is automatically satisfied in the present case because of the boundary conditions on (. Thus a necessary condition for the existence of a solution of Equation (1-3) with c = o is satisfied. Although solutions with c = o exist for symmetrical and boundary-layer flows, Tollmien shows that they exist only for a=O, i.ed for disturbances of infinite wave length. In the case of antisymmetric plane Couette flow, to which the above arguments also apply, -=0 so the only solution of Equation (1-3) with c = o is the trivial solution, =O 0 Body-force generated antisymmetric flows would appear to be unique among parallel flows in fulfilling all the above conditions. If c = o and- > 0 in L',, Sturm-Liouville theory states that only imaginary values of OY exist. (See for instance Reference 26). In the present case - < 0 throughout the domain and thus the existence of solutions of Equation (1-3) with c = o is possible, but not proved. (24) Stuart has applied similar reasoning to the somewhat analogous problem of the stability of flow near a rotating disk.

APPENDIX II EVALUATION OF INNER PRODUCTS In the application of Galerkin's method is was necessary to evaluate a large number of inner products of the form (hi1Ki1' 2 ) = /S?/,2 dy. (II-1) -I Here j2 was either the C or, function or the cosine or sine function and! was one of these functions or a derivative of one of them. These integrations were carried out by integrating by parts four times and recalling that IV = 4 C, 5/^ u'= 4 $ (11-2) ao;/= -^ _ 2 w = -xh,,// As a simple example, the evaluation of CCIy)/ ) is carried out in detail: = (c cj' l ) - (c j)+ (c S =-(C^ IS, ) IS)- C - 2(C IS ) ~= (CfCJ) + (Ca|h) )+ 3 (Cl| C-).*. (cjlys) = 4(A4-ge4)- C / 2S/) -98

-99 But, according to Reference 38, S(o,id') = (4 u \1 r, // 2 So, ( cj ^J = 8,- (A4 4_ 4-2 Cr /// I; / I) The only cases in which such a simple scheme does not suffice are those in which C2, Sin, Cco2, or 4A- appear. In these cases (38) the integration was carried out directly. Reid and Harris3 have tabulated such inner products for cases in which they involve just the G and S functions and their derivatives. The following tabulation of 56 non-zero inner products is presented in the hope that it may be of benefit to others concerned with related analyses. 1.( C"Icc,1) -1 30 (Cl< /c C, ) 45 (0 l C ) (cl I'l/S,~) 50 (c,,lelS, ) 60 (CImIYI m) 76 (C/m 21 I q ) - a 4- C (X)cC ()( — (-) ( ) A4 = (As - 4 ) (2 C/'L)C,/ f 2C// ( C/z()) = 2( -t) {)-+ C() ) C()() + C ( -) {$ (4) - 2 4 AII ( - / (4) = 2 (A4 - 4)1 { C(-C,() 54 -) Cow C 2fL) n (i4) -C ") S () - 4 (At -^ ) C ) ad ( )} a e nv in7M\

-100 8. (c l2la s ) = (A4 -A)-4 C ^) -8A4 (C'Il') ^,VI 2 ~wm /'Il I WI l 9g (Ce I/21"') t 24 (4 - C'() s' () = (A"l -M4)f {2 C'() (i)- 2C ()'(//r ) 8 4H (C l yl ^) - 24 4 -l-4 )- C ) / l l)} = (A4,- e 4 )^(C)I'I ) +36(CI 1") 10o (CJ 13I $^ ) /// -L -48(A-4 - )-4 C ())} = (A4 - ){ c () - (I) -4 j (i ) t/ (/) ^In e 4 ^-# m -4 " -r.m ii.o IYC^ 31 m) + 2 C' (~) " ( ) - 2 AI (C/' - 122 S) J-Ct(C l) +8/2A(C4K 2I) q) - 313 \4 (c v4) -48( -)'C/, iJ 12 o, (SI S/) = 13o (\ ) I) = 2" 4 4S/( 1 -() (__ c^ r+ )2 -3 2 -m Z a 2 l4. ( S 1 V ) 15 o I ( I 4) = 4,4 )-" (2 1-7"/ (2) - 2 " (-L ) Nn.M ~m 2e c5r3nl - rNH 2 ( ~ 60.JC )= 8 (C-A )-2 -() C; /) 17.0 (lC)=2(, - 4 )- L2 )) -, t/C. + Sr (2 )C) - 4 h4 ( 4-A4 )-1 4, ) I, )} 18- (S^/ lc^;) = 2 (C4- A) { ) - a ( e ( 4- /(i)u'(1)_4 4 _/-4 -'2 8]C' /

-101 + 24 ( -A4 )- (2-)c () 20. (sM.l qIC,) = ( -A { $ (,-)C$ "(4)-2S'"(f)C c"(-) + (8A4 (, | c") - 24 A (2-) S' (I)C C)} ~ a2- s -'C'c7 - -20. ely21CM14 82 (r L c) C 1 C } 25. (CsUJ) = 2(p'r, (A4 - ^^4)/ 5 (2 ) (+) 4 ( I- 4(-A4f/^ (l)c 26. (C I~4 C~ ) 4 - 280.- y 31 8 C4- ='- _ 23. (Ct^.l8mC^ ) = 4 Df)< n(S)2 1\ ^))r( 29. (^^Cl j=4 -^ A4)^)) 23. (CA-4 | 4 (2 => L + 29. (c^^l~j^^) = -4^ K^(e-l c"^ ^ ) f 30. C^ lel3) ^^) = 48C - ) )3 + 2XI( 2- )( )4 ) ( i) -3 (J2_)2) (2 )^()

-102 31. (~ "nj| yi ) = - 2 32. Jir,1 n 2_ 323. (,,;1 I} ) = - 2 33 (<2,| JC^J) =- 4X~(0^ _-~2)-2~(C [ +?-l) 34. C ( 1 )Y 3} t) )= 4 8 )<^ A ()fi /2 _-2)-3 pl + 2 2 X 2) _,2)-1 ], 35,r, (cn/i^) (E) SjIZ' 35 ( CN/ " 21+m) 2(A4 -( A'C,() -2 c$ (2()) -+A L(cfr IY1,V) + 24( 4< )-'' ("'1l $/) 36o (cI2 0) = (A a -^ )-'- 4 c) (2) ) + 24 ( - )- " / )} 370 (CJ,, ) (A4-^ {4Crl) (-C^ ^) +4 (At-4 -Y4C-(1 C)~ (/) } 39. C' )1y4 )) - (' - (c4(,c v i S,) -72 (C -l )- 2') -2 96) (C 1) - 48(A4' _4)-1 // )$ /^ )

-103 40[ (C "// e51) S)= )- 4 c ) ( + (Cr"' $NE (2 ) " 6\N m a 16 - Ctn, () fl /(2/ ~)) - 20 4 LCt 1 Sm )+ (C4.m I S 2)] -/ 240 ( 5) - 12o(CISm eC"I(C JAY21Sm 1 0 ( "',Y 240 ~~~- ch MYJI -Im 41, (Ce l ey J )m = (A4 -)-' { (c\20(Cy 41') 20(C | I /31s'") 61, ol~~~~~~l mN 42. (S"/ J = +240 (CJ,)/+2 ) 2+ 120(Cjy/SI)} (lZ 1- A: )-' {2'() cr (4)- 2-)C ()c'"(f)) - 4 [8 (" | VI Cm) + 24 (f -A ) "C (1n (/ fl)]} 43. (117Y2 C ) (4-A4 I )- {-2) f 8 ( a C) O"Y, 48A, 5r Ie j C,4 +24 (C,4-, )-1 "(~) "//( 4} (-A ) {4 (_" (c-' () - - /s ~) C"'(/) ) 44. (4S, 3 1 1 ) = -2 (So I y/ C. /) + 36 ( cS/ )] +2 - )~'(l C') ] + 3 ~ (4, I ) C, } + ~ ( -A e'I, I -/.' "'-~ 45. ( 1 IC 4 1 C ) = (, j ()- s') { )C () "- ") C$ ( -) +A4 6l6(S^l|cl) + 72(S Ic21c ) +96(yCSJ~cj] + 48 (-4g )C' (4)$ (2)}j = (g - / )4 ) - { S (^) C () - I 4 (S I G31C ~~~C~,4 L a 4 2 m a,,'/ C 46. (S I 41 C ) -72 (S I'2) c -48 (4-A4 )-' - 96} (^'I C) t (2) Ce (1) } ll~!2

-104 47. (4i 5cj) = -A4/)l) + 6 2 L - SF(1) C(L)) - 20 4 (4C)+G(S,34C - 24C0 ( IY"21c) - o2 c " c )} 48. (a C.51C,) = ( -A )-' {20( (S ey41 C ") + 120( I 31 C ) + 240 (Sh,^C/ ) + 120 (S iCm)} 49. Cjm Y)1 ) = a(A4 X{14 - C C(;)4+(#^) + 4-P,4(C, CJ c)J 50.o (Cot |lc d) (A4-St4C d/\ () 0a) Op2 [8 (crje}) cA) + /2 (cl )} 51. ^(Jyl) =2 (L - 2 L ^ ac )K0 -- ( +4X)<e(,|an,)} 52. (S I0,^Jm) = ( -,4-K4 {-1" 2 0 \) -20 (cutel y831 I)} 54~~ (ca~~, IY~~ln2 y ) fQ 1 /) ryJol~

0s n 0 ^^ - x- + l u ( e' o X II ^ + c I 1 J2 J.* L\ L\ LAr iT\.

BIBLIOGRAPHY 1. Ames, W. F. Nonlinear Partial Differential Equations in Engineering. Academic Press, 1965. 2. Batchelor, G. K. "Heat Transfer by Free Convection Across a Closed Cavity Between Vertical Boundaries at Different Temperatures." Quart. J. Appli Math., 12 (1954), 209-233. 3o Becker, M. The Principles and Applications of Variational Methods. Research Monograph Noo 27, M.I.T. Press, 1964. 4o Benard, H. "Les Tourbillons Cellulaires dans une Nappe Liquide Transportant do la Chaleur par Convection en Regime Permanent." Ann~ Chim. (Phys ), 23 (1901), 62-144o 5. Birikh, R. V., Gershuni, G. Z., and Zhukhovitskii, E. M. "On the Spectrum of Perturbations of Plane —Parallel Flows at Low Reynolds Numbers." P.M.M., 2, No. 1 (1965), 93-104. 6. Birikh, R. V. "On Small Perturbations of a Plane Parallel Flow with a Cubic Velocity Profile." P.M.M., 30, No. 2 (1967) 432-438. T, Bodoia, J. R. and Osterle, J. F. "The Development of Natural Convection Between Heated Vertical Plates." J. Heat Transfer, Trans. ASME, C84 (1962), 40-43. 8. Caddell, R. M. "An Experimental Study of Magnetohydrodynamic Flows Induced by Applied Electric and Magnetic Fields." Ph.D. Dissertation, The University of Michigan, Dept. of Mech. Eng., 1963. 9. Chandrasekhar, S. "On Characteristic Value Problems in High Order Differential Equations Which Arise in Studies on Hydrodynamic and Hydromagnetic Stability." American Math. Monthly, 1, No. 7 (1954), 32-45. 10. Chandrasekhar, S. Hydrodynamic and Hydromagnetic Stability. Oxford, 1961. 11. Collatz, L. The Numerical Treatment of Differential Equations. Springer —Verlag, 3rd Ed., 1960. 12. Damaskos, N. J. and Young, F. J. "The Stability of Thermally Induced Flow in the Coreless Induction Furnace." Int. J. Heat Mass Transfer, 8, (1965), 721-728. 13o Di Prima, R. C. "Application of the Galerkin Method to Problems in Hydrodynamic Stability." Quart. Appl. Math., 13 (1955). -106

-107 14. Eckert, E. R. G. and Carlson, W. O. "Natural Convection in an Air Layer Enclosed Between Two Vertical Plates with Different Temperatures." Int. J. Heat Mass Transfer, 2 (1961), 106-120. 15. Elder, J. W. "Laminar Free Convection in a Vertical Slot." J. Fluid Mech., 23 (1965), 77-98. 16. Elder, J. W. "Numerical Experiments with Free Convection in A Vertical Slot." J. Fluid Mech., 24 (1966), 823-843. 17. Elenbaas, W. "Heat Dissipation of Parallel Plates by Free Convection." Physica 9, No. 1 (1942), 1-28. 18 Finlayson, B. A. and Scriven, L. E. "The Method of Weighted Residuals and Its Relation to Certain Variational Principles for the Analysis of Transport Processes." Chem. Eng. Science, 20 (1965), 395-404. 19. Finlayson, B. A. and Scriven, L. E. "The Method of Weighted-Residuals —A Review." Appl. Mech. Reviews, 19, No. 9 (1966), 735-748. 20. Gershuni, G. Z. "Stability of Plane Convective Motion of a Liquid." Zhurnal Teknischeskoi fiziki, 23 (1953), 1838-1844. 21, Gershuni, G. Z. "On the Question of Stability of Plane Convective Motion of a Liquid." Zhurnal Teknischeskoi Fiziki, 25, No. 2 (1955), 351-357. 22. Gershuni, G. Z. and Zhukhovitskii, E. M. "Stability of the Stationary Convective Flow of an Electrically Conducting Liquid Between Parallel Vertical Plates in a Magnetic Field." Sov. Phys. JETP, 34, No. 3 (1958) 465-470. 23. Gill, A. E. "The Boundary-Layer Regime for Convection in a Rectangular Cavity." J. Fluid Mech., 26 (1966), 515-536. 24. Gregory, N., Stuart, J. T., and Walker, W. S. "On the Stability of Three-Dimensional Boundary Layers With Application to a Rotating Disk." Philo Trans., A, 248 (1955), 155-199. 25. Harris, D. L. and Reid, W. H. "On Orhtogonal Functions which Satisfy Four Boundary Conditions. I. Tables for Use in FourierType Expansions." Astrophys. J. Supp. Ser., 3, (1958), 429-447. 26, Ince, E. L. Ordinary Differential Equations, Dover, 1956. 27. Jakob, M. Heat Transfer, Vol. 1, John Wiley and Sons, 1949. 28. Kantorovich, L. V. and Krylov, V. I. Approximate Methods of Higher Analysis, trans. by C. D. Benster. Interscience, 3rd Ed., 1964.

-108 29. Lietzke, A. G. "Theoretical and Experimental Investigation of Heat Transfer by Laminar Natural Convection Between Parallel Plates." NACA TN - 3328 (1954). 30~ Lin, C. C. "On the Stability of Two-Dimensional Parallel Flows." Quart. Appl. Math., 3 (1946), 117-142, 218-234, 277-301. 31. Lin, C. C. The Theory of Hydrodynamic Stability. Cambridge, 1955. 32. Mordchelles-Regnier, G. and Kaplan, C. "Visualization of Natural Convection on a Plane Wall and in a Vertical Gap by Differential Interferometry. Transitional and Turbulent Regimes." 1963 Heat Transfer and Fluid Mech. Institute, 94-111. 33. Ostrach, S. "Laminar Natural Convection Flow and Heat Transfer of Fluids with and without Heat Sources in Channels with Constant Wall Temperatures." NACA TN-2863 (1952). 34o Ostrach, S. and Maslen, H. "Stability of Laminar Viscous Flows with a Body Force." International Developments in Heat Transfer, Pt. V, International Heat Transfer Conf., 1961. (Publ. ASME 1963). 35~ Ostrach, S. "Laminar Flows with Body Forces," in Theory of Laminar Flows, ed. by F. K. Moore. High Speed Aerodynamics and Jet Propulsion, Vol. IV, Princeton, 1964. 36~ Pellew, A. and Southwell, R. V. "On Maintained Convective Motion in a Fluid Heated From Below." Proc. Roy. Soc. A, 176 (1940), 312-343. 37~ Rayleigh, L. The Theory of Sound, Vol. I. Dover, 1945. 38o Reid, W. H. and Harris, D. L. "On Orthogonal Functions which Satisfy Four Boundary Conditions. II. Integrals for Use With Fourier-Type Expansions." Astrophys. J. Supp. Ser., 3 (1958), 448-452. 39. Roberts, P. H. "Characteristic Value Problems Posed by Differential Equations Arising in Hydrodynamics and Hydromagnetics." J. Math. Analo Applic., 1 (1960), 195-214. 40o Rudakov, R. N. "On Small Perturbations of Convective Motion Between Vertical Parallel Walls." P.M.M., 30, No. 2 (1967) 439-445. 41. Squire, H. B. "On the Stability for Three-Dimensional Disturbances of Viscous Fluid Flow Between Parallel Walls." Proc. Roy. Soc. A, 142 (1933), 621-628o 42. Stuart, J. T. "Hydrodynamic Stability," in Laminar Boundary Layers, edo by Lo Rosenhead. Oxford, 1963.

-10943. Taylor, G. I. "Stability of a Viscous Liquid Contained Between Two Rotating Cylinders." Phil. Trans., A, 223 (1923), 289-343. 44~ Tollmien, W. "Ein Allgemeines Kriterium der Instabilitat Laminarer Geschwindigkeitsverteilungen." Nachr. Ges. Wiss. Gottingen, 1, No. 5 (1935) 308-330. (Trans. as "General Instability Criterion of Laminar Velocity Distributions, " NACA TM 792. ) 45. Yih, C. S. Dynamics of Nonhomogeneous Fluids. Macmillan, 1965.

UNIVERSIOF MICHIGAN 3 9016 03527 1876