ICHIG NIVERSITY LIBRARIES Illllllllllll l llllll“ Ill ill 3 1293 00897 5074 This is to certify that the thesis entitled ANALYSIS OF THE EFFECT OF AXIAL CONDUCTION ON THE PERFORMANCE OF A COUNTERFLOW HEAT EXCHANGER presented by KEVIN J. DOWDING has been accepted towards fulfillment of the requirements for MASTERS MECHANICAL ENGINEERING degree in Major professor Date 3/3/ M3 0-7639 MS U is an Affirmative Action/Equal Opportunity Institution LIBRARY Mlchlgan State University PLACE IN RETURN BOX to remove this checkout from your record. TO AVOID FINES return on or before date due. DATE DUE DATE DUE DATE DUE l ‘* JE ‘1 % Uffil i MSU Is An Affirmative Action/Equal Opportunity Institution . c:\cIIcIdI-dm.m3-p.t ANALYSIS OF THE EFFECT OF AXIAL CONDUCTION ON THE PERFORMANCE OF A COUNTERFLOW HEAT EXCHAN GER By Kevin J Dowding A THESIS Submitted to Michigan State University in partial fulfillment of the requirements for the degree of MASTER OF SCIENCE Department of Mechanical Engineering 1993 “l An in. Within the “'4' exact solutio: While incluch Within the w; The efl exchange, W; 65m The re to 45% under Sinall ratio of WWII); ratio exchanger de} PICSeme of a; Pannieter Cal Ction was Cond‘lcfion, u ABSTRACT ANALYSIS OF THE EFFECT OF AXIAL CONDUCTION ON THE PERFORMANCE OF A COUNTERFLOW HEAT EXCHAN GER By Kevin J Dowding An investigation into the effect of axial (longitudinal) conduction within the wall of a counterflow heat exchanger has been completed. The exact solution was obtained for the temperature fields of the wall and fluids while including the possibility of heat flow in the direction of the fluid flow within the wall of the heat exchanger, so-called axial conduction. The effect of axial conduction on the performance of the heat exchanger was quantified by comparison to a solution that neglects this effect. The results showed that the effectiveness could be over-estimated up to 45% under certain conditions. These conditions were a small Biot number, small ratio of the heat exchanger wall length to thickness, and large flow heat capacity ratios. The severity of the effect on the performance of the heat exchanger depended on the magnitudes of these conditions. Predicting the presence of axial conduction was shown to depend on a dimensionless parameter called the Mondt number. For magnitudes less than .01 , axial con- duction was negligible; while magnitudes greater displayed a nonzero axial conduction, which grew larger as the Mondt number increased. to my parents, Roger and Mary Ann Dowding iii Ian his guidan gained in t sions cone: The Beck and 5 live. Addju‘ was comple The I Continued gt mom and F Final my Educafio aPPTCCiated collid have I] undersmndirl ACKNOWLEDGMENTS I am grateful for the opportunity to work with Professor Somerton; for his guidance and encouragement conducting this research and the experience gained in the undergraduate teaching laboratory. Also, his open, frank discus- sions concerning professorial life were enlightening and appreciated. The commitment of the remaining committee members, Professors Beck and Schock, is also recognized. Their comments and input were instruc- tive. Additional thanks to Professor Beck for his patience while this research was completed. The Department of Mechanical Engineering is acknowledged for the continued generous support. Special thanks are extended to Professors Som- erton and F oss for the involvement in their teaching laboratories. Finally, I would like to thank my family, who have endured through my education and given love and encouragement. My Mother and Dad are appreciated for their endless unassuming support, without which this point could have never been reached. While I am eternally grateful for the patience, understanding, and support of my wife Tracey. iv list of Ta List of Fig Nomenclc 1 Introdt 1.0 I 2 Dacrit 2.0 1 TAB LE OF CONTENTS List of Tables ............................................................................. ix List of Figures ............................................................................ x Nomenclature .......................................................................... xiv Chapter 1 Introduction l 1.0 Heat Exchanger Analysis .......................................................... 4 1.0.1 Basic Analysis .............................................................................. 4 1.0.2 Dimensionless Analysis ............................................................... 6 1.1 Literature Review ..................................................................... 9 2 Describing Equations and Boundary Conditions 11 2.0 Dimensional Equations ........................................................... 11 2.0.1 Wall Conduction Equation ......................................................... 11 2.0.2 Fluid Energy Balance Equations ................................................ 12 2.1 Dimensionless Equations ........................................................ 14 2.1.1 Nondimensionalization of Lengths ............................................ l4 2.2 2.3 3 Metht 3.0 3.1 3.3 2.1.2 Nondimensionalization of Temperature ...... - - - -- ...15 2.2 Dimensionless Parameters ...................................................... 16 2.3 Heat Exchanger Performance Analysis .................................. 18 3 Method of Solution 20 3.0 General Method ...................................................................... 20 3.1 General Solution of the Wall Conduction Equation .............. 22 3.1.1 General Solution for Nonhomogeneous Boundary Condition atx = L - - 23 3.1.2 General Solution for Nonhomogeneous Boundary Condition on y .......................................... 27 3.2 Solution of the Fluid Energy Balance Equations .................... 28 3.2.1 Hot Fluid Formulation ............................................................... 28 3.2.2 Hot Fluid Solution - - - - -- - ..................................... 29 3.2.3 Cold Fluid Solution .................................................................... 31 3.3 Application of Nonhomogeneous Boundary Conditions ........................................................................ 33 3.3.1 Application of Boundary Condition at x = L ........................... 34 3.3.2 Application of Boundary Condition at y = 0 ........................... 36 3.3.3 Application of Boundary Condition at y = 8 ........................... 38 3.3.4 Summary of Applying Orthogonality ........................................ 40 3.4 Computer Program AXCOND .................................................. 44 4 Raul 4.0 4.1 I 4.2 4 Results and Drscussron47 4.0 Identification of Input Parameters .......................................... 47 4.0.1 Presentation of Results -- ..................................... 50 4.1 Influence of Axial Conduction at a Constant NTU ................. 52 4.1.1 Influence of the Wall Aspect Ratio ............................................ 52 4.1.2 Influence of the Wall End Biot Numbers ................................... 63 4.1.3 Influence of the Fluid Biot Numbers ......................................... 69 4.1.4 Summary of Investigation with a Constant NT U ..................... 74 4.2 Characterization of Axial Conduction as a Function of NTU .............................................................................. 77 4.2.1 Ineffectiveness-NTU Results .................................................... 83 4.2.1.1 Heat Capacity Ratio of 1 ............................................. 83 4.2.1.2 Heat Capacity Ratio of .75 .......................................... 85 4.2.1.3 Heat Capacity Ratio of .50 .......................................... 87 4.2.1.4 Heat Capacity Ratio of .25 .......................................... 89 4.2.1.5 Discussion of Results .................................................. 91 4.3 Predicting the Presence of Axial Conduction with the Mondt number .................................................................. 96 4.4 Application of Results .......................................................... 102 4.5 Comparison to Published Results ......................................... 108 5 Conclu List of Re Appendic App Appt App< Appt Appc App: ADpe Appe 5 Conclusions and Recommendations for Future Work .. 112 List of References ................................................................... 114 Appendices ............................................................................. 116 Appendix A.1: General Solution of the Wall Conduction Equation ............................ 116 Appendix B. 1: Solution of Hot Fluid Energy Balance ................................................ 125 Appendix B.2: Solution of Cold Fluid Energy Balance ................................................ 130 Appendix B.3: Summary of Energy Balance Solutions ............................................... 136 Appendix C. 1: Application of N onhomogeneous Boundary Conditions ............................ 138 Appendix C.2: Summary of Applying Nonhomogeneous Boundary Conditions ............................................ 156 Appendix D.1: Solution for Constants ........................... 157 Appendix E.1: Derivation of Effectiveness-NTU Relationship for a Counterflow Heat Exchanger Neglecting Axial Conduction ........................................... 163 Appendix E.2: Solution of Fluid Temperatures for a Counterflow Heat Exchanger Neglecting Axial Conduction ............... 167 Appendix F: Listing of Computer Program AXCOND ............................................. 17 3 TABLE 4.3 TABLE 4.: TABLE 4.3 TABLE 4.4 TABLE 4.5 LIST OF TABLES TABLE 4.1: Input variables and dependencies for determining heat exchanger performance 49 TABLE 4.2: Relationship between the slepe of the hot and cold fluids temperature profile and heat capacity ratio 55 TABLE 4.3: Ineffectiveness as a function of the fluid Biot numbers from Figure 4.14 69 TABLE 4.4: Classification of effects seen as the NTU varied 95 TABLE 4.5: Relationship among the Mondt number, wall aspect ratio, NTU, and fluid Biot numbers 96 TABLE 4.6: Magnitude of dividing Mondt number as a function of the heat capacity ratio 99 TABLE 4.7: Operating conditions for axial conduction to exist 105 FIGURE FIGURE 1 FIGURE FIGURE FIGURE . FIGURE 1 A? '1 l noun-:2 t. i FIGURE 4 FIGURE 4. FIGURE 4.. FIGURE 4.4 FIGURE 4. FIGURE 4.61 LIST OF FIGURES FIGURE 1.1: Possible heat flow paths for wall of a heat exchanger FIGURE 1.2: Generic heat exchanger, nonspecific design exchang- ing energy between two fluids FIGURE 1.3: Thermal circuit for a heat exchanger wall neglecting any fouling FIGURE 2.1: Heat exchanger geometry studied FIGURE 2.2: Thermal circuit for heat exchanger wall FIGURE 3.1: Flowchart of solution procedure FIGURE 3.2: Flowchart of FORTRAN program AXCOND used to compute solution FIGURE 4.1: Ineffectiveness as a function of the wall aspect ratio forNTU = 7 and CR =1 FIGURE 4.2: Ineffectiveness as a function of the wall aspect ratio for NTU = 7 and CR = .75 FIGURE 4.3: Ineffectiveness as a function of the wall aspect ratio for NTU = 7 and CR = .50 FIGURE 4.4: Ineffectiveness as a function of the wall aspect ratio for NTU = 7 and CR = .25 FIGURE 4.5: Median dimensionless wall temperature as a function of dimensionless position and Biot number for NTU= 7,L =100,andCR =1 FIGURE 4.6: Ineffectiveness as a function of the wall aspect ratio and heat capacity ratio for NTU = 7 and all Biot num- bers equal .001 11 19 21 46 59 59 6O 61 61 FIGURE 4. FIGURE 4. FIGURE 4. FIGURE 4. FIGURE 4. FIGURE 4. FIGURE 4.- FIGURE 4.7: Median dimensionless wall temperature as a function of dimensionless position and heat capacity ratio for NTU = 7, L = 100, and all Biot numbers .001 FIGURE 4.8: Median dimensionless wall and fluid temperature as a function of dimensionless position and heat capacity ratio for NTU = 7, L = 100, and all Biot numbers equal .001 FIGURE 4.9: Ineffectiveness as a function of the end Biot numbers forNTU = 7, L =100,and CR =1 FIGURE 4.10: Ineffectiveness‘as a function of the end Biot numbers for NTU = 7, L = 100, and CR = .75 FIGURE 4.11: Ineffectivenessps a function of the end Biot numbers for NTU = 7, L = 100, and CR = .50 FIGURE 4.12: Ineffectiveness‘as a function of the end Biot numbers for NTU = 7, L = 100, and CR = .25 FIGURE 4.13: Ineffectiveness as a function of the end Biot numbers and heat capacity ratio for NTU = 7, L' = 100, and fluid Biot numbers equal .001 FIGURE 4.14: Ineffectiveness as a function of the fluid Biot numbers forNTU = 7,L =100,andCR =1 FIGURE 4.15: Ineffectiveness as a function of the fluid Biot numbers for NTU = 7, L = 100, and CR = .75 FIGURE 4.16: Ineffectiveness as a function of the fluid Biot numbers for NTU = 7, L = 100, and CR = .50 FIGURE 4.17: Ineffectiveness as a function of the fluid Biot numbers for NTU = 7, L = 100, and C, = .25 FIGURE 4.18: Ineffectiveness as a function of the hot fluid Biot number and heat capacity ratio for NTU = 7, L = 100, and a cold fluid Biot number equal .001 FIGURE 4.19: Ineffectiveness as a function of the NTU, wall aspect ratio, and maximum fluid Biot number for Bil“! = 0.”)1 and CR = 62 62 66 67 67 68 68 72 72 73 73 74 78 FIGUI FIGL'R FIGL'R FIGURE 4.20: Ineffectiveness as a function of the NTU and heat capacity ratio for Bi”, = 0.0001, Bi”, = 0.001, and L = 100 FIGURE 4.21: Ineffectiveness as a function of the NTU and wall aspect ratio for Bin," = 0.0001 and CR = 1 FIGURE 4.22: Ineffectiveness as a function of the NTU and wall aspect ratio for Bin“, = 0.001 and CR = 1 FIGURE 4.23: Ineffectiveness as a function of the NTU and wall aspect ratio for Bin," = 0.01 and CR = 1 FIGURE 4.24: Ineffectiveness as a function of the NTU and wall aspect ratio for Him,” = 0.1 and CR = 1 FIGURE 4.25: Ineffectiveness as a function of the NTU and wall aspect ratio for Bin,” = 0.0001 and CR = .75 FIGURE 4.26: Ineffectiveness as a function of the NTU and wall aspect ratio for Bin..." = 0.001 and CR = .75 FIGURE 4.27: Ineffectiveness as a function of the NTU and wall aspect ratio for Bi,“ = 0.01 and CR = .75 FIGURE 4.28: Ineffectiveness as a function of the NTU and wall aspect ratio for Bi = 0.1 and CR = .75 min FIGURE 4.29: Ineffectiveness as a function of the NTU and wall aspect ratio for Bin," = 0.0001 and CR = .50 FIGURE 4.30: Ineffectiveness as a function of the NTU and wall aspect ratio for Bin..." = 0.001 and CR = .50 FIGURE 4.31: Ineffectiveness as a function of the NTU and wall aspect ratio for Bi”, = 0.01 and CR = .50 FIGURE 4.32: Ineffectiveness as a function of the NTU and wall aspect ratio for Rim." = 0.1 and CR = .50 FIGURE 4.33: Ineffectiveness as a function of the NTU and wall aspect ratio for Bin,” = 0.0001 and CR = .25 81 83 83 84 84 85 85 86 86 87 87 88 88 89 F'IGUR FIGUR' FIGUR] FIGURE FIGURE FIGURE 1 FIGURE 4 FIGURE 4.34: FIGURE 4.35: FIGURE 4.36: FIGURE 4.37: FIGURE 4.38: FIGURE 4.39: FIGURE 4.40: FIGURE 4.41: FIGURE 4.42: Ineffectiveness as a function of the NTU and wall aspect ratio for Bi,“ = 0.001 and CR = .25 Ineffectiveness as a function of the NTU and wall aspect ratio for Bin,” = 0.01 and CR = .25 Ineffectiveness as a function of the NTU and wall aspect ratio for Bin," = 0.1 and CR = .25 Predicted region affected by axial conduction in terms of the NTU, Mondt number, wall aspect ratio and minimum fluid Biot number Magnitude of pertinent heat exchanger operating conditions as a function of wall thermal conductivity, which is needed for axial conduction to exist Comparison of present results to past result of Pan, Welch, and Head [15] and Rohsenow [17]. Ineffec- tiveness as a function of the NTU and wall aspect ratio for Bin," = 0.0001 and CR = 1. Comparison of present results to past result of Pan, Welch, and Head [15] and Rohsenow [17]. Ineffec- tiveness as a function of the NTU and wall aspect ratio for Bin..." = 0.001 and CR = 1 Comparison of present results to past result of Pan, Welch, and Head [15] and Rohsenow [l7]. Ineffec- tiveness as a function of the NTU and wall aspect ratio for Bin," = 0.01 and CR = 1 Comparison of present results to past result of Pan, Welch, and Head [15] and Rohsenow [17]. Ineffec- tiveness as a function of the NTU and wall aspect ratio for Bin," = 0.1 and CR = 1 89 98 104 110 110 111 111 Arabic 0‘ Arabic A An’ Bn’ Cn Bi m- 63.9 . tau-a» 1~1-‘ 23- N TU Q.~Q"u a Q‘l" E NOMENCLATURE area [m2] constants Biot number specific heat [Pi—IE] heat capacity [WI K] enthalpy [J] specific enthalpy [— lg] convective heat transfer coefficient [n ] thermal conductivity [in—K] length of heat exchanger [m] wall aspect ratio (L’ = 55') mass flow rate [s—fic— 3] number of termssec in series solutions number of transfer units perimeter [m] heat [W] heat flux [—] thermal resistance [— W] time [sec ] temperature [C or K] thermal conductance [if—IE ] heat exchanger wall windth [m] xiv Ly Greek min ”101 out Superscr- x,y Greek comm 1:: @@ Subscripts C H in L min max out Superscripts + Cartesian coordinates [m] effectiveness wall thickness [m] kronecker delta function scaled temperature [C or K] dimensionless temperature eigenvalues cold fluid wall side hot fluid wall side inlet 1: = Lorx” =1wallend minimum maximum x = 0 wall end outlet ratio wall dimensionless variable XV He: fluids. The cessin g, ht ognizable exchanger tion. There characteris form of a h exchanger. moving fluj moving in t site directio Desi g heat transfer energy need: the heat excl beat u”ansfer "met yet 11 Obtained by e; traIISfex' CODVE all IIICreaSe in {milking in an Chapter 1 Introduction Heat exchangers provide for the transfer of heat between two moving fluids. These devices are used in power generation, chemical and food pro- cessing, heating and air conditioning, and motor vehicles. It is the most rec- ognizable heat transfer device and one of the most widely used. Heat exchangers are classified based on flow arrangement and type of construc- tion. There are many types of heat exchanger designs, each type with its own characteristics that make it suitable for a particular application. The simplest form of a heat exchanger would be a double pipe (or concentric tube) heat exchanger. The construction of this type of heat exchanger consists of two moving fluids separated via a wall parallel to the fluid motion; with the fluids moving in the same direction for a parallel flow heat exchanger and in oppo- site directions in a counterflow heat exchanger. Design of heat exchangers in general, requires consideration of the heat transfer occurring between the two fluids in addition to the mechanical energy needed to overcome the frictional forces to move the fluids through the heat exchanger. These two design criteria can be generally classified as heat transfer and pressure drop. It is typically desired to achieve large heat transfer yet maintain a small pressure drop. Large heat transfer rates can be obtained by either having a larger heat transfer area or having a larger heat transfer convection coefficient. Unfortunately, both of these conditions cause an increase in the pressure dr0p, since a larger area gives more frictional area resulting in an increase in the pressure drop and the larger flow rate to increase the convection coefficient would likewise increase the pressure drop. There is a trade-off in these two design criteria; a beneficial gain in one crite- ria is usually at the expense of the other criteria and a compromise must be established. Focusing on the heat transfer, consider a simplistic heat exchanger that has two fluid streams separated by an infinitesimally thin wall. The analysis Cold Fluid 4c (1) Heat Exchanger Wall 9110‘) Hot Fluid Figure 1.1. Podble heat flow paths for the wall on heat exchanger of this heat exchanger could be accomplished by considering the energy bal- ances of the two fluids since the wall has negligible thickness. In reality, the wall will be required to have some thickness to perform its function of sepa- rating the two fluids. As the thickness of the wall becomes nonnegligible, some resistance to the flow of heat between the two fluids is added and now the wall must be considered in the analysis. In addition, an alternate path for the heat to flow is provided. Figure 1.1 shows the heat flow paths: q” is the heat convected from the hot fluid to the wall, q Ma, is the heat conducted along the wall, and qc is the heat conducted to the cold fluid from the wall. When the wall is very thin, it is unlikely that heat transferred to the wall from the fluid will flow parallel to the wall and the energy exchanged will balance locally (qH (x) = qc (x) ). This may not be the case when the wall thickness is increased; some of the heat may flow along the wall. This effect is call axial (longitudinal) conduction, and when axial conduction is nonzero(qAx,.a, = 0) the fluid streams may have a local energy imbalance (q H (x) at qc (x) ). The energy is conserved, however. Axial conduction of heat is typically neglected in the analysis of the heat transfer in a heat exchanger. This assumption may be appropriate for most cases, but under certain circumstances the effect of axial conduction can become important. These circumstances depend on more than the wall thick- ness, which was used to introduce the efi‘ect. In general, axial conduction depends on the material and physical dimensions of the wall and the proper- ties and fit be the focr quantified ble identifi heat excha The heat exchai tions and b two, while three. Chap sponding d in chapter f ties and flow rate of the fluids. Assessing the effect of axial conduction will be the focus of this thesis. Specifically, the effect of axial conduction will be quantified and the conditions of the wall and fluids for which it is nonnegligi- ble identified. Also, the adverse effect on the performance of a counterflow heat exchanger will be determined. The remaining sections of this chapter will review the analysis of a heat exchanger and conclude with a literature review. The describing equa- tions and boundary conditions for the problem will be addressed in chapter two, while the methods employed to solve the problem are shown in chapter three. Chapter four will present the results of the investigation with corre- sponding discussion. A summary and the resulting conclusions will be given in chapter five, in addition to recommendations for future work ci 116 1.0 Heat Exchanger Analysis 1.0.1 Basic Analysis Figure 1.2 shows a heat exchanger that transfers energy between two moving fluids through a wall; the geometry and construction of this heat exchanger may be considered arbitrary. A general thermal analysis of this heat exchanger will be performed. First, taking the heat exchanger as a con- trol volume and applying an overall energy balance, assuming no interactions (work or heat) with the surroundings and steady state, results in an enthalpy balance (it) HI» = Hm (1.1) which can be written in terms of the inlet and exiting conditions using spe— cific enthalpy and mass flow rate as mufiH|In+mcficlin = mflfiulout+mcficlout (1'2) and rearranged to group the fluids ) = martial - incl“) (1.3) Assuming the fluid behaves as an incompressible liquid or an ideal gas and a negligible pressure change, the change in enthalpy can be expressed as "In (511", - 511' 0 I“ 0“! din = deT (1.4) Hot Fluid «1 Cold I‘Iuld figmlJ.Gmeflehuuchng«,noMpuEededgnuchngtngenugybemnthetwoflutds 5 which can be approximated by differences(dlr ~ 712 - 71, ). Using the approxi- mate form of equation (1.4) in equation (1.3), the results for the overall energy balance are q = mHCp, H (TH, in - TH, out) = mCCp, C (TC, out — TC, in) (1'5) Introducing the heat capacity Ca Inc, (1.6) equation (1.5) can be rewritten as q = CH (TH, in ' T1101“) = Cc (Team " Tc, in) (1:7) In general, the heat transfer between the two fluid streams is a function of the following six parameters: (1.8) * q = f(e my» "'0 TH,in’ Tarn, . Teen: TILoutl ) given unknown where the first four are typically given design parameters and the last two are desired results. Thus, in order to obtain a solution more information is needed. The additional information will come from the heat transfer analysis of the wall separating the two fluids. A circuit describing the thermal communication of the two fluids is shown in Figure 1.3. This circuit shows the resistance that impedes the heat flow, neglecting any fouling on the heat exchanger wall. The total resistance of the series circuit is the sum of the individual resistances T T ” e—va w W—e C Figure 1.3. Thermal circuit for a heat exchanger wall neglecting any fouling R Also, the overall conductance of the heat exchanger wall, that will give the total amount of heat transferred if the temperature difference is known, is m = Ru+RW+Rc (1.9) UA = J— (1.10) . tot These terms are analogous to an electrical circuit, with the temperature differ- ence as the voltage potential and the heat flow as the current. The heat trans- fer over a differential element of length (11: is dq = UP [TH (x) - Tc (x) ] dx (1.11) The heat transfer over the entire length of the heat exchanger is the sum of these differentials elements over the total length L L q = jdq = UP] [7,,(x) -Tc(x)] dx (1.12) 0 0 The product UP was assumed constant and thus could be taken out of the integral, but (TH - Tc) depends on x and cannot be removed from the inte- gral. Because this integral is not known in general, we will define the mean temperature difference as L (AT)... = H [1,,(x) -Tc(x)]dx (1.13) 0 Substituting equation (1.13) into equation (1.12), the total heat transfer can be written as q = UA (A1)" (1.14) For a double pipe heat exchanger, equation (1.13) can be evaluated. The result is the log mean temperature difference. However, for geometries that are more complicated, evaluating equation (1.13) is difficult if not impossible. This suggests that another approach may be necessary to remain a general analysis. 1.0.2 Dimensionless Analysis To introduce another approach, consider the parameters that the mean temperature difference is a function of It depends on the inlet temperatures and heat capacities of the fluids and the wall conductance. The number of independent parameters can be reduced by considering a second function times the temperature difference at the inlet (A1)", = 3 (CC, C”, UA) (TH, .1.- TC, ,n) (1.16) substituting equation (1.16) into equation (1.14) the functional dependence for the heat transfer is given by ‘1 = UAg (Co CH: UA) (711.111 " Tarn) (1.17) As typical in heat transfer, scaling will be introduced to provide dimension- less parameters. Define the maximum possible heat transfer as qrnax = Cmin (TH, in - TC, in) (1°18) where Cm." = min (CH, CC) . Dividing equation (1.17) by equation (1.18) _‘1 = ”A 3 (CC, C”, 0.4) (1.19) anax Carin and introducing another function that is in terms of the scaled independent variables gives = g (C , _) (1.20) qrnax Cmin R Cmin where C . C = "'"' (1.21) R Curax Finally, noting that equation (1.20) can be written for a final function as UA C ruin This demonstrates that the performance of a heat exchanger can be expressed in terms of three dimensionless variables. The first was given previously in equation (1.21) as the ratio of the heat capacities. The second is the effective- ness i=h(CR, anax ) (1.22) q q a = —_ = ._ — — (1.23) gm Carin (TH, in - TC, in) ,. [:u.‘ r"~ - 3"? which is the 1 maximum p0 last dimensio which is a rat mum fluids 2 The pr: moderately si analysis of th. sider a double describing the in Appendix I for the case of which is the ratio of the actual heat transfer given by equation (1.7) to the maximum possible heat transfer. The number of transfer units (NTU) is the last dimensionless parameter UA min which is a ratio of the heat exchangers ability to transfer energy to the mini- mum fluid’s ability to retain energy. N TU = (1.2/I) The previously mentioned dimensionless parameters, at least for the moderately simple heat exchanger geometry, can result during the analytical analysis of the heat transfer occurring in a heat exchanger. For example, con- sider a double pipe counterflow heat exchanger arrangement. The correlation describing the performance of this heat exchanger, whose derivation is shown in Appendix E, is 1_e-NTU(1-C.) = 1 - Cge'muu - c,) (1.25) for the case of CR :1: 1. If CR = 1, equation (1.25) reduces to e ___ NTU 1+NTU (126) This is the usual form used in evaluating or presenting the performance of a heat exchanger, as a function describing the relationship between the three dimensionless parameters. These functional relations are tabulated for many different heat exchanger geometries and are the basis for easily predicting the performance of a heat exchanger design. The derivation of these functional relationships typically assumes neg- ligible axial flow of heat in the wall. The conditions for which this assump- tion is true for a counterflow heat exchanger geometry is the focus of this thesis. Furthermore, the amount that the performance is over-estimated by assuming axial conduction is negligible will be addressed. 1.1 Literature Review. The counterflow heat exchanger theory previously presented was for the idealization of no axial conduction (in the flow direction) of heat in the wall or fluids. It is not apparent under what conditions this assumption is not valid, but the consequences can be an over-estimation of the effectiveness for a given NTU by as much as 45%, as shown by Rosenhow [17]. Some of the earliest work performed in analyzing the effect of axial conduction was done by Mondt [9]. By inspecting the finite difference equa- tions for a solution that includes axial conduction, he proposed using a con- duction parameter to correlate the results. This conduction parameter was _ Kw... " 1.me and he showed that to accurately predict the convection from a surface at low Reynolds numbers requires consideration of longitudinal conduction. A (1.27) In later works, the investigation was less specific and looked to predict the influence of axial conduction on the performance of the heat exchanger. Landau and Hlinka [8] and Pan and Welch [11] developed exact solutions to predicting the performance of a counterflow heat exchanger while including the effect of axial conduction. The solutions were in the form of finite series and were algebraically involved. A general investigation was difficult, if not impossible, with the availability of the computing power at the time. However, Pan, Welch and Head [10], were able to improve on their earlier work [11] and predicted analytically the effect of axial conduction on the performance of the heat exchanger as a function of the NTU and a “modi- fied-conduction- flowrate-ratio” N ____ KM... M K mcp Lure“p Using the exact solution, the existence of axial conduction was predicted for a counterflow heat exchanger with balanced symmetric flow, that is equal mass flow rates and convective heat transfer coefficients for both fluids. Although more general than their past work, the balanced symmetric flow was still rather restrictive. (1.28) T address: incorp01 tion par: Mondt 11 He then s Mondt nr worst and mance of large Mon nant Thi: and solutic flow heat C ([8]. [10], l 90‘ needed W 011' tron to desc diffePential to solve the eaSIIy cyalu “has to obta 0f the heat e 10 The most recent work was performed by Rohsenow [17], who addressed reasons that laminar flow heat exchangers perform poorly. He incorporated a similar formulation as [8], [10], and [11], but used the conduc- tion parameter above in equation (1.27), from Mondt [9], which he called the Mondt number Kw... M0 = mep He then solved for the performance of the heat exchanger at the limits of the Mondt number, zero and infinity. The extremes of the Mondt number provide worst and best case scenarios. A Mondt number of zero gives the perfor- mance of the heat exchanger when axial conduction is negligible; while a large Mondt number gives the performance when axial conduction is domi- nant. (l .29) This investigation shall differ from previous work in the formulation and solution of the problem. All previous investigations involving a counter- flow heat exchanger have neglected temperature gradients normal to the wall ([8], [10], [l l], [17]). Although this assumption does seem reasonable it is not needed for this investigation to obtain a solution; whereas past studies required this assumption. This diflerence results in a partial differential equa- tion to describe the wall temperature for this study, instead of the ordinary differential equation obtained for past studies. However, the techniques used to solve the partial differential equation will allow for a solution that is more easily evaluated; while past solution techniques required special flow condi- tions to obtain a solution. Therefore, by considering the two-dimensionality of the heat exchanger wall the problem was more mathematically compli- cated in comparison to past studies, but allowed for a more general investiga- tion of the heat exchanger design parameters. Chapter 2 Describing Equations and Boundary Conditions 2.0 Dimensional Equations 2.0.1 Wall Conduction Equation The geometry is shown in Figure 2.1, a two dimensional wall with hot fluid entering at one end and cold fluid entering at the other. The problem was formulated as a boundary value problem on the wall assuming steady state and constant properties. The differential equation and boundary conditions corresponding to the plane geometry with the given assumptions are a’rw a’rw a—XE +5?- = O (2.1) or", 1&5; x=0 = halTWWJ) -TOI (2.2a) at}, ““53 “L = hLlTw(L, y) -TL] (2211) Cold Fluid Tom 4— 3 :2 WW 1: = 0 x = L '—'> T3 (I) Hot Fluid Figure 21. Heat exchanger geometry sturflINote, the analysis was perErmed on a per unit Width basis. 11 12 -K"3§wly=o = h” IT" 0‘) - T”(x’0)] (2.2c) BTW _K‘Q} ”5 = hen}, (x, 8) -Tc(x)] (2.211) where boundary conditions of the third kind were used to represent the most general case and still allow for boundary conditions of the first and second kind by adjusting the magnitudes of the constants (Kw, ha, 11,, I)”, he) . 2.0.2 Fluid Energy Balance Equations The describing equations for the fluids were formulated using an energy balance on a difl’erential length element of the heat exchanger. For the hot fluid the energy balance on a differential element of length Ax can be written Hlx = HI” Ax+ q"w,,PHAx (2.3) Substituting mass flow rates and the specific enthalpy and rearranging gives m£|x+ Ax "' litiil Ax Taking the limit as Ax -) 0 and assuming the fluid behaves as an ideal gas or an incompressible liquid (deT = din) gives x +quwauPH = O (2.4) . 47' H (mHCp. Ifld—x '1’ ‘1"wauPH = 0 (25) The heat flux from the wall from Newton’s law of cooling is 4"...11 = hHITH (x) - T... (x. 0)] (2.6) substituting into equation (2.5) and rearranging, results in the final describing equation for the hot fluid (IT 111’ ”+111! E mucp,[T"(x) -T,,,(x. 0)] = 0 (2.7) Prescribing a constant temperature at the inlet gives the initial condition 1,,(0) .-.- TH’ 1.. (2.3) Using a similar analysis, the describing equation for the cold fluid and corresponding initial condition are l3 ch hcpc ,. — T = O (2.9) T + .CCpc [Tw (x, 8) C (x)] TC (L) = TC, 1.. (2.10) 14 2.1 Dimensionless Equations 2.1.1 Nondirnensionallzation of lengths To remove any dimensional length dependence on the formulation of the problem the length scales were nondimensionalized using the characteris- tic lengths of the wall 3 L (2.11) + y y = 8 (2.12) Applying this nondimensionalization to the previously presented equations for the wall and fluid temperatures results in equations that only have units of temperature, which were made dimensionless after solving (to be shown later) for algebraic ease and for the physical insight available with units of temperature. The wall conduction equation and boundary conditions are 371,, L2 3’1, — — —- = 0 (2.13) 8x+2 52 ay+2 31,, + "T +BioTw(0.y ) = BioT, (2.1421) ax x’=0 or”, — +BiLTw(19y+) = BiLTL (2.1411) ax x‘ l a W . + . + ”—I +B'nTw(x :0) = Br,,T,,(x ) (2.14c) ay yO_o an, . + . + 57 +B'cTw(x 11) = Bic Tc(x ) (2.141) 2’ y’=1 where dimensionless parameters introduced in the equations are defined in the next section. Similarly, nondimensionalization of the equations for the hot and cold flu- ids gives 15 j{IEIHM(THO?) -T.,(x*.0)] = 0 dx“ 711(0) = 771.1» £0+NC[T,,(x+, 1) —Tc(x*)] = 0 dx+ Tc“) = Tc... 2.1.2 Nondimensionalization of Temperature (2.15) (2.16) (2.17) (2.18) The problem was solved as presented in equations (2.13-2.18), giving the solution with units of temperature only. The following equations were used to present the results in a dimensionless form: _ T..(x*.y*)_— Ta... 9 (x+9 +) "" w y TH,in"-TC,in T (x+)-T 9110‘” = I; -T C'in H.111 C,in T x+ -T 9C(x+) = C( ) C,in TH, in - TC, In (2.19) (2.20) (2.21) 2.2 Dim' Th alization . madcal 11 wall cond thickness As this pa tion must must be is likely the duction is The dir the wall, e 16 2.2 Dimensionless Parameters The dimensionless parameters that resulted during the nondimension- alization of lengths in the equations will provide a link between the mathe- matical model and the physical problem. The dimensionless parameter in the wall conduction equation, equation (2.13), is a ratio of the wall length to thickness L‘ 3% (222) As this parameter is increased, the change in heat conduction in the x-direc— tion must also increase. Therefore axial conduction must increase; or the wall must be isothermal in the y-direction (normal to the wall). The later is most likely the dominant effect, but there may be values of L’ for which axial con- duction is the dominant effect. The dimensionless parameters that appear in the boundary conditions on the wall, equations (2.14a-d), are Biot numbers 71,1. Bin 1! 7?: (2.23) Bi L I Elf—f (224) Bi :1 a hr”; (225) Bi C I hTC-s (2.26) These parameters are expected to have the most influence on the solution and axial conduction. The Biot numbers on the hot and cold sides of the wall, equations (2.25-26), should be more influential than the Biot numbers at the ends of the wall. Biot numbers, in general, give a ratio of the convection from a surface of the body to the conduction through the body. This provides insight into the temperature difference that will exist across the body com- pared to the temperature difference between the surface and convecting fluid. For large Biot numbers the convection dominates, resulting in a large temper- ature gradient across the body and a surface temperature that is close to the convecting fluid temperature. Conversely, for small Biot numbers the con— duction v and a larg Th of transfe ters are where the Using equ I61] In [em which relal 17 duction will dominate, giving a small temperature gradient across the body and a large difference between the surface and convecting fluid temperatures. The dimensionless parameters for the fluids are similar to the number of transfer units (NTU), but based on the pertinent fluid only. These parame- ters are h h N” = WA” -.- .55! (2.27) mHCpn C” h I: NC = .cgc = 34C (2.28) m Pc C where the area in equations (2.27-28) is A” = AC = Lw (2.29) Using equation (2.29) with the heat capacity, equations (2.27-28) can be writ- ten in terms of the corresponding Biot number Bi K L’w N1! = JG";— (2.30) Bi K L'w ”c = _Cé_ (2.31) which relates the wall parameters to the fluid parameters. 23 Heat] Solutior ture fields 1 wall. But ti heat excha] mance fi'on lation neglc terms of tht These dime eters from I exchanger’: Effectivl considered may be con imbalance 1 energy bala and the effe 18 2.3 Heat Exchanger Performance Analysis Solution of the equations presented up to this point provides the tempera- ture fields for the wall and fluids, without neglecting axial conduction in the wall. But the temperature as a function of the position along the length of the heat exchanger is not the desired result. It is desired to compare the perfor- mance from this solution with one that neglects axial conduction. The corre- lation neglecting axial conduction was presented in equations (1.25-26) in terms of the effectiveness, number of transfer unit, and heat capacity ratio. These dimensionless parameters must be related to the dimensionless param- eters from the mathematical solution to make the comparison in the heat exchanger’s performance when axial conduction is included. Effectiveness, defined in equation (1.23), will depend upon which fluid is considered when axial conduction is included because some of the energy may be conducted out the ends of the heat exchanger and result in an energy imbalance between the fluids. Thus, when including axial conduction the energy balance on the fluids are, from equation (1.5) ‘IH = CH (TH, in " TH, out) if 9c = Cc (Tc, out " Tc, in) (232) and the effectiveness must be defined as a function of the fluid qt! 8 = — (2.33) ” 9...... 5c = _qC (2.34) 9 max The definition of NTU is given in equation (1.24) and will not depend on the inclusion of axial conduction, but must be related to the dimensionless parameters used to model the physical problem. This was accomplished by evaluating the thermal circuit that describes the heat exchanger wall. The thermal circuit, neglecting any fouling, is shown in Figure 2.2, which is simi- lar to the previously presented thermal circuit (except values for the individ- ual resistances are shown for the geometry studied). The total resistance of the circuit is evaluated using equation (1.9). The conductance for the heat exchanger wall given by equation (1 . 10) is l9 O—W 'ww «MA—o 1 6 1 hHAH Kwfiw hcAc Figure 2.2. Thermal circuit for heat exchanger wall UA=-1—- = [_L.+_5 + 1 I1 (2.35) Rm hHAH KM») ”MC and NTU easily follows: HA 1 1 8 1 '1 NTU = — = + + ] (2.36) Cnu‘n Cmin[hHAH KHAN! hC‘AC Relating equation (2.36) to the mathematical formulation was accomplished through equations (2.27-28); substituting these equations into equation (2.36) gives 1 1 8 1 '1 ”T” ‘ c7; [~74 " V. * N—ccc] Noting that A, = Lw, equation (2.37) can be rewritten as a function of dimen- sionless parameters and the wall thermal conductivity using equations (2.30- 31), after some arranging it becomes (2.37) waL’ [ BiHBiC ] 7:," 31,,+ BiHBic+ Bic N TU = (2.38) The equations to be solved for the temperature fields of the wall and fluids are given in equations (2.13-2.18). Using the solution, the performance of the heat exchanger can be presented with equations (2.33), (2.34), and (2.38) in terms of heat exchanger dimensionless parameters, which can be compared to solutions neglecting axial conduction to quantify the effect of axial con- duction. Chapter 3 Method of Solution 3.0 General Method The problem presented in equations (2.13-18) is a system consisting of a second order homogenous partial differential equation (PDE) with four nonhomogeneous boundary conditions coupled with two first order nonho- mogeneous ordinary differential equations (ODE). The coupling occurs through the boundary conditions of the PDE and in the nonhomogeneous term for the ODEs. Thus, if a general solution can be obtained for the PDE without evaluating all boundary conditions, that is in terms of a product of unknown constants and functions, the solution can be substituted into the ODEs and the ODEs solved as a function of the unknown constants. The ODEs solution’s can then be substituted back into the boundary conditions of the PDE and the unknown constants evaluated. With the constants specified, the explicit solutions for all differential equations can be generated. This is the general procedure used to solve the set of differential equations and is shown schematically in Figure 3.1. This section is intended to present the methods that were applied to solve the problem. However, when the same methods were applied more than once they are not shown in detail after the first application. Also, some meth- ods contain details that are not necessary in understanding a method, and whose inclusion would overburden this section. For these reasons, the reader is referred to the Appendices (A-D) for additional details about the solution. 20 = a — anew—NEON 23:20th m \ 21 _ FEE 200 To » cm :Eonaoa .mi 55350 ~ Emu a ifiewfiom 33m — 4i T :55 cue—om 33m , H _ Eco _ 3:335 8:28 35... _ he czem H 5.88 Eves—a: .. mam—enowoatOkl . 2.2—AN noun—em n3:— - + a .1 3.525 Q dam Essa—BU .Ud 88:52 23 . gee—eav RE .«o ud> .. ESE—=3 eoquEomam .. one Awam ~30:qu .—om E U .wofiefl 4.3m $35.25 no new . L 682:8 EV): 358.50 BEE—.5 .. nous—om :85 » 1 .33 .«o ua> . dash * _ REM .woaom b €23:va 3:58me m Il— eeE a: _ FIGURE 3.1. Flowchart of solution proccdm 3.1 Genet Sino ferent matk separation of this met ferential cc lem has f0] introducin; one of the 1 matical for 22 3.1 General Solution of the Wall Conduction Equation Since the problem is linear, the PDE could be solved using many dif- ferent mathematical approaches; Fourier transforms, Green’s functions, or separation of variables. The method chosen was separation of variables. Use of this method requires that only one nonhomogeneous term exist in the dif- ferential equation or boundary conditions. As presently formulated the prob- lem has four nonhomogeneous boundary conditions equations (2.14a-d). By introducing scaling on the temperature 9,,(x*.y*) = T,.(x*.y*) -T, (3.1) one of the nonhomogeneous boundary conditions is eliminated. The mathe— matical formulation with the scaled temperature becomes ' a’e , a’e 3+‘;’+L2 83:0 (3.2) x y 30w . + ___; +8106“: (0’)) ) = O (3.3a) x . x =0 86w . + . 5;: +BzL6w(l,y ) = BtL(TL—T0) (33b) x‘=1 33:" + Rifle“, (x+, 0) = Bi,,[T,,(x+) - To] y y’=0 (33¢) ml? 0 + _ . 4- 5y— +3100», (3: ,1) _ Blc[Tc(x )-T,] (33d) y’ =1 Three nonhomogeneous boundary conditions, equations (3.3b-d), still exist. Since the problem is linear the principle of superposition may be applied. The problem is partitioned into a set of simpler problems that contain one nonho- mogeneous boundary condition, and can be solved using separation of vari- ables techniques. Then the solution to equations (3.2-3) is obtained by application of superposition. The problem in equations (3.2-3) can be repre- sented by three simpler problems, since there are three nonhomogeneous boundary conditions. Using superposition the solutions are related as fol- lows: ’thh‘v‘Tfl'Z‘i ‘ III-In). I: where 9‘. is 23 3 9,,(x".y+) = 20,.(x+.y+) (3-4) i=1 where 0, is the solution to the problem given mathematically as 3’9, ,2 3’9, +L — = 0 (3.5) ax+2 ay+2 -'—: +Bioeg (0’)”) = O (3.63) x . x =0 5—} +BiL0‘.(O,y+) = SliBiL(TL—To) (3.6b) x x’=l 803 . + . + y y0=o Bel — +Bice‘.(x*, 1) = 53gBicch(x+) —T0] (3.6d) 8y y*=1 wherei = 1,2,3 and 8.. = kroneckerdelta = { 0 “‘1' (3.7) 1' 1 i=j It can be shown by adding the three problems (i = 1, 2, 3) stated above that the originally posed problem in equations (3.2-3) is obtained. 3.1.1 General Solution for Nonhomogeneous Boundary Condition at x = L The separation of variables technique begins by assuming a product solution for 91 of the form 91 = F(x*) G (y*) (3.8) Substituting into the differential equation and boundary conditions, equations (3 .5-6). yields F, ,G+L"G , ,F = 0 (3.9) H H - GFx, +BiaGF(O) = o (3.103) {=0 where the f( Rearranging tion (3.10) 1 Where a 001 (3.12). Sing Only a func StanL This tial (limiter 6(Nation (3 Sign is ch01 Omogeno‘ 24 CF", 2 =1+13i,,(;F(1) = BiL(TL— To) (3.101.) —FG , +Bi,,m(0) = o (3.10c) y y* =0 FGy, +BiHFG(1) = 0 (3.10d) y’ =1 where the following convention for representing derivatives was used: 2 a F 3F . . (3.11) 8x+2 “‘ ‘ Rearranging equation (3.9) to group similar variables and simplifying equa- tion (3.10) gives l Fx‘x‘ = _ G”). = i 2 (3 12) [1.2 F G u ' —F, +BioF(0) = o (3.13a) x x0 = 0 017‘, . 1+BioGF(l) = BiL(TL-To) (3.131.) x = -G . +Bi,,o(0) = o (3.13c) y y* = 0 Gy, +Bi,,G(1) = 0 (3.13.1) y‘ = 1 where a constant, 1:112 (the eigenvalue), has been introduced in equation (3.12). Since the left hand side is only a function of x+ and the right hand side only a function of y+ , to obtain equality both sides must be equal to this con- stant. This key result allows for separation of the variables, reducing the par- tial differential equation to two ordinary differential equations. Before equation (3.12) can be separated a sign must be chosen for the constant. This sign is chosen to produce an eigenvalue problem for the function that has two homogenous boundary conditions to evaluate. The eigenvalue problem will have a general solution in terms of sine and cosine functions for this Carte- sian coordinate system. After choosing the sign on the eigenvalue, two ordi- nary differential equations can be written from equation (3.12) with the appropriate boundary conditions from equation (3.13) _ ‘2 2 _ Fr; L [1 F — 0 (3.14) l'w- . :5 Note that eq (3.16) has tv tion (3.13b), to either fun tion Thus, t complete so' Equat The general Applying th the two C01]: WhiCh Can t tiOn 25 -F . +BioF(0) = 0 (3.153) x {=0 G , , +1120 = O (3.16) y y -G . +BiHG (0) = 0 (3.173) y y0=o G , +Bi”G(l) = 0 (3.17b) y y*=1 Note that equation (3.14) only has one boundary condition, whereas equation (3.16) has two boundary conditions. The reason for this can be seen in equa- tion (3.13b), which depends on both functions F and G and cannot be applied to either function singularly. This is the nonhomogeneous boundary condi- tion. Thus, this boundary condition must be applied after assembling the complete solution. Equations (3.14) and (3.16) have easily obtainable general solutions. The general solution for equation (3.14) is F(x+) = Alcosth‘x”) +A,sinh(uL‘x+) (3.13) Applying the boundary condition, equation (3. 15a), gives a relationship for the two constants Bio A2 = uL' A1 (3.19) which can be substituted into equation (3.18) and rearranged giving the solu- tion F(x”) = A2 [uL‘ cosh (uL’x‘) +Biosinh (uL‘x+)] (3.20) The solution of equation (3. 16) is G (3") = A3008 (11f) +A4sin (11)”) (3.21) Applying the boundary conditions in equation (3.17a) solves for the relation- ship between the constants _ BiHA3 ' u Substituting into equation (3.21) and rearranging gives A4 (3.22) .‘V'1'I‘I‘I ' Lg Applying C any inform equation. 1 equation 51 To meet [ht satisfied fo: it has an in. solved to (11 is then With 11, g1V« Subm Vide the sol; ‘0 equation solution, A]. one constam dition, aqua rate and Was 26 60*) = A. [ucos my“) +Biysin (uy*)] (3.23) Applying the final boundary condition, equation (3.17b), does not provide any information about the constant A 4 because the constant cancels from the equation. However, it does provide information about the eigenvalue ()1). The equation simplifies to (BiH+ Bic) u” tan (tin) = 2 . . To meet the boundary condition in equation (3. 17b), equation (3.24) must be satisfied for all values of u". Equation (3.24) is a transcendental equation and it has an infinite number of solutions, hence the subscript on n. It will be solved to determine the possible values of u". The solution for the function G is then (3 .24) G(u,,.y*) = A2 [u,cos my) +Biusin (u,y*)] (3.25) with u” given by the solution of equation (3.24). Substituting equations (3.25) and (3.20) into equation (3.8) will pro- vide the solution. However, since there exist an infinite number of solutions to equation (3.25) all possible solutions will be summed to obtain the final solution. Also, the undetermined constants, A2 and A 4, will be grouped into one constant that will be determined later by applying the last boundary con- dition, equation (3.6c), which was the boundary condition that did not sepa- rate and was not evaluated. The solution for problem one is 91(x+,y+) = A" [an’ cosh (unL'f) +Biosinh (unL‘xtn *- n = 0 [uncos (unf) + BiHsin (uny+)] (3.26) and the eigenvalues are found from the solution of equation (3.24). For ease in presenting these large solutions, the following functions are defined: x1 (11» x") a put} cosh (unL‘xt) + Biosinh (nnL’xt) (3.27) Y1 (u, y‘) a u..cos (u,y*) + Bigsin (u..y*) (3.28) and the solution can be rewritten as 3.1.2 General: The so, ous section. " repeated here for problems 93(x The eigenvalue Defining the fol allows the 80111111 27 01(x*,y”) = A,x1(un,x*) Y1 (amy‘) (3.29) =0 3.1.2 General Solution for Nonhomogeneous Boundary Conditions on y The solution for i=2, 3 apply a similar procedure as used in the previ- ous section. The details of which are shown in Appendix A and will not be repeated here. Applying these procedures results in the following solutions for problems two and three: “ a a 02 (1:3)”) = Z B..[°°Sh (4f) - Cnsinh (4f N" n=0 L L [ancos (an?) + Biosin (anx+)] (3.30) 6H“+ “cm”shmH B'inha"+* w»; new (9 )+ (a )1 [ancos (Gina?) + Biosin (anx+)] (3.31) The eigenvalues are given by the positive roots of the transcendental equation (Bio + Bi L) an tan (an) = 2 . . (3.32) (an—BrLBro) Defining the following variables: X2 (“m x”) a: ancos (053+) 4» Bio sin (ans?) (3.33) all an Y2(an,y+) I cosh (—;y*)-§nsinh (7)”) (3.34) L L Y (a y“) sficosh (25y+)+8i sinh (5f) (3 35) 3 n’ Lo Lo H Lo allows the solutions to be presented in a compact form 02 (x*, y+) = 13”)!2 (an, y+) X2 (an, x”) (3.36) n = 0 03 (x*, V) = CuY3 (an, y+) X2 (an, x*) (3.37) - 0 LF 3.2 Soluti! 3.2.1 Hot F. The hot fluid, p equation in r, (H, 0) 1 from equat Try (1+: 0) = Where the 12 (3.36) and ( (2-15) and a dTH dx‘ HT With the initi; 1111's eqllation 131' 801111101]. T 28 3.2 Solution of the Fluid Energy Balances Equations 3.2.1 Hot Fluid Formulation The ordinary differential equation describing the temperature of the hot fluid, presented in equations (2.15-16), is coupled to the wall conduction equation in the differential equation itself. This is seen by the presence of Tw (x+, 0) in the differential equation. The value of this term can be obtained from equations (3.1), (3.4), (3.29), (3.36), and (3.37) as Tw(x+,0) = 6w(x+,0) +70 = 70+ ZAnX1(un,x+)Y1(un,O) n=0 + 2‘, [B,Y2(a,,.0) +C,,Y3(an,0)]X2(ocn,x+) (3.33) n = 0 where the last term in equation (3.38) is obtained from combining equations (3.36) and (3.37). Substituting this into the difi'erential equation in equation (2.15) and arranging gives d1” " , EF+NHTH = NH{T0+ 2A,,X1(un.x )Y1(ll,,.0) n=0 + 2 [BnY2(an,0) +CuY3(0t,,.0)]X2(a,,.x*) } (3.39) n=0 with the initial condition from equation (2. 16a) TH (0) = TH, ,, (3.40.1) This equation is a linear nonhomogeneous ordinary differential equation. Its solution may be written as the sum of a homogeneous solution and a particu- lar solution. The homogenous problem is given as dTH’ h x‘.’ +NHT,“ = o (3.41) and a particular problem is 41,,er The particr pler proble eral solutio where TH h dTH. pi d1+ + NH 29 dT a. dSP+NHTH,p = NH{T0+ ZA”X1(un’x+)Yl(11n,O) n=0 + 2 [3..Yz(a,,.0) +CuY3(a,,,0)]X2(an,x+) } (3.42) n=0 The particular problem in equation (3.42) was further divided into three sim- pler problems for each term on the right hand side of the equation. The gen- eral solution for equation (3.39) in terms of the simpler problem’s solutions is 3 TH(x+) = 1H,,(x+) + 2 TH,p,.(x*) (3.43) i: l where T”, n is the solution to equation (3.41) and T”, p, is the solution to (IT -—H-;p‘+NHTH.pr = NH{811T0 +821201!" X1 (umx )Y1(|.ln,0) + 53, 20 (Buy, (an, 0) + CRY3 (an, 0)]x2 (01”,?) } (3.44) (i = 1, 2, 3) and 8], is the kronecker delta function as defined in equation (3.7) 3.2.2 Hot Fluid Solution The homogeneous problem in equation (3.41) is solved giving T”, h (x*) = Dlexp (-N,,x*) (3.45) where D1 is a constant that will be evaluated after assembling the complete solution as shown in equation (3.43). The particular solution was obtained by using standard variation of parameter techniques. This method utilizes a general solution of functional form similar to the nonhomogeneous term in the problem and then solves for the unknown constants in this general solution. For the first particular solu- tion (i=1) the general solution would be a constant 7‘”, p1 (x‘) = k = constant (3.46) Substituting equation (3.46) into equation (3.44) and solving for the constant gives Thesecon . x, (1*) . sin TEPZ (1 Substituting constants, bu TH, p2 L Finally. the th TH,p3 (1 Which after su stants can be ‘ TH,p3 (J Assembling m. (3'49)! and (3 _5 dance is 30 r”, p, (x*) = r, (3.47) The second particular solution (i=2) would be of the same functional form as X, (x‘) , since all other terms in this summation are constants 1,”, (x+) = for, (x+)) = 2 bucosh (11,131?) + casinh (unL‘r‘) (3.43) n = 0 Substituting equation (3.48) into equation (3.44) and solving for the unknown constants, b, and c", equation (3.48) can be written as 0- . (BioNH— ”ZL‘Z) t T (x*) = A 2L 1 - ", cosh L Jr“) H,p2 ago n {11, [ (Ni, " 11:1. 2) (”a (BioNH - ufiL‘z) + ” " (Ni—93L”) Finally, the third particular solution is of the functional form :Isinh (unL’f) } (3.49) T“, (x+) = f(X2 (an, in = 2 ancos (an?) + dusin (an?) (3.50) n = 0 which after substituting into the differential equation and solving for the con- stants can be written as B' N 2 = 2(B.+1:c.){a.[1-‘ '0 W] ":0 (N;+a'2.) BiN +042 +NH[ o H n 2 2 NH-t-an ]sin (aha?) } (3.51) Assembling the solutions to the simpler problems, equations (3.45), (3.47), (3.49), and (3.51), using equation (3.43) the solution for the hot fluid energy balance is TH(X+) = 7 Applying m 31 T”(x*) = T0+D1exp (-N”x") " , (BioNH-nfiL‘z) , + + 2A,,{ufiL [1- (N;_u:L,2) ]cosh(unL x ) n=0 (BioNH — ufiL”) + H n inh L‘ + (Na-11:14.2) ]S (11,, x)} .- B.” 2 + 2 B +F—:C a 1-( 1° ”+a") cos(0t. x+) n=0 n L n II II (Ni, + (1%) - 2 +N BroNH-Hxn H N2“:2 H n ]sin (0%”) } (3.52) Applying the boundary condition, equation (3.403), gives the constant (BiaNH — ng‘z) ] DI = TH,fn-To- zAnugL [1- (N2 _uzLez) H u n=0 _(BgNfi+ab an 2w()[ 145» 3.2.3 Cold Fluid Solution For the cold fluid the ordinary differential equation describing the tem- perature is slightly more involved due to the wall temperature being evalu- ated at (y” = 1). This adds algebraic difficulty to the problem, but procedurally the solution is obtained exactly as the hot fluid. For this reason, only the final solution will be shown. All the details for this solution can be found in Appendix B. The solution for the cold fluid energy balance equation is ram . is. V} i M Tc(x+) Where the i; and the fOIIm 32 Tc(x+) = To+Elexp (ch+) " . NcBio + ij’2] . + A ’ unL 1+ , cosh (unL 1:”) ago I! { [ [Viv-[12L 2 N i+ ZL‘2 .. +NC[ C3" 5"” :lsinhulnL m} NC-p'uL " NCBio—az + 2 (B ’+C ’) {a [l+——-—"]cos(a x+) "=0 " " " Ni”: " N Bi -a2 +NC[-C_2°__2_"]sin(anx*) } (3.54) Nc+°‘,. where the integration constant is n=0 " _ NCBi +1121.” , 151 = exp (-NC) [(TC’ m-To) — 2A":{un1, [1+ N2 0 “25.2 ]cosh(unL ) C- n NcBiawfiL’z ,, +N , sinh(u L ) c[ Né-uib 2 " " N B' - 2 .. z (Bn'+cn') {anl:1+_ci_§2 cos(a) ”=0 Nico-a: ] " N i-oz2 +NC iii—J smut") (3.55) ”6"“: and the following constants were introduced to simplify the expressions: An’ = A,I [uncos (1.1") + Bi ”sin (un)] (3.56) Bn’ = Bn[cosh (5—2) - Cusmh (:2 J] (3.57) all all . all Cn’ = C" [—. cosh (7) + Bi Hsmh (7 )] (3.58) L L L 33 A991 Ha‘ boundary condition: condition: to be appl equation ( neous bou The difflcu an infinite : (3.37), are 33 3.3 Application of Nonhomogeneous Boundary Conditions Having obtained the solution for the temperature of the fluids, the final boundary conditions (nonhomogeneous) can now be applied. The boundary conditions that need to be applied are the three nonhomogeneous boundary conditions given in equations (3.6b-d). These boundary conditions only need to be applied to the respective 6‘. solution, (i.e. equation (3.6b) applied to 6,, equation (3.6c) to 92 and equation (3.6d) to 03 ). These three nonhomoge- neous boundary conditions are 39 ——1 +BiLe,(1.y*) = Bun-1,) (3.59) ax x’ :1 392 --—, +3510; (x+.0) = BingHU“) - To] (3.60) By y, ___ o 393 F +Bice3 (x+: 1) = Bic[Tc (x+) - To] (3.61) y y’ =1 The difficulty with equations (3.59-61) is that the relationships are in terms of an infinite series for 9." These series, from equations (3.29), (3.36), and (3.37), are 61(x*,y+) = ZAnX1(“nvx+)Y1(lln:>’+) (3.62) II 92(x*,y*) = 23nY2(°‘n’Y+)X2(°‘n:x+) (3.63) n 63(x+,y*) = ECuY3(°‘.»Y+)X2(%x+) (3.64) II This means that theoretically an infinite number of constants will need to be determined. Of course, experience with series solutions has shown that less than an infinite number of terms (and constants) need to be considered since the series converges after a finite number of terms; beyond which, the indi- vidual terms are zero or negligible in comparison to the summation up to that point and future terms decay further. However, the number of terms needed for the series to converge is not known and must remain arbitrary and assumed large. 34 The procedure used to determine the constants that will satisfy the nonhomogeneous boundary conditions involves applying orthogonality. Using the orthogonality of the eigenfunctions, the constants will be deter- mined. An appropriate mathematical statement of this property is l gwny,(u,.y*)Y1(u,,.y*)dy* = { ”(fim) '21:: (3.65) 01' 1 + + + 0 main {wxyX2(an,x )X2(am,x )dx = { N(a,,,) m=n (3.66) where wxy is a weighting constant, which for the Cartesian coordinate system is equal to one. After substituting expressions for the temperatures and deriv- atives of temperatures into the boundary condition in equations (3.59-61), the equations will be multiplied by a second eigenfunction and integrated over the boundary. Since the eigenfunctions have the property shown in equations (3.65-66), the equations will be simplified and the unknown constants can be determined. 3.3.1 Application of Boundary Condition at x = L Beginning with the solution for (31 in equation (3.26) and taking the derivative with respect to x+ and evaluating both expressions at the boundary gives 91(1,y+) = A, [unL‘ cosh (unL') +Biasinh (unL')] * =0 [u,,cos (my‘) + Biusin (11f) ] (3.67) 39 .- 5—:(l.y*) = ZAquiL’zsinhmuL') +BiouuL‘c6sh(unL‘)]* x n=0 tu,cos (u,y*) +Bigsin (u,y*)1 (3.68) Substituting these expressions into the boundary condition, equation (3.59), produces 35 A, [ufiL‘zsinh (unL') + BiounL‘ cosh (uuL‘n * - 0 [u,,cos my”) + BiHsin (u,y*)] + ML 2 A" [unL' cosh (unL’) +Biosinh (unL')] * n=0 [uncos (Hf) + BiHsin (unfll The orthogonality of the eigenfunctions is used to determine the constants. The eigenfunction for this problem is Y ("3' y+) = uncos (any) + BiHsin (unf) (3.70) and these eigenfunctions have the property shown in equation (3.65). There- fore, equation (3.69) will be multiplied by a second eigenfunction Y (um, y”) and integrated over the boundary. Because this integral is only nonzero for m=n, only this term will remain from the summations. The integrals of each term in equation (3.69) after multiplying by a second eigenfunction are shown in Appendix C. After applying orthogonality, equation (3.69) reduces to A," [uiL‘zsinh (uML') + BiaumL‘ cosh (umL‘) 1 N (pm) + AmBiL [umL‘ cosh (umL') + Biosinh (umL') 1 Mum) = _ . Bi H Bi H BIL (TL- To) [srn (pm) - T008 (um) + T] (3.71) where all summations that were present in equation (3.69) have been reduced to a single term through the application of orthogonality. Equation (3.71) can now be solved for the unknown constant . . Bin Bill 3'1, (TL‘ To) [8111(um) "' T608 (11".) + T] A - "' "' (3.72) "' - [(uiL'2+BioBiL) sinhmML’) + (Bio+BiL)umL'cosh(umL')]N(um) 36 The functional relationship for N (um) is N(u ) =1 (u2+Bi2) 1+ 8’" +Bi (373) "' 2 "' ° ( (ufi,+BiL2)) ° ' This completes the solution of problem one, equation (3.62) could be evalu- ated since A" is known from equation (3.72), noting that m is a dummy vari- able. 3.3.2 Application of Boundary Condition at y = 0 The difficulty in applying the boundary condition in equation (3.60) is recognized by noticing that the right hand side contains the temperature of the hot fluid. Using the solution for 62 from equation (3.30) and taking the derivative with respect to y+ and evaluating both at the boundary gives 92 (x+, O) = 2 B" [ancos (an?) + Biosin (anx+)] (3J4) n: 0 a62 0- an ——(x*, 0) = - z Bra-C, [ancos max”) +Biosin (anx*)] (3.75) (3y+ "=0 L In addition, the hot fluid temperature is, from equation (3.52) TH (x‘) = To + Dlexp (-NHx+) - c (BioNH - “2L. 2) t + A (FL 1- '1 cosh (an x+) "go n{ n [ (”z-“7231’ 2) ] Putting these three expressions into equation (3.60) yields 37 " or 2 3,2,5; [ancos (anxt) + Biosin (63+)] n = 0 + Bi” 2 B" [ancos (0:31?) +Biosin (anx+)] = n=0 BiH[Dl exp (—N,,x+) " . (Rwy-uh”) . + A 2I. 1- '1 cosh( L x”) n=0 .- B.” 2 + 2 B +?—:C a l-( to ”+059 cos(a x‘”) u=0 l L n n n (N3, + 6:) BiN +02 +1)!" "2” " sin(anx+) (3.77) NH-t-ot: for which orthogonality will be applied. Each term in equation (3.77) will be multiplied by the eigenfunction X (an, x*) = aucos (anal?) + Bio sin (anf) (3.78) evaluated at term m, and integrated over the boundary. The eigenfunction has the same property as discussed earlier, which is given in equation (3.66). After applying orthogonality, equation (3.77) becomes Bulk»: = -0)le - 2 An (pmn,1+pmn.2+pnm,3+pmn.4) -2(Bn+cni7)wmn(3'79) n = 0 n where the order of the terms in the two equations, equation (3.77) and equa- tion (3.79), has been maintained, except the first two terms in equation (3.77) have been combined into a single term in equation (3.79). Additional vari- ables were introduced to represent the results of the integration, the variable definitions and integration details are shown in Appendix C. The complexity added by the hot fluid being the nonhomogeneous term in the boundary condition is seen in equation (3.79), where summations are pres: cation 01 Because perature . term. Als ally two 5 1 equation, temperatu the other I 3.3.3 Appli The very similz exPression both at the afld the temp 38 are present after applying orthogonality. This was not the case for the appli- cation of orthogonality to 01, where all summations reduced to a single term. Because the summations in equation (3.77) that represent the hot fluid tem- perature do not contain eigenfunctions, they were not reduced to a single term. Also, equation (3.79) is a function of all three unknown constants (actu— ally two since A” was just obtained). Nevertheless, the dependence of this equation, which would only be a function of the constant 3,, if the hot fluid temperature were not present in the boundary condition, is now coupled to the other boundary condition, equation (3.61), through the constant C”. 3.3.3 Application of Boundary Condition at y = 5 The evaluation of the boundary condition given in equation (3.61) is very similar to the previously discussed method for 62. Starting with the expression for 63 in equation (3.31) and taking the derivative and evaluating both at the boundary gives + a all an . an 93 (I . 1) = C" —;cosh (—,)+BiHsmh (7) * n= 0 L L L [ancos (anf) + Biosin (anx+)] (3.80) EBB . an an an an —(x+, 1) = 2: Cn—,[—.sinh [.7)+ai,,cosh {—7)}- ay“ n=0 L L L L [ancos (an?) + Biosin (anx*)] (3.81) and the temperature of the cold fluid, from equation (3.54), is TC (1+) = 39 TC (x+) = To + Elexp (ch+) L, 1 + (BioNC + ufiL‘z) "l N2 _ 2L‘2 ( c “n ) + ZAn’ n80 J cosh (an' x+) (Bi N +u2L‘2) , +NC 02 C 2 ".2 sinh (unL x") (Nc‘unL ) “ (Bi N —62 + 2 (Bn’+ Cn’) {an[1+ 02 C 2 ")Jcos (unit?) ":0 (Nc+orn) Bi N —62 +NC[_°_2.C_2_'.'] sin (anx+) } (3.82) Putting these temperatures into equation (3.61) gives C n "sinh n B h n I B. . )1 "2 n ,I , ( ,) lHCOS ( J, ancos (aux ) rosrn (aux “ or or or +Bi C —:cosh(-—:)+Bi sinh(-—:')]‘ CE”, "[L L ” L [ancos (uni?) + Biosin (anx+)] = Bic [El exp (ch+) “ , (Bi N + 2L”) "I uflL 1+ 02 C “'12 (NC " 11,2,L ) ] cosh (unL' x") (Bi N + 2L”) . +NC °2c ”'12 sinh(unL x“ (Ne-“zl‘ ) (BioNC-ag) 2 (NC + 63) (BiONC- a?!) . + +NC[ (Ni-+ai) Jsrn (crux ) } J (3.83) ] cos (ornf) writ-:1- ‘ for which ’ multiplied ‘ evaluated a the same pr After applyir where the or( tion (3.83) ar (3.83) have b ilar procedure to represent ti details of this Notice : (3.61) had the Summations e t0 equation (3. After ap was obtained f Where No NOE that the sub s own since it is 40 for which orthogonality will be applied. Each term in equation (3.83) will be multiplied by the eigenfunction X(au, x‘) = ancos (uni?) + Biosin (an?) (3.84) evaluated at term m and integrated over the boundary. The eigenfunction has the same property as discussed earlier, which is given in equation (3.66). After applying orthogonality, equation (3.83) becomes Cmrbm = Quiz, 4» 2A,} (PM 1 + PM 2 + PM 3 + PM 4) + Z (Bn’ + Cn’) \PM (3.35) n I! where the order of the terms has been maintained for the two equations, equa- tion (3.83) and equation (3.85), except that the first two terms of equation (3.83) have been combined into a single term in equation (3.85). Using a sim- ilar procedure as the last boundary condition, new variables were introduced to represent the results of the integration when orthogonality was applied. All details of this step are shown in Appendix C. Notice that the presence of the cold fluid temperature in equation (3.61) had the same effect as the hot fluid in the previous boundary condition. Summations exist after applying orthogonality and equation (3.85) is coupled to equation (3.79). 3.3.4 Summary of Applying Orthogonality After applying orthogonality to equation (3.59) a closed form solution was obtained for the 01 constant, from equation (3.72) . . Bin Bin BrL(TL- To)[s1n(un) - u—cos (11") 4» ll ] An = 2 .2 c c c (3'86) [unL sinh(p.nL )+ (Bio-tBiLwnL cosh (unL )]N(un) where _ l 2 . 2 BiL . Nut”) - §[(ll,,+3‘o ) (1+ (1134'3'12) )+Blo] (3.87) Note that the subscript n was substituted for m in the form of the equations shown since it is a dummy subscript. Applica the simple a1 gel The diflerence tion (3.59) this a function (infir present two pro the series will n Secondly, they < in equations (3.4 pled equations g ranging, are Burl. +2.: (11?,l + Cmd’u'Z (gn'. These equ for the unknown exist in the equai number of terms “V0 «muons; h hence the subscr' series (m = N), t Stants. The unkn (" = 1. 7. -..N). of B, and C", no Also’ D1 and El mbsfiWIBd. Fro 41 Application of orthogonality to equations (3.60-61) did not result in the simple algebraic equation with one unknown as it did for equation (3.59). The difference is the terms on the right hand side of the equations. For equa- tion (3.59) this term is a constant, whereas for equations (3.60-61) the term is a function (infinite series). These functions are TH (x”) and TC (x*) , and they present two problems. First, they are not orthogonal functions; and therefore, the series will not be reduced to a single term after orthogonality is applied. Secondly, they contain both unknown constants, B” and C", which will result in equations (3.60-61) being coupled after applying orthogonality. The cou- pled equations given in equation (3.79) and equation (3.85), after some rear- ranging, are n+2 (313+ Cu (5;)wnm + 0) 1:301 = -2Afl (pm, 1 + pm». 2 + pmn, 3 + prim. 4) (3‘88) 1: C,,,"l — 2 (B; + Cn’) ‘PM - one, = 2A,;(PM1 + Pm M 2 + PM 3 + PM 4) (3 39) II II These equations represent a set of simultaneous equations to be solved for the unknown constants B" and C" (n = 1, 2, 3 ....... ). Because summations exist in the equations it is now required to truncate the series after a finite number of terms (n = N). This gives 2N constants to be determined from the two equations; however, both equations can be written for each eigenvalue, hence the subscript m. Since there are as many eigenvalues as terms in the series (m = N), this totals 2N equations to solve for the same number of con- stants. The unknowns in equations (3.88) and (3.89) are 3,, and C" (n = 1, 2, ....N). The primed constants in equation (3.89), which are functions of B, and C", need to be substituted for, giving a new variable for that term. Also, D1 and E, are function of the unknown constants and will need to be substituted. From equations (3.53) and (3.55), D1 and E1 are respectively .[1_ (BioNH- pita] (Nir nil-’2) D1 = TH,in—To— ZAnuyer n=0 B'N 2 ' 2 (B,+%C,)a,.[l-( 1" “69] (3.90) ":0 (Ni-Fai) Introducing (3.91) and ( reduced to Xmas"- (\{1 “WP ( ‘N. These CC Own Const [0 the left of [he were emPle 42 B, = exp (-Nc) [(13, ,,-T,) - 2A, n80 , NCBi,+n¢L‘2 , ’ L 1+ ', cosh L ) {,5 l Né-uiL 2 l m" NCBi,+ufiL’2 , +N , sinh( L ) Cl Né-uib 2 u" a. N ' _ 2 - 2 (B,'+c,') {6,[1+—C—B'°—fi]cos(a,) 2 2 ”=0 ”0+0." N Bi -62 +NC —c—2i——" sin(a,) (3.91) NC-t-or: Introducing variables to represent the large terms in equations (3.90) and (3.91) and converting B,’ and C,’ to B, and C, allows these equations to be reduced to D, = (THJu—T -2A, [3,,- 2 o,(B +C, 2"”) (3.92) 1130 E, = exp(-NC) (TQM-T0) - 2 A.’°,." Z (B,a,+c,c,) (3.93) n=0 =0 where all new variables are defined in Appendix D. Substituting for D1 and El from equations (3.92) and (3.93) into equations (3.88) and (3.89) and defining new variables to convert the primed constants gives a 3,3, + 2,119. N... - 7.0)..) + ZCnL—i' W... - 7.9)..) = - mm (TH, in - To) - 2A7! (prim, 1 + prim, 2 + pm», 3 + prim, 4-mmBn) (3'94) n 2 [9,6,- (‘l’B) M] B, + go, + 2 [9,5,- (‘PC) M] C, = Qmexp (-NC) (Tc, ,, - T,) + 231; (PM 1 + PM +,,P 3 + P — 0,6,) (3.95) n mn,4 These equations are now in the form that they were solved for the unknown constants. The equations are organized with functions of unknowns to the left of the equal sign and known quantities to the right. Standard meth- ods were employed to solve the set of simultaneous equations described by 43 equations (394-95) and these methods are addressed in Appendix D. In short, the two equations were separated into coefficient of B, and C, for each eigenvalue and equation. These terms were put into a matrix and then solved for the unknown constants using a partial pivoting scheme with the Gauss method. Once the constants were determined, all solutions could be evalu- ated. fish) of th subr tines 0U) the s were 44 3.4 Computer Program AXCOND The development of the solution to this problem has been completely analytical up to this point. A closed form solution is not possible since a set of simultaneous equations must be solved for series constants, which depend on a truncation of the general series solutions obtained for the temperature of the wall and fluids, as discussed in the previous section. With the solution framework complete a computer program needed to be developed to perform the evaluation of the analytical solution. The main functions of the program were to: 1. Determine the eigenvalues from the transcendental equations. 2. Assemble the set of simultaneous equations and solve for B, and C,. 3. Evaluate the series to calculate the temperatures of the wall and fluids. 4. Check that the boundary conditions are satisfied and that energy is conserved. A FORTRAN program was developed to perform these functions and a listing of the program is given in Appendix F. In addition, a general flowchart of the program is shown in Figure 3.2. It contains a main program that calls subroutines to perform each of the individual tasks. There are nine subrou- tines in the program: INPUTO, INPUT], ROOT, BUILD, GAUSS, PROFIL, OUTO, OUTl, and OUPUT. Each of these will be briefly discussed to sum- marize its function. In the program, variables were named to coincide with the solution shown herein, and deviations from the variables in the solution were noted. The first two subroutines INPUTO and INPUTl, as there names imply, import the data needed for a particular run. The former reads dimensional data in and computes the dimensionless parameters needed for computation; while the later reads in the dimensionless parameters, depending on the option selected. In addition, INPUT] allowed for a particular input parameter to be varied, with the first call to the subroutine reading in the data and the parameter to vary and subsequent calls updating the selected parameter for another computation. This option allowed for assessing the eflect of various parameters on the solution by observing the results as one parameter is var- ied. Subroutine ROOT is used to locate positive roots of the transcendental equation. It employs a marching scheme that searches for a change in sign of the function then backs up to locate the zero crossing, i.e. the root. A great deal of tin: the large c1 tion and m ASSt constants l obtaining 2 partial pivc Tern solutions it the derivat lated. The . ated in the and X=L) t Ewation (1 series. Usir the fluid en ditions On 1 enabled a I! two CheCks Used ‘0 Selc ity 0f the S( 45 deal of time was spent to make this subroutine sufficiently robust to handle the large changes in magnitude of the parameters in the transcendental equa- tion and not skip roots of the equation. Assembling the matrix of simultaneous equation to be solved for the constants B, and C, was done in subroutine BUILD. This was the key step in obtaining a solution. After building the matrix, the subroutine GAUSS used a partial pivoting scheme to solve for the constants. Temperatures for the wall and fluids were generated from the series solutions in the subroutine PROFIL. In addition to computing temperatures the derivatives of the temperatures with respect to the coordinates are calcu- lated. The derivatives of the analytical series solutions were taken then evalu- ated in the subroutine. Also, to calculate heat loss from the wall ends (x=0 and x=L) the integral of the heat flux over the ends was taken 1 8T q'end = IKw r 0 dy+ (3.96) and 31: Equation (3.96) was determined analytically and programmed to evaluate the series. Using the derivatives of temperature it was possible to check whether the fluid energy balances equations were being met and if the boundary con- ditions on the wall were being satisfied. The heat losses from the ends enabled a total energy balance to be performed insuring conservation. These two checks, conserving energy and satisfying the boundary conditions, were used to select the number of terms to truncate the series and insure the valid- ity of the solution. The final three subroutines generate the output file for a run. Subrou- tines OU'FO and OUTl correspond to the option listed earlier for the input. The main difference is that OUTO produces the results of one computed run while OUTl varies a particular parameter giving the results as a function of the varied parameter, but not listing the unchanged parameters again. Subrou- tine OUTPUT can write the fluid and wall temperatures as a function of posi- tion if desired. These temperature are only written if a flag is set in the input file and there are separate flags for the fluids and wall. Likewise, another flag in the input file can check that the boundary conditions are satisfied and write the results to the output file. 4.0 Ident' In the dimensionles form a link dimensionless r— the wall aspec and Biot numb The Parameters Chapter 4 Results and Discussion 4.0 Identification of Input Parameters In the formulation of the mathematical model for this problem certain dimensionless groups appeared (Chapter 2). These dimensionless parameters form a link between the mathematical model and the physical problem. The dimensionless parameters needed to solve the wall conduction equation are the wall aspect ratio L L - 5 (4.1) and Biot numbers for each surface of the wall 3' h°L 42 '0 - 7?: ( . ) B' hLL 43 'L - T: ( - ) . M5 . has BIC - —K: (4.5) The parameters needed in the fluid energy balance equations were I: Bi K L w N, = _"A” = _” W (4.6) Ca Ca h Bi K L w NC = _CAC = ___c .. (4.7) In addition to these dimensionless parameters the ambient temperatures at the wall ends (T, and TL) and the initial or inlet condition for the fluids (TH, 3, and Tc, ,,) need to be specified. Mth these eleven parameters the mathematical model can be solved for the wall and fluid temperatures as functions of posi- 47 ”I tic pd Cl] [16 de; ass ne: enc get the 11101 the deli puti part do inde corn tion tains solun ccipa mUm tions. Signir 48 tion. There are then eleven independent variables that may influence axial conduction as the problem was formulated. The investigation of all eleven parameters to determine their influence on axial conduction presents a diffi- cult task. To reduce the number of independent variables, quantification of axial conduction will be shown by evaluating the solution using effective- ness-NTU relationships. ‘ Presenting the solution in terms of the effectiveness eliminates the dependence of the solution on any of the ambient or fluid temperatures, assuming a constant ambient temperature (T, = TL); because the effective- ness is based on temperature differences and scaled by the temperature differ- ence at the fluid inlets. The NTU, which represents the physical and thermal geometry of the heat exchanger, was earlier shown in Section 2.3 as (4.8) NTU K” L’ _ ...,L [ B 1' ”Bi c Biu+ BiHBiC+ Bic] This equation introduces the need for specifying the thermal conductivity of the wall and the minimum heat capacity; whereas for the mathematical model, the thermal conductivity was accounted for in the Biot numbers and the heat capacity accounted for in NC and N H and did not need to be explicitly defined. It was decided that since the additional parameters of the wall ther- mal conductivity and minimum heat capacity need to be specified for com- puting NTU, they could also be used to compute NC and N,, from the other parameters, as shown in equations (4.6) and (4.7). This allows for the more physical parameter of heat capacity to be input. The wall thermal conductivity and minimum heat capacity do not independently affect the solution because they appear only as a ratio used to compute the NTU in equation (4.8) and the parameters NH and NC in equa- tions (4.6) and (4.7). It does not matter that the ratio in equations (4.6-7) con- tains the individual heat capacities rather than the minimum. Since the solution depends on the heat capacity ratio and the magnitudes of the heat capacity will be set by the ratio of the wall thermal conductivity and mini- mum heat capacity, the ratio will be accounted for in at least one of the equa- tions. Therefore, the individual magnitude of these variables is not significant. C. The mum heat c cific pr0Pel' dimension]: adjusted to these variat group of ter variable wil where the w. definition. 5 The Mondt 1 After Specif number is see .DUIIIber can ‘ ‘3 important 49 The appearance of the ratio of the thermal conductivity to the mini- mum heat capacity provides some additional flexibility to this study. The spe- cific properties of the wall and fluids do not need to be specified, and a dimensionless analysis can be performed. Furthermore, this ratio can be adjusted to provide the desired magnitude of NTU. To present the ratio of these variables, a new dimensionless variable will be defined to represent the group of terms that are not dimensionless in the NTU in equation (4.8). This variable will be the Mondt number, from Rohsenow [1 7] wa M, s , (4.9) CminL where the wall aspect ratio was included to be consistent with Rohsenow’s definition. Substituting the Mondt number into equation (4.8) gives (4.10) NTU = M L‘2[ BiHBiC ] 0 BiH-t- BiHBiC-t- Bic The Mondt number will be used as the adjustable parameter in this study. After specifying the fluid Biot numbers and wall aspect ratio, the Mondt number is set to provide the desired NTU. The magnitudes of the Mondt number can be looked at to correlate the conditions in which axial conduction is important TABLE 4.1. Input variables and dependency for determining heat exchanger performance Input Parameters Mathematical Model NTU Effectiveness L' L' L' L‘ Bi, Bi, Bic Bi, Bi L Bi L Bi H Bi L Bic Bic M, Bic Bi H Bi H Bi H Kw N C CR CC NH . Mo CH TH, in TH, in TC, in TC, in To To TL 1; Presen f number of ind allows for the with this soluti tiveness, will (1 of Table 4.1. An addit terms needed b vious discussio overlooked in t the series soluti ary conditions ( fled Since IllCSt applying this c1 verge, Another tion. This chec] from the hot flu 6°“ from the w 0f the TCSUlts, u Validity of the s The influ 301116011 WOUId 50 Presenting the results in an effectiveness-NTU format will reduce the number of independent variables. The input parameters that were specified to obtain a solution are shown in the first column of Table 4.1. These parameters were used to calculate the information needed for the mathematical model and NTU shown in columns two and three in Table 4. ]. This information allows for the mathematical model to be solved and a NTU to be associated with this solution. However, the solution, which will be in the form of effec- tiveness, will depend only on the seven parameters shown in the last column of Table 4. 1. An additional input that is not listed in Table 4.1 was the number of terms needed before the series is truncated, N. It was not included in the pre- vious discussion because it lacks any physical significance, but it was not overlooked in the analysis. The value of N must be sufficiently large to allow the series solutions to converge. Taking this criteria a step further, the bound- ary conditions on the wall were checked to ensure that they were being satis- fied. Since these boundary conditions involve the wall and fluid temperatures, applying this criteria is equivalent to having all three series solutions con- verge. Another check performed was an overall energy balance of the solu- tion. This checked that the solution conserved energy, i.e. the energy lost from the hot fluid was accounted for in the cold fluid and/or lost by convec- tion from the wall ends. Although magnitudes of N will not be shown in any of the results, using those two criteria the value of N was chosen and the validity of the solution assured for all results reported. The influence of N was not investigated in any detail. In general, the solution would converge based on an energy balance within 10-20 terms. However, the boundary conditions would not be met until more terms were considered; but the solution (2 - NTU) did not change fi'om the values seen at a smaller number of terms. These results suggest it may be possible to obtain simplified expressions for the solution, but this point was not pursued further. 4.0.] Presentation of Results It was anticipated that the effect of axial conduction could be demon- strated through the seven variables upon which the effectiveness depends (shown in Table 4.1). These variables were to be investigated to determine their influence on the performance of the heat exchanger. In order to present 51 these results it is useful to introduce a final variable that will describe the deg- radation in the performance of a heat exchanger. It is called the ineffective- ness, and it represents the amount that a heat exchanger’s performance is reduced due to the effect of axial conduction. It is defined by 8 “-8 iC = M (4,11) eNeg 8 "'6 i” = lit—l! (4.12) 8 Ne; where cm, is the effectiveness calculated neglecting axial conduction using equation (1 .25), which assumes one-dimensional heat flow. The ineffective- ness must be defined for both hot and cold fluids since the energy balances may difl'er due to heat lost from the wall ends. For the case that there is negli- gible heat loss from the ends of the wall, the two values for inefl'ectiveness are equal ic = iuai for 8c = e” (4.13) Note that when the ineffectiveness is zero there is no degradation in the performance of the heat exchanger due to axial conduction, indicating that axial conduction is negligible. In addition to showing when axial conduction is not negligible, the ineffectiveness gives the magnitude that the effective- ness would be over-estimated or in some special cases under-estimated, using a standard effectiveness-NTU relationship which neglects axial conduction. It is the percentage of error introduced from assuming that axial conduction is negligible. 52 4.1 Influence of Axial Conduction at a Constant N TU To begin this parametric study the NTU will be held constant as other parameters are varied. Through some initial computer runs it was seen that a larger degradation in effectiveness was seen at higher values of NTU. Hence, a NTU of seven was chosen for the study. This allows the identification of the more important parameters and will give insight into the conditions under which axial conduction should not be neglected. While maintaining a constant NTU, the Biot numbers are varied over a range from .0001 to .1, the heat capacity ratios of l, .75, .50, and .25 are investigated, and the wall aspect ratio is also varied. Variation of the Biot numbers will be done by grouping all four or grouping the ends and the hot and cold side. Although grouping all the Biot numbers may be unrealistic, it will provide insight for future more realistic cases. 4.1.1 Influence of the Wall Aspect Ratio (L‘) A natural choice to begin this parameter study is with the wall aspect ratio, which will address the physical dimensions for which axial conduction should not be assumed negligible and may need further investigation. Also, by beginning with the wall aspect ratio some bounds can be established on the other parameters that need to be investigated. This is achieved by elimi- nating ranges of the other parameters for which axial conduction only occurs under unrealistic physical dimensions. All four wall Biot numbers were set equal, and the magnitudes were varied beginning with .1 and decreasing an order of magnitude each run to the final value of .0001. The results, showing ineffectiveness as a function of the wall aspect ratio, are presented for the dif- ferent Biot numbers in Figures 4.1-4.4 with each figure depicting a different heat capacity ratio. Examination of the data, prior to plotting figures, showed a negligible difference in the effectiveness based on the cold fluid and the effectiveness based on the hot fluid (2,. = e”). For this reason, only one ineffectiveness was plotted in Figures 4.1-4.4. There was a slight deviation from the two effec- tivenesses being equal for very low magnitudes of the wall aspect ratio (L' < 40). Because this effect is not related to axial conduction and merely demonstrates the effect of the wall aspect ratio becoming very small, it was not included in Figures 4.1-4.4. .1 N .A I» .3 I! . 1' v. . I“: The 1 represents 2 nonzero. D in the effec The degrad maximum 1 numbers th( istic values The t can be expl. exchanger 1 For a plane 1 53 The results indicate that as the aspect ratio decreases, which physically represents a thick wall or short heat exchanger, the ineffectiveness becomes nonzero. Depending on the heat capacity ratio (CR), one finds a degradation in the effectiveness from 20% at small C, to nearly 50% at C, equal to one. The degradation is further enhanced as the wall Biot numbers decrease. A maximum Biot number of .1 is shown in the figures, since for larger Biot numbers the ineffectiveness only departs fi'om zero for very small and unreal- istic values of the wall aspect ratio. The trend of increasing ineffectiveness with decreasing aspect ratio can be explained by considering the thermal resistances across the heat exchanger wall in the two principal directions L RStreainwiu = K (4’14) #‘x R - 8 (4 15) - — 0 normal K 4y For a plane wall A, =-- 195 and A, = wL giving L RStreanrwiu = K W (4-16) 19 ] waL. It can then be noted that as L' becomes small, the streamwise thermal resis- tance becomes small while the normal thermal resistance becomes large. This will lead to a greater energy flow in the streamwise direction and ultimately a larger degradation of the effectiveness. Rnonnal = (4.17) In considering the influence of the Biot number on the degradation of the effectiveness, the increase in ineffectiveness with decreasing Biot number is expected. Decreasing the Biot number results in a reduction in the convec- tive heat transfer at the surface of wall and/or an enhancement of the conduc- tive heat transfer within the wall, which is shown by considering the definition of the Biot number R 5/K Bilxfl = __“fl = E]: (4.13) convection y w I - I I ‘t 1 , a .3] :--: :3 Cal exc 54 Both these trends lessen the influence of axial conduction since energy moves more freely from the surface while encountering more difficulty conducting within the heat exchanger wall. This is further seen by observing the temper- ature profile along the wall. These results are shown in Frgure 4.5. Notice that the temperature profile along the wall flattens as the Biot numbers decrease, that is, the temperature at the ends of the heat exchanger becomes closer to the temperature at the center of the heat exchanger as the Biot numbers decrease. This is the effect of axial conduction moving energy along the length of the heat exchanger wall. Coupling the effect of the Biot numbers to the wall aspect ratio follows. As the Biot numbers decrease, the alternate path of heat flow along the wall becomes more influential. At the smaller Biot numbers the wall aspect ratio has less influence because the convective resistance is so high. There- fore, the wall aspect ratio requires a larger magnitude to reduce the normal resistance or increase the streamwise resistance to the point that axial con- duction will not occur. Hence the smoothing of the curves in Figures 4.1-4.4 as the Biot numbers decreased. The final dependency that can be observed in Figures 4.1-4.4 is the influence of the heat capacity ratio. To aid in this comparison, Figure 4.6 was generated, which illustrates how the heat capacity affects the solution for a Biot number of .001. This figure does not provide any new information, it instead combines a curve from Figures 4.] through 4.4 into a single figure. In Figure 4.6 it is seen that, at a particular Biot number, the ineffectiveness increases as the heat capacity ratio increases.'lhis trend is not a physically obvious result. To better understand the effect, the dimensionless temperature profile along the heat exchanger at a wall aspect ratio of 100 was generated as a function of heat capacity ratio and is shown in Figure 4.7. In this figure the median wall temperature shows little or no variation in shape with heat capacity ratio, the curves increase a constant amount as the heat capacity ratio decreased. This can be explained through the existence of a greater dis- parity in the amount of energy per degree of temperature between the two fluid streams as the heat capacity ratio decreases. The stream with the maxi- mum heat capacity will undergo a smaller temperature change in the heat exchanger than the fluid with the minimum heat capacity. Thus, the fluid with the maximu with the mi) change, as I energy to 1’3 perature as 1 wall does 111 examinatior For n sionless terr dimensionlc capacity rat effect of ax heat exchan energy bala However, i: energy excl for CH2 Cc 55 the maximum heat capacity can cause a larger temperature change in the fluid with the minimum heat capacity while experiencing less of a temperature change, as the heat capacity ratio decreases. Therefore, there exists more energy to raise the temperature of the wall which results in a larger wall tem- perature as the heat capacity ratio decreases. The elevated temperature of the wall does not explain the decreased ineffectiveness, however; and further examination is needed. For more insight into the effect of the heat capacity ratio the dimen- sionless temperature profiles of the hot and cold fluids were added to the dimensionless wall temperatures and are shown in Figure 4.8, with a heat capacity ratio of .75 omitted for clarity. This figure clearly demonstrates the effect of axial conduction, which moves the energy along the length of the heat exchanger within the wall. The magnitude of which can be gauged by an energy balance on the fluids. For negligible axial conduction 011:?! = Ccff (4.19) However, if axial conduction is present there exist a local imbalance in the energy exchange between the fluids and equation (4.19) becomes d9 d9 H c ( dx+ ' dx” )°‘ qAxial (4.20) R for CH2 CC, which can be rearranged to give TABLE 4.2. Relationship between the slope of the hot and cold fluld’s temperature profile and the heat capacity ratio day/(1f c, ((19,!de CR x+=0 x+=0.5 x+=1 1.0 15.68 1.“) ”15.63 .75 15.53 1.28 1/11.58 .50 14.29 1.81 1n.75 .25 13.5 2.67 114.95 t4 Using excha the h< show: and a.- for we when excha shows nituc' Cold on tl the l the \ 56 d9 ”1 d3:+ (c, (dGC/dx“ ) Using equation (4.21), axial conduction at different locations along the heat exchanger can be assessed. To help quantify the differences in the slope of the hot and cold fluid temperature profiles Table 4.2 was created, which shows the ratio on the left hand side of equation (4.21) at various locations and as a function of the heat capacity ratio. This ratio should be nearly one for negligible axial conduction and considerably smaller or larger than one when axial conduction is appreciable. The data in Table 4.2 demonstrates that as the heat capacity ratio decreases, axial conduction shows a local decrease at the ends of the heat exchanger; while it increases near the center. Comparison of the magnitudes shows the dominant effect at the x+ = 1 end of the heat exchanger, while the x“ = 0 end and center show less significant change. This can also be seen in Figure 4.8, more qualitatively, by a comparison of the slopes of the tempera- ture profiles. The net effect on the heat exchanger is that axially conducted heath ineffectiveness decrease as the heat capacity ratio decreases, as dem- onstrated in Figure 4.6 where ineffectiveness decreases with heat capacity ratio. -1)... 4M“, (4.21) The difficulty with interpreting this case is that both solutions, includ- ing and neglecting axial conduction, change with the heat capacity ratio. Therefore, it is required to analyze why the solution more closely resembles the case of neglecting axial conduction as the heat capacity ratio decreases. There are two related reasons. The first is the change in the available energy on opposite sides of the wall. Second is the change in the driving potential on the hot side of the wall, which is the difference between the fluid and wall temperatures. The smaller heat capacity ratios were obtained by increasing the mag— nitude of the hot fluid heat capacity while maintaining the heat capacity of the cold fluid. Providing more energy per degree on the hot side of the wall than on the cold side and this disparity in energy across the wall becomes larger as the heat capacity ratio decreases. This situation was earlier shown to increase the wall temperature. j ' (“'91: innv‘a 'u' tion‘ obse hast pote: ture. incre cons eneq toth iBecz rado hste ontl ng it im 03133 853; thei the 1 57 The increased driving potential on the hot fluid side of the wall results from the increase in the heat capacity of this fluid, which can be shown by observing the describing equations for the hot fluid, from equation (2.15) dT , CHE? +Bi,,K,,L w[T,,(x*) - T, (1:20)] = o (4.22) The dimensionless parameter previously introduced in this equation (N H) has been separated into components to isolate the heat capacity. The driving potential balances the product of the heat capacity and slope of the tempera- ture. Therefore, when the heat capacity increases, the driving potential must increases to balance the equation since the slope of the temperature remains constant. The outcome of these two effects is that the cold fluid acquires more energy at or near its inlet, due to the hot fluid having. more energy to transfer to the cold fluid and the increased driving potential to move more energy. Because both these effects become more pronounced as the heat capacity ratio decreases, the cold fluid will depend less on axial conduction to provide its temperature rise and the ineffectiveness will decreases. Having an understanding of the overall effect of the heat capacity ratio on the performance of the heat exchanger, it is now possible to refer back to Figure 4.6 and the influence of the wall aspect ratio can be addressed. In this figure it is seen that the wall aspect ratio reduces the ineflecfiveness to zero as it increases, and the convergence to zero is more pronounced at lower heat capacity ratios. The increased convergence to zero is a consequence of the effect of axial conduction already reduced by the heat capacity ratio, and then the influence of the wall aspect ratio requires a smaller magnitude to decrease the normal resistance and eliminate axial conduction. The product of MaL‘2 is shown in the legend of Figures 4.1-4.4, from which the Mondt number can be calculated. The use of the Mondt number in this study was as a variable parameter, which could be adjusted to provide the desired NTU as the other parameters were varied. In general, as the Biot num- bers and wall aspect ratio increased; the Mondt number was decreased to maintain a constant NTU, while the Mondt number was independent of the heat capacity ratio because the minimum heat capacity remained constant. Referring to equation (4.8) for NTU, as the wall aspect ratio and Biot num- bers were inc1 by the ratio K unit width. The tre reaction to ch essence, the 1 properties to the Mondt nu increase sinc be available. by only the r 58 bers were increased, the noted trend in the Mondt number was accomplished by the ratio Kw [Cum becoming smaller because the analysis was based on a unit width. The trends seen in the Mondt number to maintain a constant NTU are a reaction to changes made in the wall aspect ratio and Biot numbers. In essence, the Mondt number is compensating the wall (Kw) and fluid (Cm) properties to match the conditions (Bi H, Bic, L‘ ). Therefore, it is fitting that the Mondt number decreases as the wall aspect ratio and fluid Biot numbers increase since more energy can be transferred for these conditions and must be available. This outcome suggest that axial conduction may be described by only the magnitude of the Mondt number. A point to be investigated later. ”a ' Tue-u 0.8-1 0.6-1 /| 0.4- 0.2- 0.0‘fi— FIGURE 4.3» 0.8— 0.6- 0.4~ /l 0.2.. 0.0J~,\ 59 0.81 MOL“ —"" .0001 140007 . BiH=Bic=BiL=Bi°= —"" .001 14007 —" .01 1407 —- .1 147 0.6- Ca=1 « NTU=7 .. tTa—T—‘W‘ .. 0.4 .\.. \.,\ \-.., . \ \~. '°\ \-. '\ \-- "\ \. 0.2- ' \ \. - \ \. - ‘\, \.. \. 0.0 l I I I 1‘- r . I UT..“_1 r r 1 10 100 1000 10000 L.* FIGURE4.1.IneffectivenessasafunctlonofthewallaspectratloforNTU = 7endCR = 1 0.8-1 M°Lu ----- .0001 140007 . BiH=Bic=BiL=Bio= —--~ -001 ”007 —-- .01 1407 _. .1 147 0.6- C,=.75 . NTU=7 0.4J ___....____,.;--\——..\ . \. ... .n\ 0.2« \ -\ \ . . \ -- ‘\ ‘\. \.\ .\.. \‘-. \" 0.0 u . Vv‘:"‘>" " \. 10 100 1000 10000 LII! “CURE“.heflecflveneuasefuncflonofthewaflaspectratioforNTU = 7 and CR = .75 0.8'1 0.6“ 0.4-1 /\ O.2~ 0.0 10 FIGURE“ Inc 0.84 0.6- 0.4- 0.81 MoLn —"" .0001 140007 - Bi =81 =31 =81 = —m .001 14007 " c L ° ——-- .01 1407 —' .1 147 (15- C'=.50 .. NTU=7 (14‘ .- ...IIO""—‘--.\.. \H. \\ '\.. \... ' \' a. .\ . \. .\ \.. \.~S-- .,\..—__'-\. .._-:.\--..__ 0.0 r I I I V I I V I r‘ ' T I ' I 10 100 1000 10000 Lil! FIGUREaneffectlveneuasafunctlon ofthewallaspectratlo for NTU = 7 and CR = .50 0.8- MOL” —°"° .0001 140007 - 81H=Bic=BiL=Bio= —'" .001 14007 _.. .01 1407 —- .1 147 0.6- c,=.25 - NTU=7 0.4-1 div 02‘ ‘;::"-..\ /e .. \ . -——-"-\,. \ \ \ \., \ \ 0.0 I T \igr $iu—‘r . \n--_I‘ l r r 1 10 100 1000 10000 LII! FIGUREMIneffectlveneuuefuneflonoflhewallupectntioforNTU = 7andCR = .25 l.( 0.E 0.6 =.5) 0.4 6w(y+ 3‘“ numb; 004a 61 I t I I —- [.0001r :.\ Blfl=Blc=Bl°=BlL= : .361 1 0.8“ "\..\.\ __ .1 -. \_ = . .\:\ ET}! 7 ‘ .\"\ L:—100 {0‘ 0-6 h' ‘\ .\ \‘H = —J n. 1 ' '\,_.\>‘ 3; Q ~-. ’ - \. \' "~- -.~.. @ 0.4 ..\\ fi . ° \. i 0.2- 'Q-.\- \. « \< 0.0 . l , r j 1 r T 0.00 0.20 0.40 0.60 0.80 1.00 x+ FIGURE 4.5. Median dime wall temperature as a function of dimensionless position and BiotnumhersforNTU = 7,L =100,ande =1 0.8- —--- Cu-1 ‘ — C.-=.75 --- c.=.50 ----' C =.25 0.6- ' —— “..\ NTU=7 0-4‘ ——--—--\.j§_ Biu=Bic=Bi°=BiL=0.001 . .....———. \\. MoLn=14007 ...—'- \. ..\ \. \. 0.2e ,,—-. .. .\ § . \ \,.\..\ 0.0 . . , .\-->‘—--+‘—’==1L—- . . , 10 100 1000 10000 LII! FIGURE 4.6. Ineflectiveness as a function of the wall aspect ratio and heat capacity ratio for NTU = 7,L = 100, and all Biot numbers equal .001 0.8 0.6- O.(y*= .5) 0.4' 0.21 0.8 0.5 62 1-0 ' l ' l ' 1 F l C.= .25 .. A ‘0. 1| ... >\ V . Q " NTU=7 " 0.2- LI=100 .. BiH=Bic=Bio=BiL=.OO1 0.00 0.20 0.14.0 0.60 0.80 1.00 FIGURE 4.7. Median dimensionless wpli temperature as a function of dimemionless position and heatcapacity ratio for NTU = 7, L = 100, and allBiot numbers equal .001 1 o o ff‘? ...... 6 j 1 I I I I ' NTU=7 "\. -\\, 1 1.7-100 \‘ \ j 0.2~ 31H=Bic=Bio=BiL=,001 \.\. _ — W0" -.5 \... ‘ ------- Hot Fluid . ""'-' Cold Fluid 0.0 1 1— ' I ' 1 v I ' fl 0.00 0.20 0.40 0.60 0.80 1.00 x41- FIGUREMMedhndhnensimhuwaflandluidgempentunuafuncflonotdhnemionku positionandheatcapadtyratioforNTU = 7,L = 100,andallBiotnumbersequal.001 4.1.1 C01] end 0011 as 1 as t nitt Bio hat tic 1 hek the ity 1 fun: tion Not two and alw the nun incr fOrx bers incr CDC; theh 63 4.1.2 Influence of the Wall End Biot Numbers (Bi 0 and Bi L) During the investigation of the wall aspect ratio it was determined that convective loses from the ends of the heat exchanger wall are negligible if the end Biot numbers are less than or equal to the fluid Biot numbers (hot and cold side Biot numbers). This dependence will be further investigated for cases that have end Biot numbers greater than the fluid Biot numbers. These cases are physically possible since the end Biot numbers have the wall length as its characteristic length and the fluid Biot numbers have the wall thickness as the characteristic length, as shown in equations (4.2)-(4.5). Thus, the mag- nitude of the end Biot numbers being greater than the magnitude of the fluid Biot numbers does not require the same dependence on the convective heat transfer coefficients, a case that would be unrealistic, because the characteris- tic lengths proportionately increase the end Biot numbers. Following a similar procedure as the previous section the NTU will be held at seven while the end Biot numbers are varied. The same magnitudes of the hot and cold side Biot number will be investigated at the four heat capac- ity ratios for a prescribed wall aspect ratio of 100. The ineffectiveness as a function of the end wall Biot numbers (Bio = Bi L) are shown for these condi- tions in Figures 4.94.12, with a different heat capacity ratio in each figure. Notice for these figures that each curve has a starting point at approximately a zero end Biot number, then is double valued for larger end Biot numbers. The two values at each end Biot number are the inefl‘ectiveness based on the hot and cold fluids. For a specified end Biot number, the cold ineffectiveness is always greater than or equal to the hot ineffectiveness in absolute magnitude. Observing Figure 4.9, which shows the ineffectiveness as a function of the end Biot number for a heat capacity of one, the effect of the end Biot number is seen. As the end Biot numbers increase i H decreases while iC increases, for all fluid Biot numbers investigated. This divergence of the per- formance based on the hot and cold fluids is expected as the end Biot num- bers increase because the heat lost by convection from the ends similarly increases resulting in an energy imbalance between the fluids and a differ- ence in the ineffectiveness based on the hot or cold fluid. Recalling that an ineffectiveness of zero implies axial conduction can be neglected, it is seen for smaller fluid Biot numbers that the performance of the heat exchanger is actually enhanced, as noted by the inefl’ectiveness based 64 on the hot fluid being negative. That is, the effectiveness calculated based on the hot fluid is greater that the effectiveness calculated neglecting axial con- duction. However, this enhancement of the performance based on the hot fluid is at the cost of the performance based on the cold fluid, as seen by the hot ineffectiveness being negative while the cold ineffectiveness approaches one. This outcome exists because the heat that is acquired from the hot fluid is lost from the wall ends and is not transferred to the cold fluid. The appro- priateness of this result depends on the use of the heat exchanger. For remov- ing heat from the hot fluid this outcome is beneficial, but for adding energy to the cold fluid the results are discouraging. In systems which have one stream as a waste stream, such as an automobile cooling system, this could be a use- ful result. The fluid Biot numbers affect both the starting point of the curve and the amount that the ineffectiveness calculated based on the hot fluid diverges from the ineffectiveness based on the cold fluid. Figure 4.9 shows that the curves begin at a higher ineffectiveness and diverge a greater amount between the two ineffectivenesses as the fluid Biot numbers decrease. Fur- ther, the divergence of the two ineffectivenesses as the end Biot number increase occurs at a faster rate at smaller magnitudes of the fluid Biot number. There are three relevant observations to be made as the fluid Biot numbers decrease 1) the starting ineffectiveness increases, 2) the magnitude of the divergence of the two ineffectivenesses increases, and 3) the rate of the diver- gence with respect to the end Biot number increases. These points will be addressed in subsequent paragraphs. The higher starting ineffectiveness of the curves at is a consequence of the effect of axial conduction increasing as the fluid Biot numbers decrease. The basis of this outcome was addressed in the previous section. The magnitude of the divergence of the ineffectivenesses increasing as the fluid Biot numbers decrease is also due to the effect of axial conduction increasing as the fluid Biot numbers decrease. Since energy is more likely to proceed along the heat exchanger as the fluid Biot numbers decrease, it fol- lows that increasing end Biot numbers will result in more energy being removed. The additional energy removed results in an increase in the magni- tude of the divergence as the fluid Biot numbers decrease, which correspond- ingly grows with the end Biot numbers. 65 Finally, addressing the rate of the divergence with respect to the end Biot numbers, which increases with the decreasing fluid Biot numbers, is again linked to the increased effect of axial conduction at lower fluid Biot numbers. With axial conduction being more prevalent at lower fluid Biot numbers it is expected that the sensitivity to the end Biot numbers will be greater at the lower fluid Biot numbers. With axial conduction more preva- lent, a change in the end Biot number will have more effect than when axial conduction is less and the effect of the end Biot number is diminished. This sensitivity to the end Biot number is seen by the change in ineffectiveness as the end Biot numbers change, which explains the change in rate of diver- gence as the fluid Biot numbers decrease. The influence of the heat capacity ratio requires a comparison of infor- mation from Figures 4.9-4.12. To aid in visualizing this information Figure 4.13 combines the curve for a fluid Biot number of .001 from each of the four figures, onto a single figure. In this figure the ineffectiveness is seen to increase as the heat capacity ratio increases, which was covered in an initial discussion of the heat capacity ratio, where this trend was also seen. How- ever, as the heat capacity ratio varied; the end Biot number showed more influence on the hot ineffectiveness. At a heat capacity ratio of one the curve for ineffectiveness based on the hot and cold fluids is symmetric. As the heat capacity ratio decreased, the curve for the ineffectiveness becomes asymmetric, sloping more on the lower leg. This leg of the curve represents the ineffectiveness based on the hot fluid. The change with heat capacity ratio represents the performance of the heat exchanger actually improving, as seen by the curving sloping more negative. The effect is due to the increased driving force, or larger difference between the hot fluid and wall temperature, on the hot side resulting in more energy available to the wall. This additional energy is then lost from the wall ends through convection. For this case the Mondt number depends only on the hot and cold Biot numbers, and therefore was constant with respect to the end Biot numbers. For this reason little can be learned about the Mondt number by varying the end Biot numbers. 10 . a e T . . . ././".’._——-.——w_— J 08- /-/ - 0.61 /. . ’ .____ ....— -'"" "7" -1 0 -° ‘* 04- ,,/-'/ _ NTU=7 ,g'z ' ../ C'=1 ‘ .\ ‘ L‘=1OO 02- - '\-\. ...... - \, ......... \ -'->-<:: """"""" 3159‘. MO _I' _\.. ......... .— 1%.“... ...... 1 oo_— ...... :KV- oooooooo s—esIs-ssss-8218=::::=:;.=_:.3 -t —.... '1 '0147 . .\‘\_\. -""‘ .01 .1407 ‘°“-~---..—-———— 4 —-- .001 1.4007 —° .0001 14.0007 —0.2 r l l I ' l 0.0 02 0.4 06 08 81° or BIL 22%“ ”i Ineffectiveness as a function ofthe end Biot numbers for NTU = 7, L' = 100, R = 67 a———--—-'——° .. 1 --"' . NTU=7 .‘r . c,=.7s O . 0'2'\ \"\.. ‘ L*=1oo . \ ...—...*OOO—.IC—I.I‘I..”'C _ 8=::::..-J-Brmin) Mawimwim) = M0(Bim=mm,n) .2 (4.26a) MoL . . k NTU (Blmax=Blmin) J while the NTU is w." . . \ W (Blmax=Blmin) NTU(Bimu>Bim-n) = NTU(Bim=Bi,,,,.,,) 2 (4.26b) MOL. . . K-N—TIT (Blmax>Blm‘-n)) where M0 (Bimx=Bim.n) in equation (4.26a) and NTU (Bimeimin) in equa- tion (4.26b) are read from Figure 4.37 and the ratio on the right hand side is obtained from Table 4.5 at the corresponding maximum and minimum fluid Biot numbers. Equation (4.26a) corrects the Mondt number to correlate to the NTU in Figure 4.37 for unequal fluid Biot number and equation (4.26b) cor- rects the NTU to correlate to the Mondt number in Figure 4.37 for unequal fluid Biot numbers. Only one of the variables (NTU or M) needs to be cor- rected; the other can be read from the figure directly. A final restriction to the dividing line in Figure 4.37 is that it applies for a heat capacity ratio of one. Because this heat capacity ratio is the most adversely affected, these results will be correct in predicting when axial con- duction is negligible for smaller heat capacity ratios. However, the results may indicate axial conduction is nonzero when in fact it is negligible at smaller heat capacity ratios. To alleviate this restriction, the dividing line can TDABIJE4J£.thpflhmkaofdhfldhmghimndtnunfiuusurlhuunhniofflue heateapaeflyraflo C M 1.0 .010 75 .015 .50 .030 25 .070 100 be shifted for smaller heat capacity ratios, where the effect of axial conduc- tion is less, if all other conditions are the same. Table 4.6 gives the approxi- mate Mondt numbers that divide axial conduction as a function of the heat capacity ratio. Because the line does not shift very much, the additional lines were not added to Figure 4.37. The Mondt number varies with the wall aspect ratio, NTU and fluid Biot numbers, as seen in equation (4.24). The Mondt number can be inter- preted as the parameter that matches the properties of the wall and fluids (Kw and Cm) to the operating conditions (Binm, Him, L', NTU). Figure 4.37 in conjunction with the results of Section 4.2.1 shows that as the Mondt number approaches zero the effect of axial conduction is negligible, while at larger values axial conduction is present. The smaller Mondt number results in less axial conduction because to become smaller the product of the streamwise resistance and minimum heat capacity must get larger, which is shown using equations (4.25) and (4.16) Kw 1 M0=CWL'=C R min min streamwise (4.27) A larger minimum heat capacity means the fluid retains its energy better, and the energy that is transferred to the wall is hindered from conducting axially by the large resistance in that direction. However, for larger Mondt numbers the heat capacity is small, and the fluid gives up its energy easily. This energy can then be conducted down the wall because of the smaller streamwise resis- tance. Physically the Mondt number represents the capabilities of the fluids and wall, the fluid’s ability to carry energy and the walls ability to move energy in the axial direction. Because the Mondt number varies with fluid Biot number, wall aspect ratio and NTU it will reflect any changes seen in these variables, which explains why it is possible to base the presence of axial conduction on this parameter The results in Figure 4.37 are intended to serve as a design guide to predict if axial conduction will exist. A certain amount of caution must be exercised in applying the results in Figure 4.37 and Table 4.6, however. The line presented on the figure for a heat capacity ratio of one was estimated by observing the ineffectiveness-NTU plots presented earlier, as were the magni- 101 tudes listed in Table 4.6. The criteria was an ineffectiveness less than 1% was negligible. For the regions near the dividing line, it is recommended to check the magnitudes of the ineffectiveness on the appropriate figure in Section 4.2.1. The corrections allowing Figure 4.37 to be applied for unequal fluid Biot numbers in equations (4.26a) and (4.26b) will be in error for small val- ues of the NTU, wall aspect ratio, and fluid Biot numbers if the relationship between the heat capacity and fluid Biot numbers associates the opposite extremes (Cm, - Bin“). This is due to the influence of the secondary effects when the thermodynamic energy becomes large. Even though the Mondt number is getting smaller the ineffectiveness increases, and this outcome reaches larger values of the NTU as the heat capacity ratio decreases. It is rec- ommended not to use the corrections for this association of the fluid Biot numbers and heat capacity ratio. 102 4.4 Application of the Results Up to this point, all operating conditions of the heat exchanger have been in terms of the Biot numbers, wall aspect ratio, heat capacity ratio, and Mondt number. Also the material properties of the wall or fluids have not been an issue since they were absorbed into dimensionless parameters. In this section the dimensional variables will be presented to show specific operating conditions (,i.e. flow rates, fluid types, and wall material) when the previ- ously presented results will apply. The results of this section could be used to establish experiments to verify the analytical solution. The parameters that will be addressed in dimensional terms are the wall aspect ratio L L=8 minimum and maximum fluid Biot numbers where heat capacity ratio where and the Mondt number Bi in = min (Bic, Bi”) Bi a = max (Bio Bi”) Cm." = min (Co CH) C max = max(Co CH) min (4.28) (4.29) (4.30) (4.31) (4.32) (4.33) (4.34) (4.35) (4.36) 103 which were all defined earlier, but are shown again for discussion purposes. The NTU is then NTU = M L‘2[ “mm“ ] (4.37) 0 Bimx + BimmBim‘ + Bin,“ In the previous section it was shown for a Mondt number less than .01 axial conduction was negligible. Consequently, the Mondt number must be larger to see axial conduction, and this is the parameter to begin the investiga- tion of the dimensional quantities. Figure 4.38 gives an indication of the con- ditions required on the wall and fluids (CminL‘ Iw) as a function of the wall thermal conductivity, that are needed to provide the required Mondt number. Also, for the given wall thermal conductivity the product of the convective heat transfer coefficient and wall thickness, calculated from equations (4.29- 32), are given as a function of the minimum fluid Biot number by the addi- tional labelled axes. The results in Figure 4.38 give insight into the magnitudes and trends of the basic heat exchanger parameters required for axial conduction to exist. The thermal conductivity that intuitively would improve the performance if it were larger actually increases the possibility of axial conduction because the value of (CmmL'lw) becomes larger as the thermal conductivity increases, as does the product 1:5. The larger these two groups of parameters are the more likely the chance they will be seen in a typical application. Noting that typical magnitudes of the wall aspect ratio and minimum heat capacity are at mini- mum 102 - 103, the required Mondt number will surely not be reached unless a moderate thermal conductivity exist. This result gives an indication of the conditions necessary for axial conduction to exist; it is not an ordinary occur- rence. To experience axial conduction when the wall aspect ratio is relatively large (L’ > 1000), the heat capacity must be small or the wall width must be large. These conditions require a fluid with a low specific heat (such as a gas), a creeping flow to provide the small mass flow rate, or a large heat transfer area with a correspondingly large wall width. The wall thickness is present on both axis and produces opposite out- comes for each. As the wall thickness becomes smaller, the wall aspect ratio gets larger and the possibility of axial conduction decreases because the 104 10‘ 10’ 10‘ 10’ CmL/dw 102 10‘ 10° 1o-t i 10-' 10° 10' 102 10 KW 164 . ....qu . ....qu ....qu . '11‘73- h6(Bim=.OOO1) 1(5-4 . --.-.wd . ”Hid-2 . .84 . ””"tba h6(Bim=.001) 154 . H.164 . ””164 . .T....1.bo . ......{10‘ h6(Bim=.01) 16-2 - ””164 . ”Hub . ......{b‘ . ......PO hd(Bim=.1) FIGURE 4.38. Magnitude of pertinent heat exchanger operating conditions as a function of the wallthermalconductivity,whichisneededforaxialconductiontoexist 105 required minimum heat capacity is very snmll. However, the convective heat transfer coefficient on the horizontal axis becomes larger and more likely, which is due to the scaling scheme on the fluid Biot numbers while, the former effect of the wall thickness is due to its physical role. From the previous discussion it is clear that to obtain the Mondt num- ber necessary for axial conduction to exist requires one or more of the follow- ing conditions: 1. High wall thermal conductivity 2. Small ratio of wall length to thickness 3. Larger heat transfer area and corresponding wall width 4. Low specific heat fluids 5. Small mass flow rates TABLE 4.7. Operating conditions for axial conduction to exist niL/w (kg/sec) h (W/mzK) “it B' ' Wall .01 .1 1 .0001 .001 .01 .1 Material Air 1.5E—3 1.5E-4 1.5E-5 1.5 15 150 1500 Swel Water 3.75E-4 3.7555 3.75E-6 Air 4.&-2 4.&-3 4.&-4 40 400 4000 40000 Copper Water l.&-2 LOB-3 LOB-4 = .01 Air 1515-2 1513-3 1513-4 .15 1.5 15 150 Steel Water 3.75E-3 3.75134 3.7SE-5 Air 4.&-1 4013-2 4.&-3 4 40 400 4000 Copper Water 1.&-l 1.&-2 l.&-3 = .1 Air 1.5E-1 1.5E-2 1.55-3 .015 .15 1.5 15 Steel Water 35755-2 3.75E-3 3.75B-4 Air 4.0E+0 4.&-1 408-2 .4 4 40 400 Copper Water 1.0500 1.&-1 1.&-2 106 Consider the conditions needed for axial conduction in a wall material of steel (Kw ~ 15) or copper (Kw ~ 400) and fluids of water (Cp ~ 4000) or air (Cp ~ 1000). The possible combinations of wall material and fluid type are shown in Table 4.7 at various Mondt numbers and wall thicknesses. The results shown in Table 4.7 are the mass flow rate combined with the wall length and width and the convective heat transfer coefficient as a function of the fluid Biot number. The fluids, wall materials, and wall thicknesses were chosen to give a perspective over the range of the parameters that may be encounter for typical applications. Even though a wall thickness of .1 is quite extreme, it gives an indication of the magnitudes required for other parame- ters to become reasonable. It was not possible to reduce to only the mass flow rate because the case would be too specialized. But noting that Llw is a minimum of one, and most likely larger than one; the magnitude of the mass flow rate is even smaller than the values shown for niL/w in Table 4.7. It is apparent from the data presented in Table 4.7 that small mass flow rates and corresponding convection coefficients are required for axial con- duction to exist. The required mass flow rate decreases as the Mondt number increases while the convection coefficient decreases as the minimum fluid Biot number decreases. Both trends demonstrate the requirement for axial conduction to become larger since large Mondt numbers and small fluid Biot numbers increase the effect of axial conduction on the performance of the heat exchanger. Recalling that the results predicting axial conduction were based solely on the Mondt number, the results presented in Figure 4.38 and Table 4.7 do not particularly depend on the magnitude of the convective heat transfer coef- ficient. In essence, the magnitude of the convective heat transfer coefficient is arbitrary as long as the requirement on the Mondt number is met. But for a required Mondt number the magnitude of the convective heat transfer coeffi- cient and wall dimensions will set the NTU. To obtain a particular NTU with a prescribed Mondt number and convective heat transfer coefficient may require unrealistic wall dimensions. Likewise, if the Mondt number and wall dimensions are prescribed, unrealistic convective heat transfer coefficients may be required. The magnitudes of the convective heat transfer coefficients are intended to represent the range that axial conduction was shown to exist 107 in physical situations. The range of the convective heat transfer coefficient, however, may not be physically realistic for the entire range of Mondt num- bers shown. 108 4.5 Comparison to Published Results All previous studies of the effect of axial conduction have been consid- erable less general than the present work. For this reason, there is only one specific case that can be compared to previously published results. This is the case of balanced symmetric flow in which the mass flow rates, heat capaci- ties, and convective heat transfer coefficient between the fluid and wall are all equal for the two fluids. Therefore, the heat capacity ratio is one (CR = 1) and the fluid Biot numbers will be equal (Bi = Rim“). min For this case, Pan, Welch, and Head [15] derived a closed form expres- sion for the effectiveness as a function of three variables NL = ”—611 = NTU (4.33) Kw“ hA WK 2 . ”K = (—'E':) (a) = (_C—w) 3' (4.39) 1 2 2 p =(1+IVK) (4.40) which are related to the variables of the present study as shown. Note that the subscripts on heat capacity and Biot numbers have been dropped since the conditions are the same on both sides of the wall. The effectiveness is then given ([15] equation (34)) as _ NL+ (NR/2p) tanh (pNL) - NL+ (NK/Zp) tanh (pNL) + (NR/2) +1 Using equation (4.41) the effectiveness (and ineffectiveness) was cal- culated and is shown in Figures 4.39-4.42 as a function of NTU at various minimum Biot numbers and wall aspect ratios. The results show excellent agreement with the present work. The results (i - N TU ) differ only beyond the third decimal place, which can be attributed to computational round-ofl’. 8 (4.41) The same restrictive case was investigated by Rohsenow [17]. How- ever, he proposed the effectiveness as a function of NTU and the Mondt num- ber (Mo), which was earlier shown to be KVAW wa C M: = . CL 0 L (4.42) 109 and solved for the effectiveness at the extremes of the Mondt number, ([17] equations (32) and (34)) M0 = o NTU 8 " TIN??? “'4'” Mo = on e = %(1-e'm”) (4.44) Obviously, with M, = 0 the ineffectiveness is zero, since this is the case of neglecting axial conduction by allowing the wall thermal conductivity to shrink to zero; whereas, the wall conductivity becomes infinite for M0 = co. Both of these cases are shown also in Figures 4.39-4.42, and they provide bounds on the present results. 110 — Bi_-ai__ J Bi*=.0001 X Pan, Welch, Head [15] 0'6 81:81 = 0001 + Rohsenow [17] MO-.. ‘ ° ' x Rohsenow L17] tin-=0 CI=1 L*= 100 + + 0.4“ + . + b + 500 4- 0.2- + 1000 0.0 ‘ ‘ 0 2 4 6 8 10 NTU FIGURE 4.39. Comparison of present results to past result of Pan, Welch, and Head [15] and Rohsenowblaxn ffectivene: as a function ofthe NTU and wall aspect ratio for — Bin-Bi“ Bi*=.001 3! Pan, Welch. Head [15] 0-5" Bi =3; =.0001 + Rohsenow [17] Mos—o- ‘ ° X Rohsenow [lzLMQ=0 4 C.-.-1 ... + + + + 0.4" + L*= 100 FIGURE 4.40. C of present results to past result of Pan, Welch, and Head [.15] and Rohsenow [17]. Inc ectiveness as a function of the NTU and wall aspect ratio for B t = O. 001 111 . — Bi_=Bi_. Bi_.=.01 3r Pan. Welch, Head [15] 0'60-4 BlL=Bio=.0001 "l’ Rohsenow [17] M0309 X Rohsenow [17] MO=0 q Cn=1 + 1’ + + 0.40-« + 1 , + a + ‘I’ 0.20q + L*= 100 0.00 NTU FIGURE4.41. ComparisonofpresentresnltstopastresultofPan,Welch,andHead [15]and mmwuplnsffectivenessmafunctlonoftheNTUandwallaspectratiofor Bin“ =0.01 R — a 0.701 — ai_=ai__ Bug-7.1 * Pon, Welch, Head [15] Bi =5; =.0001 + Rohsenow [17] M08- " ° X Rohsenow [17] Mn=0 050-] C.-1 + "' + + + + + a 0.30- 1’ 9’ + 0.10- + + 'I- 100 W— #4 x x ’0-10 I l l l l 0 2 4 6 8 10 NTU FIGURE“2.C0mparbonofpresentruultstopastresultofPan,Welch,and Head [15] and {Isabenowfl‘ln .IneffectivenessasafunctionoftheN’I‘Uandwallaspectratiofor Bi R- ruin" Chapter 5 Conclusions and Recommendations for Future Work An exact solution for the analysis of the heat transfer occurring in a counterflow heat exchanger including the efiect of axial conduction has been presented. In addition to the solution methodology a computer program to evaluate the numerical aspects of the exact solution was given. The program allowed for the investigation of the influential parameters in the solution and verification of the solution’s validity. The role of all describing parameters was investigated to determine their importance in affecting axial conduction. The effect of axial conduction was quantified by the amount that the perfor- mance of the heat exchanger was degraded in comparison to the case of neglecting axial conduction. The following results and conclusions are drawn from this investiga- tion: 1. Axial conduction has less effect on the performance of the heat exchanger at larger wall aspect ratios. 2. Axial conduction has less effect on the performance of the heat exchanger at larger fluid Biot numbers. 3. Axial conduction is reduced at smaller heat capacity ratios. 4. The effect of the end Biot numbers on the performance of the heat exchanger is negligible for magnitudes less than or equal to the magnitude of the fluid Biot numbers. 5. Axial conduction is negligible for magnitudes of the Mondt number less than .01. 6. High wall thermal conductivity promotes axial conduction. 7. Axial conduction is more likely for large heat transfer areas. 8. Applications with low heat capacity fluids promote axial conduction. The most obvious continuation of this work would be to experimen- tally verify the results presented herein. Although the efl’ect requires some 112 113 atypical conditions, it was shown possible for common fluids. The applicabil- ity of this solution to an experimental set-up, which would typically be con- centric tubes, will require the set-up to be approximated by the model, a plane dividing wall. Therefore, it must be possible to “cut and unroll” the concen- tric tubes and represent them as parallel plane walls without incurring an appreciable difference in the areas on opposite sides of the wall. In addition, the issue of identifying the type of boundary condition to model the physical problem at the ends of the wall will need to be addressed. The computer program AXCOND could be coupled with existing ther- mal design software to predict the performance of a counterflow heat exchanger. Use of this program will be computationally expensive but will only need to be used when certain conditions exist, otherwise the traditional effectiveness-NTU relationships can be used. Incorporating the program into existing software will require conversion of the main program into a subrou- tine and consideration of the number of terms needed to produce a converged solution. Also, if the required analysis is a sizing problem, solving for the heat exchanger dimensions needed for a particular operating condition, an iterative scheme needs to be implemented. The issue of the number of terms needed to obtain an accurate solu- tions could also be further investigated. The outcome of this additional work may lead to possible simplifications if only a small number of terms are needed. Furthermore, simplified equations may result under certain condi- tions, which may allow for approximate solutions. In general, the complex solution presented may be studied for possible simpan cases. The solution presented could also be modified rather easily to deter- mine the performance of a parallel flow heat exchanger. The formulation of the problem for the heat exchanger wall is exactly the same. The formulation of the fluid energy balance equations will differ by a sign and the location of the specified initial condition. The first difference will infiltrate the applica- tion of orthogonality to apply the nonhomogeneous boundary conditions on the wall and change one of the fluids temperature solution, but only by sign change(s). Whereas the initial condition will change the form of the integra- tion constant on one of the fluid temperature solution. The conversion may be tedious, yet would require minimal computational efl’ort or programming changes since the outline of the solution and program are given. LIST OF REFERENCES [1] [2] [3] [4] [5] [6] [7] [8] [9] [10] [11] [12] [13] 114 Beyer. W. H. Wallabies. New York: CRC Press Inc., 1987 Campo, A., and J. C. Auguste, “Axial Conduction in Laminar Pipe Flows with Nonlinear Wall Heat Fluxes.”11;te;n_ational_193gal_gf_lieat_ MW; 1987, volume 21, p. 1957-1607 Chiou, J. P., “The Advancement of Compact Heat Exchanger Theory considering the Effects of Longitudinal Heat Conduction and Flow Nonunifonnity.” Compag Heat Exghangers, Chiou, J. P., “Thermal Perforrnanoe Deterioration in Crossflow Heat Exchangers Due to Longitudinal Heat Conduction in the Wall.” ASME paper no. 76-Wa/HT-8 Chiou, J. P., “The Effect of Longitudinal Heat Conduction on Cross- flow Heat Exchangers.” f H Tr f r May 1978, volume 100, p.346-351 Edwards D. K., V. E. Denny, and A. F. Mills, W, New York: Hemisphere Publishing Corporation, 1979 Eisensmith, S. P., W, September 1990 Hewitt. G. F.. W. New York: Hemisphere Publishing Corporation, 1990 Incropera. F. P.. and D. P. Dertt. W New York: John Wiley and Sons, 1985 Kays W. M., and A. L. London, W, New York: McGraw Hill Book Company, 1984 Landau, H. G., and J. W. Hlinka, “Steady State Temperature Distribu- tion in a Counterflow Heat Exchanger Including Longitudinal Conduc- tion in the Wall.” ASME paper no 60-Wa/HT-236 Mondt, I. R., Effects of Longitudinal Thermal Conduction in the Solid on Apparent Convection Behavior, with Data for Plate-Fin Surfaces.” WM; 1960-61. P 614-621 Nyhoff. L.. and S. Leestma. W New York: Macmillan Publishing Company, 1988 [14] [15] [16] [17] [18] [19] [20] [21] 115 Ozisik M. N., W90, New York: John Wiley and Sons, 1980 Pan, C. L., N. E. Welch, and R. R. Head, “A New Concept to the Gen- eral Understanding of the Effect of Longitudinal Conduction for Multi- stream Counter Flow Heat Exchanger.” Proceedings of the 3rd International Heat Transfer Conference, 1966, volume I, p159-169 Pan, C. L., and N. E. Welch, “Exact Analytical Wall and Fluid Temper- ature Field for a Counterflow Heat Exchanger with the Effect of Longi- tudinal Heat Conduction.” ASME paper no. 65-HT-63, 8th National Heat Transfer Conference, August 1965 Rohsenow, W. M ., “Why Laminar Flow Heat Exchangers Can Perform Poorly.” W, New York: McGraw Hill Publishing Corpo- ration, 1981, p. 1057-1072 Shah, R. K., “A Correlation for Longitudinal Heat Conduction Effects in Periodic-Flow Heat Exchangers.” Journal of Engineering For Power, 1875, volume 97, p. 453-454 Somerton, C. W., Notes on Heat Exchanger Experiment, ME413 Course Notes, Michigan State University, 1990 Somerton, C. W., Personal Communication on Axial Conduction, 1990 Wark, K., Thermodynamics, New York: McGraw-Hill Book Company, 1983 APPENDICES APPENDIX A A.1 General Solution of the Wall Conduction Equation Beginning with the partial differential equation for the wall and bound- ary conditions that were derived in Chapter 2 flq-L'l'z 326” = O ax+2 ay+2 80 ”—7 +Bi00w(0,y+) = O x e x =0 W +BiL0w(l,y+) = BiL(TL—To) 3):“ x. =1 30w --—+' +3139... (1+a0) = BiHITH(x+) -T0] 8y y’=0 aew , + . . 5? ”teem .1) = B‘CITCO‘ >-Tol y‘=1 (A.l) (A.2a) (A.2b) (A.2c) (A.2d) The problem is split into three simpler problems each with only one nonho- mogeneous boundary condition. Problem 1 329‘ +L"‘2 9:31- = 0 ax+2 ay+2 301 ' -—+ +Bi091(0,y+) = 0 81 {=0 1 $7 +BiLe,(1,y*) = Brim-T.) x’=l 801 . + 73—; +Bz,,el(x ,0) = 0 y ,.=0 116 (A3) (A.4a) (A.4b) (A.4c) 117 +3.09, (x+, 1) = o Problem 2 3’9, 3’9, 2 .. ax+2+Lit ay+2 - O + BiaG2 (0, y+) = 0 x’=0 +BiL02(l,y+) = o x’=1 892 ay" 4495102 (x*,0) = BiH[TH(x+) —To] y‘ =0 892 +3369, (x", 1) = 0 y'=l Problem 3 3’9, ”*2 air, = o ax+2 ay+2 303 +Bi003(0,y+) = o +BiL03(1,y+) = o +Bi303 (x+,0) = 0 y’=0 __1 a)!“ «t-Bice3 (x‘, 1) = Bicch(x+) -To] y‘ =1 (A.4d) (A5) (A.6a) (A.6b) (A.6c) (A.6d) (A.7) (A.8a) (A.8b) (A.8c) (A.8d) It can be shown by adding the three problems in Eq. (A.3-8) that the original problem is obtained in Eq. (A.1-2). The solutions are related as follows: 118 0“, (x1, y”) = 01(x+, y“) + 02 (x‘, y”) + 03 (x+, y+) (A.9) A.l.l Solution for Nonhomogeneous Boundary Condition at x“ = 1 (Problem 1) Assuming a product solution for 01 of the form 01 = F(x”) G(y*) (A.10) Substituting into the partial differential equation and boundary conditions, Eq. (A.3-4), gives F , .G+L’ZG . .F = 0 (A.11) x X y y -GF . +BioGF(0) = 0 (11.123) x x1 =0 orx. e =1+BioGF(l) = BiL(TL- To) (A.12b) -FG . +Bi,,FG(0) = 0 (A.12c) y ’9 = o my. +Bi,,ro(1) = 0 (1112.1) r’ = l with the following convention for representing derivatives used: air ax” - Fm. (A.13) Rearranging Eq. (A. 1 l) to group similar variables and simplifying Eq. (A. 12) gives _L Ff” ._ _GY'Y _ :l: 2 (A14) Leg 1; - G - 11 ° —F, +Bi,r(0) = o (A.15a) x x‘ =0 3 = -G . +Bi,,G(0) = 0 (A.15c) y ’930 G . +Bi,,G(l) = 0 (A.15d) 1' f =1 119 where a constant, 1:111, was introduced in Eq. (A. 14). With the left hand side only a function of x“ and the right hand side only a function of y+ for equality both sides must be equal to a constant, hence the introduction of the constant. This key result allows for separation of the variables, reducing the single par- tial difl‘erential equation to two ordinary differential equations. Before Eq. (A. 14) can be separated however, a sign must be chosen for the constant. This sign is chosen to produce an eigenvalue problem for the function that has two homogenous boundary conditions to evaluate. The eigenvalue problem will have a general solution in terms of sine and cosine functions for this cartesian coordinate system. After choosing the sign, two ordinary differential equa- tions can be written from Eq. (A.14) with the appropriate boundary condi- tions from Eq. (A.15) F . . -L'2u2F = 0 (A.16) x I -r. +Bi,r(0) = 0 (11m) x x’ 80 G . . +1126 = 0 (A.18) Y y -G, +BiHG(0) = o (A.19a) y y9=o Gy. . +Bi,,G(1) = o (A.19b) y =1 Note that Eq. (A.16) only has one boundary condition, whereas Eq. (A.18) has two boundary conditions. The reason for this can be seen in Eq. (A.15b), which depends on both functions F and G and cannot be applied to either function singularly. Thus, this boundary condition must be applied after assembling the complete solution. Eq. (A.16) and Eq. (A.18) have easily obtainable general solutions. The general solution for Eq. (A. 16) is r(x+) = Alcosth’xt) +A2sinh(uL‘x+) (.4120) Applying the boundary condition, Eq. (A.17a), gives a relationship for the two constants 120 A2= "HA, (A.21) which can be substituted into Eq. (A20) and rearranged giving the solution F(x”) = A, [pL‘ cosh (uL‘xt) +Biosinh (pL‘x*)] (A22) The solution of Eq. (A.18) is G (yl) = A3608 my") +448in (11f) (A.2a) Applying the boundary conditions in Eq. (A.19a) solves for the relationship between the constants BiHA3 A = .24 4 u (A ) Substituting into Eq. (A23) and rearranging gives 60*) = A. [ucos (11)”) +Biusin (uy“)] (A25) Applying the final boundary condition, Eq. (A.19b), does not provide any information about the constant A 4 because the constant cancels, it does pro- vide information about the constant )1 though, the equation simplifies to (11" + -BtHBrC) To meet this boundary condition Eq. (A.26) must be satisfied for all values of u". Eq. (A.26) is called a transcendental equation, it has an infinite number of solutions, hence the subscript on u, and will be solved to determine the possi- ble values of u". The solution for the function G is then (A.26) 601,. f) = 112111.003 01,3”) + Biusin (u,y*)] (11.27) and u" is given by the solution of Eq. (A.26). Putting Eq. (A27) and Eq. (A.22) into Eq. (A.10) will give the solu- tion. However, since there exist many solutions to Eq. (A.27) all possible solutions will be summed to obtain the final solution. Also, the undetermined constants, A2 and A ‘, will be grouped into one constant that will be deter- mined later by evaluating the last boundary condition, Eq. (A.4b), which was the boundary condition that did not separate and was not evaluated. The solu- tion for problem one is 121 (11 (x+, m = .1,,[unL’ cosh (unz‘xt) + Biosinh (an‘x+)] * n = 0 [mos My“) + Bigsin (u,y*)] (A28) and the eigenvalues are found from the solution of Eq. (A.26). A.1.2 Solution for Nonhomogeneous Boundary Condition at y+ = 0 (Problem 2) Using similar methods problem two can be solved. Defining an assumed product solution 02 = F(x*)G(y*) (1129) and substituting into the partial differential equation and boundary conditions and choosing a sign to make the proper eigenvalue problem allows the partial differential equation to be separated into two ordinary difierential equations F . . +0121r = 0 (A30) x x "Fe +BioF(0) = 0 (A.3la) x’ = 0 F. +BiLF(1) = o (A.31b) x x‘ =1 G . . -a20 = 0 (A32) I y G . +Bi,,o(1) = 0 (11.33...) Y y. =1 The solutions are of the form, respectively F(x*) = Blcos (ax‘) +stin (of) (11.34) + (I + - (l + G(y ) = B3cosh(-Ey )+B4s1nh (Ey ) (A.35) Applying the boundary conditions in Eq. (A.3la-b) to Eq. (A.33) gives F (16‘) = Br [ancos (an?) + Biosin (anx*)] (A36) where the values of an are obtained from the solution of the transcendental equation 122 (Bi, + 31,) an (ah-aim» Applying the boundary condition given in Eq. (A.33a) to Eq. (A.35) gives tan (01") = (A.37) 0 (y‘) = B.[cosh (1%)“) -Csinh (fl-fl] (1.33) where C is a constant, given by fimflb +Biccosh (3;) c = L L L (A39) (1 a a -;cosh (-—;) +Bi sinh (-;-) L L C L As before, putting the two solutions together and noting that the solution to Eq. (A36) is the sum of all possible solution, and lumping undetermined constants together allows the solution to problem two to be written as " or or 02 (x+. y“) = 2 Bn[cosh (4f) - 9.8m (_fy‘t )] e n = o L L [ancos (of) + Biosin (anx*)] (A.40) and the values of an are obtained fi'om the roots of Eq. (A.37). A.1.3 Solution for Nonhomogeneous Boundary Condition at y+ = 8 (Problem 3) Finally, applying the method again to solve problem three. Defining an assumed product solution 03 = F(x") G (y”) (A.41) and substituting into the differential equation and boundary conditions and choosing a sign to make the pr0per eigenvalue problem allows the partial dif- ferential equation to be separated into two ordinary differential equations F . .+a2F = 0 (A42) 1 x - r1. +Bi,r(0) = o (A.43a) x’ = 0 17‘. + BiLF(1) = o (A.43b) x’ = l 123 G . . -—azG = 0 (A44) H G . +Bi,,o (0) = 0 (11.45.) y y. g 0 Realizing that Eq. (A.42-43) is exactly the same problem solved in the previ- ous section, Eq. (A.30-3l), the solution is F(x+) = Cl [ancos (aan) + Biosin (anx+)] (A46) where the values of an are obtained from the solution of the transcendental equation in Eq. (A.37). The solution to Eq. (A4445) is G(y ) = C3[:—cosh(— —y )-1-Bil,,sirrh(i1 y )] (A47) Putting the two solutions together and noting that theL solution to Eq. (A.46) is the sum of all possible solution, and lumping undetermined constants together allows the solution to problem two to be written as or e (x+,y CnLl:-§cos 1):!“ )+Bi sinh(—,y+)]* 3 )=n§0n Ly H L [ancos (an?) +Biosin (anx*)] (A.48) and the eigenvalues are given by the roots of Eq. (A.37) A.1.4 Summary of the Solutions The solutions for the three problems without applying the final (nonho- mogeneous) boundary conditions are: 01(x+,y+) = 2An[ttnL'cosh(unL.x+) +Biasinh(unL‘x+)]* n=0 [u,cos (11,.y“) + Bigsin (u,y*)] (A.28) 92 (1+, ”.20" a "[cosh (L, y )- ; "WGWH' [ancos max“) +Biosin (anx*)] (A39) 124 " or or a 0 x1. + = Cl:—1-'cosh(-';' *)+Bi sinh(—§ +)]* 3( y ) "go a L L y H L y [ancos (up?) + Bio sin (anx*)] (A47) where the eigenvalues u" and a" are given by the positive roots of the tran- scendental equations (Bi +Bi )u taunt.) = 2" , C, " (A.26) (nu-BrHBrc) (Bi +Bi )(X mum") = ° L " (A.36) (a: + —BiLBi0) The remaining boundary conditions to be applied to determine the con- stants A", B", and C", are Eq. (A.4b), Eq. (A.6c), and Eq. (A.8d) and they depend on the fluid temperatures for B, and C", therefore all will be applied after the fluid temperatures are determined in Appendix B. The solution for the wall conduction equation is given by Eq. (A.9) as the sum of the three solutions and the unsealed temperature (needed to solve for the fluid temperatures) is T,.(X*.y*) = 9,,(x*.y+) +T, = 0101.)”) +92(x*.y+) +93(x+.y*) +7}. (448) APPENDIX B B.l Solution of Hot Fluid Energy Balance Taking the describing differential equation for the hot fluid that was derived in Chapter 2 and rearranging gives d1” + 1,, (0) = 1,“, (B.2a) The expression for the unsealed wall temperature (from Appendix A) is 1,, (x*.0) = 10+ 2.1.11" [unL'coshan’xU +Biosinh (an36)] 71:0 " or + 2 (Bn+L—:Cn)[ancos(anx*) +Biosin (anx*)] (B3) n-O Putting this into the differential, Eq. (B1), and arranging gives dTH - a s + . . t 4' Ewing, = NH{T0+ 211,41”an cosh (an x )+Brosrnh(|.tnL x )1 n=0 " or + 2 (Bn+i-:C.)[ancos(anx*) +Biasin(anx+)] } (3,4) n=0 To solve this problem it was first split into two problems, a homogeneous problem dTH, h dx+ “var“ = o (13.5) and a particular problem (1TH, p d1:+ +~,,Tu.. = N..{T.+ 2'. .u.tu,.L‘cosh(u,L’x‘) +Bi.sinh(u.L'x*)l n=0 .. an + - - + + 20(Bn+17cn)[ancos (crux )+Brosm(ornx )] } (8.6) 125 126 The particular problem in Eq. (B6) was further divided into three problems for each term on the right hand side of the equation. The general solution for Eq. (B. l) in terms of the simpler problem’s solutions is 3 7“,,(x+) = Tfl’h(x+) + z Tmp‘(x+) (3.7) i=1 Where 7'”, p, is the solution to the problem dTH, p, + +NHTH,pi = NH{5ttT, + 52, 2 Ann” [u'L' cosh (unL'f) + Biosinh (11,5151 u=0 +53, 20(Bn+;§cn)[ancos(aux+) +Biosin(anx+)] } (13.3) (i = l, 2, 3) and 5], is the kronecker delta function as defined in Eq. (3.7). Superposition of the problems in Eq. (BS) and Eq. (B.5) will give the origi- nally posed problem in Eq. (B.4). B.1.1 Hot fluid Homogenous Solution The homogenous problem in Eq. (B5) is solved giving TH, ,, (x+) = chxp (N”x+) (13.9) where D1 is a constant to determined after assembling the complete solution as shown in Eq. (B.7). B.1.2 Hot fluid Particular Solution 1 Solution of the particular problem in Eq. (B8) was obtained using standard variation of parameter techniques. This method substitutes a general solution of functional form similar to the nonhomogeneous term in the prob- lem then solves for the unknown constants in this general solution. For the first particular solution (i=1 ) the general solution would be a constants, TH,” (x1) = f(To) = k = constant (B.10) 127 Substituting this into the difierential equation Eq. (B.8) yields NH]: = NHTO (3.11) and the solution is T”, p, = To (3.12) B.1.3 Hot fluid Particular Solution 2 The functional form of the second particular solution (i =2) is T“, = 2 bacoshan'x") +cnsinh(|.tnL.x+) (B.13) n=0 Substituting into the differential equation, Eq. (B.8), (note the summations are not shown) gives buan‘ sinh (nave) + cnunL. cosh (anti) NH [bucosh (unL‘xt) + cnsmh (unL'f) 1 = NHAnun [unL' cosh (uuL‘xt) +Biasinh (an’x+)] (3.14) For this particular solution to be correct equality must hold, which is accom- plished through the constants b" and en. Grouping all constants for the sinh and cash terms separately into two equations produces bnunL' +1)!ch = NHAnttnBio (3.15) for the sinh terms and cnunL‘ 1'”an = NflAnufiL’ (B.16) for the cash terms. Solving these two equations for the unknown constants . (BiN - 2L‘z) bn=AnufiL 1- "2” 5”,: (3.17) (NH-11314 ) (33,117,, - 33L”) ] c =ANHu , (3.18) " " "[ (Ni-uiLz) and substituting into Eq. (B.13) gives the second particular solution 128 " . (3.- N wit”) . . Tmnw) = 2.1,{332 [1- (1:219:32) :Icosh (11,,L x ) H— n (31,111,, - ufiL‘z) (Ni-ufiL") +Nulln[ ]sinh (unL’f) } (3.19) B.1.4 Hot fluid Particular Solution 3 The functional form of the third and final particular solution (i=3) is 1,,” (x+) = 2 aucos (anxt) + dnsin (a,x*) (3.20) n = 0 Substituting into the differential equation —a,,ct,sin (0.3”) +dnctucos (an?) +NH[a,cos (our?) +d,,sin (anx+)] or = NH[B,+ i—an) [ancos (uni?) + Bi,sin (anx+)] (B21) then separating into two equations, the coefficients are a dud" + NHa, = NH (B, + Egg)», (3.22) for the cos terms and all . - and" + N ”d" = N H (B, + Z;- C")Br, (B23) for the sin terms. Solving for the unknown constants 0. (Bi N + a2) a" = (B,+-£C,)an[1- "2” 2" :l (3.24) L (N ”-1- or") or Bi ”3* a2 d = (B +—-f-C )N ° " (3.25) " " L " ”[ N§,+ a}, and substituting into Eq. (B.20) gives the particular solution 129 " ct (Bi NH-t-a’) 711.33 0+) = 2 (Bn+i-§Cn){au[l- (13244;; ]cos (anx+) H n +N 3i,N,,+a3, H N2+a2 H n :lsin (any) } (3.26) Assembling the solutions to the simpler problems, Eq. (B.9), Eq. (B.12), Eq. (B.19), and Eq. (B.26), using the relationship given in Eq. (B.7), the solution for the hot fluid energy balance is TH (15*) = T, + Dlexp (-NHx+) “ , (Bi,NH-|.t:L'2) , + + 2A,,{nfiL [1- (Nil-113132) :lcosh (3,1. x ) (BLNH-uiL’z) . N inh L + + ”u"[ (Na-33L”) ls ('1’ “l “ B'N 2 + 2 (Bn+a—2Cn){au[l-( r, H+a")]cos(anx+) (N3, + 01:) 3i,1v,,+ a: + N” 2 sin (our) (3.27) NH-t- or: Applying the boundary condition, Eq. (B.2a), the constant is °' , (Bi,NH-112L‘2) D=T.—T- AuzL 1- ", l H.173 o 2 n n [ (Ni-11:1! 2) n=0 - 0‘. (3i,N,,+ a3) - 3 —c 1 - .23 n§0( "+1: a)0‘n[ (Nificri) ] (B ) 130 B.2 Solution of Cold fluid Energy Balance The differential equation that was derived in chapter 2 for the cold fluid, after rearranging is are , -;: +NCTC = NCTW (x ,1) (B29) TC(1) = TC, , (3.30.1) The wall temperature solution (from appendix A) is 7,06, 1) = 70+ 2 A,[n,L’cosh(n,L‘x+) +3i,sinh(nnL‘x+)] * n=0 [uncos (11,.) + Biysin 05.)] +2:~[cosn(:—r)-c~m(%)]=~ [ancos (mg?) + Bi,sin (anx+)] ”Canha" B’sinha"* + — -— + — .2. "[L‘ M (L) '” (L' )] [ancos (an?) +Bi,sin (0131)] (B31) Defining the following constants: A,’ = A, [uncos (11») + BiHsin (11") ] (B32) 0: or 3,; = 3,,[cosh {—2) — gen). (4')] (3.33) L L or a or. C,’ = Cn[—: cosh (4) + Bi ”sinh (—f )] (B34) L L L simplifies the expression for T, (x”, 1) to Tw(x*, 1) = 70+ 2A,’[ttnL.cosh(ttnL‘x*) +3iosinh(n,L“x+)] n=0 + 2 (Bn'+C,,') [ancos (crux‘) +Bi,sin (anx+)] (3,35) n=0 131 Putting this into the differential Eq. (B29) and arranging gives dTC . , O 3 + o . . + _quvcrc = NC{T0+ 211nm}. cosh(|.tnL x )+Bzosmh (unL x )1 n=0 + 2 (3,; + C,’) [ancos (any) + Biosin ((13+)) } (3.36) n = 0 To solve this problem it was first split into two simpler problems, a homoge- neous problem ch (hf +NCTC’,, = 0 (3.37) and a particular problem dTC’ p + + NCTC’ p = NC {To + Z A.’ [unL' cosh (11.13:?) + BioSinh muff“ n=0 + 2 (Bn’+ Cn’) [ancos (anf) +Biosin (anx+)] } (B38) 3 = 0 The particular problem in Eq. (B3 8) was further divided into three problems for each term on the right hand side of the equation. The general solution for Eq. (B29) in terms of the simpler problem’s solutions is 3 raw) = Twat) + 2 Tam-(f) (3.39) i=1 Where Tc, p, is the solution to the problem given as ”Q“ NT -N 51‘ 17+ cc,p:- c{uo + 52: X A,’ [unL' cosh (unL'f) + Biosinh (an‘r‘n n=0 +5”; (3,,’+c,,') [ancos(anx+) +Biosin(ocnx+)] } (3.40) n=0 (i = 1, 2, 3) and 5;.- is the kronecker delta function as defined in Eq. (3.7) 132 32.] Cold fluid Homogenous Solution The solution to the homogeneous problem in Eq. (B.37) is T0,). (x+) = Elexp (NC?) (B.41) and E1 is a constant to be determined, after assembling the complete solution. Superposition of the problems in Eq. (B40) and Eq. (B.37) will give the orig- inally posed problem in Eq. (B.36). B.2.2 Cold fluid Particular Solution 1 Assuming the first particular solution ( i=1 ) is a constant and substitut- ing into the differential equation Eq. (B.40)and solving for the magnitude of the constant gives Tm,1 (x*) = To (13.42) B.2.3 Cold fluid Particular Solution 2 The form of the second particular solution (i =2) is Tamer) = 2 bucoshan'f) +cnsinh(unL’x+) (3.43) n=0 substituting into the differential equation Eq. (B.40) (with summations not shown) gives bnan‘ sinh (uuL'JF) + enunL‘ cosh (puL‘xt) - Nc[bncosh (fluff) + cnsinh (unL'x+)] = —NCA,,'[unL‘cosh (unL'f) +Biasinh (unL.x+)] (3.44) and the coefficients of the sinh terms and cash terms are, respectively bnunL’ -Nccu = -NCA,’Bia (13.45) canal: - NCbn = 'NCAn’unL. (13°46) Solving for the unknown constants 133 . N 1' MPH2 b" =31;an [1+ Cf" 2".2 ] (3.47) NC—unL N Bi +1121." c" = An’NC[ C2 ° 2",2 J (B.48) NC-uflL and substituting into Eq. (B.43), the second particular solution is " , NCBi +1125" , Tcfloc”) = 2A,;{an [1+ N2:112£‘2 ]cosh(|.1nL x+)+ C n=0 1! [NCBio + 33L” C big—11%” ]sinh(unL x*) } (3.49) II B.2.4 Cold fluid Particular Solution 3 The last particular solution is of functional form Tc’p3(x+) = Z ancos(anx*) +dnsin(anx+) (3.50) n = 0 The resulting equation after substituting into the differential equation (with- out summations shown) is (-anan) sin (Olaf) + dnancos (uni?) - NC [ancos (an?) + dnsin (anx+)] = —NC (Bn’ + Cn’) [ancos (an?) + Biosin (anx+)] (B31) Separating into coefficient equations for the sin and cos terms - an“, - Ned, = -NC (Bn’ + Cn’) Bio (352) dud" - Ncan = -NC (Bn’ 4» Cn’) an (3.53) and solving for the unknown constants (13.54) N '- 2 an: (Bn’+Cn’)an[l+ CBIO an] 2 2 ”CHI, 134 d (B ’ C ’)N "emf“: (355) = + —— . n n u C N2- + a: and substituting into Eq. (B50) gives the final particular solution + ” I I NOEL-0 - a: + TC,p3 (x ) = "go (8,. +C,. ) {an[l + $441: ]cos (aux ) NCBio - (xi + NC ‘7‘“?- sin (any) (B.56) NC + an Assembling the simpler problem’s solutions, Eq. (B.41) Eq. (B .42), Eq. (B.49), and Eq. (B56), using Eq. (B.39) the'solution for the temperature of the cold fluid can be written as TC(x+) = T0+Elexp (New?) '° , NCBio-t-p’L'z] , + A’ uL 1+ '1 cosh(uLx+) ago 8{ n I: [Viv-[1:14 2 II NCBio+p§L‘2 , +N . sinh( Lx") C[ ”is-nil- ” p” r N Bi -a2 + 2 (3,;+C,,’) an l+—CE°—-2—" cos(anx+) "=0 Newt-an NCBio - a: . + +NC[ Ni. +01: :lsmuxnx ) (3.57) By applying the boundary condition in Eq. (B.30a) the unknown constant E1 can be determined 135 " . NCBio-t-uiL'z . E, = exp(-NC) (Tm—To) - 2A,; uuL 1+ 2 2 .2 cosh(|.tnL ) n=0 :3 N Bi +11%” , +Nc[ 02° 2",2 ]sinh(unL)} NC-uflL " N ' -a2 - 2 (Bn’+C,,’){an[l+ 2131:“: ":lcosmn) C n It: 0 Nam-a3] . H + NC[ sm (an) (3.58) 2 2 Nc+an 136 B.3 Summary Energy Balance Solutions This section does not provide any new information it instead summa- rizes the solutions obtained in this appendix. The equation numbers will aid in locating the details of these solutions. B.3.1 Hot Fluid Energy Balance Solution Solution for the hot fluid temperature is TH(x+) = To-t-Dlexp (-NHx+) 0- . (BioNH - ”3L. 2) :l e + A 32L 1- , cosh(|.le+) .20 "{ " [ (33,—ng 2) " (BioNH - ufiL‘Z) (33,-1133‘2) .- B. N 2 + 2(Bn+g;cn){an[1-( I" "+a")]cosm,,x*) (N3, + 01%) +NHP..[ ]Sinh (lint?) } BioNH + (:2 "Jsin (anf) } (3.27) +N "[ ”31+“: where (BioNH - 33L”) ] D1=T,“-T— ApZL‘[1 , "' ° .21. " " (Na-3:1. 2) " a. (BioNH+ afi) - 3 —c 1 - . Ell " + L‘ ")a”[ (N§,+ 0:) ] (B 28) 137 B.3.2 Cold Fluid Energy Balance Solution The solution for the cold fluid temperature is Tc(x+) = T0+Elexp (NC-x”) " , N i+ ZL‘2 , +ZA,'{IJ,,L [1+ CB” l1" Jcoshut'L x‘) N26 - 113,1." n=0 NcBio + 1131.” + NC N2 _ 21:2 C ”I! Jsinh (H133) } " NCBio-afi + (B ’+C ’) a 1+ cos (a x+) E." ~{~[ ~33] ~ New), - a: + No 2 sin (of) (3.57) NC + a: where New, + 33L” N3 - #3132 E1 = exp (-NC) [(TC,in - To) - Z An’{uul‘. [1 + :ICOSh flint.) n=0 NCBio+11§L‘2 +NC 2 2 '2 NC-unL ]sinh(u,L‘) } - E B +C, a + COS a 11:0 II n n ”Cl “)2; n N Bi -a2 +NC C2" 2" sin(0tn) (3.53) Nc+an APPENDIX C C.1 Application of Nonhomogeneous Boundary Conditions In Appendix A the solution of the wall conduction equation describing the temperature distribution in the wall was solved to the point of applying the nonhomogeneous boundary condition for the three simpler problems that the original problem was split into, because two of the three boundary condi- tions depended on the temperature of the fluids. With the solutions for the fluid temperatures complete (Appendix B) the final boundary conditions can be applied. C.1.1 Application of Orthogonality at x” = 1 The solution for 61 up to the point of applying the nonhomogeneous boundary condition, and the boundary condition (Appendix A Eq. (A28) and Eq. (A.4b)) are 91(x+,y*) = 2An[uuL'cosh(uuL'x+) +Biosinh(|.tnL'x+)]* n80 [u,.cos day“) + Bil-(Sin (u,y”)] (c1) 36 57: ”31.91 (1, y+) = BiL(TL- TO) (0.22)) x’ =1 Taking the derivative of Eq. (Cl) and evaluating the expression and it deriv- ative at the boundary gives 91(1,y+) = ZAuan'coshwnL‘)+3iosinh(nnL‘)]* u=0 [1130s (11,?) + Biysin (u,y*)] (C3) 391 -+ _ 2.2. t . t c * —(1..V ) - A, [11,14 smh(unL )+Bloll,,L 00511 (11,1! )] ... ax n=0 [uncos My") + Bigsin (u,y*)l (CA) 138 139 Substituting these expressions into the boundary condition, Eq. (C23), pro— duces A, [ufiL‘zsinh (unL') + BiouuL' cosh (unL’) 1 * - 0 [woos (11.3”) + Biysin (unyin + 31L 2 A, [unL' cosh (unL') + Biosinh (unL')] * n=0 [11,608 (113*) + Biysin (11M) ] = BiL(TL- To) (05) The orthogonality of the eigenfunctions was used to solve for the constants, A". The eigenfunction for this problem is Y (nu. y*) = new (113*) + Bigsin (unyi) (as) and these eigenfunctions have the property that 0 natm l {wwi’muwhY(u,,..y*)dy+ = { N01,.) n=m ((2.7) where wxy is a weighting constant and for the Cartesian coordinate system it is equal to one. Therefore, Eq. (C.5) will be multiplied by a second eigen- function Y (11”, y*) and integrated over the boundary. Because this integral is only nonzero for m=n that is the only term that will remain from the summa- tion. The integrals for each of the three terms in Eq. (C5) are 1 a. I { 2 A. [ufiL’zsinh (unL’ ) + 31011,!) cosh (unL’ )] * 0 n=0 [(1,608 (fly) + Bigsin (u,y*)] [uncos (umf) + BiHsin (umy+)] }dy+ __- A,” [uiL’zsinh (me‘) + BioumL' cosh (11mm) No.1,") ((2.3) 140 l - j {31L 2 A, [unL’ cosh (nut) 4.3103311) (11_L‘)] * 0 n=0 [uncos (11f) +BiHsin (11.3”)1 [umcoswmyh +191},,sin(tl,,.y*)]aly+ } = AmBiL [me‘ cosh (umL‘) + Biosinh (umL') 1 N (11,") (C9) 1 Biz. (TL- T,)j[u,,,cos (umf) +Biusin (twin ally+ = 0 _ . Bi H Bi H 31L (TL - To) [3m (um) - E—cos (um) + T] (C.10) After applying orthogonality Eq. (C5) is A," [piL'zsinh (umL') + BioumL' cosh (umL' ) 1 N01") + AmBiL [umL' cosh (umL' ) + Biosinh (me‘) ] N (ll...) = . Bin 31“,, Eli (TL- To) [sin (11,") - Tcos (um) + r] (C.11) which can be solved for the unknown constant Am B ' B' BiL(TL— To) [sin (11 )- ficos (u )+ it] III ”In In “In A _ (C.12) - [ (11:1.‘2 + 31,311) sinh (umL') + (310 + 35) umL‘ cosh (me‘ ) ] N (11,") The functional relationship for N01,") is N(|.1) = l[(uerBiz) 1+ Bi" +Bi (€13) "‘ 2 "' ° ( (1133832)) °] ' Having this information everything is completely known for evaluat- ing 91. It would be possible to substitute Eq. (C.13) into Eq. (C.12) and this result into Eq. (CI) to obtain a closed form solution for 01, but because the expression would be algebraically cumbersome this will not be done. 141 C.l.2 Application of Orthogonality at f = 0 The equation and boundary condition to be applied for determination of the unknown constants (Appendix A Eq. (A.40) and Eq. (A.6c)) are " a a 02(x+.y+) = ZBn[cosh(-7"y+)-Cnsinh (4y+]]* n=0 L L [ancos (an?) +Biosin (anx+)] (C.14) “Rifle2 (x‘, 0) = BiH[T,,(x*) - To] (c.1531) y‘ = 0 The difficulty in applying this boundary condition is recognized by noting the right hand side contains the temperature of the hot fluid. The needed temper- ature expressions for evaluating Eq. (C. 15a) are: 92 (x+, 0) = B,I [ancos (anf) + Biosin (anx+)] (C.16) 11:0 392 " an — (x”, O) = - 2 337C. [ancos (anf) +Biosin (anx+)] (C.17) 3)” n = o L TH 0:“) = To + Dlexp (-N,,x+) ~ . (aw-3.1L”) .. + 24,,{113L [1- (N; #13132) ]cosh(p.an) Putting these expressions into Eq. (C.15a) 142 " a Z 3,,L—fgn [aflcos (mg?) + Biosin (up?) 1 n=0 + Bi” 2 B" [ancos (an?) + Biosin (anx+)] = n=0 " , (BiN - 2L”) , + 2A,,{113L [1— ”2” 5",: cosh(ttuL x+) "=0 (NH-1131' ) :lsinh (11,113+) } a. B. N 2 + 2 (Bn+a_:cn){an[l-( 1° "+a")]c08(a,,x1) 2 2 (NH+°‘..) Bi N +012 +NH[ 1:2: 2"]sin(anx+) }] (C.19) H a). for which orthogonality will be applied. Each term in Eq. (C. 19) will be mul- tiplied by the eigenfunction X (an,x*) = ances (an?) +Biosin (czar) (C20) evaluated at term m and integrated over the boundary. The eigenfunction has the same property as discussed earlier, which is given in Eq. (C.7). Each of the integral will be listed separately and evaluated starting with the first term in Eq. (C.19) and proceeding to the last term. 1 a. a J { 2 Bn-écnIancosmnf) +Biosin(aux+)] * 0 n=0 L [amcos(amx*) +Biasin(amx*)] }ch+ = Bm%CmN(am) (C21) 143 1 to j (31,, 2 3,, [ancos (an?) +Biosin (anxtn * 0 n=0 [amcos (an?) +Biosin(amx*)] }dx+ = BiHBmNmm) (C22) 1 31,4), 1' {exp (-1~/,,,x+ ) [amcos (amxt) +Biosin (amx+)] }dx* = 0 BiHDl {(N +B' )a + -N”[ (N +B° )a cos(a )+ l e "' l (”31+a'2n) H a m H 0 m m (Uzi-Nam.) sin (am) ] } (0.23) l on - 2 .2 t (BIONH-unl‘ ) t {Bi 3,1131. 1- , coshutnL x+)* 1 ”E. 1 (Nib-nil- ”) 1 [amcos (an?) +Biasin(amx*)] }dx* = . . t2 2 An31flam 2 , _ Name - 1131. (113132 + cf") " N2, - 11:1.” ] {sinh (unL’) unL' cos (am) n=0 + cosh (unL‘) amsin (am) } + i 11,133,131, “2 , 1_NHBio-u:L (nfit‘zwg) " Ng-ufiL‘z n=0 2 ] {sinh (11,13) u,L‘ sin (a...) -cosh (unL‘) amcos (am) + am } (C34) 144 1 c- . 2 02 (BloNH-Ll'L ) B ° A N ]sinh(uuL'x+)* 0 13:0 [aacos (an?) +Biosin (anxtn }dx+ = i AnBiHam NHBio - 33L” (11,2,L’2 4» 01,2") l4, H N}, - ufiL" n30 ] {cosh (unL’) unL' cos (am) + sinh (uuL')amsin (am) ‘11,}; } " AnBiHBio NHBio- p 2‘1. ) n H + 2 N},-: "=0 (uiL‘2+a: 22,,]{008h(u L) p. L sin(am ) -sinh (uuL’) amcos (am) } (C25) The remaining integrals depend on whether n=m. (BiaNH-I- 01:) + I“; iunzo (B +_ L' 1C")a"[l - (N},+ 01%) Jcos (“"x )* [amcos (an?) +Biosin (amx+)] }dx” = For n=m " 3' N 2 Biflz (B"+E;C")anl:l- ( '02H+:n):l* "=0 L (NH+an) (1m 1 2 [Tm+43in(amff)+ Bi sin (a )] (C26) Fornstm Bi” i (Bu-(- %Cn)an[l _ (BiaNu+ (1%)] {am|:8in(a’"-a”) + sin (“n+“n)] n=0 (Nift- (1:) 2 (am-0%,.) 2 (am+an) + Bio[l - COS (am— an) 1- cos (am+ an) 2(am-an) + 2(am+an) ]} (€27) 145 (31°,N,,+01:) . + j{31,,2(3,, +—c )NH[ (Nffiafi) ]sm(anx )2 n=0 [amcos (an?) +Biosin (amx+)] }dx+ = For n=m . u an (BioNH+atzr) Bio Bio . 1 . Bl” 2 (Bu '1' PC")NH[ (Nil... (13) ]|:T - 33:83] (2am) '1' 581112 (am)] ((3.28) n=0 Fornaem . " on,l (BioNH-t-ai) . sin(am-an) sin(am+an) B'”.§.(B"*PC")"”[ mama) 1W1 21.1,-..) " 2(a..+a,.)] 1- cos (“111 - an) 1- cos (“111 + at") +“[ mam-01") + 2(am+an) ] (C29) Instead of substituting these large expression into the Eq. (C.19) some new variables are defined. For Eq. (C23) let “{(N +31H)01 +eN”[- (N,,+31°)01m cos(a ) (OMB N2: + (ag-NHBio) sin (01,") ] } ((2.30) In Eq. (C.24-25) define: 3' N 3' - 2L’2 '21“... ) 2L'[l- H 10 ll, ]{sinh (unL') unL'cos (am) 5 p "'"1‘ (ngL‘2+o1; N§,-11,2.L“2 +cosh(unL')amsin(am) } (c.31) BiHBio 2L! 1 NHBio_u:L.2 {. ( L.) L. . ( ) p a . - , smh 11 [l sm 0: "'2 (nthag) 33-1192 " " "' -cosh (unL') amcos (“111) + am } (C32) 146 31,,01," [NHBio -11§L‘2 p "‘2" (1121.32.32) 11,, " Ni-uiL" 1km (""L ) ”"1“ com") + sinh (unL' ) a," sin (01m)-11nL‘ } (C33) 3533'} ”3310 " "3L” { sh ( L.) L" - (a ) a co sm pmn,4 (“£12.2‘1'az.) n H [Vii-11:14.2 “I! "a m -sinh (nut) amcos (am) } (C34) and for Eq. (C.26-29), for m=n _ . (BioNH+a:) an: 1 . Bio . 2 wmn-Btuau 1" (N§+a:) [7+zm(2am)+2am81n (am)] (Bi.N11+a§> 31, 31,, . 1 . 2 +NH (N;l+a:) [T ”Ia—MSU! (2am) +5811! (am):| (C35) and for man: -3° 1 (31',N,,+01§) sin(am-au) sin (amwn) “'"m’ ”“1 ’ (was) "'1 2111",-..) “ 2(a..+a.>] . 1-cos(am-an) 1-cos(am+an) +B‘°[ 2(a111'0‘11) + 2(“111+°‘11) ] (BioNH + 01%) . sin (“11 - on") sin (am + “11) +NH[ (N21411:) '0': 2(am-an) - 2(am+an) :| +am[1-cos (“111-“11) l-cos (“111+a11) nan-01,.) + 2(am-t-an) ] (c.2115) Using the newly defined variables given in Eq. (C.30-36) in the integrals of Eq. (C.21-29); Eq. (C.19) after applying orthogonality, can be written as 147 “111 . ” BMFCmNfllm) + BmB'HNuxm) = (0le + Z A” (plan. 1 + pm”. 2 + pm": 3 + pm", 4) n=0 " a + 2 (Egg-[gum (C37) n=0 where the orthogonality relationship is N01,, ) 1 ((114-32) +——-Bic +B° (C38) = _ 1” l . [ ( .2...) .1 Defining one last variable a A,” = —[ZT-CM+BiH]N (“111) ((3.39) and rearranging the equation gives 3317‘". +1.20"? +C,, gm)wm+o) D1: “zAu(P11111,1+P11111,2+P11111,3+Pm11,4) n=0 (C.40) where the sign changes were introduced (in A," and Eq. (C.40)) so that the variables are exactly as programmed in the computer program used to evalu- ate the problem (in Appendix F). 148 C.l.3 Application of Orthogonality at y” = 1 The equation and boundary condition to be applied for determination of the unknown constants (Appendix A Eq. (A.48) and Eq. (A.8d)) are 93(X+.)’+) = ZCn[%cosh(%gy*)+3iflsinh(%y+)]* n=0 [ancos (anf) + Biosin (anx+)] (C41) 803 a— +BiC93 (x+, 1) = BiC[TC(x+) - To] ((2.423) y . y = l The evaluation of this boundary condition is very similar to the previously discussed method for 62, except the algebra is a little more involved because the boundary condition is evaluated at y“ = 1 instead of O. The needed tem- perature expressions for evaluating Eq. (C.42a) are u a" all . an 9303.1) = 2 Cn[—,cosh (7) +BiHsmh (7)] * n g 0 L L L [ancos (anf) + Bio sin (anx*)] (C43) 89 " a a a a —3-(x+.l) = ZCn-§[—fsinh(—:)+Biucosh(—§)]* ay+ ":0 L L L L [ancos (anf) + Biosin (anx+)] (C.44) 149 Tc (x‘) = To + Elexp (NC?) (Ni-nib”) n=0 .0 B. N 2Lt2 + 2A.’{uuL‘[l+( 1" Cw" )]cosh(unL'x+) (BiN + 2L”) . +NC ”2C 5"” :|sinh(unL m} (NC-“"1! ) " B'N — 2 n = 0 BioNC — a: + ”c 2 sin (uni) (c.45) NC + a: Putting these expressions into Eq. (C.42a) gives “ an an a" an 2 Cr; [—, sinh (7) + BiHcosh (—, II [ancos (anx+) + Biosin (anx*)] n: o L L L L +Bic i C,[%cosh (;§)+Biusinh (2.3)]:- n=0 [ancos (Qua?) +Biosin (a,x")] = Bic [El exp (ch) " , (BioNc+ ng‘z) . + A ’ L l + . cosh L x+ (BioNC+u:L'2) . +N , sinh( L 1:”) C[ (NE-uiL 2) "" “ (BiONC-ag) B ’ C ’ l + n§0( n+ n){an[ + (Né+a:) +N (Rifle-ab ' ( +) (C46) (1 . C[ (Ni-mi) is” "x ] :l cos (053+) 150 for which orthogonality will be applied. Each term in Eq. (C.46) will be mul- tiplied by the eigenfunction X(an,x*) = aucos (ans?) +Biosin(anx*) (C.47) evaluated at term m and integrated over the boundary. The eigenfunction has the same property as discussed earlier, which is given in Eq. (C.7). Each of the integral will be listed separately and evaluated starting with the first term in Eq. (C46) and proceeding to the last term. {Ca"a"inha” 3 ha I B. ( )1... I 2 nL.[L.S (U) tycos (21)] ann”c08(ax) lsin ax [amcos (an?) +Biosin (amx+)] }dx+ = am am am an Cm—. [—;— sinh(— T)+ BIHCOSh(-— ?)]N (am ) (C.48) L L L L an an 3,131.0 20 C "E 0031‘ (L—f )+B‘33inh (f: [I [ancos (an?) +Biosin (anx+)] * [amcos (amf) +Biasin (amx+)] }dx" = am am am 31'ch [— cosh ( ) + BiHsinh(— )flN(aM ) (C.49) L L L 1 area, 1 {exp (ch+) [amcos(amx+) +Biosin(0tmx+)] }dx+ = 0 BicEi N2 :{wi -Nc)am +eN‘[(NC-Bi )am cos(am )+ (0:,+Nc3i,)sin(am) ] } (C.50) 151 l 0- ° 2 ‘2 . (B'oNc'l'llnL ) {Bi An’u L 1+ , i C 2 " [ (Ni-nil- 2) n=0 ]cosh (unL'f) * [amcos (amf) + Biosin (amx*)] }dx* = " Au’BiCamttnL' (BioNcwfiL'z) 2 -2 1+ 2 .2 (nib mi.) (NC-nil. ) n=0 ] {sinh (nnL‘) an" cos (am) + cosh (unL') amsin (am) } + i An’BicBiounL’ 1+ (BioNC+u:L'2) with 0:3,.) (N3 wil-'2) n=0 ] {sinh ‘ 2(am+a,)] + (1,,[1 - cos (am-0t") 1— cos (am+an) 2(am-au) + 2(am+an) ] (C63) Using the newly defined variables in the expressions for the integrals, Eq. (C.46) after applying orthogonality is C — [—,sinh (1)4- BiHcosh—;-]N(am) L L L ML. am am am + BiCCm[—;- cosh (T) + BiHsinh—;:|N (am) = L L L QME1+ Z An’ (PM 1 + PM 2 + PM 3 + PM .) + 2 (B; + 0,5) ‘1’” (C64) n=0 n=0 Defining a final variable 052.. . . . a... . . a... (but 2 [(14.2 +BICBIHJSlnh (F)... (BIH+BIC) COSh (FHNUIH') (C55) and substituting into Eq. (C64) and rearranging gives 155 Cmom - 2 (B; + Cn’) ‘l'm — QmEl = 221,319”, 1 + PM 2 + PM 3 + P ) (C66) I! l mu, 4 156 C.2 Summary of Applying Nonhomogen. Boundary Conditions The final boundary conditions were applied and simplified using the orthogonality of the eigenfunctions. This procedure produced a close form solution for the constant A, in 61 since the nonhomogeneous term in the boundary condition was a constant , 3 Bi” Bi” BIL(TL- To) [sm (“111) - -u—cos (“111) + T] Am = 2 t2 . o . o a ((112) [(umL +Bi0BiL)s1nh(umL )+ (B10+BiL)umL cosh(ttmL )]N(um) where _ _1_ 2 - 2 Bi!— - However, for the constants B3, in 63 and Cu in 03 the boundary condi- tions to be applied have nonhomogeneous terms that are functions containing the unknown constants B, and C3,. These functions are contained in the solu- tions for the fluid temperatures. Thus, in applying the boundary conditions and orthogonality the resulting equation is simplified, yet still contains sum- mations and is a function of both unknown constants B" and C3,. The equa- tions to be solved to determine these constants are B", 2 “+203? +C,—2:1,”)wm+m 1),: " 2 A”(91.111.12“ Pm, 2 + PM 3 + pm, 4) (C40) n=0 Cmtbm - Z (B3'+ Cn’)‘1’,..,. 13:0 (25,: 2A,’(PM3+PM +Pm m3+P3M3) (C.41) These two equations represent a set of simultaneous equations that need to be solved for the unknown constants B" and C,I (n= 1,2 ....... N) . The terms on the left hand side contain unknown constants and the right hand side contains all known information for two equations. The solution of these equations is given and discussed in Appendix D. APPENDIX D D.l Solution for Constants B, and c, After evaluating the final nonhomogeneous boundary conditions and applying orthogonality two equations were derived to solve for B, and C,, (Appendix C Eq. (C40) and Eq. (C.66), which are B 1,, + 2 0(8, +C, :Z-MJVMR'i-O) DI: -2A,(pm1+pm,,2+pm,,3+p,,,4XD.l) u=0 C,,o, - 2 (B,’+ c,') ‘l’m- 9,51 = 2 A,’(P,,,,1+P,,,,2+Pm,’3 +P,,4) (D2) n=0 n=0 where all the greek variables are defined in Appendix C and the expressions for E1 and D1 are (from Appendix B Eq. (B.25) and Eq. (B.52)) o 1 (BiaNH-flila.2) " N2_ 2L'2 H "a D1: Tam-To- 2"»sz n=0 .- B. N 2 -2 (B,+%C,)a,[l- ( t°2H+a")] (D.3) n=0 2 NH+a, " , , (Bi NC+ pZL’Z) , E1 = exp (-NC) [(Tmn" T,)- 2 A, {unL [1+ 02 2 ",2 cosh (unL ) n=0 Nc‘llnL (BiN +p2L’2) . , +Nc[ ”2C 2",, ]smh(u,L)} NC-u'nL “ B'N — 2 -2 (B,'+C,'){a,[1+ 1" c a"]cos(a,)+ 2 2 ,go Neal-oz, N (BiaNC—afi) . ( ) 4 [ (Pew: ism" i“ Defining the following variables: 157 158 (Bi N” - fl") 2 t - 0 n B, a (u,L ) [I N}, __ ”it, (D.5) (BioNH + (23,) Y..'°‘..[1' ”3”“; (13.6) ll _ . (BiN wil-‘2) . Date Ncfl‘lnl‘ [1+ £2.qu .2 ]c08h(}.lnL) (BioNC+ ufiL‘z) +NC 2 2 ‘2 NC- ”EL 4 a, a, (BiaNC - afi) ‘ an?" C[cosh (—,)-Cnsinh (7)] {ual} + 2 2 cos (an) L L NC+ a, _ (BiONC - a3) ‘ + NC 2 2 NC + a, Jsinh (u,L‘) } (v.7) sin ((1,) } (D8) _N a, a, a, (BiaNc- uni) ruse C[_.cosh (—;)+Biusinh (7)] a,[l + ]cos (an) L L L N3: + a: (BioNC- a3) . + NC 2 2 sm (a,) (D9) NC + a, and substituting into Eq. (D.3-4) gives a Q a D1 = (TH, ,., - To) - 2 A,p, — 2 (B, + c,—f)y,) (D10) 1: = o n = o L E, = exp (-NC) (TC, ,,— :20 A, ',u - 20 (B, a, + C,1:(D.11) where Eq. (D.8-9) contain the necessary constants to convert B,’ and C,’ to B, and C,. Defining two last variables, to convert all unknown constants to a form without a superscript prime (VB) M, = ‘Pm [cosh (g) — Cnsinh (E; )] (D.12) 159 a a a ($10,, = W,,[—3 cosh (4)+Bi,sinh (4)] (0.13) L L L Introducing Eq. (D.10-13) into Eq. (D.1-2) and rearranging gives " " a 2 W.,-m.) + 2 C,,_ W.,-m.) = n=0 n=0 -mm(THJfl—T0) - 2Al(pmn,l+pmn,2+pmn,3+pmn,4-mmfln) (1114) n=0 2 [annou- (YB) mu] 3,, + Cm¢m + 2 [flm1n_ (WC) mu] C, ___ n=0 Qm‘xp(‘ c) (TQM-To) + zAn’u’nm.1+Pmm2+Pmm3+Pmm4-Qmun) (13°15) n = 0 These two equation represent a set of simultaneous equations to be solved for the constants B,and C,. To accomplish this the series must be truncated at a value n=N giving 2N constants to solve for (B, and Ci 2' = 1, 2 ..... N).There are only two equations, however these equations can be written for each eigen— value; hence the subscript m denotes the eigenvalue considered. Therefore, since there are m=N eigenvalues and two equations for each eigenvalue the total number of equations are 2N to solve for the same amount of constants. For illustrative purposes consider a matrix representation of the equa- dons [ME] = [f] (D.16) where [K] is the coefficient matrix of order (2Nx2N) and [f] is the forcing vector (2le) and the unknown constant vectors given as It‘xli‘et :11:- 'FDFJ- 160 81 Cl B2 C2 LB] = LC] = - (D.17) _BN. -CNJ Further, the coefficients matrix can be quartered and the forcing vector halved K K [K] = " 1 " 2 (D.18) L. .J F) [El = (D.19) [2] Each quarter of the coefficient matrix is a square matrix of order (MW) and half the forcing vector is of length (NxI). Referring to Eq. (D. 16) the simulta- neous equations can be represented as [§:§:]E’;~]=[§;] This shows where the terms that need to be put into the coefficient matrix and forcing vector will come from. The rows in the coefficient mauix in Eq. (D20) represent the coefficients of B, and C, for each equation. Row one corresponds to Eq. (D.14) and row two to Eq. (D.15). Similarly, the entries in the forcing vector correspond to Eq. (D. 14) for row one and Eq. (D.15) for row two. Presenting the equations in this form, Eq. (D20), makes it possible to write out some of the terms of the matrices to show the developing pattern and indicate where these terms originate. The individual matrices and forcing vectors that form Eq. (D.18) and Eq. (D.19) are listed below. 161 (A, My" -11co,) (Wu-72ml) .. (Wm-won (V21 "' 7102) (A2 + V22 " 7202) - - (WZN — YNmz) [K1] = (0.21) (VNI - 710)”) (VNZ — 720)”) . . (2%” + VNN -' yNtoN) " T (anl a2 a ((V11'71m1)27) (le’yzml)f) '°((W1N-7le)zg) ((v vana‘) ((v vanaz) (( )a") 21‘12'7 22‘227 -- W'Y‘D'T [K2] = L L 2N N 2 L (( - (0)11) (( - (0)93) a” WM 71 N L‘ Wm 72 N L' -- ((WNN'YNlezr) (D22) (alol-(WB)11) (Gig—(1%)”) . . (91°N‘(‘FB)1N) (Q,a,-(\PP)21) (Q,a,-(\PB)22) . . (azoN-(WB) 2N) [.Kg] = (Owl-(WM) (flNaz-(WB),2) .. (Q,a,-(\PB)NN) (D23) 162 P ( CC and the heat capacity ratio is used C,, = CC/ C,,, Eq (E35) and Eq. (E36) become (7.11.17: ' Tc, in) (CR '1) C1 = ER-exp [-UP/CC(CR-1)L] (E37) _ CRTC’ 1.. - 73,1119” [-UP (CR - 1) L] ‘ - CR-exp[-UPICC(CR-1)L] Substituting Eq (E35) and (E36) into the solutions for the hot and cold flu- ids, Eq. (E31) and (E32), gives (E38) 169 CCTC,iu “Tu,tnCH¢xP[’UP(Z~l—"C—)L] 1,,(x) = C + C -C x -UP(— —)L] C- He pl: C1” -CC (TH, in - TC, in) (CC - CH) (1 - C,,/CC) {Cc— C,,exp[-UP(EI; - —1— L]} exp[- UP(-— -2}- )x](E (13.39) C ) CH C C CCTC ,n— TH inCHexp[- UP(— - —)L] Tc(x) = CH CC — CC- Cuexp[-UP(— -—)L] CH Cc ”min" Tarn) (CC' CH) exp ”PFC-l- --)x (13.40) (l-CCJCH){CC-Cuexp[-UP(£;—EIE)L]} [ Cu Cc ] For the case of CH> CC Eq. (E37) and (E38) can be used in place of Eq. (E35) and (E36) giving the following expression for the hot and cold fluid temperatures: + - - .— THO‘ ) ' CR-exp[-UP/CC(CR-1)L] (TH,in-TC,in) CR CR- exp [—UPICC(CR- 1)L] ”‘1’ [‘UP (CR ’1)“ “‘14” T (x+) = CRTC,“in-TH,in:xp[-UPI_CC(CR-1)L]_ C (TH, in - TC, in) CR - exp I-UIVCE (CR ‘1)L] exp l-UP (CR -1)” 05.42) Nondirnensionalizing the length scales by the total length of the heat exchanger +-§ x -L (13.43) and introducing the number of transfer units (NTU) 170 UPL C min N TU = (E44) the final equations for the hot and cold temperatures become - CRTC, in "' TH, mexp ["NTU(CR - 1)] H- CR-exp[—NTU(CR-1)L] _ (TH, in "' TC, in) CR CR-exp [-NTU(CR- 1)] “P [‘NTU (CR - 1)x*] (545) T = CRTC, in " Tn, inexp [-NTU (CR - 1)] C CR-expl-NTU(CR— 1)] (Tu, m' TQM) CR - exp [-NTU(CR -1)] exp [-NTU(CR - l)x+] (E46) Note that Eq. (E45) and (E46) apply for C,, > CC, which could easily be con- verted if the opposite were true. E.2.2 Equal Heat Capacities (C,, = 0,) If the heat capacity ratio of the two fluids are equal the solution of the differential equation, in Eq. (E10), for the temperature difference between the fluids is a constants. Which upon substituting into Eq. (E8) and (E9) give the following differen- tial equations and initial conditions for the temperatures of the hot and cold fluids, respectively. dTH -UP 7; — WC, (5.48) T” (0) = T", in (E.483) dTC -UP TC (L) = To». (E49a) The general solution of these differential equations is 171 -UP TH(x) = C —Clx-t-C2 (E50) 1! Tc(x) = Z-Cl—J-I—iClx-t-C3 (13.51) c where C2 and C3 are constants, giving a total of three constants to be deter- mined. The initial conditions in Eq. (E48a) and (E49a) provide two condi- tions and the third condition comes from Eq. (E25). Applying these three conditions gives the following set of equations to be solved for the unknown constants. T”, M = C2 (5.52) TC, in = -g—:+:C_ L+C3 (5.53) The solution for C2 is obvious, solving fore the other constants gives (Ta, 1.. - Tc 1..) C1 " (1 +_ UPL (355) Cc UPL TC,’ in + TH, in CC C3 = (1+___UPL (E56) CC which can be rewritten using Eq. (E44) as (TH, in - TC, in) CI " 1+ NTU (EN) Tc, in + T”. “N TU C3 ' 1 + NTU (E'SB) Substituting Eq. (E57) and (E58) into the general solutions, Eq. (E50) and (E51), the expressions for temperature of the hot and cold fluids are _ UP (TH, in - TC, in) CH(1+NTU) 711(1) = x + T”, in ([559) 1‘72 UP(THiu-TCIn) (TCin+TH,inNTU) Tc(x) = — ’ ' x+ ’ CC(1+NTU) 1+NTU If length dimensions are scaled according to Eq. (B.43), expressions for tem- perature become (E60) NTU (Tu, :. - Tarn) (1+ NTU) NTU (TH, in ' Tc, in) x” + (Tc, in "' TH, inNTU) (1+NTU) 1+NTU TH (x) = - x” + T”, in (E.61) TC (1:) = (5.62) APPENDIX F Program AXCOND BOUNDRY CONDITIONS OF THE THIRD KIND ORTHORGONALITY APPLIED WITH COSINE AND SINE FUNCTION ASSUMING To = Tl (CONSTANT AMB TEMP) LAST MODIFIED 4127/92 ************ implicit double precision (a-h,o-z) character para*10,runfil*60 integer size,size2 double precision len,Nc,Nh,Kw,Lstr,mu,NTU parameter (size=200,size2=50) dimension coeff(size,size+l),BC(size),A(size), + Tcold(size2),Thot(size2),Twall(size2,size2), + mu(size),alph(size), + chold(size2),dThot(size),dexxO(size2),dexx1(size2), + deyy0(size2),deyy1(size2) common delta,len,Lstr,Nc,Nh,Bic,Bih,Bio,Bil,Tcin,Thin,To,'l‘l it ** INPUT DATA *1- ” delta - wall thickness (m) ** len - total lenth (m) """ Cc - heat capacity cold stream (WIK) ** Ch - heat capacity hot stream (WIK) *" hc - convection coeff cold side (W/m**2 K) ** hh - convection coefl“ hot side (Wlm**2 K) *"' Kw - thermal conductivity of wall (W/m K) ** ho - convection coefl‘ to ambient at x=0 (W/m* *2 K) ** hl - convection coefl‘ to ambient at x=L (Wlm**2 K) ** To - Temp of ambient at x=0 (C) ** Tl - Temp of ambient at x=L (C) *"' Tcin - temp of cold fluid at x=L ** Thin - temp of hot fluid at x=0 ** Lx(L) - scale length for x-direction 173 I” 13* it it it it it *4! *II it II"! it *1! I” ** iii IMF ill ** It'll It'll C *I-I'i-‘I-‘I'i‘I-‘I-‘I-‘I-‘I‘I'i-‘I'I-M *1- 174 Ly(delta)- scale length for y-direction CALCULATED DATA Lstr -Lley Nc -thch/Cc Nh -thth/Ch Bio - ho Lx I Kw all variables dimensionless Bil -h1Lx/Kw Bic -thy/Kw Bih -thy/Kw read the input data in and define looping method iopt - 0 -input raw data (dimensional) 1 -input Non-dimensional data iloop = 0 irun = 1 continue READ IN NEXT RUN FILE NAME AND BEGIN CALCULATIONS if(irun.eq.l)runfil='[dowding.thesis.rundata]Bih.one' if(irun.eq.2)runfll='[dowding.thesis.rundata]Bih.low' if(irun.eq.3)runfil='[dowding.thesis.rundata]Bio.one' if(irun.eq.4)runfll='[dowding.thesis.rundata]Bil.one' if(irun.eq.5)runfil='[dowding.thesis.rundata]BihBic.one' if(irun.eq.6)runfll='[dowding.thesis.rundata]BioBil.one' if(irun.eq.7)runfll=’[dowding.thesis.rundata]Lstr.one' if(irun.eq.l)runfil='[dowding.thesis.rundata]crl.low' if(irun.eq.2)runfll='[dowding.thesis.rundata]cr75.low' if(irun.eq.3)runfil=‘[dowding.thesis.rundata]cr50.low' if(irun.eq.4)runfil='[dowding.thesis.rundata]cr25.low’ if(irun.eq.l2)runfll=’[dowding.thesis.rundata]nterm.one' if(irun.eq.5)go to 300 print *,'running data set ',irun open(20,status=‘unknown',file=runfil) read(20,'(i10)')iopt print *,'enter ioption O-single run l-vary parameter 2-stop' I'ead("'.'(il)')iOPt 175 if(iopt.eq.2)go to 300 100 continue *‘I'i'l‘l' {ii-1" **'I-* if (iopteq.0)then call input0(Cc,Ch,Kw,wid,nterrn,ibound,ifluid,iwall,ider, deltax,deltay,istop) elseif(iopt.eq. l)then iloop=iloop + 1 call inputl (Cc,Ch,Kw,wid,ntenn,ibound,ifluid,iwall,ider, deltax,deltay,istop,iloop,para,ipara,del,parmax) endif Nc = Bic*Kw*I.su*widICc Nh = Bih*Kw*Lstr*wid/Ch print *,' read in input' if(istop.eq.l)go to 200 SOLVE FOR THE EIGENVALUES print *,'going to calculate eigenvalues‘ SOLVE FOR ALPHA AND MU call root(alph,size,Bio,Bil,nterm) call root(mu,size,Bic,Bih,nterm) if(alph(l).gt.3.l41592654)then write(10,"')alph(l),' eigen troubles' go to 300 endif if(mu(l).gt.3.l41592654)then write(10,*)mu(l),' eigen troubles' go to 200 endif print *,' calculated eigenvalues' BUILD MATRDK OF UNKNOWN COEFFCIENTS call build(mu,alph,coefl',nterm,size,A) print *,' built matrix' SOLVE THE COEFF MATRIX FOR THE NEEDED CON STANT S call gauss(coefl‘,size,nterrn,BC) print *,' solved matrix ' CALCULATE THE FLUID AND WALL TEMPERATURES call profll(BC.A,size,size2,nterrn,deltax,deltay,Thot, 176 + Tcold.Twall,chold,dThot,mu.alph,d'I‘dxx0,dexx1, + deyyO,deyyl,qx0,qxl) print ‘3' solved for temperatures ' CALCULATE THE EFFECI'IVNESS AND NUMBER OF TRANSFER UNITS WITH AND WITHOUT AXIAL CONDUCTION INCLUDED *§**§ n = nint(lldeltax) + 1 m = nint(lldeltay) + 1 if(Cc.le.Ch)then Cmin = Cc Cr = Cc/Ch else Cmin = Ch Cr = Cthc endif NTU = (Bic*Bih"Kw*Lsu'*wid)/(Cmin*(Bih+Bic*Bih+Bic)) effc = Cc*(Tcold(1) - Tcold(n))/(Cmin*(’l'hot(l)-Tcold(n))) effh = Ch*(Thot(l) - 'Ihot(n))l(Cmin*('Ihot(1)-Tcold(n))) * write(*,*)'NTU ',N'IU if(Cr.eq.1)then eff = NTU/(1.0 + NTU) else eff = ( 1-exp(-NTU*(l-Cr)))/(1-Cr*exp(-NTU*(l-Cr))) endif print *,' calculated Effect and NTU ' CALCULATE THE HEAT FLUXED AND INSURE THAT ENERGY IS CONSERVED ***** qhot = Ch "' (Thot(l) - Thot(n)) qcold = Cc * (Tcold(l) - Tcold(n)) qu = -wid"‘Kw/Lstr‘qx0 qxl = -wid*Kw/Lstr'"qxl *write(10,*)'heat fluxes' *write(10,*)'qx0 =',qx0 *write(10,*)'qxl =',qxl *write(10,*)'qhot =',qhot I"write(10,"')'qcold =',qcold **§*** 4: 177 qtot = (qhot - qcold + qu - qxl)/(rnin(qhot,qcold)) print ‘, ' calculated heat fluxes ' GENERATE THE OUTPUT FILE OF THE TEMPERATURE DISTRIBUTIONS OF THE FLUIDS AND WALL AND CHECK THAT THE BOUNDRY CONDITIONS ARE MET if(iopt.eq.0)then call output(deltax,deltay,size2,dexxO,dexx1,deyyO, + deyyl,Thot,Tcold,Twall,iwall,ifluid,ibound, + idenponplxpoymly) call out0(NTU,efl',efl'c,eflh,nterm,Kw,Cc,Ch,qtot,p0x,p1x, + P0y.ply) else call output(deltax,deltay,size2,dexxO,dexx1,deyy0, + deyyl,ThoLTcold,Twall,iwall,ifluid,ibound, + ider,p0x,p1x,p0y,ply) call out1(NTU ,efl'pfl‘c,eflh,nterm,Kw,Cc,Ch,para,qtot,p0x,p1x, + p0y,ply,iloop,ipara) endif print *, ' Generated output files ' if(iopt.eq.0)close(10) goto 100 200 continue if (iopLeq. l)then close(10) close(40) irun=irun+l ismp=0 iloop=0 go to 50 endif 300 continue a: t a. Stop end DATE 3105/92 SUBROUTINE TO GENERATE THE OUTPUT FILE FOR OPTION 0 178 ‘I' - SINGLE DATA FILE OPTION VARYING SOME PARAMETER * subroutine out0(NTU,efl,effc,eflh,nterm,Kw,Cc,Ch,qtot, + Pokplwomly) implicit double precision(a-h,o-z) double Precision NTU,Nc,Nh,len,Lstr,Kw common deltaJen,I..str,Nc,Nh,Bic,Bih,Bio,Bil,Tcin,'I‘hin,To,Tl 1"} WRITE THE NONDIM DATA AND INLET PARAMETERS USED IN CALCULATION write(10,"')" write(10,*)' Nondimensional Data' write(lO,'(a,f12.5)')' L“ =',Lstr write(lO,'(a,f12.5)')' Hot side(Bih) =',Bih write(10,'(a,f12.5)')' Cold side(Bic) =',Bic write(10,'(a,f 12.5)')' Wall End(Bio) =',Bio write(10,'(a.f12.5)')' Wall End(Bil) =',Bil write(10,‘(a,f12.5)')' 'ntu" hot (Nh) =',Nh write(10,'(a,f12.5)')' 'nm" cold(Nc) =',Nc write(10.*)' ' write(10,*)' Fluid and Wall conditions ' write(10,'(a,fl 2.5,a)')' Wall Conductivity(Kw) =',Kw,’ W/mAZK' write(10,'(a,f12.5,a)')' Heat Cap cold (Cc) =',Ch,‘ W/K' write(10,'(a,fl 2.5,a)')' Heat Cap hot (Ch) =',Cc,‘ W/K' write(10,'(a,f12.5,a)')' Inlet Temp cold (Tcin) =',Tcin,‘ C' write(10,'(a,f12.5,a)')' Inlet Temp hot (Thin) =',Thin,’ C' write(10,*)' ' write(10,*)'Ambient Temperatures ' write(10,'(a.f 12.5)')' At x=L (Tl) =',Tl write(10,'(a,f12.5)')' At x=0 (To) =',To write(10,"‘)" write(10,'(a,iZ,a)')' Summations terminated at ',nterm,' terms' write(10,"')" write(lO,'(l7x,a,5x,a,4x,a,4x,a)')'I-lot',’Cold', + 'Negl', 'Eng Boundry Condition err Bi' write(l0,'(5x,a,8x,a,6x,a,5x,a,5x,a)')'NTU','Efl‘,'Efl‘,'Efl‘, + 'Bal Bio Bil Bih Bic' write(10,'(2x,t7.4,2x,3f9.4,lx,e8.2,4f7.4)')NTU,effh, + efic,efl’,qtot,p0x,plx,p0y,ply return end * DATE 3/05/92 i-I-‘I'fii * + 4.. + 179 SUBROUTINE TO GENERATE THE OUTPUT FILE FOR OPTION 1 - OPTION VARYING SOME PARAMETER subroutine outl(NTU,efl‘,efl'c,efl‘h,nterm,Kw,Cc,Ch,para,qtot, pOx.p1x.pOy.ply.iloop.ipara) implicit double precision(a-h,o-z) double Precision NTU,N'IUlam.Nc,Nh,len,Lstr,Kw character para" 10 common delta,len,Lstr,Nc,Nh,Bic,Bih,Bio,Bil,Tcin,Thin,To,Tl WRITE THE NONDIM DATA AND INLET PARAMETERS USED IN CALCULATION if(iloop.eq.1)then write(10,*)' ' write(lO,"')' Nondimensional Data' write(10,'(a,f12.5)')' L“ =',Lstr write(10,'(a,f12.5)')' Hot side(Bih) =',Bih write(10,'(a,f12.5)')' Cold side(Bic) =',Bic write(10,'(a,f 12.5)')' Wall End(Bio) =',Bio write(10,'(a,f12.5)')' Wall End(Bil) =',Bil write( lO,'(a.f 12.5)')' "ntu" hot (Nh) =',Nh write(10,'(a,f12.5)')' ”ntu" cold(Nc) =',Nc write(10.*)' ' write(10,*)' Fluid and Wall conditions ' write(10,'(a,f12.5,a)')' Wall Conductivity(Kw) =‘,Kw,‘ W/m"2K' write(10,'(a,f12.5,a)')' Heat Cap cold (Cc) =',Cc,’ W/K' write(10,'(a,fl 2.5,a)')' Heat Cap hot (Ch) =',Ch,‘ W/K' write(10,'(a,fl 2.5,a)')' Inlet Temp cold (Tcin) =',Tcin,’ C' write(10,'(a,f12.5,a)')' Inlet Temp hot (Thin) =',Thin,‘ C‘ write(lO,‘)' ' write(lO,‘)’Ambient Temperatures ' write(10,'(a,f12.5)')' At x=L (Tl) =',Tl write(10,'(a,f 12.5)')' At x=O (To) =',To write(lO,*)' ' write(10,'(a,i3,a)')' Summations terminated at ',nterm,' terms' write(10,*)' ' write(10.'(a.a)')' varying parameter 29m write(10,*)' ' write(lO,'(3x,a9,l6x,a,5x,a,4x,a,4x,a)')para,'Hot','Cold', 'Negl', 'Eng Boundry Condition' write(lO,'(l9x,a,6x,a,5x,a,5x,a,5x,a)')'NTU', 'Efl‘, 'Efl‘,'Efi‘,'Bal max error' endif 180 BCerr = max(p0x,plx,p0y,p1y) uneff = (efl‘ - (efl'h+eft‘c)l2)leff Cmin = min(Ch,Cc) Cmax = max(Ch,Cc) NTUlam a Kw*Lstr/Cmin if(iparaeq.l)then write(10,'(3x,f7.2,5x,f8.3,1x,3f8.5,lx,e8.2,3x.f7.5)') + Isu',NTU,eflh,efl'c,eff,qtot,BCerr write(40,'(lx,f10.3,lx,f4.2,2f10.6,2x,4f9.5)')Lstr,Cmin/Cmax, + Cmin,NTU,eflh,efl'c,eff,unefl‘ return elseif(ipara.eq.2)then write(10,'(3x,f7.4,5x,f8.3,lx,3f8.5,1x,eB.2,3x,f7.5)')Bih,NTU, + eflh,effc,efl',qtot,BCerr write(40,'(1x,1f 10.6,lx,f4.2,2f10.6,2x,4f9.5)')Bih,Cmin/Cmax, + Cmin,NTU,eflh,effc,eff,uneff return elseif(ipara.eq.3)then write(10,'(3x,f7.4,5x.f8.3,lx,3f8.5,1x38.2,3x.f7.5)')Bic,NTU, + efl'h,efl'c,eff,qtot,BCerr write(40,'(lx,f 10.6,1x,f4.2,2f10.6,2x,4f9.5)')Bic,Cmin/Cmax, + Cmin,NTU,effh,efl'c,efl‘,uneff return elseif(ipara.eq.4)then write(10,'(3x,f7.4,5x,f8.3,lx,3f8.5,lx,e8.2,3x,f7.5)')Bio,NTU, + effh,effc,efl‘,qtot,BCerr write(40,'(lx,f10.6,lx.f4.2,2f10.6,2x,4f9.5)')Bio,Cmin/Cmax, + Cmin,NTU,eflh,efl'c,eff,uneff return elseif(ipara.eq.5)then write(10,'(3x,f7.4,5x,f8.3,1x,3f8.5,1x,e8.2,3x,t7.5)')Bil,NTU, + eflh,effc,efl',qtot,BCcrr write(40,‘(1x,f10.6,1x,f4.2,2fl0.6,2x,4f9.5)')Bil,Cmm/Cmax, + Cmin,NTU,eflh,efl’c,eff,uneff return elseif(ipara.eq.6.and.Ch.lt.Cc)then write(10,'(3x,f7.4,5x,f8.3,lx,3f8.4,lx,e8.2,3x.t7.5)')Cthc,NTU, + efl‘h,efl'c,eff,qtot,BCerr write(40,'(1x,f4.2,1x,2fl0.6,2x,4f9.5)')Cmin/Cmax,Cmin,NTU,eflh, + efl'c,eff,unefl‘ return elseif(ipara.eq.6.and.Ccle.Ch)then write(10,'(3x,f7.4,5x,f&3,1x,3f8.5,lx,e8.2,3x.f7.5)')Cc/Ch,NTU, **** + + + + + + + + + + + 181 eflh,efl'c,eff,qtot,BCerr write(40,'(1x,f4.2,1x,2flO.6,2x,4f9.5)')Cmin/Cmax,Cmin,NTU,effh, efl'c,eff,unefl‘ return elseif(ipara.eq.16.and.Cc.le.Ch)then write(10,'(f7.4.f8.2,f8.3,lx,3f8.5,lx,e8.2,3x,f7.5)')Cc/Ch,LsU,NTU, eflh,efl‘c,eff,qtot,BCerr return elseif(ipara.eq.l6.and.Ch.lt.Cc)then write(10,'(t7.4,fB.2,f8.3,lx,3f8.5,lx,e8.2,3x,fl.5)')Ch/Cc,lsu,NTU. efl‘h,effc,efl‘,qtot,BCerr return elseif(ipara.eq.23)then write(10,'(3x,t7.5,5x,f&3,lx,3f8.5,1x,e8.2,3x,f7.5)')Bih,NTU, eflh,effc,efl',qtot,BCerr write(40,'(lx,2f8.5,lx,f4.2,lx,f8.4,lx,f$.2,4f8.5)')Bic,Bih, Cmin/Cmax,Cmin,NTU,eflh,effc,efl‘,unefl‘ return elseif(ipara.eq.45)then write(10,'(lxfl.5,lx,f7.5,f8.3,1x,3f8.5,1x,e8.2,3x,t7.5)')Bio, Bil,NTU,efl'h,efl'c,efl',qtot,BCerr write(40,'(1x,2f8.5,lx,f4.2,1x,f8.4,lx,f5.2,4f8.5)')Bio,Bil, Cmin/Cmax,Cmin,NTU,eflh,efl‘c,efl',unefl’ return elseif(ipara.eq.8)then write(lO,’(3x,iS,7x,f8.3,lx,3f8.5,1x,e8.2,3x,f7.4)')nterm,NTU,eflh, efl‘c,eff,qtot,BCerr write(40,'(lx,3i,f10.6,2x,4f9.5)')N,NTU,efl'h, efl'c,eff,uneff return endif end DATE 2/17/92 SUBROUTINE TO OPEN FILE AND READ IN RAW DATA subroutine input0(Cc,Ch,Kw,wid,nterm ,ibound,ifluid,iwall,ider, deltax,deltay,istop) implicit double precision (a-h,o-z) double precision len,Lx,Ly,Kw,Nc,Nh,Lsu' 182 character outfll‘60.infil*60 common deltalenlsu',Nc.Nh,Bic,Bih,Bio,Bil,Tcin,Thin,To,Tl print ‘,Enter the input file' read(*.'(a)')infil print ‘,Enter the output file' read(*.'(a)')outfil CHECK IF EOF FLAG TO TERMINATE THE PROGRAM HAS BEEN READ IN if(outfll.eq.'eof')then close(20) istop = 1 return else istop=0 endif READ IN THE RAW DATA open(10,status='new',file=outfil) 0pen(40,status='new',flle=outfil) open(50,status='new',filwoutfil) open(30,status='unlrnown',file=infll) read(30,'(l7x,i20.0)')l..str read(30,'(l7x,f20.0)')Bih read(30,'(l7x,f20.0)')Bic read(30,'(l7x,f20.0)')Bio read(30,'(l7x,i20.0)')Bil read(30,'(l7x,f20.0)')Kw read(30,’(l7x,i20.0)')Ch read(30,'(l7x,f20.0)')Cc read(30,'(l7x,f20.0)')'lhin read(30,'(l7x,f20.0)')Tcin read(30,'(l7x,f20.0)')To read(30,'(l7x,f20.0)')'1'l read(30,'(l7xf20.0)')wid read(30,'(l7x,i20)')nterm read(30,'(l7x,i20)')ibound read(30,'(l7x,iZO)')ifluid mad(30,'(l7x520)')iwall read(30,'(l7x,iZO)')ider read(30,‘(l7x,i20.0)')deltax read(30,'(l7x,i20.0)')deltay close(30) *****§ 1. at a: 183 return end DATA 3105/92 SUBROUTINE TO OPEN FILE AND READ IN RAW DATA subroutine inputl(Cc,Ch,Kw,wid,nterm,ibound,ifluid,iwall,ider, + deltax,deltay,istop,iloop,para,ipara,del,parmax) implicit double precision (a-h,o-z) double precision len,Lstr,Nc,Nh,Kw character outfil*60,para*10,infil*60 common deltaJeanu;Nc,Nh,Bic,Bih,Bio,Bil,Tcin,Thin,To,Tl if(iloop.eq. l)then READINPARAMETERTOVARYJNCREMENTANDMAXVALUETHENTHERAW DATA FILE NAME AND PERTINENT DATA print *,'enter the para to vary' l'ead("‘.'(i10)')il>fim print *,'enter the del' read(*,'(f20.0)')del print*,'enter the max value' read(*,'(f20.0)')parmax print ‘,‘enter the input file name' read(*.'(a)')infil print *,'enter the output file name' read(*.'(a)')outfil open(40,status='new',file=outfll) Open(10,status='new',file=outfil) open(30,status='unhiown',file=infil) read(30,'(l7x,f20.0)')Lstr read(30,'(l7x,f20.0)')Bih read(30,'(l7x,i20.0)')Bic read(30,'(l7x,i20.0)')Bio read(30,'(l7x,f20.0)')Bil read(30,'(l7x,f20.0)')Kw read(30,'(l7xfl0.0)')Ch read(30,'(l7x,i20.0)')Cc *i't‘l' 184 read(30,'(l7x,i20.0)')Thin read(30,’(l7x,i20.0)')Tcin read(30,'(l7x,i20.0)')To read(30,'(l7x,i20.0)')Tl read(30,'(l7x,i20.0)')wid read(30,'(l7x,iZO)’)ntenn read(30,'(l7x,i20)')ibound read(30,'(l7x,i20)')ifiuid read(30,'(l7x,im)')iwall read(30,'(l7x,i20)')ider read(30,'(l7x,f20.0)')deltax read(30,'(l7x,f20.0)')deltay close(20) close(30) OTHERWISE INCREMENT THE PARAMETER THAT IS BEING VARIED AND SEND DATA BACK TO PROGRAM endif if(iparaeq.l)then para = 'Lstr' if (iloop.eq. l)then Cmin = min(Cc,Ch) Cr = Cchh ratio = Lstr/Cmin else Lsu' = Lstr + del write(*,"')'Lstr ',Lstr Cmin = Lstr/ratio if(Cc.le.Ch)then Cc = Cmin Ch = Cc/Cr else Ch = Cmin Cc = Cr " Ch endif endif istop = 0 if(I.str:gt.parmax)istop=l return elseif(ipara.eq.2)then para = 'Bih' istop = 0 if(iloop.eq. l)then Cmin = min(Cc,Ch) Cr = Cchh 185 ratioi = Lstr/Cmin“(Bih"Bic*Kw*wid/ (Bih + Bic‘Bih + Bic)) else Bih = Bih * del ratio = Lstr *(Bih*Bic*Kw*widl (Bih + Bic*Bih + Bic)) Cmin = ratio/ratioi if(Cc.le.Ch)then Cc = Cmin Ch = Cc/Cr else Ch = Cmin Cc = Cr " Ch endif endif if(Bih.gt.parmax)istop=l return elseif(ipara.eq.3)then para = 'Bic' istop = 0 if(iloop.eq. l)then Cmin = min(Cc,Ch) Cr = Cc/Ch ratioi = Lstr/Cmin*(Bih*Bic*Kw*widl (Bih + Bic*Bih + Bic)) else Bic = Bic * del ratio = Lstr *(Bih*Bic*Kw*widl (Bih + Bic‘Bih + Bic)) Cmin = ratio/ratioi if(Cc.le.Ch)then Cc = Cmin Ch = Cc/Cr else Ch = Cmin Cc = Cr * Ch endif endif if(Bic.gt.parmax)istop=l return elseif(ipara.eq.4)then para = 'Bio' kmp=0 if(iloop.gt.1)Bio = Bio * del if(Bio.gt.parmax)istop=l return elseif(ipara.eq.5)then ********************* 186 para = 'Bil' istop = 0 if(iloop.gt.l)Bil = Bil * del if(Bil.gt.parmax)istop=l return elseif(ipara.eq.6)then para = ‘Cr' istop = 0 Cr = Cthc if(iloop.gt.l)Ch = Ch 4» del if(iloop.gt.l)Cc = Cthr if(Cc.gt.parmax)istop=l return elseif(ipara.eq. l6)then para = ‘Cr & L*' istop = 0 Cmax = max(Cc,Ch) Cmin = min(Cc,Ch) if(iloop.eq.1)then ratio = Lsu'lCmin deLsu' = del * (parmax-Lstr) delCr = del "' ratio/(Cmax - Cmin) write(*,*)'delCr delLstr =',delCr,deLstr else if(Ch.le.Cc)then Lsu' = Lsu' + deLstr Ch = Ch + delCr else Lstr = Lsu' + deLstr Cc = Cc + delCr endif endif if(Lsu'.gt.parmax)istop = 1 return elseif(ipara.eq.23)then para = 'Bih & Bic' istop = 0 if(iloop.eq. l)then Cmin = min(Cc,Ch) Cr = Cc/Ch ratioi = Lstr/Cmin*(Bih*Bic*Kw*wid/ (Bih + Bic‘Bih + Bic)) else Bih = Bih + del Bic = Bic + del ratio = Lstr *(Bih*Bic*Kw*wid/ 187 + (Bih + Bic‘Bih + Bic)) Cmin = ratio/ratioi if(Cc.le.Ch)then Cc = Cmin Ch = Cc/Cr else Ch = Cmin Cc = Cr "' Ch endif endif if(Bih.gt.parmax)istop = 1 return elseif(ipara.eq.45)then para = 'Bio & Bil' istop = 0 if (iloop.eq. l)then Bioi=Bi0 Bili=Bil else Bio = Bioi "‘ del Bil = Bili "‘ del Bioi = Bio Bili = Bil endif if(Bio.gt.parmax)istop = 1 return elseif(ipara.eq.8)then para = 'N' kmp=0 if(iloop.gt.l)nterm = nterm "' del if(nterrn.gt.parmax)istop=l return endif end date 2109/92 SUBROUTINE TO BUILD THE MATRIX OF UNKNOWN COEFF *i-ii-I'l’l- subroutine build(mu,alph,coefl',nterm,size,A) implicit double precision (a-h,o-z) integer size {ii-I- 188 double precision 1en,Nc,Nh,mu,mulsq,rnul,Lstr,Ncsq,Nhsq,lam,musq common deltalen,Lsu',Nc,Nh,Bic,Bih,Bio.Bil,Tcin,Thin,To,T1 dimension coefl‘(size,size+l),mu(size),alph(size) dimension A(size) GENERATION OF THE MATRIX, ORDER WILL BE (2*NTERM X 2*NTERM +1) SET CONSTANTS AND INITIALIZE THE PARAMETERS pi = 1141592654 Ncsq = Nc * Nc Nhsq = Nh "' Nh rhsl = 0.0 rhs2 = 0.0 rhsB = 0.0 rhs4 = 0.0 exch = dexp(-Nc) BEGIN COMPUTTIN MATRIX -i CORRESPONDS TO ROW NUMBER do 10 i=1,nterm mul = mu(i)*Lstr mulsq = mul "' mul musq = mu(i) * mu(i) alpl = alph(i) I Lstr alphsq = 81111110)“ alph(i) alplsq = alpl * alpl hsin = sinh(alpl) hcos = cosh(alpl) sine = sin(alph(i)) cosin = cos(alph(i)) zeta = (alpl*hsin + Bic‘hcos)l(Bic*hsin + alpl*hcos) orth =((alphsq+Bio"'*2.0)*(1.+Bill(alphsq+Bil**2.))+Bio)l2.0 orthmu =((musq+Bih"2.0)*(l.0+Bicl(musq+Bic**2.0))+Bih)/2.0 lam = -(alpl*zeta + Bih) "' orth if(Nh.lt.100)then 0mg = Bih/(Nhsq+alphsq)*((Nh+Bio)*alph(i)+dexp(-Nh)* ((alphsq-Nh‘BioY’sine - (Bio+Nh)*alph(i)*cosin)) else omg = Bil1l(Nhsq+alphsq)*(Nh+Bio)*alph(i) endif 189 cphi =((alplsq+Bic*Bi11)*hsin+(Bic+Bih)"'alpl*hcos)*orth 1} ABSORBED EXP(-NC) INTO COMG FROM LOADING OF THE MATRIX if(Nc.lt.100)then comg = Bic/(Ncsq+alphsq)‘(exch*(Bio-Nc)*alph(i)+ + ((alphsq+Bio*Nc)*sine + (Nc-Bio)*alph(i)*cosin)) else comg = Bic/(Ncsq+alphsq)*((alphsq+Bio*Nc)*sine + + (Nc-Bio)*alph(i)"'cosin) endif I- j CORRESPONDS TO COLUMN NUMBER do 20 j=1,nterm alphi = alph(i) alpl = alph(i)/Lstr alphsq = alph(i) " alph(i) musq = mu(i) "' mUG) mul = mu(j)*Lstr mulsq = mul*mul * hsinmu = sinh(mul) "' hcosmu = cosh(mul) hsinal = sinh(alpl) hcosal = cosh(alpl) sinal = sin(alph(i)) cosal = cos(alph(i)) sinmu = sin(mu(i)) cosmu = cos(mu(i)) chl = (Nh*Bio-mulsq)/(Nhsq-mulsq) ch2 = (Nh*Bio+alphsq)/(Nhsq+alphsq) ccl = (Nc‘Bio+mulsq)/(Ncsq-mulsq) * cc2 = (Nc*Bio-alphsq)l(Ncsq+alphsq) * CONSTANTS 01.62.03.04 AND cps HAVE A COSH FACTORED OUT OF THEN WHICH "' CANCELS WITH THE COSH FACTORED OUT OF THE CONST AG) WHEN MULT :- TOGET'HER WITH rho AND crho IN RHS CONSTANTS AS COMPARED TO SOLUTION c1 = tanh(mul)*mul*cos(alphi) + + alphi*sin(alphi) c2 = tanh(mul)*mul*sin(alphi) + - alphi‘cos(alphi) + alphi ++ 190 c3 = mul*cos(alphi) + tanh(mul)*alphi*sin(alphi) - mul c4 = mul‘sin(alphi) - tanh(mul)*alphi*cos(alphi) denoml = mulsq + alphi*alphi zeta = (alpl‘hsinal + Bic*hcosal)l (Bic‘hsinal + alpl*hcosal) rhol = Bih*alphildenom1*musq*Lstr*(l-chl)*cl rh02 = Bih‘Bio/denoml*musq*Lsu'*(1-chl)*c2 rh03 = Bih‘Nh/denoml*mu(i)*alphi*chl*c3 rho4 = Bih*Nh*Bio/denom l *mu(i)*chl *c4 crhol =Bic‘alphildenoml*mul*(l+ccl)*cl crh02 =Bic*Bio/denom1*mul*(l+ccl)*c2 crho3 =Bic*Nc/denoml*alphi*ccl*c3 crho4 = Bic*Nc*Bio/denom1*ccl*c4 beta = (l-ch1)*musq*Lstr 8am = (l-ch2)‘alph(i) eps = (l+ccl)*mul+ Nc*ccl*tanh(mul) sig = (hcosal-zeta*hsinal)*((1+cc2)*alph(j)*cosal + Nc*cc2*sinal) tau =(alpl*hcosal+Bih*hsinal)*((1+cc2)*alph(j)*cosal + Nc*cc2*sinal) ammn = alphi - alph(i) ampn = alphi + alph(i) anmm = alph(i) - alphi anpm = alph(i) + alphi if(i.eq.j)then psi = (.5+1/(4*alphi)*sin(2*alphi))*(l-ch2)*alphsq + ((Bio*(l-ch2) + Nh*ch2)l2) *(sin(alphi))**2 + (.5-ll(4"'alphi)*sin(2*alphi))*Nh*Bio*ch2 cpsi = (.54-ll(4*alphi)*sin(2*alphi))*(1+cc2)*alphsq + ((Bio*(1+cc2) + Nc*cc2y2) *(sin(alphi))*"'2 + (.5-1/(4‘alphi)*sin(2*alphi))*Nc‘Bio‘ch *‘I-i-I'i- 191 else sinlm = sin(ammn)/(2*ammn) sinlp = sin(ampn)l(2"'ampn) cos2m = cos(ammn)l(2*ammn) cos2p = cos(ampn)/(2*ampn) cos3m = cos(anmm)l(2*anmm) cos3p = cos(anpm)l(2*anpm) argl = sinlm + sinlp arg2 = -cos2m - cosZp + l/(2*ammn) + 1/(2*ampn) arg3 = -cos3m - cos3p + 1/(2*anmm) + ll(2*anpm) arg4 = sinlm - sinlp psi =arg1‘alph(j)*alphi*(1-ch2) + arg2‘alph(i)*(l-ch2)*Bio + arg3*ch2*Nh*alphi + arg4’ch2*Nh*Bio cpsi =argl’alph(i)*alphi*(l+cc2) + arg2*a1ph(i)"'(l+cc2)"‘Bio + arg3‘cc2‘Nc‘alphi + arg4*cc2*Nc*Bio endif psiB = psi * Bih psiC = psi * Bih * alpl cpsiB = cpsi "' Bic * (hcosal - zeta "‘ hsinal) cpsic = cpsi * Bic * (alpl‘hcosal + Bih*hsinal) LOAD THE MATRIX, LOADING IS DONE IN QUARTERS(i.e.TOP LEFT, top right, bottom left. ..... ) but concurrently lST EQUATION LOADING - TOP LEFT AND TOP RIGHT coefl‘(i,j) = -omg*ga1n + psiB coeff(iJ+nterm) = -omg*gam"'alpl + psiC 2ND EQUATION LOADING - BOTTOM LEFT AND BOTTOM RIGHT coeff(i+nter1n,j) = comg‘sig - cpsiB coefl'(i+nterrn,j+ntenn) = comg*tau - cpsiC if(i.eq.j)then coeff(i,j) = coefl'(i,i) + lam coeff(i+nterrn,j+ntenn) = coefl‘(i+nterm,j+ntenn) + cphi 192 endif CONSTANT AG) HAS A COSH FACTURED OUT OF THE DENOMINATOR THAT IS CANCELLED BY THE COSH FROM C1..C4 WIEN MULTIPLIED BY rho AND crho IN TIE RHS CONSTANTS CALCULATIONS *‘I'O‘I'fi if((Tl-To).ne.0)then A0) = (Bil"(Tl-To)*(sin(mu(i))—Bih/mu(i)*cos(mu(i)) + +Bihlmu(j))) l(orthmu"' + ((mulsq+Bil*Bio)*tanh(mul)+(Bil+Bio)*mul)) Aprm = A(j)*(mu(i)*cosmu + Bih’sinmu) else A0) = 0.0 Aprm = 0.0 endif rhsl = (rhol+rho2+rho3+rho4)*A(i) + rhsl rth = (crhol-+crho2+crh03+crho4)*Aprm + rhs2 rhs3 beta*A(j) + rhs3 rhs4 eps‘Aprm + rhs4 LOADING THE RHS CONST ANTS - DONE FOR EQUATION 1 THEN EQ 2 **'I-* if(i.eq.nterm)then coefl‘(i,2"'nterm+1) = -rhsl + omg’rhs3 - omg*(Thin-To) coeff(i+ntenn,2*ntenn+l)=rhs2 - comg*rhs4 + + comg'CTcin-To) endif 20 continue rhsl = 0.0 rhs2 = 0.0 rhs3 = 0.0 rhs4 = 0.0 10 continue return end DATE 51 12192 *I'i'i' *I'I'i Mia-{area- 193 SUBROUTINE TO CALCULATE EIGENVALUES FROM THE TRANSCENDENTAL EQUATIONS. USING A SCHEME THAT MARCHES OUT UNTIL TIE FUNCTION CHANGES SIGN THEN BACKSTEPS TO TIE ROOT subroutine root(eigen,size,B1,B2,ntenn) implicit double precision(a-h,o-z) integer size double precision Nc,Nh,len,Lsu' common deltalenJAu',Nc,Nh,Bic,Bih,Bio,Bil,Tcin,Thin,To,Tl dimension eigen(size) func(x) = tan(x) - (x*(B l+B2))/(x**2-B1*B2) irep = 0 delin = .5d+0 ndel = 2 pi = 4.0d+0*atan(1.0d+0) eps = 1d-8 ep = 1d-10 del = delin/ndel K = 0 arg =sqrt(Bl*B2) idenom = 0 CHECK IF FUNCTION CHANGES SIGN DUE TO CHANGE OF SIGN IN DENOM- INATOR OF THE FUNCTION. (X