fir“ WW‘ I" ‘ “n I h' ‘ IQ l-( '1 L“ "‘ _- 1- 3‘ ‘ . " a WWW WWhQQmW ' QWNQ AQHW"QLJW$§753E Q“ "M"! [Q QIQIQKQQQIII 11091111“ I Q'QIWQQQWWWQ 11- 1' EM Q. 511$" TQQ'Q’Q m-QQIIQ’Q" I1. 2.1.; l‘ Q“Mmd;)3*;jlt?n1nu% . Q: Q wa‘I-Lki‘ . ' '1 QQ’QQQQQ‘ ‘ )"Q' Q5".QQ‘1I(‘QQ;IQ:;1L \: . "I m“ ”11m“ ‘W‘ II IQ'W 1‘ .'WQIEII; . ' “Wm -1 ‘IH'WQQ' 1W QW'I'QQQ " 21:14; I III ‘I'III ‘ I W; ' h" In .. . .' Q Q“ I ‘ IQ Q I ' '1'“! I1, ‘I "Q Q ' Q1!” .1 . _-r.-_;;.“:_J :QQQ Q Q Q» . Q * :5 1Q] . W ' W {1}, ' ' 4.1 i '23"?! "I 3.4% «t1 "12.1"" ‘ . . :IL’Q . .1 1- .:1 gig-65 gQfiQi. 5 V 4} 'QQ 1 kW " ‘Qf'ffi' W ‘Q’gUI'II AQQ‘J'Q Q" . . QM 1 ‘ 21.111” Jan-1% .4 Q ,l 11.89 . . b I 1{ 0b.“). 5.! .51“! . 1d,!” Q. I . 13.1. ' 1:21- ;: .u‘ 9 . Q . Q‘QQQQ. ' U‘ ‘ J ~ .' "Q" . Q 'JQV? g? —. '. . . 1 -_-_ .f - ., I —~.., _. . -‘-—v .A .‘ mum 1 . . . fix. 3': .'.: “' - -'— 7 .. . , .Q . , -1... ’t-D ‘ , _ .-.. _. . _ .‘.,“ .. . 1'. -. -. ‘5 115215;“? . 1. 1’1" IKE-Y; ‘1 1‘ '. “1' =1 r ----- :11. I M‘65 q __.... .. ,» .'m'. . v: u . .- 3-, ‘4' . > I;. 1... -..... ' J: . *l ::::' I '4 .I ,1"? '. I Z .‘5' ‘ : .IQH‘Qg‘QQI‘ 3Q}; 1‘5: 'fi'l‘ . 1" . ,. _ . - % .‘- iu' 5‘ In} .01.”, VIY.I V J .H V '1. ;. ..§;. .718 $111 -. .-. '1‘?- 1 .QF‘f: l‘ 33:. W ‘I"5{:'IQI’QS\I.}' 121%!“ 'l 7:31;; 15:: “(1% Q'QQ-Hvtrfikg' "13' 3:;11 IQQW 33:0." Q '1.'Q "2'3: Q 1‘1. 1. .‘ ' "'3: ' .lgll '. 'Q'ffig’: A ... ‘31 fi‘lt‘! I'M 419311“szde IQQ’I IQ .. ,s'i‘w 'Z'QJ'III-Q IQ" IIIIm ,1“... - Big".- ~': . I II m 1'" ! ..IH. IIQH:.~=.'.I:IIQM'13111 , 1. a...“ ,"WIQQQ . Q WWWQMW' W'WVHngmWflWQqQQQ thQQQBQ’Q ; ' Q1111] QQfll'Q‘Q' Q11 ' ‘QQQQ "QQ"I QQQ ~I' ""1‘Q""WQQQQ;'WQ' .'IWQ ‘- "' "Q Q : . ' .I.u~... .3 I‘ .'II QIIIQUIIIIQ II-"IQL'QIIIIQI ’~QQQQJQ'IQHQLIQL‘l-IIFQQHII. IQ! ' PP THESIS LIBS-EAR 0 q, .. m , U .1' _ ' '. ‘ i ‘0‘ -s‘ , .‘ .17 "“) El". 1' This is to certify that the thesis entitled "An Analytical Study of the Convection Stage of a Cryomicroscope System" presented by Hasmukhbhai Patel has been accepted towards fulfillment of the requirements for Master's dame“! Mechanical Engineering Q71“ Qf [a £79762?! L1 (\/ Majéy‘rofessor Date U/H (/8/ 0-7639 MSU LIBRARIES “ ‘- RETURNING MATERIALS: Place in book drop to remove this checkout from your record. FINES will be charged if book is returned after the date stamped below. sail: 3: 20011 AN ANALYTICAL STUDY OF THE CONVECTION STAGE OF A CRYOMICROSCOPE SYSTEM BY HASMUKHBHAI K. PATEL A THESIS Submitted to Michigan State University in partial fulfillment of the requirements for the degree of MASTER OF SCIENCE Department of Mechanical Engineering 1981 ABSTRACT AN ANALYTICAL STUDY OF THE CONVECTION STAGE OF A CRYOMICROSCOPE SYSTEM BY Hasmukhbhai K. Patel The present study describes a numerical analysis of the hydrodynamic and thermal characteristics of a heat transfer system designed to produce thermal control of a small sample placed on a light microscope. A simplified model of an actual system has been developed. The model consists of a laminar flow in a rec- tangular duct of aspect ratio 0.3. The velocity field in the hydrodynamically fully deve10ped region is solved. This solution is used in the thermally developing region where the duct experiences uniform heat flux at the top (heater) surface and is insulated on the other three sides. The numerical solutions have been obtained using the finite difference technique with the Gauss-Siedel method (using successive over relaxation) and the Alter- native Direction Implicit Method. The temperature field on the surface of the heater (the site of the biological specimen) is presented and it reveals unacceptably severe temperature gradients in both the axial and lateral directions (280°C/CM and 320°C/CM Hasmukhbhai K. Patel respectively). These results are in qualitative agreement with experimental data obtained from the actual system. Results for the velocity field indicate that the solution yields accurate results based on mean velocity criteria (Um n) but may be inaccurate by as much as 20% ea based on local velocity gradient criteria (Cf). Thus the results for the temperature field on the heater surface may be in error by a comparable magnitude even though an energy balance based on mean temperatures is accurate to less than 1%. Modifications in the computer program required for model improvement and further recommendations for thermal design studies are discussed. Dedicated to my dear Brother and my friend Chunibhai and Hirokiyoda without whose support and encouragement this work would never have been completed. ii ACKNOWLEDGMENT The author is deeply indebted to Professor J. McGrath for his efforts, guidance, encouragement and patience. The author avails this opportunity to thank all the faculty members of the Mechanical Engineering Depart- ment, Michigan State University, East Lansing. The author wishes to express appreciation to Professor J. V. Beck, Professor J. F. Foss and the Honor- able Chairman Professor John Brighton for their timely as- sistance on many occasions in spite of their heavy engage- ments. Also the author wishes to express appreciation to Mr. Kim Hull for his technical assistance and help, and Mr. Narendra Dahotre for his help and encouragement. iii TABLE OF CONTENTS LIST OF TABLES . . . . . . . . . . Vi LIST OF FIGURES . . . . . . . . . . Vii NOMENCLATURE . . . . . . . . . .Viii Chapter 1. INTRODUCTION . . . . . . . . 1 1.1 Background . . . . . . . 1.2 The Actual Cryomicroscope System . . 2 1.3 Fluid Mechanics and Heat Transfer Aspects of the Actual System Compared to the Model System . . . . . 1.4 Previous Solutions Available . . . 1.5 Specification of Model Duct and Heater Based Upon a Proposed Experimental System . . . . . . . . 1.6 The Model Duct/Heater System . . . 2. HYDRODYNAMIC CONSIDERATIONS . . . . 14 2.1 Non-Dimensionalized and Finite Differ— ence Momentum Equations With Specific Boundary Conditions . . . . . 14 2.2 Hydrodynamic Solution of the Problem . 17 2.3 Conclusions and Comments on the Hydro— dynamic Solution . . . . . . 27 3. THERMAL CONSIDERATIONS . . . . . 33 3.1 Dimensional Approach for the Thermal Solution of the Problem . . . . 33 3.2 Thermal Problem Solution by the Finite Difference Approach . . . . . 40 3.3 Conclusions and Comments on the Thermal Solution . . . . . . . . 53 iv Chapter 4. DISCUSSION AND CONCLUSIONS . 4.1 Discussion and Conclusions 4.2 Suggestions for Future Work APPENDIX A: Methods of Iteration APPENDIX B: ADI Method . . APPENDIX C: Tridiagonal Matrix Algorithm APPENDIX D: Program ‘HEAT' . APPENDIX E: Convection Program APPENDIX F: Non-Dimensionalization of Momentum Equation for Fully Developed Flow . REFERENCES . . . . . 63 63 64 67 74 77 79 80 86 90 LIST OF TABLES Table . . . . + + + 1 Non-DimenSIonal Veloc1t1es U (y ,z ) as a Function of Non-Dimensional Position . 2 Velocities U(y,z) in (m/sec) . . . . 3 Temperatures at Different Locations in the Rectangular Duct at Various Distances from the Leading Edge of the Heater . . . . . 4 Temperature of the Heater Centerline as a Function of Distance from the Leading Edge of the Heater . . . . . . . . 5 Results by Cess and Shaffer Method . vi 50 61 Figure 10 11 12 13 14 15 16 LIST OF FIGURES Schematic Diagram of Cryomicroscopic System Temperature Gradients in Original Cryomicro— scopic system . . . . . . . Developing and Developed (Hydrodynamically) Regions . . . . . . . . Reynolds Number and Mass Flow Rate as a Function of Total Pressure Drop . . . Finite Difference Control Volume . . . Non-Dimensional Hydrodynamic Boundary Condi- tions on Rectangle Cross Section Quarter . Non-Dimensional Velocity and Hydrodynamic Boundary Conditions on One Quarter of the Rectangular Cross Section . . . . Finite Difference Approximation of he Non- Dimensional Velocity Variation in Y DIRN at 2+ = 3. 333 . . . . . . . Finite Difference Approximation of the Vari- ation of Non-Dimensional Velocity in Z+ Direction at Y+ = 1.0 . . . . . Thermal Boundary Conditions on Surfaces of Rectangle Duct . . . . . . . Finite Difference Cross Section for Thermal Solution . . . . . . . . . Velocity Assumption in Calculating Velocity on Boundary for Thermal Solution . . . Non-Dimensional Temperature Parameter R(I,J) Used in the A D I Method for Solving the Energy Equation . . . . . . . The Product of Velocity and Temperature (UT) at Different Nodes in the Duct at the Trailing Edge of the Window (x = 5 CM) . . . Predicted Variation of Temperature at Window Heater Surface in 'Z' Direction at a Given 'X' Variation of Thermocouple Temperature with Distance from Leading Edge . . . . vii 10 12 16 19 23 30 31 34 38 41 43 52 54 55 P NOMENCLATURE matrix containing coefficients of the unknowns area of cross section of the duct half the height of the duct coefficient of matrix A at 1th row and jth column matrix containing constant of equations half the width of the duct surfaces specific heat of compressed air (kJ/kgok) Cf fully developed L“ l." WU- th th Ap Apl Ap2 friction factor in fully developed region of the duct friction factor in developing region of the duct hydraulic diameter of the duct specific gravity of compressed air row number column number thermal conductivity of compressed air (w/mok) length of the duct (m) thermal length required for full development of temperature profile (m) dimensionless thermal length required for full development of temperature profile (m) mass flow rate of the fluid (kg/sec) total pressure drop in the duct pressure drop intflmadeveloping region of the duct pressure drop:hathe fully developed region of the duct viii q" number of parts made from the height of the duct for finite difference method heat flux input (w/mz) total heat added to the fluid number of parts made from the width of the duct for finite difference method constant dependent on velocity, molecular thermal diffusivity at the position (I,J) temperature (0C) mean bulk fluid temperature (0C) bulk temperature of fluid at the exit of the window region ( C) bulk temperature of fluid at the exit of the window region (0C) velocity at given location (m/sec) non-dimensional velocity mean velocity of the compressed air in the duct (m/sec) constant depends on pressure drop for non- dimensionalization of velocity (m/sec) maximum velocity at the center of the duct (m/sec) non—dimensional mean velocity Optimum overrelaxation factor in S. O. R. method the hydrodynamic length required for full develop- ment of the velocity profile matrix containing unknowns x coordinate of the cross section non-dimensional x coordinate of the cross section y coordinate of the cross section non-dimensional y coordinate of the cross section 2 coordinate of the duct length non-dimensional z coordinate of the cross section ix 0(8) Greek Symbols molecular thermal diffusivity (mz/sec) fluid density (kg/m3) dynamic viscosity coefficient (pa ° sec/m2 sec) aspect ratio of the duct specific radius of the matrix of coefficients CHAPTER I INTRODUCTION 1 . 1 BACKGROUND The cryomicroscope is a research tool which allows a researcher to directly observe freezing and thawing pro- cesses. The better systems of this type are computer- controlled and allow precise and virtually independent control over temperature and its rate of change. The desirable characteristics of such a system include the capability to subject a specimen to a uniform temperature process T(x,y,z,t) = T(t) over a rather wide dynamic range of temperature rates of change. These char— acteristics make it possible to identify particular events of interest with specific temperatures and/or temperature histories. For small samples, this has been possible and a single temperature measurement provides a reasonable basis for feedback control. The present work considers a heat transfer system which has been used to study freezing damage in the micro- circulation of a hamster cheek pouch--a macroscopic speci- men that was viewed through the cryomicroscope. There is experimental evidence showing that temperature gradients of approximately lOOOC/cm may exist in the immediate vicinity of the specimen for this type of system. The research presented here is a first step at im- proving the thermal design of such systems. A simplified model of the actual system is developed and numerical solu- tions of the velocity and temperature fields are presented. The computer model developed here and modifications of it 2 can be used to optimize the thermal design of such sys- tems in the sense of reducing temperature gradients at the site of the specimen. The remainder of the introduction describes the actual system more completely and presents the simplifi- cations made in modelling the system. 1.2 THE ACTUAL CRYOMICROSCOPE SYSTEM The cryomicroscope, developed by Diller and Cra- valho at the M.I.T. Cryogenic Engineering Laboratory (1), was the primary investigative tool used by Thomas Hrjcaj (2) for his experimental investigation of the freezing damage in the microcirculation of a hamster cheek pouch. The cryomicroscope consists of a Zeiss Universal light micro- SCOPe (see Figure 1) that incorporates a specially designed convection stage. This convection stage is placed in the optical path of the microscope between the objective and a long working distance (7 mm) condenser. The convection stage consists of a closed rectangular chan- nel (0.118)(0.375)<2.69 inch) fitted with circular quartz windows (0.9" diameter, 0.010" thick) on the top and on the bottom of the channel. These windows are located at that point along the channel which intersects the micro- scope objective. In this case, biological specimens are placed on the top heater window for cooling. The cheek- pouch of a sedated hamster could be everted and mounted flatly over the top window. For the convection stage used by Hrjcaj (2) in his research, pre-cooled gaseous SMBmNm UHmOUmomUHEOMmU m0 ZflmwdHO UHBflmeUm .H mmeHm momaom emquulllwxnmmmmwwv e Kmogm zoEom>zoo w\\\\\\L _\\\\\\\_ af/medm BZflmmemmmm I mmwz¢muxm Bdmm u\\\\\\N\h _\\\\\\\\L - \.\\\\\t DMBHmOmmn m0m¢> m>HBUMUmOl|lV g mamaooozmmme BHZD domBZOU 00A42< nitrogen was convected through the closed channel. The gas stream supplied the necessary cooling capacity to freeze the cheek pouch. The top quartz window of this system was manufac— tured with a thin transparent resistant film on its bottom surface. This film was therefore between the biological sample and the refrigerant fluid so that dissipation of electrical energy in the film provided thermal energy when desired to offset cooling by the refrigerant. An analog thermal control system measured one tem- perature in the specimen and used this single measurement to control "the" specimen temperature by varying the heat generation in the film heater. As shown in Figure 2, the temperature gradients in the original system could be severe. A quantitatively accurate model of this system is desired. 1.3 FLUID MECHANICS AND HEAT TRANSFER ASPECTS OF THE ACTUAL SYSTEM COMPARED TO THE MODEL SYSTEM The actual cryomicroscope heat transfer system used by Hrjcaj is characterized by several important complexi- ties that were beyond the scope of this initial study. The following factors are included in the actual system: 1) flow separation at the entrance of the rec- tangular duct due to sudden expansion from an upstream circular tube 2) turbulent flow 3) both hydrodynamically and thermally developing flow EWEmNm UHmOOmomUHZONmU AflZHmeo ZH wBZmHQAme mmbeflmmnHEMB .N mmDGHm TZHL MUZ/NBmHD Nm\wH Nm\mH Nm\m Nm\m o P — — _ u EMMEm BEWOHMhflM I n N _ -3- \ y \ If 5: \ OO 0 .. 1 95.383 \ eve/p \\ 3 momg \ \e 90 \ . m \ oo \ .vv. 1 \ WOO \ «(0 ... fl \ \ «we \ ,ve . m \ CA .9 o .... okay .. .w \ oo O \\900 11 ( \ \ .v? \ WW. I... O m \ \ e . e L. 6 4) heat dissipation over circular boundary re- gions of the rectangular duct system 5) ill-defined natural convection boundary con- ditions over the remaining three duct surfaces Since the present research focuses primarily on the thermal aspects of the problem, a simplified flow field was employed. Hence the flow is taken to be: (1) laminar; (2) fully developed hydrodynamically; (3) steady flow. In addition, the thermal model primarily considers the development of the thermal boundary layer in the rec- tangular channel in the heater window region in order to predict quantitatively the thermal gradients on this sur- face. The following assumptions are made for the heat transfer calculations: 1) uniform heat flux from the heater window into the refrigerant fluid 2) adiabatic walls on the other three duct sur— faces 3) constant fluid properties. 1.4 PREVIOUS SOLUTIONS AVAILABLE Hornbeck (12) provides a solution to the heat transfer problem described above. However, this solution does not in- clude the non-dimensional axial length, aspect ratio, or Reynolds number of interest here. Details of the Hornbeck technique are given later. 7 1.5 SPECIFICATION OF MODEL DUCT AND HEATER BASED UPON A PROPOSED EXPERIMENTAL SYSTEM Since the model to be analyzed is not completely comparable to the actual system, the conditions to be mod— elled are defined in terms of a proposed experimental system which could be used to verify the numerical results ob- tained from the current simplified model. The conceptual experimental system serves the useful purpose of providing a specific illustration of the fluid mechanics and heat transfer concepts involved and does so on a scale that would be easily realizable in an experimental laboratory. The model fluid assumed here is compressed air since it would be readily available in most laboratories. 1.6 THE MODEL DUCT/HEATER SYSTEM The model duct is taken to have the same cross- sectional dimensions as the original Hrycaj system (1.0 cm width by 0.3 cm height). The model heating length is arbitrarily chosen to be 5.0 cm compared to an actual length of 2.9 cm. This increased length allows predic- tion of temperatures for longer heaters of interest. The simplified model assumes a completely devel- oped hydrodynamic boundary layer in the window region with the laminar flow assumption. In order to assure the fully developed condition, the behavior of the flow was examined for a range of Reynolds numbers less than 2300, and a developing length was defined. This criteria places the model heater far enough downstream of the duct 8 entrance such that full development of the hydrodynamic boundary layer occurs. It has been assumed that a uniform velocity distribution would be present at the entrance of the duct. The Reynolds number will be based upon the hydraul- ic diameter since this approach is found to describe flow through non-circular ducts very well (8,9). The hydraulic diameter is defined as: 4 x cross-sectional area Dh = perimeter of the cross-section In the present case for a duct of dimensions (0.3 cm x 1.0 cm): D = 0.4615 cm h The criteria for the criticalhydrodynamicdevelopment inside the duct is given by Kays and Crawford (8) in terms of the hydraulic diameter approach: NI?) DID X— D; _ Since the fluid in the duct is assumed to be air, the following prOperty values have been used: 1.1766 kg/m3 9 = u = 18.53 x 10'6 kg/m-sec K = 26.14 x 10-3 w/mok cp= 1.005 kJ/kgok These values are taken at 23°C and are assumed constant. A Reynolds number of Re = 1827 is arbitrarily DH chosen to have laminar flow and the corresponding develop- ment length is X = 0.4215 m. Thus the model heater will be 9 placed downstream at X = 0.475 m and the solution will be for ReDH = 1827 with a total duct length of 1.0 m (see Figure 3). Although the total pressure drop in the Hrycaj system was estimated as 10-20 psig, the pressure drOp required to maintain laminar flow will be considerably less than this. The pressure drop mechanism and expressions quantifying the pressure drop are different in the undeveloped and developed regions. Knowing the criti- cal development length, we consider the two pressure drops separately and calculate the total pressure drop for the entire tube. Let AP = AP + AP 1 2 = total pressure drop where: AP 1 - pressure drop in the developing region AP2 = pressure drop in the fully developed region Shah and London (9) present graphical data for C 'Re as a function of the duct aspect ratio and: f DH ._17.51 _ FEE; (1.1) Cf fully developed where the aspect ratio is 0.3 here. Shah and London (9) discuss the developing nature of the flow in which the friction factor is a function of axial position in the tube. The same authors point out that the more useful quantity is the apparent average friction factor over the development length C This f app' quantity is a function of aspect ratio where the aspect ratio is defined as: Aspect ratio = (%%) 10 mZOmem ANAQmD 02¢ wszOAm>MQ .m mmDUHm 11 (2a = height of channel; 2b = width of channel) For the present case, the aspect ratio has a value: Aspect ratio == 0.3 For this value, the apparent friction factor is: Cf app - ReDH = 22.25 (9) (1.2) Therefore: sz }( AP1 = 4 Cf app (76;) 51} (1'3) 2 = 0V L-x AP2 4 Cf fully developed (29c)(l%i) (1.4) where L = length of duct = 1.0 m Results of the calculation of AP for various Reynolds numbers are given in Figure 4. These values are considerably lower than the estimated maximum figures for the flow rate: 1.3 x 10-3 kg/sec, and pressure drop: 9.0 x 3 10 NT/mz, for the turbulent flow in the actual system. To complete the fluid mechanics description of this model duct flow, we note that the mean velocity of the fluid calculated for ReDH = 1827 is V'= 6.237 m/sec. The pressure drop per unit length in the fully developed region will be: dP 0V2 1 _ 3 d? - 4(Cf fully developed) fag'fig” 190'012 N/m (1'5) for this case. This information will be used in the solution of the energy equation. +4) Mass flow rate (kg/sec) ——> (x10 ) (m 12 2460 d 2050 d 1640 Q C) C! ' 820 _ 410 l I O 0 60 120 180 240 300 360 420 P (N/Mz) FIGURE 4. REYNOLDS NUMBER AND MASS FLOW RATE AS A FUNCTION OF TOTAL PRESSURE DROP. 13 Finally, we consider the nature of the thermal boundary layer development in the model system developed here. The non-dimensional thermal length is defined by Shah and London (9) as follows: + _ Lth th DHPe L where Lth = duct length from the point of duct heating to the point at which the temperature profile becomes fully developed. No value of L+ was found for the particular case th of interest here. However, an estimate of the behavior is obtained by considering laminar flow in a duct with aspect ratio = 0.3 when one wall is maintained at a con- stant heat flux and the other three walls are maintained at a constant temperature (9). The thermal development length Lzh obtained ushxrthis criteria is 41.5 cm, which is much greater than the heater length of interest. Thus the thermal boundary layer will be deve10ping along the heater window surface. CHAPTER 2 HYDRODYNAMIC CONSIDERATIONS 2.1 NON-DIMENSIONALIZED (FINITE DIFFERENCE) MOMENTUM EQUATION WITH SPECIFIC BOUNDARY CONDITIONS 2.1.1 NON-DIMENSIONALIZED MOMENTUM EQUATION The velocity field in the fully deve10ped laminar flow (directly beneath the heater surface) is obtained by solving the Navier-Stokes momentum equation in the x-direction (see Figure 3). The solution is developed in non-dimensional form in order to be general and the deri- vation of the non-dimensional equations from the dimen- sional form is given in Appendix F. The non-dimensional momentum equation is given as: a2U++3—2[.'J.«:=—l 3y+2 az+2 where U+=9—; U * 1 dp 2 U — F dx)a + y = y/a 2+ = z/a and 2a is the height of the channel shown in Figure 3. The boundary conditions applied to this equation are: i) U+ = 0 at y+ = 0 ii) 0+ = 0 at 2+ = 0 iii) §;;-= 0 at 2+ = b/a 2.1.2 FINITE DIFFERENCE APPROACH A finite difference numerical method of solving 14 15 this equation has been applied using the finite difference technique. The finite difference nodal representation of the control volume is given in Figure 5. In the case where + _ + _ + _ + _ + = + (Ay )n - (Az )e — (Ay )S — (Ay )w — Ay Az we get 2 + U+ U+ U+ U+ 8, 3 32" ((Jl—ahn- -—U—:> S)/Ay =(————P-- -P——U-S-)/Ay 3y 3y+ 3y (Ay )n (Ay )S therefore: + + + 3y 2 (Ay+)2 Similarly: + + - 20 + a2U+ = UW p + UE az+2 (Az+)2 Substituting these expressions into the governing equa- tion and simplifying, we get: + + + + + _ _ + 2 Un - 4Up + US + Ue + Uw — (Ay ) (2.1) 2.1.3 METHODS OF SOLUTION We can solve equation (2.1) by the following methods methods: i) Jacobi Method ii) Gauss-Seidel Method iii) Gauss-Seidel with S.O.R. Method These three methods are discussed in Appendix (A). The methods are iterative methods which are used when the problem is two-dimensional and there are large numbers of equations to solve simultaneously. It is very difficult 16 N o l---'"n--_I I l I | Wf J ! 43 W e | I l I I___._._.___! O S FIGURE 5. FINITE DIFFERENCE CONTROL VOLUME. 17 to handle two-dimensional equations with large numbers, so iterative methods are used with an initial estimated solution:hiorder to get new values. The rate of con- vergence of such methods depends on the type of iterative method used. Such iterative methods may not converge at all due to aspects discussed in Appendix (A). Other iterative methods,such as the Alternative Direction Implicit (ADI) Method, can also be used. This method has been used for a solution of the thermal problem (discussed in Appendix (B)). Patankar (3) suggests using the Gauss- Seidel technique with the line-by-line method or the successive overrelaxation method. We will discuss the methods used here and the conditions for convergence will be considered in Appendices (A) and (B). 2.2 HYDRODYNAMIC SOLUTION OF THE PROBLEM 2.2.1 NODAL CONVERSION OF THE FINITE DIFFERENCE EQUATION INTO NODAL FORM In the previous chapter, we have developed an appro- priate non-dimensional finite difference equation for the equation of motion in the 'X' direction: + + .+ + + _ _ + 2 on - 4Up + Us + 08 + Uw — (Ay ) + + + + . . . . s’ Ue' Uw' Up are non-d1mens1onal veloc1t1es where U+, U n defined previously and Ay+ is the non-dimensional distance between nodes. Due to symmetry about the center line, we only need to take into account one quarter of the rectangular 18 cross section and determine the velocities at other points within the rectangular cross section. Non-dimensional hy- drodynamic boundary conditions are shown in Figure 6 and for our solutions, we define four rows and eleven columns so that Ay+ = Az+ = 0.333 These node arrangements on one quarter of the rectangular cross section in terms of non-dimensional distances are shown in Figure 6. In this case, equation (2.1) from the previous section will be reduced to 0+ - 40* + 0* + U+ + U+ = -0.1111 (2.2) n p s e w Now defining i and j as follows: I + I 0 row number in z d1rection i j = column number in y+direction we can write this equation in terms of (i,j) as follows: 4U+(i,j) =tfki,j+1)-+U+(i,j-1) +tfiki+1,j) + U+(i-l,j) + 0.1111 with the following boundary conditions i) U+(l,j) = 0, j = 1 to 11 ii) U+(i,l) = o, i = 1 to 4 iii) at i = 4: au+ = 0 so U+(5,j) - U+(3,j) 0 + ' + By 2(Ay ) Thus U+(5,j) = U+(3,j), j=l to 11 iv) at j = 11: 19 mmfimM¢OZDOm UHEdzwmomn—wm QZOHmZmZHn—IZOZ .w mmeHh x 26 u +n V. II + 5 8 C Amm.m.ov 1* 1o ’1. l 4a 'H H + ++ 1.1+: H1. 1 1r 11 H \\\ o I Dru] (om \\\ + xxx‘xxx \\\\k\\\$\\\ xxx \\\%\\\r\\\\\\\u Amm.m.: 8.: \\\\ 0 ll .12. (o (O 20 + +. +. 30+ = 0, so U (1,12) - 0 (1,10) 32 2VD mmHeHoogm> .m mamas bmmv.o mm¢.o mv.o mmv.o mhv.o mv.o mmv.o mmm.o mam.o mmH.o o o.H mmv.o wm¢.o mmv.o Hm¢.o mm¢.o Hv.o 5mm.o mvm.o vmm.o hhH.0 o v.0 th.o vhN.o mnm.o hm.o mom.o mmN.o va.o NNN.o vma.o mHH.o o mmm.o o o o o o o o o o o o m.o mmm.m o.m hm.m mmm.N o.N hm.H mmm.H oo.H bow.o mmm.o o N +% + ZOHBHmOm AdZOHmZmSHQIZOZ m0 ZOHBUZDW fl md A+Ne+>v+D mMHBHUOQm> AdZOHmzmSHQIZOZ .H Hflmdfi 22 using the formulae discussed in Appendix A. The Jacobi, Gauss—Seidel, and S.O.R. methods are discussed in Appendix A along with the respective criteria of convergence. Using equations (a) and (b) in Appendix A: 6(B) = S (cos n/p + cos n/9) for a unit square mesh of dimension 'h' and a total rec- tangular cross section of dimensions ph by qh, we get in our case, ph=l qh = 3.333 where h = Ay+ = Az+ = 0.333 p = 3 and q = 10 Substituting the values of p and q above, we get 6(B) = % (cos n/3 + cos n/lO) = k (0.5 + 0.951) = 0.725 Substituting the calculated value of 6(B) as discussed in the appendix: 2 w = = 1.185 b 1 + /1 - (0.725)2 The computer program used to solve the hydro- dynamic problem, using the successive overrelaxation method along with the Gauss-Seidel method, is given in Appendix D under the name 'HEAT.‘ 2.2.3 NON-DIMENSIONAL VELOCITIES AT DIFFERENT NON-DIMENSIONAL POSITIONS OBTAINED BY THE PROGRAM 'HEAT' The results obtained from the computer program after fifty iterations are given in Figure 7. 23 ZOHaumm mmomU M§ADGZ¢BUHM mmB m0 mmfimdDO MZO ZO mZOHBHOZOU MfidQZDom UHEdZNQOfiDflE 924 NBHUOQfl> AdZOHmZmSHDIZOZ .n WMDOHM x mm.m u +6 m. Amm.m.ov £3 24 2.2.4 CALCULATION OF VELOCITY AT DIFFERENT POSITIONS IN THE DUCT To convert the non—dimensional velocities to dimensional velocities at different positions in the y and 2 directions, we have to first determine the value of (u*) which is used to non-dimensionalize the x-direction momentum equation. The quantity (u*) is defined in the previous chapter as follows: it: 3.92.2 u u(dx)d (m/seC) where a = 0.15 cm. We have calculated the pressure drop dp . . (3;) 1n the prev1ous chapter. Calculation of the velocity at any position in the fully developed region can now be made from the non- dimensional velocity results obtained from the program: (m/sec) Substituting the value of g3 and a from the previous chapter, .519 _N_ = 3; 190.012 m3 and a 0.15 cm u* = 1 x 190.012 (N) x (0.45 x 10’2)2(m2) U(lg see) 33 m . l 190 012 o 15 0 15 10"“) ( / ) ‘l8.53x10"5x . x . x . x 10 sec 23.07 m/sec 25 Thus u* = 23.07 m/sec for the fully developed flow. Calculating the velocities now using the appropriate non- dimensional velocity result and substituting into the following, u(y.2) _r__. u = u+(y.2) where u+(y,z) is the non-dimensional velocity obtained at a given position. For example, we are interested in evaluating the velocity at the center of the duct which is the maximum velocity, €¥ = 0.4937 (from Figure 7) ucenter = 23.07 x 0.4937 = 11.39 m/sec 2 . 2 . 5 CALCULATION OF FRICTION FACTOR The friction factor is now calculated from the results in order to compare with the assumption that Cf = %¥4§l-for fully developed flow. e DH We defined the following earlier: * _ l d 2 u — E(3)2961 (2.3) dp 0v2 1 26 Substituting (2.4) into (2.3), we get 2 2 * _ P_V__ 3.. uDh so, C = 12 f 2 D(%?) a2 2 UDh C = 5— f 20V+VDh 8.2 where V+ = um/u* and 11m = the mean velocity = V and considering the definition of the Reynolds number, 2 h 2v+ (Re) (a2) D Cf = Now to calculate V+ from the results obtained, the trap- ezoidal rule is applied which is in final form below. A . C + * 1k . S1nce um : iL- I udAc and U = u/u where u is constant c 0 we derive: V+ = i%- [%((U+) in corners of the duct) + c %(sum of all U+ on the surface of the duct excluding corners) + (sum of all U+ inside the duct):] 4. In our case, V = 0.346 and the computer program calculated friction factor is: Cf==l4/Re We believe that the dif— DH' ference between this value and Cf = 17. Sl/ReDH is due to the finite difference approximations in the solution. The ap- proximation in the integral above is another sourcecnferror. No calculations were performed to isolate the difference. This friction factor calculated from the non-dimensional velocity results is less than the Shah and London value by approx- imately 20%. 27 2.3 CONCLUSIONS AND COMMENTS ON THE HYDRODYNAMIC SOLU- TION 2.3.1 CONCLUSIONS ON APPROACH USED TO SOLVE THE PROBLEM The hydrodynamic solution for a rectangular duct flow is achieved numerically by the finite difference numerical approach. We observe from the output of the computer program (which is given under the name 'HEAT' in Appendix D) that after forty-eight iterations using the S.O.R. method with the Gauss-Seidel method, convergence is achieved. The first attempt to solve the problem was using the Jacobi method. A convergent solution was not achieved. This is due to constraints for convergence which are stated in Appendix A. These constraints state that the solution converges only when the coefficient of diagonal element ex- ceeds the sum of the off-diagonal elements. Another method, which can be used for solution of the finite difference equation, is the Alternative Direc- tion Implicit Method (ADI), which is discussed in Appendix B. Taking large numbers of odd and even steps, in the 'X' direction, after a certain number of steps the solution will be constant for all steps in the 'X' direction, which will reflect the fully developed velocity profile. Successive overrelaxation is also used for this type of problem and it is a very good alternative to the ADI method, using an optimum overrelaxation factor (w). In the present study, the successive overrelaxation method 28 is used to solve the hydrodynamic finite difference equation. 2.3.2 CONCLUSIONS AND COMMENTS ON THE SOLUTION ACHIEVED BY THE S.O.R. METHOD (1) The solution converged after forty-eight iterations. (2) From Figure 7, we can conclude that at the centerline in the z direction (i.e., z = 0.5 cm), the velocity at y+ = 1.0 (y = 0.15 cm) is a maximum and the ratio: Ucenter:(23.07)x(computer non-dimensional) U 6.237 velocity at center mean so Ucenter = 23.07 x 0.4937 = 1 826 U 6.237 ’ mean h U - ( *)( + t ) d * - 23 07 / w ere center — u u cen er an u - . m sec. Shah and London (9) and Knudson (10) also gave experi- mental results for rectangular cross sections for the same problem as follows: IJ _ _ y n _ z m U — (1 (6) )(1 ‘4’ ) (2.5) max where m = 1.7 + 0.5 (0L*)-1'4 and n = 2 for 0* E 1/3 = 2 + 0.3(o*-l/3) for 0* 1 1/3 *= . =E=0015= where a aspect ratio b THE? 0.3 Shah and London (9) also gave the final result of the ratio of maximum velocity to mean velocity as follows: 29 U max _ m+l n+1 U - (——In >(-—n) <2-6> m _ -1.4 _ so, m - 1.7 + 0.5 (0.3) — 5.395. Substituting in equation (2.6), we get Umax _ 5.395+1 2+1 _ 6.395 x 3 U ‘ ( 5.395 )( 2 ) ‘ 5.395 x 2 = 1'77 Comparing the results of the present computer model with the Shah and London value yields a three percent differ- ence : 1.826 - 1.77 1.77 = 3% This agreement indicates good agreement between the present computer results and previous work. However, these values may agree without good agreement at the wall (and therefore Cf factors) since Um and Um are less ax sensitive to the approximations made at the boundary than the actual boundary values themselves would be. (3) If we consider Figures 8 and 9, the wall gradient of velocity in the 2+ direction is smaller than the gradient of velocity in the y+ direction. This is 30 ’l‘ 1]? Y+ *——3.333——‘* Z 1.o-—— — —— __ _$ CL __ 0.667- 0.333d 0 0.1 0.2 0.3' 0.4 0.5 NON-DIMEN§IONAL VELOCITY ——————9 U (Y ,3.333) FIGURE 8. FINITE DIFFERENCE APPROXIMATION OF THE NON-DIMENSIONAL VELOCITY VARIATION IN Y+ DIRN AT Z+ = 3.333 31 Y+T' 1.0 Y+' 5L K——3. 333—) 2+ NON-DIM. VELOCITY, U+(1.0,Z+) I 0.333 0.61 1.0 1.3331.61 2.0 2.333 2.61 3.0 3.333 z+-———> FIGURE 9. FINITE DIFFERENCE APPROXIMATION OF THE VARIATION OF NON-DIMENSIONAL VELOCITY IN z+ DIRECTION AT y+ = 1.0 32 due to the fact that the viscous stress in the z+ direction is less than the viscous stress in the y.'. direction (i.e., Bu Bu EVA—2" CHAPTER 3 THERMAL CONSIDERATIONS 3.1 DIMENSIONAL APPROACH FOR THE THERMAL SOLUTION OF THE PROBLEM Although a non-dimensional solution to the energy equation could be developed as was the case for the velo- city field, only a dimensional solution is presented at this time for the specific model system presented here. 3.1.1 INTRODUCTION As described in Chapter 2, the present problem is assumed to consist of boundary conditions on each side of the cross section as follows: (1) three sides of the rectangular cross section are insulated, i.e. q" = 0 (2) on one side of the rectangular cross section, q" = constant (axially and laterally) The boundary conditions are shown in Figure 10. The sole thermocouple used in the cryomicroscope system is placed on the center of the window heater typically and this is shown in Figure 10, where the position is denoted by "B." The energy equation for hydrodynamically fully developed flow and thermally developing flow is derived by the control volume approach and is discussed by Keys and Crawford (8): 2 2 2 3T 3 T 3 T 3 T U-—— = a ( + ———-+ ) (3.1) 3 8x2 By 322 82T where a 5—7 is the axial conduction term. X 33 34 BUDQ mqwz¢+ 4 Addaid E #1. 0 L 3&4: 8.3 WI) {'0 39 32T + BZTI = u (Tx+Ax - Tc) 822 x+Ax 3y2 x a AX and in finite difference form, 0 o o 0 TN 2Tp-l-TE + TN -2Tp -+TS = U(I,J) (Tp-Tp ) (Az)2 (Ay)2 a Ax In our program (U(I.J) (yz) a Ax ) = R(I,J) = constant at given position So our finite difference equation for even steps will be ° = (R(I,J) + 2) + T P O o TN+TE+TN +TS +(R(I,J) 2)Tp (3.3) For odd steps, and similarly as described in even steps, the finite dif— ference form of the energy equation will be as follows: T °+T °+T +T +(R(I,J)-2)TPO= (R(I,J)+2)TP w E N S (3.4) The thermal entry length problem was solved by digital computer using the finite difference equations by the ADI method. The only difficulty which was encountered was that on the sides of the rectangular cross section, the velocity is actually zero but the finite difference -v—-'-( 40 equations must account for the convection properly near the wall. Therefore, due to the linearity considered in the finite difference equation while formulating the problem, the velocity at the wall at a given point is assumed to be: I UP+UE+US+US U(se) = 4 (see Figure 12) A s1m11ar procedure 15 used for U(ne)’ U(nw)’ and U(sw)’ After calculating all the velocities above, UP is calculated: U = U(se) +U(SW) +U(nW) +U(ne) P 4 This linear approach is assumed because throughout the study, (Ay = A2) is considered. 3.2 THERMAL PROBLEM SOLUTION BY THE FINITE DIFFERENCE APPROACH 3 . 2 . 1 INTRODUCTION We have discussed the assumptions and the dimen- sional approach to the thermal problem of interest in the previous section. We have applied the Alternative Direc- tion Implicit (ADI) method for this three-dimensional problem which is discussed in Appendix B, along with the TDMA method which is discussed in Appendix C. The governing equation is formulated in finite difference form along with the boundary conditions in the previous chapter. The computer program which is designed to solve this problem is given in Appendix E under the name "Convection Program." 41 N (UN) 'w'-n (U ) w 4 :2 E (0W) Use S'- O O s (Us) Us, FIGURE 12. VELOCITY ASSUMPTION IN CALCULATING VELOCITY ON BOUNDARY FOR THERMAL SOLUTION. 42 In solving the problem by the convection program, the R(I,J) values which are defined in the previous chapter U(I.J) (Ay)2 a A): Figure 13 and used in the program. The computer results as R(I,J) = ( ) are calculated and shown in are shown in Table 3 at different y and 2 positions as a function of the x direction. The thermocouple is placed along the centerline of the window. The predicted variation of temperature at the centerline location of the thermo- couple is shown in Table 4 with respect to the x direction (i.e., distance from the leading edge). 3.2.3 CALCULATION OF MEAN FLUID TEMPERATURE AT EXIT OF WINDOW REGION The mean temperature of the fluid at x = 5 cm is calculated in order to perform an energy balance on the system for examining the validity of the calculated solution. The trapezoidal rule in two dimensions is used as given below. We can use the Simpson rule also by modifying it for two dimensions according to Hornbeck (12). It is quite accurate to use the Simpson rule for hydrodynamically and thermally developing flow in the rectangular duct (12). The mean temperature is defined in the general case by Rays and Crawford (8) which is as follows: _ 1 TM — AcUm I UTdAc A C 1 b a TM = A U I I UT-dY-dz C m oneAom mom comes: H O a mme zH ammo In.va mmemz¢m¢m amoedmmmzme qeonmzmzHouzoz .me mmDon xx 20 H x 44 oo.mm oo.mm oo.mm oo.mm oo.m~ oo.mm om.o oo.mm oo.mm oo.mm oo.m~ oo.m~ oo.m~ mm.o oo.m~ oo.mm oo.mm oo.m~ oo.mm oo.m~ om.o oo.mm oo.mm oo.mm oo.mm oo.mm oo.m~ mH.o oo.mm oo.mm oo.mm oo.mm oo.m~ oo.m~ H.o oo.mm oo.mm oo.m~ oo.m~ oo.m~ oo.m~ mo.o oo.mm oo.mm oo.m~ oo.mm oo.mm oo.mm o m.o v.0 m.o N.o H.o o mwmmww 80 o u x .mmfidmm HEB m0 MOON UZHQdmA mmB 20mm mm02¢8 ImHD mDOHm¢> Bfl RUDD m¢QDw24Bumm mmB ZH mZOHBGUOQ BzmmmthD Ed mflMDBfiMMQZmB .m mqm<9 ..ZOHBUM>ZOU. dewomm MMBDmZOU wm DMZHdBmO mBADmmm N.N.m 45 em.eme eH.mmH me.eme ee.mee we.mma -.mm~ om.o me.~m em.~m om.mm mm.em me.~e me.me m~.o mm.m~ em.mm om.m~ me.m~ Ha.e~ ee.mm ~.o o.mm oo.mm oo.mm Ho.m~ mo.mm m~.e~ mH.o o.m~ oo.m~ oo.m~ o.m~ oo.m~ ee.mm oa.o e.m~ oo.m~ oo.m~ o.mm oo.m~ efl.m~ mo.o o.m~ o.m~ o.mm e.m~ o.m~ H.m~ o m.o e.o m.o m.o H.o 0 So» EON/y =6an poscflucoo .m mqm om.o mN.mN mm.m~ Nm.mm mv.mm mw.vN mm.hm mH.o No.mm No.mm No.mm mo.mm vN.mm wo.hm 0H.o oo.mm oo.m~ oo.mm Ho.mm ho.mm mm.vm mo.o oo.mm oo.m~ oo.mm o.mN oa.mm mm.mm o xsow m.o v.o m.o N.o H.o o hWW// Eo m nx conceucoo .m mqmde 48 mm.mmm ww.oom vm.mmm mm.m>m mm.hom om.mmm om.o mm.mm m¢.mm hm.mm mm.¢m mv.mma om.mma mm.o mH.Hm mm.am om.am Nw.mm Hw.mv Hm.hm om.o mo.mm H>.mm m>.mm bo.vm hh.mm vo.wv mH.o mo.mm mo.mm mo.mm mo.mm >m.mm mv.m~ oa.o oo.mm oo.mm oo.m~ ~o.mm m~.m~ Hm.m~ mo.o oo.mm oo.mm oo.mm ~o.mm m~.mm om.em o 50% m.o ¢.o m.o N.o H.o o EON 80 v u x toseeucoo .m memes 49 sea.msm mm.ms~ sa.~m~ em.mm~ mm.omm sm.msm om.o mm.aoa mm.moa oa.moa m~.eHH me.meH ms.mo~ mm.o mo.mm ~m.mm ~m.om mm.mm em.mm ma.aoa om.o mm.em mm.em Hm.e~ so.mm os.m~ mo.am ma.o mmH.mN mH.mm 6H.mm m~.mm He.e~ mm.~m oa.o Ho.mm Ho.mm mo.mm eo.mm ~m.m~ ms.o~ mo.o Ho.m~ Ho.mm Ho.m~ mo.m~ mm.m~ www.mm o m.o e.o m.o m.o H.o o ewwx =8me eosceucoo .m mqmee 50 TABLE 4. TEMPERATURE OF THE HEATER CENTERLINE AS A FUNCTION OF DISTANCE FROM THE LEADING EDGE OF THE HEATER. x cm TherfigggupIg (0C) 0 23.00 0.2 58.10 0.6 116.12 1.0 154.57 1.4 181.19 1.8 200.63 2.2 215.63 2.6 227.85 3.0 238.28 3.4 247.52 3.8 255.94 4.2 263.76 4.6 271.12 5.0 278.12 51 According to the trapezoidal rule for one direction (ref. 18) for: a A (JU-Tdy) =71 (f0+2f1+2f + +2f O 2 n-l + fn) f f f are 0’ l’ 2' "' n values of (UT) at y = 0, Ay, 2Ay, 3Ay, ... nAy. Now where f = f(i,j) and in our case f integrating again with respect to the z direction: _ l y, z Tm‘A‘EIE‘z" —2— [:f(0,0)+2f(0,Az)...+2f(0,(m—l)Az)+f(0,mAz) +2(f(Ay,0)+2(Ay,Az) ... +2f(Ay,(m-l)Az)+f(Ay,mAz) +2(f(n-1)Ay,0)+2f(n-1)Ay,Az) ... +f(n-l)Ay,mAz)) +f(nAy,O)+2f(nAy,A2) + f(mAy,mAz):) The product (UT) at different locations is given in Figure 14. Because of the symmetry about the center line in the z direction of our problem, using the trapezoidal rule, where C: II V = 6.237 m/sec 0.3 x 1 x 10"4 m2 11’ ll Ay = A2 = 0.05 m we can calculate 'Tm' from the half cross section as follows: = ZAy-Az m A U [%(UT(7IO)+UT(1:0))+%(Sum of product C m 52 :8 m n as 38sz was m0 mwflm OZHAH/wmm. mmB HR 8030 mam. ZH mMQOZ BzgmmmHQ Ed ABDV mmbfigmmEMB 92¢ NBHUOQMNV .mO BUDoomm mmB .QH mmeHm smemzzwm W so m .o 0— .m0 mHXAN cw...» 2...; .23» rd.» «.3 we...» 58... $23 3.: to. .5 t... i :3... r so), 5.: $5.42.. 0 . .mwte v» +r~ C 3 W 3.3». «t. as. ofl¢~3 flaws, 11w «:5 5a 53 (UT) at all nodes on boundary y=0, z=0, y=0.3 cm) + (sum of product (UT) at all nodes inside the cross section[] Calculating the mean temperature from Figure 14 using the above approach, we get: _ o Tm — 48.1 C This value will be compared to that computed from a first law energy balance. The lateral variation (in the z direction) of temperature on the window at given distances from the leading edge of the heater is shown in Figure 15. In Figure 16, the variation of the centerline temperature at different distances from the leading edge of the heater is shown. 3.3 CONCLUSIONS AND COMMENTS ON THE THERMAL SOLUTION 3.3.1 GENERAL COMMENTS ON THE RESULTS The solution to the energy equation with the specified thermal boundary conditions is discussed and pre— sented in the previous chapter. This solution was obtained by using the finite difference method with the ADI approach. Figure 15 shows the variation of temperature at different lateral positions on the window surface at given distances from the leading edge. Figure 16 shows the centerline tem- perature variation (of the thermocouple) from the leading edge of the window. The following conclusions are derived from the results obtained from the computer program: T (°C) 54 y=(L3CM i=7 WINDOW HEATER SURFACE (0'0) i=0 390 z = 0.5 CM 2 360 330 300 TX x = 5 CM 270 1x = 4 CM 4 . 240- x = 3 CM 210- x = 2 CM 180- % O 120- 90- 6o- 30' x = 0 CM 0 0.2 0.3 0.4 0.5 Z (CM) 5 FIGURE 15. PREDICTED VARIATION OF TEMPERATURE AT WINDOW HEATER SURFACE IN 'Z' DIRECTION AT A GIVEN 'X' 55 V 280 260 240 220 qfifig 200 180 160 140 120 100 80 60 4. mean TEMPERATURE AT TYPICAL THERMOCOUPLE POSITION, B, OC€> 20 0 0.4 0.8 1.2 1.6 2.0 2.4 2.8 3.2 3.6 4.0 4.4 4.8 5.2 X (CM) ————> FIGURE 16. VARIATION OF THERMOCOUPLE TEMPERATURE WITH DISTANCE FROM LEADING EDGE 56 1) On the window surface, the temperature varies both laterally and axially. The lateral variation is shown in Figure 15. The lateral temperature gradient (g?) depends on the position in the rectangular cross-section. To illustrate the magnitude of these gradients we calculate the approximate gradient at x = 1 cm: 3T é T(0.3,0.l) - T(0.3,0.0) Az _ 195.68 — 252.22 _ -56.54 _ _ 0 — 0.1 — T—i— - 565.4 C/cm Examining the results of Table 3 indicates that the largest lateral thermal gradients appear closest to the beginning of the heated length. This region warrants a closer examination in future work. The maximum temperature gradients in the x direc- tion are also a function of distance from the leading edge of the heater because of the developing temperature profile. The maximum axial temperature gradient also occurs at the leading edge of the heater and is approximately: 8T ; T(0.4,0.3,0.5) - T(0.0,0.3,0.5) .3? - —0.4 max é 58.100-223.00 = 175.50C/Cm The temperature gradients in the x direction (i.e., axial direction) decrease with distance from the leading edge of the heater. 2) The temperature at the leading edge is 230C 57 while temperature at the trailing edge of the window is different for different locations on the window surface. However, the maximum temperature is always at the corners of the duct near the heater window. 3) The temperature at the centerline increases rapidly within one centimeter from the leading edge, then the gradient decreases considerably. As shown in Figure 16 the mean temperature of the fluid increases linearly, a consequence of the first law for uniform heat flux. The 3Tmean 3x the length of the heater (5 cm) indicates that the flow is fact that the slopes of and HTS/8x do not match over still thermally developing at x = 5 cm which confirms our earlier estimate predicting that this would be the case. 4) The temperature along the line '8' where the thermocouple is typically placed is the minimum temperature on the surface of the window at a given distance from the leading edge of the window. 5) We have neglected the axial conduction term while solving numerically the thermal problem, which is justified on an order of magnitude basis. The Peclet number in our case is equal to 1827 x 0.72 (=Re'Pr) = 1315 which is much higher than five (5). The axial conduction term can be neglected for cases having Peclet numbers more than '5' as discussed in Keys and Crawford (9). 3.3.2 ENERGY BALANCE CHECK OF THE RESULTS Checking our solution with the help of the energy balance by the first law of thermodynamics which 58 states that for the steady state energy in = energy stored + energy out i.e., IECPAT = Q (3.5) where AT = sz - Tbl sz = mean bulk fluid temperature at trailing edge of window surface Tbl = mean bulk fluid temperature at the leading edge of the window. (kg/m3) (m2) In our case, fi1== p x A x V (m/sec) (kg/sec) = 1.1766 (kg/m3) x 0.3 x 10"4 (m2) x 6.237 (m/sec) 2.207 x 10'4 (kg/sec) cp = 1.005 x 103 joule/kgok Q = 5 watts Substituting into equation (3.5), we get 5 (w) P -4 3 w.sec 2.207 x 10 (kg/sec) x 1.005 x 10 (E65E-) _ _ = 0 AT — sz Tbl 22.70 c T = 22.70 + 23.0 = 45.70°c b2 In the previous chapter, the solution technique for mean temperature was discussed. We can use either the trape- zoidal rule or Simpson rule in two dimensions to calculate the mean temperature at the trailing edge (i.e., x = 5 cm). Using an approximate method, i.e., trapezoidal method, the mean bulk temperature at the trailing edge 59 is calculated in the previous chapter and it is equal to 48.1°c, i.e., = 48.1°c. T52 Comparing the mean bulk temperature obtained by the energy balance and by the finite difference method, they are quite close and the error of the finite differ- ence solution is: 48.1 - 45.70 _ Tb2(finite diff.)'.Tb2(energy bal.) 45.70 ‘ Tb2(energy balance) 5.3% This error is due in part to the velocity approximation on the surfaces of the duct and the approximate solution tech- nique (i.e., alternative direction implicit) to solve the thermal problem with specified boundary conditions. The trapezoidal method of solving the mean temperature at the trailing edge of the window is also approximate. Though the error in the velocity solution is multiplied by the error in the temperature solution, i.e., error in (UT) = (error in u)(error in T) the error as a whole is 5.3% which is acceptable for our present purposes. 6) This problem for a variety of boundary conditions is also solved by V. Preingerova and P. H. G. allen (14) by adimen- sionless approach for both thermally and hydrodynamically develOping flow. However, the problem was solved for 60 quite a small magnitude of length in developing and de- veloped (hydrodynamically) region. 3.3.3 COMMENT ON PARALLEL PLATE APPROXIMATIONS Comparing the actual rectangular duct solution with a parallel plate approximation solution which is discussed below, we can examine the conclusion given by Hornbeck (12) which states that if the aspect ratio is less than '5', a parallel plate approximation is rejected in order to avoid misleading results for a rectangular duct flow. There are two methods of solving the parallel plate problem which are considered here and discussed by Kays and Crawford (9), who give an hydraulic diameter approach tabular form to solve the flow behavior through the parallel plate. Cess and Shaffer (16) used the analytical approach in which the parallel plates are heated equally axially but at different levels of heat flux on each plate. For our approach Um should be calculated for Reynolds number = 1827 according to the parallel plate approximation, 4 x a (m) x Um (m/sec) Re = 1827 = 2 V(m /sec) where Dn = 4a (for parallel plate) Um = 4.8 m/seC' Now defining the Peclet number: 4Uma (m/sec)m a mZ/sec 61 where a = 22.106 x 10-6 mz/sec (in our case) Then Pe = 1315 Now using the graph which is derived analytically for ]_ x Twi - Tb 155-5), as a function of' (q,a) for q2 = 0 (i.e., ( one side heated and one side insulated), we get the following results: TABLE 5. RESULTS BY CESS AND SHAFFER METHOD. x cm Twi-Tb(ot) Tb(oc) Twi(oC) 0 0 23.0 23.0 0.5 116.68 28.37 144.85 1.0 146.32 33.74 180.06 2.0 188.2 44.78 232.68 3.0 220.36 55.25 275.60 4.0 245.6 65.95 311.55 5.0 263.98 74.54 338.50 Kays and Crawford also obtain results for the same problem using the hydraulic diameter approach and obtain very similar results. Comparing our solution for the flow in the rectangular duct which is given for point '8' in Table 4 and the parallel plate solution by Cess anui Shaffer which is given in Table 5, we conclude that the results are quite different for a given aspect ratio if we were to simplify the problem using a parallel plate assumption 62 and the solution will be misleading in all respects. So for an aspect ratio less than five, the effect of corners are significant and cannot be neglected. CHAPTER 4 DISCUSSION AND CONCLUSIONS 63 4.1 DISCUSSION AND CONCLUSIONS A preliminary quantitative understanding of the temperature gradients on the heater surface of a convec- tive cryomicroscope heat transfer system has been realized by developing a simplified model of the actual complex situation. The primary aim of the present study was to con- sider the effects of the develOping thermal boundary layer with respect to the thermal gradients on the window heater surface. The numerical solution to this simplified problem, to the best of our knowledge, represents the first solu- tion to a thermally developing,laminar boundary layer in a rectangular duct of aspect ratio 0.3 with constant heat flux on one wall and with the other three walls insulated. The solution is in qualitative agreement with experimental results obtained by Hrycaj on the actual complex system indicating very large thermal gradients in the axial di— rection. (See Figure 2.) Furthermore, the present model predicts severe thermal gradient problems in the lateral direction as well as the axial direction. (See Figures 15 and 16.) The velocity field solution has been checked by calculating the friction factor, Cf, in the hydrodynam- ically fully developed region and the ratio (Umax/Umean) and comparing these values with those available in Shah and London(9). We find that there is excellent agreement with the values of (U /U ) (+_w_(b,_ z a, .x,_ X0 1 at. l o_ I 1 11 3—1 3 aijxj(n)), i=l(l)m "Ma i The above equation can be simplified and written as fol- lows: 1"]. m i ii j=1 3 3 j=i+1 3 3 The factor w' in above equations is called ac- celeration parameter or relaxation factor. It generally lies in the range 1 < w < 2. The determination of the Optimum value of 'w' for maximum rate of convergence will be discussed in the next tOpic. The value w=l gives Gauss Seidel Iteration. This method is also suggested by Robert Hornbeck (ref. 5) at National Aeronautical Research Center for developing flow hydrodynamically in rectangle duct which is quite complicated problem compared to fully developed flow (hydronamically) in solving large number of two- dimensional equations simultaneously. 73 Patankar (ref. 3) suggests on page 66 that to solve the two—dimensional conduction problem of similar formulation as above formulations, S.O.R. method with line-by—line method which is described in the same book by the author. TO CALCULATE THE OPTIMUM RELAXATION FACTOR 'w' Smith (ref. 4) gives the following criteria to decide optimum overrelaxation factor for S.O.R. method (page 252): 0(8) = % (cos % + cos 3) (A) for square mesh of side 'h' and for rectangular cross section of sides 'ph' and 'qh' using five point differ- ence approximation and wb = optimum relaxation factor = 2 (B) l + /Ql-02(B)) where p(B) is the spectral radius of the Jacobi method iteration matrix. This condition is discussed by Smith in detail on pp. 243-50. Roache in his book Computational Fluid Dynamics also suggested above approach by Smith for calculating suc- cessive overrelaxation factor (w) which is given on pp. 118 in the same book. 74 APPENDIX B: ADI METHOD The equation we used in thermal solution of the problem in general is three-dimensional equation and it is as follows: N 2 2 M + 3% = R(I,J) g}; 82 By For the lst, 3rd, ... (odd) steps in x direction, 2 2 3—%-| + §_% | = R(I,J) ¢x+AxA - ¢x 32 X By x+Ax AX where R(I,J) is constant and it is either dependent on position or not according to the problem to be solved. In finite difference form, taking Figure 5 into account, o o o 0 ¢ - 2¢ + ¢ ¢ - 2¢ + ¢ ¢ -¢ w 92 e + N E 2 S = R(I,J) —E-——S p (62) (6y) X In our case if Az = Ay, then the above equation is re- sulted in the following final form for odd steps in x direction: 2 _ = R(I,J)(Ay) _ o + 4n 2¢p + ¢S Ax (¢p 4p ) o o o w - 2dip + ¢e ¢ (A) For the 2nd, 4th (even) steps in x direction: ¢x+Ax -¢x) Ax (:6 2 1% + | = R(I,J) ( Bz x+Ax 3y x N and simplifying into finite difference form, 75 ¢w-2¢p+¢e (<82)2 + o o o ¢n 2q’p +¢s = R(I,J) (¢ _ o) 2 6x P P (6y) and if Sz = Sy, we get the final finite difference equa- tion for even space steps in x direction as follows: _ o _ o o = 2 ¢w 2¢p + 0e + ¢n 2¢p + 63 R(I.i;(4y) (¢p‘¢p°) (B) In above equations (A) and (B), we used 0 superscript which denotes the present time and no superscript means future time for the unknown property (¢)- This is called Alterative direction implicit method due to the fact that for odd steps to solve the equations, previous step iteration value which are available are used for the present step as the known value at points 'w' and 'e' and the ¢ value for point p, n and s are unknown value so solving one dimensional problem with TDMA method which is discussed in Appendix C, we get the unknown value of 0 at 'p', 'n' and '5' points. Similarly, for even steps to convert the equation into one dimen- sional, previous step 0 values at 'n', 's' and 'p' are used and ¢ values at w,p and s which are unknown at given step are found by solving one-dimensional prob- lems with the help of TDMA method. Another alternating direction method is dis- cussed by Yarenko (ref. 7) which is called method of fractional steps. 76 Step 1,3 . . . (odd) 2 a _ 29¢ 3:; _ R(I,J) .5; Steps 2, 4 . . . (even) 2 1%: R(I,J) 39. By 3X But this method cannot be used for the problem concern- ing developing nature of '0' in 'y' and 'z' direction and also in 'x' direction. This method can be used for solving fully developed velocity profile or fully de- veloped temperature profile. 77 APPENDIX C: TRIDIAGONAL MATRIX ALGORITHM This method is used to solve one-dimensional (Patankar, ref. 3) discretization equations. This is also called Gaussian elimination method. The designa- tion of TDMA refers to the fact that when the matrix of the coefficients of these equations is written, all the non-zero coefficients align themselves along three diag- onals of the matrix. If we do note equations in general, Bi¢i = Ci¢i+1 + Ai¢i-1 + Di (C1) for i = l,2,3,... Thus the quantity '¢i' is related to the neighboring '0' quantities ¢i+l and ¢i-l' To account for the special form of the boundary point equations, let us set A1 = 0 and CN = 0 (C2) so that the 00 and ¢n+l will not have any meaningful role to play. These conditions imply that '01' is known in terms of '02'. The equation for i=2 is a relation between 01, ¢2, and ¢3. Since ¢1 can be expressed in terms of ¢2, ¢2 can be expressed in terms of ¢3. This process of substitution can be continued until ¢n is expressed in terms of ¢n+l' But because (¢n+l) has no meaningful existence, we actually obtain the numerical value of '¢n' at this stage. This enables us to begin the 'back substitution' process in which '¢n-l. is ob- tained from ¢n’ ¢n-2 from ¢n-l' ¢2 from ¢3 and ¢1 from 78 '02.' This is called TDNA. Let us describe forward substitution process as follows: ¢i = Pi¢i+1 + Qi ¢i-1 = Pi-1¢i + Qi-l (C°2°a) Bi¢i = Ci¢i+1 + C1 (Pi-l¢i + 91-1) + Di where P = C1 Q _ Di+AiQi-l (C 3) i B. - A. P. ' i ’ B.-A.P. ° 1 1 1-1 1 1 1-1 for i=1, Ci Di Pi '-'-' E— and Qi = fi— (C.4) l 1 At the other end of Pi' Qi sequence, we note that Cn=0 so P = 0. n So we can write SUMMARY OF TDNA (1) Calculate Pi and Qi from equation (C.4) (2) Use the recurrence relations (C.3) to obtain Pi and Qi for i=2,3,...n (3) Set 9n = Qn (4) Use equation ¢i-I= Pi-l¢i + Q1-1 for i=n-l, n-2, ... 3,2,1 to obtain ¢n-l’ ¢n-2' ... ¢3I¢21 ¢lo This technique is very powerful means to solve one dimen- sional equations and it requires computer storage and computer time proportional only to 'N'. 11 12 10 14 13 30 35 50 17 16 15 79 APPENDIX D: PROGRAM 'HEAT' DIMENSION U(50,5,12),A(50,5,12) CALL SRCH$$(2,'0UTPUT',6,1,N1,N2) DO 10 K=l.50 DO 11 J=1,ll U(K,1,J)=0.0 DO 12 I=l,4 U(K,I,1)=0.0 CONTINUE D0 13 J=l,ll DO 14 I=l,4 U(1,I,J)=0.0 CONTINUE CONTINUE WRITE(5,30) FORMAT(' K I J U') DO 15 K=1,49 D0 16 I=2,4 U(K,I,12)=U(K,I,10) D0 17 J=2,1l FORMAT(I3,SX,13,SX,I3) U(K,5,J)=U(K,3,J) A(K,I,J)=U(K+l,I-1,J)+U(K,I+l,J)+U(K,I,J+l)+U(K+l,I,J-l)+0. $1111 U(K+l,I,J)=—OO.185*U(K,I,J)+l.185*A(K,I,J)/4.0 WRITE(5,50) K,I,J,U(K+1,I,J) FORMAT(I3,5X,13,5X,I3,5X,E15.6) CONTINUE CONTINUE CONTINUE CALL SRCH$$(4,'OUTPUT',6,1,N1,N2) CALL EXIT END 80 APPENDIX E: CONVECTION PROGRAM PROGRAM THERMAL(INPUT,OUTPUT,TAPE2=INPUT,TAPE3=0UTPUT) COMMON A(7,21,26),8(7,21,26),C(7,21,26),D(7,21,26),T(7,21,26),R( $7,21),J1 READ(2 :* ) ((R(1 :3) 93:1,11), 1:154) D0 5 J=1.11 D0 10 I=1,4 10 R(8-I,J)=R(I,J) DO 15 I=1,7 15 R(I,22-J)=R(I,J) 5 CONTINUE WRITE(3,78) ((I,J,R(I,J),I=l,7),J=l,21) 78 FORMAT(212,E15.6) DO 20 I=1,7 DO 25 J=1,21 25 T(I,J,1)=23.0 20 CONTINUE WRITE(3,30) 30 FORMAT<20X, * CONVECTION PROGRAM *) J1=1 lOO J1=Jl+l CALL EVEN CALL TDMAEV CALL PRINTT IF(J1.EQ.26)GOTO 110 J1=Jl+1 CALL ODD CALL TDMAODD CALL PRINTT IF(J1.LT.26)GOTO 100 110 STOP END 50 45 33 37 36 48 81 SUBROUTINE EVEN COMMON A(7,21,26),B(7,21,26),C(7,21,26),D(7,21,26),T(7,21, $26),R(7,21),J1 DO 45 I=2,5 DO 50 J=1,21 B(I,J,Jl)=R(I,J)+2.0 A(I,J,Jl)=-l.0 C(I,J,Jl)=-l.0 CONTINUE CONTINUE DO 33 J=1,21 B(6,J,J1)=R(6,J)+2.0-2.0/(R(7,J)+2.0) B(1,J,J1)=R(1,J)+2.0 A(1,J,J1)=0.0 C(6,J,J1)=0.0 C(1,J,Jl)=-2.0 A(6,J,Jl)=-l.0 DO 36 I=1,7 DO 37 J=2,20 D(I,J,J1)=T(I,J—l,Jl-l)+T(I,J+l,J1-1)+(R(I,J)-2.0)*T(I,J,J $1-l) D(I,1,Jl)=2.0*T(I,2,Jl-l)+R(I,l)-2.0)*T(I,l,Jl-1) D(I,21,Jl)=2.*T(I,20,J1-l)+(R(I,21)-2.0)*T(I,21,J1-l) DO 48 J=1,21 D(6,J,Jl)=D(6,J,J1)+(D(7,J,Jl)+382.55)/(R(7,J)+2.0) RETURN END 23 82 SUBROUTINE TDMAEV COMMON A(7,21,26),B(7,21,26),C(7,21,26),D(7,21,26),T(7,21, $26),R(7,21),J1 DO 23 J=1,21 C(1,J,Jl)=C(1,J,Jl)/B(l,J,Jl) D(l,J,Jl)=D(l,J,J1)/B(1,J,Jl) D0 3 I=2,6 E=l./(B(I,J,Jl)-A(I,J,Jl)*C(I-1,J,Jl)) D(I,J,Jl)=E*(D(I,J,J1)-A(I,J,Jl)*D(I—l,J,Jl)) C(I,J,Jl)=E*C(I,J,Jl) T(6,J,Jl)=D(6,J,Jl) DO 4 I=2,6 MI=7-I T(MI,J,J1)=D(MI,J,Jl)-C(MI,J,J1)*T(MI+1,J,J1) T(7,J,Jl)=(2.*T(6,J,Jl)+382.55+D(7,J,J1))/(R(7,J)+2.0) CONTINUE RETURN END 11 12 14 13 21 83 SUBROUTINE PRINTT COMMON A(7,21,26),B(7,21,26),C(7,21,26),D(7,21,26),T(7,21, $26),R(7,21),J1 WRITE(3,11)J1 FORMAT(2X,2HK=,I2) D0 21 I=1,7 WRITE(3,12)I FORMAT(2X,2HMs,12) D0 13 J=1,21,6 WRITE(3,14)T(I,J,J1),T(I,J+1,J1),T(I,J+2,J1),T(I,J+3,JI),T $(I,J+4,Jl),T(I,J+S,J1) FORMAT(5X,6E15.7) CONTINUE CONTINUE RETURN END 17 11 84 SUBROUTINE ODD COMMON A(7,21,26),B(7,21,26),C(7,21,26),D(7,21,26),T(7,21, $26),R(7,21),J1 Do 7 I=1,7 DO 8 J=2,l9 B(I,J,Jl)=R(I,J)+2.0 A(I,J,Jl)=-l.0 C(I,J,J1)=-1.0 CONTINUE B(I,1,J1)=R(I,1)+2.0 B(I,20,Jl)=R(I,20)+2.0-2.0/(R(I,21)+2.0) A(I,20,J1)=-1.0 A(I,l,Jl)=-0.0 C(I,l,Jl)=-2.0 C(I,20,Jl)=0.0 DO 17 J=1,21 D0 9 I=2,6 D(I,J,Jl)=T(I-l,J,Jl-1)+T(I+1,J,J1-l)+(R(I,J)-2.0)*T(I,J,Jl-l) D(1,J,Jl)=2.*T(2,J,Jl-l)+(Rl,J)—2.)*T(1,J,Jl-1) D(7,J,J1)=2.*T(6,J,J1-l)+(R(7,J)-2.0)*T(7,J,Jl-1)+382.55 DO 11 I=1,7 D(I,20,J1)=D(I,20,J1)+D(I,21,J1)/(R(I,21)+2.) RETURN END 66 67 65 85 SUBROUTINE TDMAODD COMMON A(7,21,26),B(7,21,26),C(7,21,26),D(7,21,26),T(7,21, $26),R(7,21),J1 DO 65 I=1,7 C(I,l,J1)=C(i,l,Jl)/B(I,1,Jl) D(I,l,Jl)=D(I,1,Jl)/B(I,1,Jl) Do 66 J=2,20 E=l.0/(B(I,J,Jl)-A(I,J,Jl)*C(I,J-1,Jl)) C(I,J,Jl)=E*C(I,J,J1) D(I,J,Jl)=E*(D(I,J,Jl)-A(I,J,J1)*D(I,J—1,Jl)) T(I,20,Jl)=D(I,20,Jl) DO 67 J=2,2o NJ=21 -J T(I,NJ,J1)=D(I,NJ,J1)-C(I,NJ,J1)*T(I,NJ+1,J1) T(I,21,J1)=(2.*T(I,20,J1)+D(I,21,JI))/(R(I,21)+2.0) CONTINUE RETURN END 86 APPENDIX F: NONfDIMENSIONALIZATION OF MOMENTUM EQUATION FOR FULLY DEVELOPED FLOW We are concerned with the solution of the convec- tion problem through a rectangular duct having dimensions stated in the previous chapter. The flow through the rectangular duct is partially developing hydrodynamically and partially developed after a certain length of the tube. These aspects have been discussed in the previous chapter. We are interested in examining the thermal conditions on the quartz window heater of the convec- tion cryomicroscope system. In the present study, the actual complex system has been simplified in order to develop a preliminary understanding of the problem. The window is situated at a distance beyond the critical length for the flow being fully developed so the flow in the region of our interest is hydrodynamically developed. In order to achieve the temperature profile along the surface of the heater window, we must determine the hydro- dynamic characteristics of the channel flow velocity pro- file. After obtaining the velocity profile, we can use the hydrodynamic velocity information in the energy equa- tion which will be discussed in the next chapter. To achieve the hydrodynamic solution, we begin by stating the Navier Stokes equations which are: x direction: uau Bu Bu - l ___+ _..+ ___=—.— 8x V 8y 82 p y direction: qu 3v 3v _ 1 TF+V87+W55“ 6 z direction: C. I + < I II I ‘OIH 87 2 2 2 3P 3 u 8 u 3 u ——-+ u(——§-+ ——— + ——— (Fl) 3x 8x Byz 322 2 2 9.2+6(__3‘2’+__3‘2’+——3‘2’) (F2) 3y 3x 3y 32 8 32 32 82 _E + v(——;-+ ——§-+ -—§9 (F3) 3 3x 3y 8z In the above equations, we have made the following assumptions: (1) the flow is steady (2% 0 where ¢ = velocity in given direction) (ii) three dimensional flow (iii) incompressible flow (0 = CODSt-) (iv) constant property flow We are interested in finding the fully developed velo- city profile, so an O O 0 v = 0, w = 0 and §§'= 0 in each direction in (F1, F2, F3) So 2 2 2.§+i_121_l22 By 32 )1 ax . 3 X dir : SE = 0 BY . 3 z dir : AB = 0 az So from equations F4, F5 .28 = constant ( . n . , x dir momentum equation becomes (F4) (F5) (F6) and F6, we get U = f(y,z) but u 75 f(x)) 88 So p = f(x) only. Rewriting the equation (4) we get 2 2 3 + 3 = By 32 C1 (F7) N .I. tlH can. 9):: The x direction momentum equation is non-dimensionalized by defining U+ = 9+ where U _ :1_dp 2 l n 2 = U* — u ax a (_N/MZX N——/M3 x m ) (m/sec) + = Z. T y a (m) + — E. E z — a (m) where a = width of the rectangular cross section as shown in Figure 3. The x direction equation is non- dimensionalized in order to generalize the solution which we obtain for various aspect ratios and various pressure drop magnitudes. The non-dimensionalized solu- tion can be used to get various information about the effects of pressure drops and dimensions of the duct on the velocity magnitude at given locations in the duct. So converting the coordinate system Y and Z into non-dimensionalized coordinates y+ and 2+, the width of the duct equals unity (the length of rectangular duct equals (b/a). Now converting equation F7 into a non- dimensionalized form, 2: 82U+ + (U*) 32U+ % QE 2) 3‘y+z g3. 8y+2 - dx a 89 2 + 2 + 2 3 U 3 U _ l_dp a 3'72” ___—ay+2 ‘ u UHF" (F3) 0 I * _ 1 d 2 0 Substituting the value of U — - 313% a on the right hand side of the above equation (8), we get 2 + 2 8 U + 8 U ay+2 32+2 — (F9) REFERENCES (l) (2) (3) (4) (5) (6) (7) (8) (9) (10) (ll) (12) (13) REFERENCES Diller, K.R., and Cravalho, E.G., "A Cryomicroscope for the study of Freezing and Thawing Processes in Biological Cells,"Cryobiology,7,(1970):l9l. Hrycaj, Thomas, "The Effects of Freezing on the Micro- circulation of the Hamster Cheek Pouch" (M.S. thesis, Massachusetts Institute of Technology, Patankar, Suhas V., Numerical Heat Transfer and Fluid Flow, McGraw-Hill, New York, 1980. Smith, G.D., "Numerical Solution of Partial Differential Equations)"Eiei:s-9£§ezsass-!s&heél an ed- Hornbeck,'Robert W., "Numerical Marching Techniques for Fluid Flows with Heat Transfer" National Aero- nautics and Space Administration, NASA SP-297, (1973):241. Roache: Patric J-v QQTPBESEEQEEE-ElEiQ_PY£§TiS§r Hermosa Publishers,:AJbuguergue, N.M.,1972. Yonenko, Springerkerlag; Shovery, M.M., and Gaane, R. R-: 399£E§l_9§-EQTBi-EPY§iS§v 23: (1977>=242- Kays and Crawford: QQEYSSEEY§_§§§E_§P§-E§§§_IEEP§ESEr 2nd Edition McGraw-Hill, 1980. Shah, R.K., and London,A.L., "Laminar Flow Forced Convection in Ducts)" in ASYEPS§§_iP-§§§E_T£§E§§S£r Academic Press, New York 1978; Knudsen: James: and Katzv El229_PYP§TES§-229-§§§E-T£§E§ESEv McGraw-Hill, New York, 1958. Roshsenow, Warren M., and Choi, Henry, HeatLM§§§_§nd Mgmentum_Tran§f§§, Prentice-Hall; Englewood Cliffs, N.J.:l96l. Hornbeck, Robert, "Marching Solutions for Fluid Flows .with Heat Transfer" (NASA Research Center), (1973), p.241. Preingerova, Vera, "Laminar Entry Length Heat Transfer in Ducts of Rectangular Cross Section with Boundary Conditions of the Second Kind," in Heat_Tr§n§fer Source Book (Minsk: Fifth All Union Conference, 1976), p. 74. 90 (14) (15) (l6) (17) (18) 91 Preingerova, Vera, and Allen, P.H., "Laminar Flow Entry Length Heat Transfer with Varying Physical Properties in Simple and Complex Duct Geometries" (paper presented at the Heat Transfer Conference, Tokyo, 1974) (NC 5-4). Sastri,V.M., and Chandrapatla, Ashok, "Laminar Flow and Heat Transfer to a Non-Newtonian Fluid in an Entrance Region of a Square Duct with Presented Constant Axial Wall Heat Flux," in Numerical Heat Transfer, l,(l979):243. Cess and Shaffer, "Heat Transfer to Laminar Flow Between Parallel Plates with Prescribed Wall Heat Flux," Montgomery,S.R., and Wibulsatas, P., "Laminar Flow Heat Transfer for Simultaneously Developing Velocity and Temperature Profiles in Ducts of Research Journal, 18, (1967): 247. -Lenert, Louis W., Advanced Technical Mathematics, Merril, Columbus, Ohio, (1970): 327. MICHIGAN STATE UNIV. LIBRARIES IIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIII 31293009962568