|UFDC Home||myUFDC Home | Help|
This item has the following downloads:
PRES SURE DROP AND HEAT TRANSFER INT INVERTED FILM BOILINTG HYDROGEN
A DISSERTATION PRESENTED TO THE GRADUATE SCHOOL
OF THE UNIVERSITY OF FLORIDA IN PARTIAL FULFILLMENT
OF THE REQUIREMENTS FOR THE DEGREE OF
DOCTOR OF PHILOSOPHY
UNIVERSITY OF FLORIDA
The effort put forth over the last four and a half years to complete this Ph.D. is dedicated to my
children, Nicholas and Connor. This is one component in my continuing efforts to be a good
father and role model for them. Life is much more interesting and rewarding when you remain
I would like to thank Dr. Samim Anghaie for agreeing to work with me on this effort that
started four and half years ago. I understand that working with a long-distance student is
difficult all the more reason I appreciate his patient support to achieve this goal. I thank my
wife, Melynda, who supported my efforts by giving me time to study. I express my gratitude for
having a great and supportive family; John and Alice Pasch, brother Jack, and sisters Alison and
Lorelei. I also express my gratitude to Robert Hendricks for giving freely of his memories of
these experiments in which he was centrally involved. His efforts, then and now, provide the
engineering community with unique information.
TABLE OF CONTENTS
ACKNOWLEDGMENT S ........._..... ...............4.._._. ......
LIST OF TABLES ........._..... ...............7..____ ......
LI ST OF FIGURE S .............. ...............8.....
AB S TRAC T ............._. .......... ..............._ 16...
1 INTRODUCTION AND STATE OF THE ART ................ ...............18........... ..
Introducti on ................. ...............18.................
M otivation .............. ...............18....
Obj ectives ................. ...............19.......... .....
Pressure Drop............... ...............20..
Heat Transfer .............. ...............23....
2 MODELLING APPROACHES FOR TWO-PHASE FLOW .............. .....................3
Angular Simplifications ................. ...............38.................
Basic M odels .............. ...............38....
Flow Regime Analysis............... ...............40
3 TEST DATA DESCRIPTION AND EVALUATION AND MODEL DEVELOPMENT....44
Description of Experiments ............ ......_.. ...............44....
Experimental Setup .............. ...............44....
Experimental Conditions ....._.. ............... ......._.. ..........4
Heat Leaks ..... ._ ................ ......._.. ..........4
In strmentati on ................. ...............46........... ....
Data Validation ................. ............ ...............49.......
Comparison with Similar Data ................. ...............49................
End Effects ......... .. .. .................... ...............5
Hydrogen States: Parahydrogen and Orthohydrogen ........._.._.. ....._.._ ........._.....55
Model Development .............. ...............58....
Nature of Data .........._.... ........____ ...............60.....
Magnitude of Radiation Heating ........._._._. ....___ ...............63...
Conservation Equations ........._._.. ..... .___ ...............64.....
Entrance Lengths ........._._.. ..... .___ ...............66.....
Boundary Conditions............... ...............6
Closure Conditions .............. ...............70....
Vapor super heat............... ...............7
Liquid energy fl ow ........._... ...... ..... ...............7 1...
W all friction .............. ...............74....
Model Implementation .............. ...............75....
4 ANALYSIS AND VALIDATION OF MOMENTUM MODEL RESULTS ........................90
Data Referencing .............. ...............90....
Data Refinement ................. ...............90.................
Om itted Data .............. ...............90....
Problematic Data .............. ...............90....
Data Representation............... .............9
Problematic Runs............... ...............93..
Vapor Super heat ................. ...............95.......... .....
M odel Results ................. ........... ...............96......
Validation of Model Results ................. ...............97................
Range of Validity ................. ...............98........... ....
5 EVALUATION AND CORRELATION OF DATA AND CORRELATION
AS SES SSMENT ................. ...............106................
Data Correlation................. .............10
Low Pressure Slip Correlation ................. ........... ...............111 ....
Low Pressure Slip Correlation Assessment ................. ...............113........... ...
High Pressure Slip Correlation ................. .......... ...............115.....
High Pressure Slip Correlation Assessment ................. ......... ......... ............1
Accuracy of the Slip Correlations............... .............11
Validation of the Slip Correlations ................. ...............118........... ...
Observations ................. ...............119......... ......
6 HEAT TRANSFER ANALYSIS .............. ...............168....
Data Omission ................. ...............168.
The Nature of IFB Heat Transfer ................. ...............168..............
The General HTC Profile .................. .. ........ ......... .........16
An Interpretation of Controlling Effects in IFB Heat Transfer ................. .................1 70
Assessment of Various Models ................ .......................... ....................172
7 CONCLUSIONS AND RECOMMENDATIONS ................ ...............................180
General Conclusions ................. .. .......... ............. ............18
Pressure Drop Conclusions and Recommendations ................ .......... ...............181
Heat Transfer Conclusions and Recommendations............... ............18
Recommendations for Future Efforts .............. ...............182....
LIST OF REFERENCES ................. ...............184................
BIOGRAPHICAL SKETCH ................. ...............190......... ......
LIST OF TABLES
3-1. Table of experimental conditions. ............. ...............77.....
3-2. Comparison of Core et al. and Hendricks et al. heat transfer coefficients ............................78
3-3. Comparison of heat transfer coefficients for Wright and Walters data and TN 765.............79
3-4. Summary of test conditions for maj or hydrogen heat transfer studies ................ ...............79
3-5. Result of parametric sensitivity study of end axial heat conduction. ................ ............... ..79
3-6. Tube wall axial heat transfer analysis............... ...............79
4-1. List of tube numbers, dimensions, and runs executed with the tubes. ................ ...............99
4-2. Statistical analysis of pressure data ................ ...............100..............
5-1. Accuracies of some common slip correlations. .............. ...............121....
5-2. Comparison of pressure drop prediction accuracies ................. .............................122
6-1. Comparison of predictive accuracy of various IFB models. ............. ......................7
LIST OF FIGURES
2-1. Various flow regimes for IFB. ............ ...... ._ ...............43
2-2. Flow regime map generated by Takenaka for IFB (1989). ............. .....................4
3-1. TN 765 experimental setup............... ...............80.
3-2. TN 3095 experimental setup............... ...............81.
3-3. 1961 data test section. .............. ...............82....
3-4. TN 3095 test section ................. ...............83........... ..
3-5. TN 3095 instrumentation. .............. ...............84....
3-6. Nodal distribution and heat generation distribution used to model end effects ....................85
3-7. Radial metal temperature profiles as a function of metal thermal conductivity...................85
3-8. Radial metal temperature profiles as a function of metal thickness. .............. ................86
3-9. Effect of specified parameters on tube end wall axial heat transfer. .................. ...............86
3-10. Difference in wall to liquid temperature for all data considered .................... ...............8
3-11. Wall to liquid hydrogen temperature differences for four runs ................. ..........___.....87
3-12. Theoretical liquid core temperature profile at the exit of the heated test section. ...............88
3-13. Flow diagram for momentum and energy analysis of data. ............. ......................8
4-1. Sample of 1961 data wall temperatures............... .............10
4-2. Tube 3 exhibits a consistent reduction in wall temperature at 34 cm................. ...............101
4-3. Comparison of runs 7 and 8 pressure profiles .....__.....___ ..........__ ........10
4-4. Run 14 energy and momentum balances .....__.....___ ..........._ ...........0
4-5. Results of modifying the coefficient in Burmeister' s equation............. .. .........___....103
4-6. Culled data momentum and energy balance results from model............. ... .........___...103
4-7. Calculated void fraction from model for the culled data set. ............. .....................10
4-8. Velocity slip ratio vs quality from model for the culled data set. ................... ...............10
4-9. Void fraction vs. equilibrium quality for three runs of Ottosen' s experiments ................... 105
5-1. Vapor velocity vs. superficial velocity. ............. ...............122....
5-2. Comparison of model slip and slip predicted from correlations. ............. .....................12
5-3. Predicted versus measured pressure gradients for all data used in correlating slip.............123
5-4. Model and prediction results for run 1. ............. .....................124
5-5. Model and prediction results for run 2. ............. ...............125....
5-6. Model and prediction results for run 3. ................ ....__ ....__ ........__.........126
5-7. Model and prediction results for run 4. ............. ...............127....
5-8. Model and prediction results for run 5. ............. ...............128....
5-9. Model and prediction results for run 6. ............. ...............129....
5-10. Model and prediction results for run 7. ............. ...............130....
5-11. Model and prediction results for run 9. ............. .....................131
5-12. Model and prediction results for run 10. ............. ...............132....
5-13. Model and prediction results for run 11. ............. .....................133
5-14. Model and prediction results for run 12. ............. ...............134....
5-15. Model and prediction results for run 13. ............. ...............135....
5-16. Model and prediction results for run 15. ............. ...............136....
5-17. Model and prediction results for run 16. ............. ...............137....
5-18. Model and prediction results for run 17. ............. ...............138....
5-19. Model and prediction results for run 18. ............. ...............139....
5-20. Model and prediction results for run 19. ............. ...............140....
5-21. Model and prediction results for run 20. ............. ...............141....
5-22. Model and prediction results for run 21. ............. ...............142....
5-23. Model and prediction results for run 33 .......__......_.__.... ...._.. ......._.._.......143
5-24. Model and prediction results for run 34. ............. ...............144....
5-25. Model and prediction results for run 35. ............. ...............145....
5-26. Model and prediction results for run 37. ............. ...............146....
5-27. Model and prediction results for run 38. ............. ...............147....
5-28. Model and prediction results for run 39. ............. ...............148....
5-29. Model and prediction results for run 40. ............. ...............149....
5-30. Model and prediction results for run 41. ............. ...............150....
5-31. Model and prediction results for run 42. ............. .....................151
5-32. Model and prediction results for run 43 ........_......_._._... .....___......_.........152
5-33. Model and prediction results for run 45. ............. ...............153....
5-34. Model and prediction results for run 46. ............. ...............154....
5-35. Model and prediction results for run 47. ............. ...............155....
5-36. Model and prediction results for run 48. ............. ...............156....
5-37. Model and prediction results for run 49. ............. ...............157....
5-38. Model and prediction results for run 50. ............. ...............158....
5-39. Model and prediction results for run 51. ............. .....................159
5-40. Model and prediction results for run 23 ........_......_._._... .....___......_.........160
5-41. Model and prediction results for run 24. ............. ...............161....
5-42. Model and prediction results for run 25. ............. ...............162....
5-43. Model and prediction results for run 27. ........._. ...... .___ .....___......_.........163
5-44. Model and prediction results for run 32. ............. ...............164....
5-45. Model and prediction results for run 36. ............. ...............165....
5-46. Model and prediction results for run 44. ............. ...............166....
5-47. Model and prediction results for run 8. ............. ...............167....
6-1. Variation of the HTC as a function of quality in IFB flow. ............_. .. ...___............176
6-2. Variation of HTC versus equilibrium quality in the IAFB flow regime. ............................177
6-3. Variation of HTC versus mass quality for runs 39-42. ................ ............................177
6-4. Variation of HTC versus mass quality for runs 44-47............... ...............178.
6-5. Variation of Dittus-Boelter vapor properties with pressure and temperature. ....................178
6-6. Comparison of predicted HTC using the TN 3095 correlation with the experimental .......179
6-7. Comparison of predicted HTC using the modified equilibrium bulk Dittus-Boelter. .........179
AIAFB agitated inverted annular film boiling
As surface area
b y-intercept of line
Bo boiling number
C conversion constants for ortho-para conversion
CHF critical heat flux
Co Colburn number
Co drift flux model distribution parameter
cp specific heat at constant pressure
cy specific heat at constant density
DFB dispersed film boiling
f friction factor
F Chen's enhancement factor
fi low pressure slip correlating parameter
f2 high pressure slip correlating parameter
Fr Froude number
G mass flux
Go reference mass flux
Gr Grasshof number
h mass-specific enthalpy
h, or HTC heat transfer coefficient
IAFB inverted annular film boiling
IFB inverted film boiling
ISFB inverted slug film boiling
j superficial velocity
k thermal conductivity
K conversion factor for ortho-para conversion
LOCA loss of coolant accident
m slope of line
n number density
Nu Nusselt number
Pr Prandtl number
q" heat flux
qo" reference heat flux
Q heat flow rate
r radial direction, radial distance
Re Reynolds number
s velocity slip
S Chen's suppression factor
mass flow rate
oc void fraction
P volumetric quality
AT temperature differential
X Lockheed-Martinelli parameter
4 friction multiplier
a surface tension, Stefan-Boltzmann constant
z shear stress
u specific volume
c cross section
crit critical condition
f film conditions
i inlet, interface
I liquid phase
10 all fluid flowing as liquid
m mean conditions
mac macroscopic, in Chen's correlation
mic microscopic, in Chen's correlation
s saturated conditions
tt turbulent-turbulent liquid-vapor phases
v vapor phase
Abstract of Dissertation Presented to the Graduate School
of the University of Florida in Partial Fulfillment of the
Requirements for the Degree of Doctor of Philosophy
PRES SURE DROP AND HEAT TRANSFER INT INVERTED FILM BOILINTG HYDROGEN
Chair: Samim Anghaie
Major: Nuclear Engineering Sciences
Two-phase boiling hydrogen pressure drop and heat transfer is studied in the context of
high velocity upflow in a constant, high heat flux, steady state, internal pipe flow environment.
These data were generated by NASA in the early and mid 1960s in support of the manned space
flight programs. Measurements taken were local pressure, temperature, and voltage drop.
System measurements included mass flow rate, and test section inlet and discharge pressure and
This effort establishes the nature of the flow as inverted film boiling, which has been
studied to some degree. In this structure, the wall temperatures are too hot to allow liquid to
remain at the surface. Therefore, a vapor film is established at the wall throughout the flow. The
approach of this analysis is to reverse-engineer the data to determine mass quality, void fraction,
and velocity slip. This is accomplished by applying a one-dimensional, fiye-equation model,
with pressure gradient being the one combined equation for the liquid and vapor phases. Other
maj or assumptions are that all of the vapor is at the mean film temperature, and the liquid core
experiences no sensible heating.
The resulting velocity slips are correlated for high and low pressure conditions, with the
cutoff established at 600 kPa. Good agreement is achieved between the pressures predicted
using the slip correlations and the measured pressures. Results are in general significantly better
than those from the homogeneous equilibrium model.
Various established heat transfer coefficient models are also applied to these data. It is
shown that pre-critical heat flux models fail absolutely to predict the heat transfer coefficient. It
is further shown that film boiling models that focus on buoyancy fail as well. While all forced
convection film boiling models are within a reasonable range of the data, recommendations for
appropriate models are made.
The range of pipe inlet conditions are 188 kPa to 1265 kPa, mass fluxes from 327 kg/m2-S
to 3444 kg/m2-S, and heat fluxes from 294 kW/m2 to 2093 kW/m2. Two heated test section
lengths are 30.5 cm. and 61.0 cm. long, and five different diameters range from 0.48 cm. to 1.29
INTRODUCTION AND STATE OF THE ART
This dissertation investigates the state of understanding of and prediction capabilities for
boiling hydrogen, and the needs for improving the current condition. It presents an engineering-
based approach to improve on the prediction capabilities for pressure drop and heat transfer.
Accurate predictions of pressure drop in and heat transfer from a pipe to hydrogen during
forced convective two-phase flow benefit engineers throughout the life of a product. During the
design phase, good pressure drop and heat transfer models will help the engineer reduce the
uncertainty in the design parameters. During the product test and development phase, good
models will help the engineer to correctly interpret test data, therefore allowing him to determine
where modifications are necessary. During the use of the product, problems inevitably arise that
require the engineer to search for the root cause. This investigation requires a good
understanding of how the product will react under off-nominal operating conditions. Accurate,
mechanistic models allow the engineer to perform this investigation with confidence that the
thermal-hydraulics related results of the investigation are valid.
The rocket industry uses liquid hydrogen as a fuel. Heat transfer to two-phase flowing
hydrogen routinely occurs during three phases of rocket operation; fuel tanking, rocket engine
conditioning, and possibly during rocket firing. Nuclear Thermal Propulsion (NTP) systems are
powered by high temperature nuclear reactors that are used to heat up hydrogen propellant to
temperatures in excess of 3,000 K. Hydrogen is the only viable propellant for the NTP systems
because of its low molecular weight that generates the highest specific impulse (Isp) at the
maximum operating temperatures of these reactors. Hydrogen is pumped at cryogenic
temperatures and relatively high pressures to cool the rocket nozzle before entering the reactor
core. Heat removal in the rocket nozzle and reactor core areas transform subcooled liquid
hydrogen to superheated hydrogen gas. The evolution of hydrogen flow in the system involves
two-phase flow and heat transfer under subcooled, saturated, and superheated thermodynamic
conditions. In addition to the rockets, a nascent industry that may require modeling of this sort is
the hydrogen-fueled car industry.
This research effort includes a number of obj ectives. First, it is necessary to conduct a
literature search to determine the best battery of two-phase hydrogen tests to analyze. Using the
data from this test series, the next obj ective is to evaluate the quality of these data. The primary
obj ective is to improve the accuracy of predicted pressure drop of and heat transfer to two-phase
hydrogen in a forced convection, highly heated, internal pipe flow environment.
Mechani stically-based models are preferred, but correlations that provide improvements to
pressure drop and heat transfer predictions are considered acceptable. Since very high wall to
bulk temperature ratios can reasonably be expected with liquid hydrogen flowing in a heated
pipe, the effect of radial temperature variation will necessarily be included. This goal will
include the generation of void fraction, quality, and slip information that must be evaluated
against data. It is an obj ective to develop a predictive model for one or more of these parameters
so the pressure drop can be predicted. An important criterion of success is to reproduce the
pressure drop data with minimal error using the predictive model.
Additionally, it is an obj ective to either improve on the accuracy of current heat transfer
models, or at least review the current understanding of this subj ect and recommend models to use
for two-phase hydrogen.
Pressure change for a vaporizing fluid is comprised of three contributing effects:
momentum decrease due to increasing fluid velocity as it vaporizes, friction between the fluid
and the wall, and pressure change due to a change in height of the fluid. From Collier (1981),
these three terms for a homogeneous flow are
a -G (1.1)
: = -p (1.3)
In these equations, p is pressure, a refers to acceleration, F refers to friction, z refers to
elevation, G is mass flux, u is velocity, frP is the two-phase friction factor, p is density, D is tube
diameter, and g is gravity. The frictional term is often determined by calculating what the
frictional loss would be if the entire flow were liquid, then applying a two-phase frictional
multiplier. This method was developed by Lockhart, Martinelli, Nelson, and others at the
University of California in the 1940's (Martinelli et al., 1944, 1946; Martinelli and Nelson, 1948;
Lockhart and Martinelli, 1949). The form of their equation is
F F ;(1.4)
The two-phase frictional multiplier, $102, iS modeled as a function of flow quality and system
pressure. Determining this multiplier is the goal of much research, particularly since the research
of Martinelli and coworkers was limited to atmospheric pressure.
The frictional multiplier can take four different forms. Two of them develop from using
only the liquid or vapor mass that is present in the flow and are represented by a single letter
subscript "l" or "v" to indicate liquid or vapor. The other two develop from using the entire flow
as either liquid or vapor, and have a two-letter subscript "lo" or "vo" to signify that the entire
flow is liquid or vapor. Traditionally, liquid conditions are used in evaporating systems, and
vapor conditions are used for condensing systems. Relations can be developed between these
various frictional multipliers.
A correlating parameter developed by Martinelli and his coworkers is Xut, which they
determined to have the following form
X,, =1- x p'(1.5)
In this equation, x is quality and C1 is viscosity. The sub scripts 1 and v refer to liquid and vapor
phases, respectively. This parameter, in various similar forms, has been used to evaluate the
frictional multipliers. The exact forms of the frictional multiplier models depend on the nature of
the flow of the liquid and vapor phases turbulent or laminar. Thus, there are four combinations
of turbulent/laminar flows.
The currently preferred model for predicting the frictional multiplier was developed by
Chisholm (1973). His model uses a property index, composed of the property terms in equation
(1.5) above without the quality term. Mass flux is also a factor in his model.
B is a function of the saturated liquid and vapor densities and the flow regime combination of the
two phases (i.e., turbulent-turbulent, turbulent-laminar, laminar-turbulent, or laminar-laminar),
and x is quality.
Hendricks et al. (1961) derived a version of Xut for the peculiar case of inverted annular
flow. The primary difference in the derivation is the position of the phases in the flow.
Martinelli et al. (1944) assumed that liquid was adj acent to the wall and vapor was at the tube
core, or that the flow is homogeneous. For convective hydrogen, the flow is usually better
described as separated and the phase adj acent to the wall is vapor, not liquid. However,
Hendricks et al. (1961) determined that the correlating parameter, presented in equation (1.5),
was the same in both cases.
Papadimitriou and Skorek (1991) processed data from two of Hendricks' tests with their
one-dimensional thermohydraulic model THESEUS. They observed that the pressure drop due
to viscous shear forces is about 100 times smaller than that caused by momentum change. The
Chisholm (1973) method was used to model the two-phase friction multiplier.
John Rogers at Los Alamos Laboratories contributed significantly to the understanding of
parahydrogen flow friction characteristics in the 1960's (1963, 1968). His efforts included
extending Martinelli's friction multiplier quantification work beyond one atmosphere. His
results were based upon theoretically determining the values of vapor void fraction and its
derivative with respect to pressure at one atmosphere and at the critical point pressure with
quality as a parameter, then interpolating the curves of void fraction verses pressure between
these boundaries for the specified qualities. The empirical equation he developed for turbulent-
turbulent flow as a result of his work is
2 1+ 8 8187 0.1324Gpent p)+ 0.03966(p cnt -P)3 1.7
E = 1.896x 2. 646x2 +1.695x3 (1.8)
In this equation, the subscript crit refers to the critical condition. Note that pressure is in
atmospheres, x is quality, and the correlation gives the multiplier for the pressure drop for the
liquid only in the tube, not for all the flow considered as liquid. Comparison of predicted versus
experimental pressure drop with one set of parahydrogen data at various pressures indicated
good agreement, with the error generally smaller at lower system pressures.
More recent work on separated two-phase flow pressure drop and heat transfer with vapor
core was performed by Fu and Klausner (1997). In their work, conservation of mass,
momentum, and energy laws are applied with closure relationships for vapor-liquid interface
friction, liquid film turbulent viscosity, turbulent Prandtl number, and liquid droplet entrainment
rate. Their results compared with 12 data sets of upflow and downflow were good. Although
this theory assumes a liquid film and vapor core, the general procedure may prove useful with
inverted annular flow after making the appropriate modifications to the various correlations.
John Chen published a correlation in 1966 that was based on the superposition of heat
transfer caused by forced convective flow and by bubble generation. These terms are referenced
with subscripts mac and mic for macroscopic and microscopic effects.
h = h,,, + h,,, (1.9)
h,,, = hF(g,,) (1.10)
h...._2 0 01i5k O079IU 9 045hl 04P49P 2
her = .012 51 2p 02 02 T Ts 024 IP,(T,)-g oyS(F, Ret) (1.11)
hi is the heat transfer coefficient associated with single phase liquid flowing alone in the pipe. In
this equation, k is thermal conductivity, c, is specific heat, o is the surface tension, hiv is the heat
of vaporization, T is the temperature, S is a boiling suppression terms, F is an enhancement term,
and Re is the Reynolds number. The sub scripts w and s refer to wall and saturation conditions.
hIOIZ = 0.02 ePl0 4 k;` (1.12)
Rez (1 ) (1.13)
In the above equation, Pr is the Prandtl number. The model incorporated heat transfer data from
water, methane, pentane, and cyclohexane in the form of two factors, F and S, that were applied
to the two different heat transfer components. His model proved to be very successful.
Modifications to the original model have been proposed. Collier provided curve fits for the
factors F and S as a function of Xtt.
Shah (1984) developed a correlation for saturated flow boiling for both vertical and
horizontal tubes as a function of the Colburn, boiling, and Froude numbers, represented as Co,
Bo, and Fr.
h = hzf (Co, Bo, Frze ) (1. 14)
C=1x p (1.15)
Bo = (1.16)
Frze =~ (1.17)
In the above equations, q" is the heat flux.
Schrock and Grossman (1959) reviewed vertical, upward flowing boiling heat transfer data
for water with the following result;
h = hC,I Bo+C2 1 06 (1.18)
where C1 and C2 are COnstants with values of 7390 and 0.00015, respectively.
Gungor and Winterton (1986, 1987) developed the following for vertical, convective flow
S0 75 041 I(.9
h,,ea = hi 1+ 3000Bo0 86 1(.9
Bj orge, Hall, and Rohsenow (1982) developed a correlation for vertical, internal, upward
forced flow boiling for qualities above 0.05. Note that this correlation is a superposition of heat
fluxes as opposed to Chen's superposition of heat transfer coefficients.
h = fr'(1.20)
qto = 4f rc 4fa 1- (1.21)
qj, = F, Pr; T, T (1.22)
F, = 015 -+032 (.3
Xtt Xtt (.3
C2 = f (Prt, Rel) (1.24)
q,= u#hi / 9,1,:p';,'~/8(~ 5/8j 1/ (1.25)
a #lh:h! Pt Pv)9i a :T::
(TK T,)1 = Th, U1) (1.26)
In the above equation, u is the specific volume.
Kandlikar (1990) developed the correlation below for vertical and horizontal flow boiling
heat transfer in tubes;
h,,,a = h, CICoc"(25Frzec+ +CzRoc F, (1.27)
where the constants C1 through Cs can each take on two different values depending on the
Colburn number. The value of the constant FK depends upon the fluid being modeled.
Hendricks et al. (1961) performed experiments with hydrogen flowing inside a highly
heated tube. Nusselt numbers were determined from measurements. The deviation from these
measurements that the calculated Nusselt numbers generate approaches 80% at large values of
Martinelli parameters, and roughly 40% at low values. As a result, the researchers found it
necessary to curve fit the Nusselt number ratio as a function of the Martinelli parameter. This
technique significantly improved the predictive accuracy, with most experimental Nusselt
numbers lying within +15% of the curve fit.
The model for the Nusselt number, Nu, that Hendricks et al. published for their forced
convective heat transfer for flowing hydrogen was as follows:
Mr ic 'c~ (1.28)
exf 0.611 +1.93X,,
Miat = 0.023 Reo.s Pr0.4 (1.29)
Re = p ,aD(1.30)
Pr,,, = ?r lx(1.31)
The result of this method can be seen in Figure 1-1.
Hendricks et al. (1966) developed a similar equation to correlate the combined subcritical
data from TN 3095 (1966) and TN 765 with somewhat worse results due to data scatter. It is
critical to note that this correlation excludes those data for which the thermodynamic equilibrium
quality indicates subcooled flow. This excludes possibly up to one-third of the 612 points in the
data set! The authors remarked that this equation should describe subcritical convective fi1m-
boiling data up to pressures near the critical pressure when non-equilibrium characteristics are
small. The correlation based on the remaining data is
exp, f + .5(1.32)
Miat, ,,,, 0.7 + 2.4Xr,x,
In this equation, subscript f refers to properties evaluated at the average of wall and bulk
temperatures. Sub script fm refers to mean film conditions, e.g., using the density defined above
with subscript f~m.
These authors also developed a correlation based on a pseudo quality with similar results in
accuracy. This correlation included some of the subcooled data, but far from all of it. Their
assessment of this correlation was that it covered the liquid-hydrogen data for convective fi1m
boiling from a slightly subcooled state through two-phase and well into the superheat region.
It should be noted that all models presented above perform very poorly on the data
addressed in this dissertation. The exceptions, of course, are the models from TN 765 and TN
3095. Chen' s (1966) model performed the best of all the others, while those of Shah (1984),
Schrock and Grossman (1959), and Gungor and Winterton (1986, 1987) predict convection
coefficients that are hundreds of times too high. This is due in large part to the form of the
Reynolds number used by Hendricks et al.
Heat transfer coefficient models have been developed that focus specifically on the flow
structure of the data in this dissertation inverted film boiling (IFB). In general, the forms fall
into two categories: those that attempt to capture the heat transfer mechanics of a highly
convective flow, and those that focus on the effects of buoyancy. The convective models
generally expand on the basic Dittus-Boelter model, while the buoyancy models usually take the
form of the Bromley model (1950), which was developed for laminar film boiling, and is
analogous to film condensation theory. These low velocity models are generally used to model
Loss Of Coolant Accidents (LOCA' s) in the nuclear industry.
Bromley's model (1950) is an extension of theory developed by Nusselt (1916) for laminar
film condensation on a horizontal tube. His heat transfer coefficient model for laminar film
boiling from a horizontal tube is
h=g 0.6 (1.33)
where hfg' is the effective latent heat of vaporization accounting or vapor superheat. Numerous
film boiling models expand on this basic form.
Bromley et al. (1953) extended his own model to include forced convection. For low
velocity flows, he determined the following:
hg g 0.6 (1.34)
co= DpD~, ATk
where the subscript 'co' refers to convection only excluding radiation heat transfer. AT refers
to the temperature delta between the wall and the centerline. For higher velocity flows, he
proposed the following;
h = 2.~IkPhf1 (1.35)
Here, the enthalpy of vaporization is defined as
h = hg 1+ h '8T (1.36)
Berenson (1961) modified the Bromley model by incorporating the hydrodynamic
instabilities predicted by Taylor instability theory. He published the following result;
h = 0.425~~T [5 2p~:j (1.37)
The vapor properties are evaluated at the mean film temperature, liquid properties at saturation
temperature, and 0.425 is used as a coefficient instead of the 0.62 in equation 1.34 above to
account for enthalpy of vaporization to superheated conditions.
The analogy to liquid film condensation has been extended to the assumption of turbulent
flow in the vapor film. Wallis and Collier (1968) presented conclusions from this theory and
=() 0.056ReU [Pr Grr *]i (1 .38)
where the modified Grasshoff number is defined as
Gr* = : p (1.39)
An obvious characteristic of the heat transfer coefficient models presented thus far is their
inclusion of buoyancy effects. The models that focus on highly convective flows ignore
buoyancy effects. In these models, heat transfer is quantified within the framework of the
traditional Dittus-Boelter forced convection concept.
Dougall and Rohsenow (1963) developed the following model for dispersed flow and
inverted annular film boiling (IAFB) of Freon 1 13:
h = 0.023 ""Reo Pr~O4 (1.40)
w x w (1 -x)P
Re, = p,~,Dh (1.41)
The velocity term applied here is the throughput velocity. This effort focused on low quality
mass flows. In this equation, w is the mass flow rate, and the subscript s refers to saturation.
A subsequent research program that focused on higher mass qualities was completed by
Laverty and Rohsenow (1964, 1967). Their IFB nitrogen studies included visual analysis of the
flow structure. Through theory, they determined that a significant amount of superheat was
present in the vapor. As a result, they determined that it is impossible to obtain a simple
expression for the overall heat transfer coefficient, although they did present a model for their
data, presented below. Instead, they presented arguments based on the Dittus-Boelter model to
set the upper bound and approximate value of the heat transfer coefficient. Their published
model is as follows:
h = 0.023 pv bD jPr04 kv (1.42
In this equation, the sub script b refers to bulk conditions.
Forslund and Rohsenow (1968) also used nitrogen to improve the analysis of Laverty and
Rohsenow (1964). Improvements focused on droplet breakup due to vapor acceleration,
modified drag coefficients on accelerating droplets, and a Leidenfrost heat transfer from the wall
to the droplets at lower qualities. Test conditions covered the quality range from saturation at the
inlet to 35% to 315% at the exit. They focused on estimating the magnitude of departure from
thermal equilibrium and droplet size. They concluded that vapor superheating was significant -
up to 50% in vapor quality. The heat transfer model they proposed, presented below, attempted
to modify the Reynolds number to reflect conditions in the vapor:
h=~ l~ i 0.1 G r"x(1- x~"1P ', (1.43)
Kays (1980) presented an analysis for heat transfer between parallel plates. This model has
been used by Hammouda (1996, 1997) in his modeling of IFB nitrogen. The Kays model is
below. Note that the length dimension is the film thickness, 6.
h= =r 03 + 0.0028 Pr o 645 Re-l (1.44)
Bailey (1972) presented a buoyancy-based heat transfer model as follows:
h = (1.45)
D= vk,3g, (T -T)( T,)hr ) 2
Takenaka at Kobe University in Japan is associated with a number of IFB studies from the
late 1980's and early 1990's. In general, his working fluids are R-113 and nitrogen flowing
upward inside a vertical heated tube. Heat fluxes and mass velocities are generally an order of
magnitude or more smaller than those addressed in this dissertation. His work is unique in that it
is the only research found in the literature search that produced a flow regime map for IFB.
Takenaka et al. (1989, 1990) found that heat transfer coefficients, as a function of equilibrium
quality, did not vary with heat flux or inlet subcooling, but segregated consistently with mass
flux. As a result, their IFB flow regime map uses mass flux and equilibrium quality as
coordinates. As equilibrium quality increased, higher mass fluxes produced higher heat transfer
coefficients at the same quality. They found that the Nusselt Number predicted using the
Dougall-Rohsenow (1963) model were reasonably close to their data. Takenaka also worked
with Fujii (Fujii et al. (2005)) to investigate pressure drop in IFB. Because the mass velocities
are very low, the pressure drops measured in the nitrogen flow are in general much smaller than
those exhibited in the data of this dissertation. They found that the pressure drop characteristics
correspond well with the heat transfer characteristic map.
Hammouda et al. (1996, 1997) investigated the effects of mass flux, inlet subcooling, and
system pressure on the heat transfer coefficient using R-12, R-22, and R134a as the working
fluids. The characteristic shape of the heat transfer coefficient as a function of equilibrium
quality is consistent with those in Takenaka' s experiments. The effect of mass flux is the same,
but varying the inlet subcooling measurably segregated Hammouda' s data while Takenaka noted
no such influence. Different results are also noted in the effect of heat flux on the heat transfer
coefficient. While Hammouda' s data show that higher heat flux increases the heat transfer
coefficient, Takenaka' s data shows very little, if any, effect. Hammouda also observed that
higher system pressure increases the heat transfer coefficient a parametric effect that Takenaka
Ishii has been involved with a number of experiments that focused on the flow regime
characteristics and transition criteria of post-critical heat flux (IFB) flows. Ishii and De Jarlais
(1986) investigated the basic hydrodynamics of this flow regime. The mechanisms that
disintegrate the liquid core were investigated, as well as the formation and entrainment of
droplets in the vapor annulus. The experimental portion of this work involved adiabatic two-
phase flow, resulting in a flow regime transition criterion based on the Weber number. Ishii and
De Jarlais (1987) presented experimental data for an idealized IFB flow generated by inj ecting a
liquid inside a vapor annulus in up-flow using Freon 113. Fluid heating was incorporated into
the test setup. Visual observations revealed the nature of the flow structure to include smooth
IAFB, agitated inverted annular film boiling (AIAFB), followed by inverted slug film boiling
(ISFB) and dispersed film boiling (DFB). Obot and Ishii (1988) extended this work with the
same fluids and test setup. More extensive results of flow regime transition are presented. Ishii
and Denten (1990) continued this work to investigate the effects of bubbles present before post-
critical heat flux is attained on the IFB flow regimes and their transitions. Three regimes were
observed; rough wavy, agitated, and dispersed ligament-droplet. They found that the flow
pattern in IFB depends upon the nature of the pre-CHF flow. A general flow regime transition
criterion between the agitated and dispersed droplet regimes is given based on conditions at
dryout. This correlation includes void fraction at this point as an important parameter. Babelli et
al. (1994) used the same experimental apparatus to continue the research. He concluded that the
most significant flow regime is the agitated regime, since the large interfacial surface generated
in this regime probably correlates with high momentum and heat transfer. A correlation for the
axial extent of this flow regime was proposed, again dependent upon the void fraction at the
point that CHF occurs. It should be noted that all of the work performed by Ishii and his
associates was performed for the purpose of better understanding nuclear reactor LOCA. As
such, the flow velocities are quite low compared with the data in this dissertation.
Per Ottosen (1980) published the first known results from the use of y-ray absorption to
measure void fraction in low Reynolds number IFB nitrogen. He observed the transition from
IAFB to DFB at void fractions between 80-90%. These void fractions were typically attained by
the point at which equilibrium quality was 20%. Given that superheat will be present, this
equilibrium quality probably relates to a lower actual quality. Since his work was in support of
understanding LOCA' s and reflooding, his fluid velocities were low. Also, the work was
executed at a constant temperature condition instead of a constant heat flux condition, as is more
often the case. Nonetheless, trends in heat transfer coefficients as a function of mass flux are
Experiments using hydrogen as the working fluid are rare. This is primarily because of the
dangerous nature of the fluid. Hendricks (personal communication, 2005) relates that, in the
series of experiments during 1961 and 1966, the building in which they worked was evacuated of
people, and emergency personnel were notified of each experiment. It is determined through the
literature search that the only published hydrogen experiments performed in the United States
that present heat transfer data occurred in support of the manned space missions in the 1960's.
Published results from hydrogen experiments in the Soviet Union and Europe, though they likely
occurred, have not been found.
Core et al. (1959) performed experiments with hydrogen similar to those in TN 3095, but
with much fewer measuring points of pressure and temperature. Twenty-seven heat transfer tests
with liquid hydrogen were completed in the series. Since only test section inlet and exit
conditions were measured, the heat transfer coefficients calculated from these measured data are
overall average coefficients for the entire tube. The authors did not present a theoretical
correlation for the heat transfer coefficient. Their primary goal was to evaluate the utility of
hydrogen as a regenerative rocket nozzle coolant. Nonetheless, the data from this study may be
considered as complementary to the data of TN 3095, and therefore useful.
Wright and Walters (1959) found that stable film boiling of hydrogen could occur for wall
to bulk temperature differences as low as about 22 K to 28 K. Also, peak heat transfer
coefficients were about 10% of the magnitudes of those in nucleate boiling. Their film boiling
heat transfer coefficients were almost constant over the range of wall to bulk temperature
differences of 22 K to 167 K.
Papadimitriou (1991) presented results of a simulated rocket engine two-phase hydrogen
chilldown process using a modified form of Dougall and Rohsenow model in the computer
program THESEUS. The modification is a temperature correction,
applied as a multiplier on the model for the heat transfer coefficient. It was stated that this better
accounts for the real film conditions at high wall temperatures.
Many of the above forced convection models are based on the classic Dittus-Boelter model.
Variations on this standard model are implemented by using properties and flow conditions
calculated in specified ways. For example, the properties used in the Reynolds number could
represent bulk calculated values for the two phases, the vapor saturated condition, or superheated
vapor conditions. Below is the standard model for later reference:
h = (1.47)
In this model, the coefficient and exponents can be adjusted to fit the data. Common values are
0.023 for the coefficient and 0.8 and 0.4 for the exponents m and n. Unless stated otherwise,
these are the values used in this research effort.
The literature search has found a number of experiments that are peripherally related to the
data in the NASA data. However, the data addressed in this dissertation are rare or even unique
in several ways. First, the working fluid is hydrogen. As stated above, there are only three other
published reports of experiments with hydrogen in a convective, IFB condition. None of these
three experiments operated at the high mass flux levels of the NASA data. Finally, and most
importantly, the extent of measured parameters makes these data extremely valuable. These
measurements provide the means to theoretically analyze the pressure drop and heat transfer
characteristics of hydrogen, and to validate any proposed model or correlation.
exp~f.1 061 e1.95 XH
-02 .04 .06i.08.1 .2 4 ,6 .8 I 2
MARTINELLI PARAMETER, Xal
Figure 1-1. Ratio of experimental to calculated Nusselt number for the 1961 data.
mr pi B
MODELLINTG APPROACHES FOR TWO-PHASE FLOW
Two-phase fluid flowing in a pipe can have characteristics that vary in the axial, radial, and
azimuthal directions. Axial dependencies of properties and flow structure can result from
entrance effects, wall friction, turbulence, and heat addition. This dependency is typically not
ignored, since it is changes in conditions in the axial direction that interest engineers typically.
Radial dependencies can result from these same sources. Since it is a great simplification
to ignore this dependency, this is commonly done. Corrections can be applied to models that
explicitly ignore the radial dependency. For example, the effect of a radial temperature gradient
on fluid properties can be accounted for by multiplying the Dittus-Boelter Nusselt number by a
ratio of wall-to-centerline temperatures, usually raised to an exponent. Another example is the
drift flux model, the purpose of which is to account for radial variations is fluid density and
velocity. These two examples speak to the duel importance of neglecting radial dependencies in
the formal conservation equations while simultaneously including radial effects through semi-
Finally, azimuthal dependencies are usually important only in horizontal flow, where
gravity strongly segregates the liquid and vapor phases due to the large difference in densities.
In vertical flow with uniform heating, this dimension is typically confidently neglected.
There are four basic approaches that can be used to model the thermal hydraulics of two-
phase flow. Each method explicitly defines the number of independent conservation equations
used. The number of closure relations that link the corresponding conservation equations
increases as the number of conservation equations increase, so that the number of conservation
equations minus the number of closure relations will always equal three. While complexity
increases as the number of conservation equations increase, the variety of information obtained
about the flow also increases. This does not necessarily mean that predictions for pressure drop
and heat transfer will be better for a six-equation model compared with a three-equation model.
It simply means that more predicted information will be generated. The reliability of these
predictions will depend directly on the validity of the closure relations and assumptions used to
develop the overall modeling approach.
The most sophisticated model is called the two-fluid model, in which there are separate
mass, momentum, and energy equations for each of the phases. Closure relations must link the
corresponding equations for each phase; mass, momentum, and energy transfer rate terms are
defined at the phasic boundaries. The mass transfer term is relatively simple and is directly
related to the change in quality. The momentum and energy transfer terms are more complicated
at they depend on the momentum and energy associated with the newly vaporized fluid. They
also depend upon the interfacial shear and heat transfer rates two terms for which data are
difficult to obtain. These terms are usually developed in terms of theory and assumptions, or a
combination of theory and experimental findings.
In addition to the interfacial closure relations, there must also be relations for momentum
and energy transfer at the fluid-wall boundary. These conditions are usually determined with
more confidence because the experimental data that have been generated to understand these
conditions are more complete. Research in single-phase flow, which has been extensively and
reliably performed, often applies. For instance, the wall friction is a term that is of fundamental
importance to engineers, and therefore has been studied since the beginning of the science of
fluid mechanics. Heat transfer also is of fundamental importance. To simplify the analysis of
data from an experiment, the wall-to-fluid heat transfer boundary is usually established as one of
two conditions constant heat flux or constant wall temperature. With either, the wall heat
transfer boundary condition is well defined.
The next simpler model includes five equations. In this, the developer can choose which
conservation equation to simplify, but usually selects either momentum or energy. In two-phase
flows, it is commonly accepted that the pressures of both phases are the same. Therefore, the
momentum equations for the two phases are usually reduced to one. This is accomplished by
equating the interfacial momentum transfer terms, since they must be the same. That is, the
momentum that one phase loses at the interfacial boundary is gained by the other. This approach
has the great advantage of eliminating the interfacial shear stress term. Alternatively, the
developer may choose to equate the energy transfer terms in a similar fashion. This eliminates
the need to determine the rate of sensible heating of the liquid phase.
A further simplification is made by reducing the number of independent conservation
equations to four. In this case, there is sometimes a specific piece of information required, such
as velocity slip.
The simplest approach is the three-equation model, also called the homogeneous
equilibrium model (HEM). In this, equations of mass, momentum, and energy conservation use
properties that represent the mass-weighted values of the vapor and liquid phases. There is no
information regarding the separate velocities of the phases. Equilibrium quality is used, which
neglects liquid subcooling or superheating, and vapor superheating. In spite of its simplicity, the
HEM is often cited as a standard against which the results of other models are compared.
Flow Regime Analysis
When the more complicated models are used, it is frequently necessary to determine the
structure of the two phases relative to each other. The various structures in two-phase flow have
been distilled down to a few flow regimes. A heated two-phase flow progresses through these
flow regimes as it increases in quality. The specific set of flow regimes may be different for
different conditions. For example, flow through a horizontal pipe can experience separated flow
with the heavier species at the bottom of the pipe, and the lighter species at the top a flow
structure not developed in vertical tubes. Flow through vertical tubes can also progress through a
different set of flow regimes, depending upon the amount of applied heat. Low heat loads will
result if pre- Critical Heat Flux (CHF) conditions. The vapor phase is generated at the wall and
migrates to the center of the tube. Liquid is always on the surface of the tube wall until dryout
occurs at high qualities. After this, the liquid is dispersed as droplets in a continuous vapor
matrix. High heat loads can produce post-CHF conditions, or IFB, at very low qualities. In this
situation, the wall is too hot for liquid to remain. Vapor stays along the wall of the tube
throughout the increase of quality. The progression of flow regimes in IFB are IAFB, AIAFB,
and DFB. These flow regimes are presented in Eigure 2-1. If the mass flux is low, then ISFB can
occur after IAFB. Note that this figure, taken from Takenaka (1989), does not include the 'B'
for boiling in the regime nomenclature that this dissertation includes. IAFB is characterized by a
relatively smooth interface between the vapor and liquid. The liquid flows through an annulus of
vapor. The interfacial area is easy to determine assuming the void fraction is known. AIAFB is
characterized by a rough interface. The liquid core is still whole, or in separate, parallel liquid
filaments, but is rough such that determining its surface area is no longer a straight forward
calculation using void fraction. The area for heat and momentum transfer likely increases
relative to IAFB even though the amount of liquid is decreasing. Finally, in DFB, the liquid core
completely breaks up into drops and is carried along in the continuous vapor matrix. This flow
structure is very similar to pre-CHF dispersed flow.
Because the physics of the flow is strongly dependent on the flow regime, it is common to
base closure conditions and other modeling decisions on the local flow regime. Of course, this
requires that the various transitions between regimes be predictable. As pointed out in chapter
one, Ishii has put in significant effort to develop predictive models. His more recent work is
with heated Freon 113 in relatively low velocity conditions. Observations are that the void
fraction at the point of dryout has a significant impact on the flow regime transition correlation.
The correlation is as follows (Babelli et al. 1994):
= 55 (2.1)
D cr 0.854]-
In this relation, L is the length at which the flow regime transitions from IAFB to DFB, D is
diameter, Clf is the fluid viscosity, jJ is the volumetric flux, o is the surface tension, and oes is the
Takenaka (1989, 1990) generated a flow regime map for IFB, as shown in Figure 2-2
where coordinates are equilibrium quality and total mass flux. Note that inlet velocity is used on
the ordinate instead of mass flux, but his final map actually used mass flux. For his test
conditions, this map predicted the IFB regimes he viewed.
- 0.1 0
0.3 0.4 0.5
Figure 2-1. Various flow regimes for IFB (Takenaka, 1989). The ISFB regime on the left is
associated with low mass flow rates.
Ifsub =10 n
Figure 2-2. Flow regime map generated by Takenaka for IFB (1989). Flow regimes are IAFB in
region (a), AIAFB in region (b), DFB in region (c), and ISFB in region (d)
TEST DATA DESCRIPTION AND EVALUATION AND MODEL DEVELOPlVENT
Description of Experiments
As referred to earlier, the data used to validate the model were generated at NASA Glenn
Research Center (formerly Lewis Research Center) and published in two separate technical
notes, NASA TN 765 and NASA TN 3095, in 1961 and 1966, respectively. These data will be
referred to collectively as the NASA data to distinguish it from other hydrogen experiments, or
as the 1961 and 1966 data when the data from the individual reports are discussed. The
experiments were performed in support of rocket engine modeling for the US manned space
The experimental setup for the 1961 experiments is presented in figure 3-1. Hydrogen was
stored in a large tank and pressurized by gaseous hydrogen to force it through the system. Piping
from the tank to the test section and the test section were enclosed in a vacuum environment to
eliminate convection heat transfer to the piping and working fluid. The vacuum container was a
stainless steel cylinder 38. 1 cm in diameter. Heat was generated inside the tube metal by
applying a voltage across its length. The power supply for heat generation was external to the
vacuumed environment. Therefore, the leads for the voltage supply, along with instrumentation
leads, were passed through the wall of the vacuum chamber. The voltage was applied to the
heated test section through copper flanges brazed to the tube. It was found that unevenly brazed
joints distributed the power unequally circumferentially in the tube. Therefore, multiple
connections to the buss bar were made and the brazed joint was X-ray inspected. After passing
through the heated test section, the hydrogen was completely vaporized and then exhausted
through the roof of the facility into the atmosphere. All system flow conditions were remotely
controlled. The system pressure and flow rate were set by valves upstream and downstream of
the test section.
The setup of the 1966 experiments is similar to that of the 1961 setup. Figure 3-2 presents
the configuration. More useful information is given in the 1966 report that will be repeated here.
It is pointed out that the liquid hydrogen storage tank is enclosed in a vacuum to mitigate
heating. This in turn is contained within a liquid nitrogen radiation shield to mitigate conversion
of parahydrogen to orthohydrogen. Finally, this is contained within a foam insulated container.
The liquid hydrogen was forced through the flow system using gaseous normal hydrogen as a
pressurant. Just upstream and downstream of the test section were mixing chambers of high
turbulence in which the fluid bulk pressure and temperature were measured. Mixing the fluid in
the mixing chambers and having an entrance length to the test section were found to be important
since there could be some thermal stratification of the liquid as it is transferred from the storage
tank to the test section, with warmer liquid adj acent to the wall and colder liquid in the center of
Five different tube diameters were used in the NASA experiments, ranging from 0.48 cm
to 1.29 cm inside diameter, and all were vertical with hydrogen flowing upwards. The heated
test section length in the 1961 and 1966 experiments are 30.5 cm and 61.0 cm long, respectively.
Straight, unheated approach lengths were included in all test sections; approximately 12.7 cm for
the 1961 tests, and 30.5 cm for the 1966 tests. Approach sections and test sections were
contained within the vacuum environment. Figure 3-3 and 3-4 present the test sections for the
1961 and 1966 data, respectively.
Heat fluxes and mass flow velocities are very high, and tube diameters are similar to those
used in regeneratively cooled rocket engine nozzles and other rocket engine piping. The
experimental conditions of these data reflect the nature of hydrogen flowing in a rocket engine.
Table 3-1 presents a summary of test conditions.
Heat was generated within the test sections by applying a voltage across the length of the
section. Care was taken to ensure a uniform weld of the copper flange around the tube so that
current would flow uniformly down the tube. As Hendricks (personal communication, 2006)
stated the problem,
The most damaging effect [on uniform heat generation] was the braze j oint between the
tube and the copper flange. Erratic j points distributed power unequally into the tube and the
current paths in turn did not heat the tubes properly.
Paths for undesired heat transfer into or out of the system that have not already been
addressed were either analyzed or otherwise considered by the authors of the 1961 data and
determined to be insignificant.
All test sections were instrumented for local static pressure, tube outside wall
temperatures, and local voltage drops. Accuracy was of paramount importance (Hendricks,
personal communication, 2006). When initial results of tube wall temperatures ran counter to
anything previously experienced or expected, double and triple instrumentation redundancy was
implemented to determine the source of the "error". Data published in the reports represent
those deemed most accurate of the redundant measurements. Figure 3-5 illustrates some
specifics of thermocouple and pressure tap installations. The 1961 report gives no information
about instrumentation measurement accuracies. The 1966 report gives information on this
subj ect, and in general, the accuracies of instrumentation and measurements in the 1961 data are
consistent (Hendricks, personal communication, 2006).
All test sections had 12 thermocouples along the outer surface of the heated lengths, plus
inlet and exit temperatures in the mixing chambers. Thermocouples were either copper-
constantan, which were silver soldered, or Chromel-Alumel, which were welded in place.
Connections to the tube outer wall were made with great care to avoid affecting the test
conditions or measurements. Leads from the thermocouples were 30 gage wire. Circumferential
thermocouple placements were intended to determine the circumferential uniformity of power
distribution in the tube and as checks for accuracy. The cold junction was atmospheric boiling
nitrogen in the 1961 data, while the 1966 data used either liquid nitrogen or ice.
Thermocouple accuracy was determined by the recording system accuracy, standard
calibration, lead wire and junction temperature gradients. The mixing chamber fluid temperature
measurements were estimated to have less than 1% probable error. Multiple thermometers in the
mixing chambers agreed to within 1.1 K at the inlet and 5.6 K. No percent accuracy is given
for the tube surface thermocouple measurements. Tube surface temperatures were checked by
comparing multiple thermocouple readings attached by different techniques. The readings
usually agreed to within 15.6 K.
The 1961 data had five static pressure taps spaced along the length of the test section and
one at each of the inlet and exit mixing chambers. These pressure measurements were not
differential relative to a datum. The other four tubes from the 1966 experiments had three static
pressure taps spaced along the test section, and one at each of the inlet and exit mixing chambers.
These pressure measurements were differential relative to the pressure reading just upstream of
the test section inlet. No pressure taps on any of the five test sections were located at the same
axial location as a wall thermocouple. To complete the pressure data set, smooth curves were
hand-fitted through the measurements. From these curves, pressure values were interpolated at
the locations corresponding to the 12 thermocouple measurements. Commercial transducers
with a maximum of 1% full-scale nonlinearity were used. Readings from these transducers were
confined to half of the full scale. Therefore, errors from the pressure readings were estimated to
be 2%. Unfortunately, the range of the transducers is not given, and efforts to discover this
information have been unsuccessful. The differences in local static pressure measurements were
found to agree with differential pressure measurements to within 20%. This was reported to
correspond to an absolute static pressure measurement uncertainty of 1%.
Mass flow rates were measured both upstream and downstream of the test section. A
venturi was placed upstream of the test section and a sharp-edged flow orifice was placed
downstream of the heat exchanger. A second venturi, primarily used for flow control, was also
used for mass flow measurements. Measurements from these were compared for accuracy, and
all agreed to within 3%.
Local values of voltage drops were measured by eight voltage taps along the length of the
heated test section to assist in determining local power generation. Two sets of voltmeters and
ammeters that had independent shunts or taps were used. These incremental measurements of
power input were summed and compared with the overall power input measured by voltage and
ammeter taps at the bottom and top of the test section. Agreement between these two methods
was good (Hendricks, personal communication, 2006). Accuracies for these measurements are
stated to be +1%.
The values of the eight voltmeter measurements were not included in either publication.
However, as will be explained later, these measured local voltage drops appear to have been used
to determine local heat transfer coefficients, and in this sense, the local voltage drops are
The literature search has revealed five maj or experimental efforts investigating the heat
transfer characteristics of convective internal pipe flow boiling hydrogen. Two of these are the
1961 and 1966 NASA reports that are the focus of this dissertation. The other three were also
performed during the early stages of the U.S. manned space program. These studies were
scrutinized for possible use to validate the NASA data set.
Comparison with Similar Data
Core et al. (1959) performed experiments with hydrogen similar to those in the 1966 data,
but with much fewer measuring points of pressure and temperature. Twenty-seven heat transfer
tests with liquid hydrogen flowing through an electrically heated stainless steel test section, 6.35
cm long and 0.213 cm inside diameter, were completed in the series. Each test comprised a
number of different steady state conditions, isolating the effect of changing inlet pressure, mass
flux, or heat flux. As a result, there are a total of 164 steady state conditions, with two points of
heat transfer coefficient measurements each, in the set. Only the inlet pressure was measured, so
a pressure loss analysis cannot be compared with data. The authors did not present a theoretical
correlation for the heat transfer coefficient. Their primary goal was to evaluate the utility of
hydrogen as a regenerative rocket nozzle coolant. This source stands out as the only one that
presents wall superheats that are likely to represent transition boiling conditions. While most
experimental results indicate that transition boiling occurs between wall superheats of 5 K and 20
K, the data in this experiment show some superheats between these values. Therefore, these data
may represent results from transition boiling. Table 3-2 presents comparisons of heat transfer
coefficients averaged from the two points of measurement on the test section, compared with
runs with similar conditions from the NASA data set. The Core et al. data set includes calculated
equilibrium qualities based on pressure and enthalpy. Negative equilibrium qualities were set to
zero. Therefore, inlet subcooling is not known. The two calculated heat transfer coefficients for
each run in the Core et al. data are averaged and compared with the average heat transfer
coefficient for runs with similar conditions over the same equilibrium quality range in the NASA
data set. Sets of compared runs are separated by bold lines in the table. The first runs listed in
each comparison is from the Core et al. set, while the second listed run is from the NASA data
set. The RMS difference between these comparisons is 46.2%.
Wright and Walters (1959) experimented with liquid and vapor hydrogen flowing in a 15.2
cm long and 0.635 cm inside diameter heated tube. Most of their 35 steady state liquid hydrogen
experiments were pre-CHF, with 11 runs showing wall-to-bulk temperature differences
consistent with IFB. In fact, their data show a marked gap in wall-to-bulk temperature
differences between 2.8 K and 22.2 K. Temperature differences between these values were not
obtained. This gap is consistent with a transition in flow regime from pre-CHF and CHF
conditions to IFB. They concluded that stable fi1m boiling could occur for wall to bulk
temperature differences as low as about 22 K. Test section pressure measurements were not
obtained. There are three runs from their data set with conditions similar to several runs in the
1961 data. Table 3-3 presents the test conditions and average heat transfer coefficient over the
tube length. Note that the average heat transfer coefficient listed for the 1961 data represent an
average of points two through six. This omits the first point that is affected by inlet conditions
and concludes at approximately 15 cm into the test section. The heat transfer coefficients from
the two different test series agree well.
Lewis et al. (1962) experimented with boiling hydrogen and nitrogen flowing upward in a
type 304 stainless steel, electrically heated vertical tube 41.0 cm long and 1.41 cm inside
diameter. Critical heat fluxes corresponding to transition to IFB were determined over a range of
flow rates, heat fluxes, and qualities. They noted that the maximum CHF increased with
increasing mass flux and decreased as the point of transition occurred farther into the tube.
These findings are consistent with the interpretation of runs 22, 26, 29, and 30 from the NASA
data in figure 3-11 that will be discussed later. The mass flow rates in these experiments were so
low that no measurable pressure drops were observed. Wall superheats were similar to those
observed in the NASA data, with a maximum wall superheat of 500 K. Since mass fluxes and
heat fluxes are an order of magnitude lower than in the NASA data set, there are no test
conditions that are similar enough to warrant a comparison.
Table 3-4 summarizes the test conditions of the three forced convection heated tube flow
boiling hydrogen experiments discussed above and the NASA data. From the data in these three
experiments and other hydrogen experiments in geometries other than internal tube flow, it can
be said that transition boiling occurs between 5 K and 20 K. Review of tables 3-2 and 3-3 show
that the data from the NASA experiments are reasonably consistent with results from other,
similar works. From this comparison, it is determined that the NASA data are, in general, valid.
From the 1961 data, it is obvious that axial heat conduction occurs in the tube wall. Using
the finite difference heat transfer theory presented by Incropera and DeWitt (2002), a Fortran
program was generated to model the end axial heat conduction effects for the purpose of
determining the data that are affected and should therefore be omitted from the analysis. It was
assumed that curvature effects on axial conduction were negligible. Therefore, a two
dimensional infinite plate with axial and radial heat conduction was used to approximate the tube
geometry. The middle of the length of the plate corresponds to the beginning of the heated test
section. Left of this position is the unheated approach section, while right of this point is the
section in which heat is generated by electrical current.
To ensure that the imposed boundary conditions did not affect the solution, lengths of 50
wall thicknesses were generated on either side of the midpoint, for a total length-to-thickness
ratio of 100. It was found that the number of radial nodal points were not crucial to generating
acceptable results, so a minimum number of five nodal points were selected in the y direction,
with nodes one and five at the tube inner and outer walls, respectively. For the length-to-
thickness ratio of 100, this required 401 nodal points in the x direction. Figure 3-6 presents the
nodal structure and applied power distribution. Note that the distribution in the x direction is too
close to discriminate separate nodes, and the power generation is typical.
The applicable energy equation is
~+ + = 0 (3.1)
Dx" 2 2 k
In this geometry, x is the axis parallel with the flow, and y is the radial direction. Also, q"' is
the heat generation rate per unit volume. The variation in thermal conductivity as a function of
temperature will have only a very small impact on the results provided a representative
temperature is used to select the constant thermal conductivity. The thermal conductivity can
therefore be assumed constant in the analysis.
The four boundary conditions applied to this problem are:
1. T(x + -oo,y) = 25 K (a representative liquid hydrogen temperature in approach section)
2. 8T/8x (x + +oo,y) = 0 (adiabatic boundary far into heated test section)
3. -k8T/8y = h(Tw-Tb) at (x,y=0) (conduction = convection at wall/liquid interface)
4. 8T/8y = 0 at (x,y=Y) (adiabatic surface at tube outer wall)
In the heated section, boundary condition three assumes that the axial heat transfer is much less
than the radial heat conduction at the wall-liquid interface. To use this boundary condition, an
estimate of heat transfer coefficient that supports the purpose of the particular scenario at hand is
For each problem, the following four parameters must be specified; wall thickness, wall
thermal conductivity, heat generation rate, and fluid-to-wall heat transfer coefficient. The heat
transfer coefficients used in this analysis come from those values listed in the 1961 data set at the
first point, which is 1.4 mm above the heated section inlet. The algorithm was iterated until the
maximum difference in temperature in adj acent iterations was less than 1.0E-6 K.
The computer model was validated through five observations. First, the boundary
conditions at the left and right hand sides of the tube are satisfied, as is the boundary condition
corresponding to the outside of the tube wall. Second, it is logical that the point of largest
temperature slope should occur at the point that heat generation starts. Every scenario has
satisfied this requirement. Third, the effect of varying the parameters listed above affect the
results in a reasonable way. For example, increasing metal thermal conductivity causes the
effect of heat conduction to be felt deeper into both sides of the point of heat generation. Fourth,
magnitude of predicted inner and outer wall temperatures are reasonably close to those published
in the 1966 report (1961 report did not publish outer wall temperatures). Two runs, seven and
11, were selected at random for comparison purposes. For run seven, the inner and outer wall
temperatures are 231 K and 269 K, respectively, while the model calculated 207 K and 238 K.
Run 11 inner and outer wall temperatures are 461 K and 482 K, with model predictions of 412 K
and 431 K. Finally, the difference in tube inner and outer wall temperatures in the heated portion
reasonably agree with published data. Again, using runs seven and 11, the published differences
are 38 K and 19 K, while the model results are 31 K and 19 K. These differences are deemed to
be well within the uncertainty in the four parameters and errors associated with the model
assumptions for the intended purposes of this analysis.
Figures 3-7 and 3-8 present inner and outer wall temperatures for the scenarios in which
thermal conductivity and wall thickness are parameters. The two dimensional effects are
noticeable in the right hand portion of the tube.
To evaluate the effect of each parameter listed above on the tube end wall temperatures,
high and low values of each parameter were run, with all other parameters set to nominal values.
The length from the heated section inlet to the point at which 95% of the Einal temperature is
achieved was determined for each run and compared. Large differences in lengths by which
95% of the Einal temperature is achieved indicate a significant parametric effect on tube end axial
heat transfer. Figure 3-9 presents the results of the computer model. Since the difference in
outer and inner wall temperatures is small, only the outer wall temperatures are presented for
each scenario for clarity. Table 3-5 shows the distances in thicknesses from the heat section inlet
at which 95% of the Einal temperatures are achieved. This analysis suggests that the effect of end
axial heat conduction in the tube metal increases with increasing thickness and thermal
conductivity, and decreasing heat transfer coefficient. It is approximately independent of heat
flux. For a given test section, wall thickness and thermal conductivity are determined. The
remaining variable that changes the end effect for a given test section is the heat transfer
To determine the maximum distance into the heated test section that experimental results
might be affected, the worst-case heat transfer coefficient of 1000 W/m2-K was used for all test
sections. This is half the lowest heat transfer coefficient in the entire data set, and should
represent a worst-case scenario in the unheated section where liquid hydrogen is flowing next to
the tube wall. That is, liquid hydrogen will have a significantly higher heat transfer coefficient
than will vapor hydrogen. Table 3-6 presents the model results and suggests that all test section
data more than 0.8 cm from the heated test section boundaries are adequately unaffected to be
used in the analysis. As a result of this analysis, all 12 points in the 1966 report will be used
since the end points in these runs are far more than 1 cm from the ends. However, points 1 and
12 in the 1961 data are theoretically affected, and the data of wall temperatures strongly supports
this conclusion. Therefore, these 40 points will be excluded from the heat transfer and pressure
drop analyses, leaving 572 points for consideration. All other data in the 1961 report are
predicted to be adequately unaffected and will be used.
Hydrogen States: Parahydrogen and Orthohydrogen
Hydrogen is naturally found as a molecule composed of two atoms of hydrogen, j oined by
a covalent bond. The proton at the nucleus of each atom has a spin associated with it giving rise
to four possible combinations of spin pairs between the two protons of a hydrogen molecule, H2-
Three of these combinations of nuclear spins are symmetric, resulting in orthohydrogen (ortho),
while the fourth combination is antisymmetric, resulting in parahydrogen (para). This two-state
nature or hydrogen is significant for several reasons. The heat of formation released during the
transition from ortho to para, coupled with the unstable nature of ortho at low temperatures, can
cause significant boil-off of stored hydrogen if ortho constitutes a large fraction of the liquid.
Ortho conversion to para is an exothermic process, with the emission of 703 kJ/kg of heat at 20
K, which is significantly more than the latent heat of vaporization of 443 kJ/kg. Secondly, the
thermal properties of specific heats and thermal conductivities of the two forms are known to be
significantly different at cryogenic conditions, causing the need to consider the issue of the
ortho-para makeup of the test fluid throughout the test section.
The relative equilibrium abundance of each form varies with temperature. At room
temperature, the ratio is 3 parts ortho to 1 part para, reflecting the number of spin combinations
available to each form. This state of hydrogen is called normal hydrogen. The ratio changes to a
larger proportion of para as the fluid is refrigerated. At 20.4 K, the ratio is 0.002 parts ortho to
0.998 parts para, at equilibrium. Note that time is needed to allow for equilibration, which can
be hastened in the presence of a catalyst.
There are four processes in which one form of hydrogen can transition into the other;
collisional, spontaneous, adsorption, and radiative. The collisional process can be further
segregated to homogeneous and heterogeneous processes. Through the homogeneous collisional
transition, an ortho molecule acts as a paramagnetic medium through which spin exchange
occurs either with another ortho molecule or a para molecule (Iverson, 2003). The
heterogeneous collisional transition requires a catalyst, such as a tank or pipe wall, that is
propitious to the transition of one form to another. This method involves the interaction between
the magnetic Hield generated by a magnetic material and the magnetic Hield associated with the
nuclear spin of the H2 HUClOUS. The interaction causes a reversal of spin in one of the nuclei,
which effectively changes the form from one to the other. In both of these collisional processes,
the transition from ortho to para is exothermic in the form of increased kinetic energy of the
participating molecules. Natterer et al. (1997) describe a method of catalyzing the transition of
ortho to para by flowing hydrogen through a tube that is charged with charcoal. Without a
catalyst, the conversion from ortho to para liquid hydrogen has a time constant on the order of
180 hours (Scott, 1959). Milenko et al. (1997) measured natural ortho-para conversion rates
within a wide region of hydrogen fluid states, including fiye different liquid temperature states.
Their Eindings indicate a conversion time constant near 12 hours.
The spontaneous transition of ortho to para produces a photon, and is therefore also an
exothermic process. Ehrlich (1991) sites theoretical results showing that the time constant for an
isolated ortho molecule to transition is on the order of 1011 years.
Chemical adsorption of the hydrogen on the metal can lead to conversion of hydrogen.
Ptushinskil (2004) addressed the physics of this process. The adsorption process is composed of
physisorption and chemisorption, which denote different levels of interaction between the
hydrogen molecule and the metallic surface. These two levels are separated by a repulsive
barrier of variable magnitude. As of yet, no theory for the time constant of transition between
the para and ortho states for this process have been found.
The fourth method considered here requires radiation bombardment of the hydrogen. In
this process, H2 mOlecules dissociate due to the bombardment. The subsequent hydrogen atoms
can recombine with each other generating, on average, the equilibrium ratio of para and ortho
forms associated with the system temperature (Kasai, 2003). Since the hydrogen storage tank
used in the NASA tests was surrounded by a radiation shield, this process is not expected to
contribute significantly to the production of ortho.
Iverson (2003) presents a method to quantify the dynamic equilibrium density of para and
ortho in a mixture of liquid hydrogen with collisions and irradiation present. He uses the
following set of equations to quantify the concentration of ortho and para, considering
homogeneous and catalyzed transitions:
= Kponlno -Ko Y1 pon p Copno (3.2)
S= Ko1pnpno -KopIIZ ponp'~ Copno," (3.3)
subj ect to the conservation equation,
no(t)+ n,(t)= NH2 (3.4)
In the above equations, no, n,, and NH2 are the densities of ortho, para, and all H2 mOleCUleS,
respectively. K,o and Kop are COHVersion factors for homogeneous conversions from para to
ortho, and from ortho to para, respectively. C,o and Cop are COnversion constants for catalyzed
conversions in a similar sense. Both the homogeneous and catalyzed conversion constants are
strong functions of system temperature and pressure. Milenko et al. (1997) provides information
about the values of the constants.
Hendricks et al. (1961) analytically quantified the various means of transition from para to
ortho and visa versa and chose to neglect the effects based on the results. For the analyses in the
NASA reports, 100% para was assumed. It is stated, though, that neglecting the presence of
ortho may introduce error into some of the heat balance calculations. An accurate quantification
of the ortho-para makeup was extremely important in the NASA analyses (Hendricks, personal
communication, 2006). While the parahydrogen flowed through the heated test section, there
was also concern about the transition from para to ortho as the fluid was heated.
To test for this possibility, one test section was gold plated and then used. This experiment
is based on the fact that any heterogeneously catalyzed transitions from para to ortho that occur
with a stainless steel test section should be eliminated by the gold plating. Since the transition
from para to ortho is endothermic, a stainless steel tube should show lower wall temperatures
than the gold-plated tube under the same test conditions. However, the opposite effect was
observed, which was attributed to experimental error. Their assessment was that the residence
time of the hydrogen molecules in the test section was not long enough to generate significant
ortho concentration from the para population as it was heated and flowed in the tube to warrant
adjusting the properties from the assumed 100% para makeup. 100% parahydrogen is assumed
in the current analysis.
Inverted annular film boiling of hydrogen is modeled in this analysis as a separated flow of
vapor and liquid. The liquid flows as a homogeneous core through an annulus of homogeneous
vapor. In this geometry, the vapor interfaces with both the wall and the liquid core, while the
liquid interfaces only with the inner boundary of the vapor annulus. All of the heat from the wall
is assumed to be absorbed by the vapor through convection. Radiation of energy to the vapor or
directly to the liquid is assumed, and has been shown, to be negligible. Additionally, momentum
loss through friction at the wall is largely a function of vapor conditions. This approach is
consistent with the experimental observations of Kawaji and Banerjee (1983, 1987). In their IFB
quench front experiments with water flowing upward in a highly heated quartz tube, bubbles
were seldom observed in the liquid core. They concluded that nearly all the vapor generated at
the liquid-vapor interface flowed upward in the vapor film. They also found no evidence that the
liquid column rewetted the tube wall.
Local static pressures, tube wall temperatures, and voltage drops were recorded. This is
enough information for only a three equation model, also known as a homogeneous equation
model (HEM), with mixture mass, momentum, and energy conservation equations. An extensive
literature search has not uncovered data-based models for vapor superheat or vapor slip in the
flow structure of this analysis. It is likely that these profiles will be unique relative to pre-CHF
flows, so that information on vapor superheat and slip from pre-CHF will not apply.
It was desired to obtain void fractions from the hydrogen data. To obtain useful void
fraction data, it was determined that a no-slip condition was not acceptable, since the slip ratio
directly affects the void fraction. In addition, a reasonable value for vapor velocity was desired
to allow for a reasonable estimate of frictional losses. Also, since void fraction and slip are
related to density, it was determined that the vapor superheat needed to be quantified. Without
information regarding superheat, vapor velocity slip, or applicable information regarding void
fraction profiles, theory and assumptions must be applied if more information is to be obtained
from these hydrogen data than what a HEM can provide.
The desired information can be obtained with a one-dimensional, five-equation model,
with separate vapor and liquid mass and energy flows, but with one momentum equation. This
assumes that the local pressure is the same for both fluids, which is commonly accepted.
Completing this model requires closure conditions for two of the following three quantities;
vapor mass-specific energy flow, vapor slip, and liquid mass-specific energy flow. Since wall
temperatures are part of the data set, it was determined that a closure condition for the vapor
energy flow, through quantifying vapor superheat, could be reasonably determined. Neither the
liquid heating nor the vapor slip is well understood. It was determined to model the liquid
energy state. It was determined that modeling the interfacial momentum effects was not
necessary for the obj ectives of this analysis. Including such effects would lead to a two fluid
Nature of Data
Consideration of figure 3-10 of tube inner wall temperatures minus liquid hydrogen
temperatures leads to the expectation that the vast maj ority of data is IFB. The vast maj ority of
data show very large temperature differences between the inner wall and the liquid hydrogen
temperature. These large temperature differences can only be sustained in an IFB flow structure.
The four runs presented in figure 3-11 exhibit at least one point with relatively low temperature
difference. It is likely that these points correspond to pre-CHF conditions, or possibly transition
The trend of the temperature differences for these runs in figure 3-11 supports this
explanation. For example, run 30 has an extremely high mass flux of 3406 kg/m2-Sec and a very
low heat flux of 310 kW/m2. These operating conditions are most likely to delay the onset of
CHF and the transition to IFB, and this is what is indicated by the data. It is not until
approximately 40 cm into the heated section that the temperature difference increases greatly.
Run 29 has a lower mass flux of 2669 kg/m2-Sec and approximately the same heat flux and
would thus theoretically depart from nucleate boiling at a lower elevation than run 30. This is
indeed what the data show, with run 30 temperature difference increasing significantly starting
after the 16 cm point. Run 22 has a similar mass flux as run 30 at 3444 kg/m2-Sec, but a much
higher heat flux of 1128 kW/m2. One would expect this run also to transition to IFB at a lower
elevation than run 30. While the temperature differences for run 22 at low elevations is higher
than run 30's, it appears that the IFB structure is not stable until after the 24 cm point earlier
than the run 30 transition. Since runs 22, 26, 29, and 30 show that pre-CHF conditions exist at
least at some points, and the model generated to analyze these data assumes IFB conditions,
these four runs will be excluded from the analysis.
There are other experimental findings that support this conclusion to omit them. Previous
research with hydrogen indicates the magnitude of wall to bulk superheat that hydrogen will
allow before departing from nucleate boiling. Walters (1960) reported a maximum wall
superheat from his single-tube forced hydrogen flow heat transfer experiments of about 2.8K.
Sherley (1963) experimented with free-convection hydrogen heated by a small flat heating
surface and reported wall superheats as high as 6.1K. Class et al. (1959) experimented with free-
convection hydrogen on various surface conditions, heating surface orientations, and pressures.
For a very thin layer of silicone grease applied to the test section, wall superheats of about 16.7K
were reported. Graham et al. (1965) presented test results from parahydrogen pool boiling that
showed wall superheats of up to 5.6K at a system pressure of 290 kPa before departure from
nucleate boiling. Kozlov and Nozdrin (1992) measured heat fluxes and wall superheats at DNB
during pool boiling of hydrogen for steel, aluminum alloy, and copper at low pressures. They
found that wall superheats at DNB varied significantly between the three metals, as did the wall
superheats during return to nucleate boiling from film boiling when they reduced the heat flux.
At one atmosphere on steel, the wall superheat was on the order of 16 K. All of these studies
support the previously stated assumption that the vast maj ority of data from TN 3095 represent
post critical heat flux conditions.
Carey (1992) states that the variables affecting critical heat flux are tube diameter, system
pressure, and mass flux. The fourth controlling variable depends on whether the bulk flow is
subcooled or saturated. For saturated flow, Carey sites the critical quality, while for subcooled
flow it is the difference between saturation and bulk temperatures. Collier (1981) also lists
length to diameter ratio as an important parameter.
Chun et al. (2000) developed a new theoretical model for predicting CHF for low quality
flows of water and refrigerants in round tubes. Chun states that there is general agreement that
for highly subcooled flow, the liquid sublayer dryout approach performs well, while for low
subcooling the bubble crowding model performs better. No one model works well in all
conditions, though. Chun attempts to improve this situation by proposing that the controlling
factor in CHF is the evaporation of the superheated liquid layer along the tube wall.
Recent research into this issue has been performed by Celata et al. (1994, 1996, 1998,
2001) in Italy. While most of his research is focused on highly subcooled CHF of water, the
general concepts will probably prove relevant to liquid hydrogen. While Carey (1992) lists three
postulated mechanisms for CHF at low quality dryout under a growing bubble, vapor
crowding, and dryout under a vapor slug Celata states that the liquid sublayer dryout theory
predicts the CHF under a wide range of subcooled conditions.
Magnitude of Radiation Heating
Heat is transported from the tube inner wall to the hydrogen primarily through convection.
However, the large temperature differences experienced in the test series raises the concern that
radiative heat transfer from the wall to the vapor and/or liquid hydrogen may be significant.
While the exact analysis of radiation heating is complex, a simplified analysis of the worst-case
scenario will reveal that radiative heating is at least three orders of magnitude less than
Sparrow (1964) presented a thorough theoretical analysis of the effect of radiation heating
from a tube wall to a vapor/liquid flow in film boiling. His work generated a quantitative
criterion by which the relative significance of surface-to-liquid radiation can be determined. A
more recent paper by Liao et al. (2005), which presents an excerpt of his Ph.D. work, addressed
this complicated problem by modeling the liquid core flow as a long inner tube at the center of a
long outer tube. The equation for radiation heating he applied to this geometry is
q" d =(3.5)
1f 1-E r,
The emissivities, s, that will lead to the largest radiative heating are 1 for both hydrogen and
wall. In this equation, r is radius, and o is Stefan-Boltzmann constant. The radiative heat flux
then reduces to
q", =o-(T -T,(3.6)
The highest wall temperature from the data is 560 K, and the fluid temperature is roughly 25 K.
Using these values to represent the upper limit of radiative heating, the magnitude is 5.6 kW/m2.
The lowest heat flux in the data set is 294 kW/m2, SeVeral orders of magnitude larger.
Additionally, this lowest heat flux does not correspond to the highest wall temperature of 560 K
used in this analysis, but instead has a much lower wall temperature of 178 K. In summary,
there is no run in this data set that has a radiative heating contribution of more than 2% of the
total applied heat flux, and in fact is certainly much less than 2%. The impact of ignoring
radiative heating of hydrogen is therefore justified.
Most of the experimental runs have subcooled liquid entering the heated test section. The
amount of subcooling is appreciable, up to 7 K in some runs, and cannot be ignored in the energy
balance. The velocities attained in some of the experiments required that the stagnation
enthalpies of the two fluids be used in the energy balance instead of the static enthalpies. Thus,
the momentum and energy equations are coupled and must be solved simultaneously.
A one-dimensional model of this system was developed to calculate mass, momentum, and
energy balances. It is assumed that the pressure is constant across the flow cross-section, and
while separate velocities of the two phases are determined, the bulk velocity for each phase is
used. Additionally, bulk thermodynamic properties are assumed.
The conservation of mass equation is simply
w = w1 + w,, (3.7)
The liquid momentum equation is
d~p~,Au,) d(PA,;, )
dz dz +r, 2nr,d: gp, A,~,d: (3.8)
where Ti and ri are the vapor liquid interface shear stress and radial location, and A, is the flow
area. The corresponding equation for the vapor phase is
d (puca~, z= d(PAUca )d: r, 2nr~ d: r,xi~d: gpg Ac,,d: (3.9)
The velocity and area terms in these momentum equations can be replaced by use of the
u = Gx(3.10)
G (1 x)
P 10 a)
Ac, = aAe (3.12)
Ac,, = (1 a)A, (3.13)
In these equations, oc is the vapor void fraction.
During the expansion of the derivatives, the vapor density was allowed to vary as a
function of z. Doing so facilitates investigating the effect of vapor superheating and its axial
variation on the pressure profile. The liquid density axial variation was also allowed to vary.
Also during the expansion, certain derivatives were replaced with equivalent expansions that
used terms more amenable to the analysis.
In a one-dimensional analysis such as this, these separate momentum equations are
combined by equating the interfacial interactions of the two phases. The result is seen in
a S -2( -x)2x Ba (1 x)2 x x S,,S
dP I8z LPt(1 ) p,,a dx Pt (1- a)2 p,,a p~a T P,
fi ,,,z1 X2 T z1+ z, + giP r 4j ,(1-a)+ p,,a
x-G r S p,apj (1 -x)z Sp P )
p a SPTP1 \1-a P
This equation is similar to that commonly presented in two-phase flow textbooks, but with
Jacobian expansions useful for this analysis. The following relation for the wall shear stress was
As previously stated, the velocities attained in some experiments were high enough that
they should be included in the energy balance. Radiation heating of the liquid is ignored based
on the previous analysis of liquid heating by radiation. As a result, conservation of energy is
modeled as follows:
L!W~1I Gx / +1 G (1 x)`' 3.6
Q = wxh, w(- ) t* < x 3.6
2 p, a 2 p 1-a
where h is the enthalpy.
In the application of this equation, the total energy flow rate of the flow is determined to be
the total energy of the flow at the first point of measurements, point 1, plus the cumulative
energy added through heating:
()= w h+ u + q"A, (3.17)
As is the cumulative tube inner surface area up to a particular point of calculation.
There are three types of entrance lengths considered here; hydrodynamic, thermal in the
fluid, and thermal in the tube metal. Although all test sections included straight entrance
approach sections approximately 12.5 cm and 30.5 cm long in the 1961 and 1966 data,
respectively, to develop the velocity and thermal profiles, this concern is obviated by the nearly
instantaneous and violent change in flow structure from single-phase liquid to IFB. Stated
another way, the history of the flow up to the start of heating is not important. Instead of
modeling liquid velocity and temperature profiles across the radius of the tube and their effects
on heat transfer and pressure drop, these processes are controlled by the conditions in the vapor,
the inception of which occurs at the heated section inlet, and in which the radial dimension is
The developing hydrodynamic and thermal profiles in the vapor from the test section inlet
onwards must still be considered. Hsu and Westwater (1960) used law-of-the-wall theory to
determine that the vapor in the annulus transitions from laminar to turbulent at a Re = 100.
Some rather arbitrary assumptions were employed in their theoretical analysis. Somewhat
marginally applicable computations from Rohsenow et al. (1956) for condensation on a vertical
plate were used to justify this transition Reynolds number. Regardless, this transition Reynolds
number appears to be commonly quoted and used to determine transition from laminar to
turbulent flow of the vapor in IFB. Note that for typical values of vapor density and viscosity,
and for typical velocities at the test section entrance, the vapor annulus dimension that produces
a Re of 100 is 0.001 cm an extremely small thickness. This film thickness is achieved at a void
fraction for the smallest tubes in the NASA data set, which will give the largest required void
fraction, of 0.008. From this, it is reasonable to assume that the vapor is always turbulent.
Additionally, it is hard to conceive of the vapor flowing in a laminar fashion after its violent
generation at the heated test section entrance.
As previously discussed, in the tube metal at the boundary between the heated test section
and the entrance piping, there will be a significant axial gradient in metal temperature. This will
lead to axial heat conduction, which in turn will affect the local heat flux and temperature.
Instead of the approximately constant heat flux established within the tube far from the
boundaries of the heated test section, the local heat flux can be significantly reduced. Measured
wall temperatures from the 0.795 cm diameter tube support this conclusion. It is important to
note that, while there is axial heat transfer in the metal, at any particular station near the inlet, all
of the energy that is calculated to be transmitted to the flow up to that point will indeed be
transmitted to the flow. Thus, the calculated total energy input to the flow up to a given point
will not be in error. At the test section exit, this is not the case. Heat flows up and out of the test
section at the exit. Thus, the flow will not receive all of the heat input until some point after the
heated section exit.
The first point at which enough information is given to determine the thermodynamic state
of the flow is the first point listed in the tables of measurements for each run. For the 1961 and
1966 reports, this point is at 0.14 cm and 6.35 cm up from the test section inlet, respectively. If
the flow at this point is subcooled, then the published quality is zero, and the published
temperature and pressure is used to determine the thermodynamic state. If a positive quality is
listed, then the published pressure and quality is used.
Quality and void fraction are determined from the momentum and energy balances. The
balances calculate changes in static pressure and total energy. Therefore, the quality and void
fraction of the first point in each run must be determined in a method other than using these
balances. It was determined to initialize the quality and void to zero at the inlet. It was assumed
that quality and void increased monotonically at each successive point.
Implementing this boundary condition required knowledge of the thermodynamic state of
the fluid at the test section inlet. This information is not given directly. However, the energy
state of the fluid at the inlet can be found by subtracting the energy added from the inlet to the
first measured point from the energy of the flow at the first point. Assuming the energy
associated with the local velocity to be negligible relative to the enthalpy of the flow, this energy
level per unit mass is used as the bulk enthalpy of the flow. The pressure at the inlet is
determined using the same technique the authors used to determine the pressure profie Sit a
smooth curve through the measured points. The cubic least squares fit of the pressure profiles, as
previously described, were use to extrapolate backwards to calculate the test section inlet
pressure. Thus, pressure and bulk enthalpy are determined for the inlet. From this, the
thermodynamic and kinetic state of the liquid and vapor is determined. The inlet was then
defined to be point one for each run, and the number of points used in the analysis of each run
increased from 12 to 13.
The momentum equation requires positive qualities. However, as stated previously, many
runs had subcooled inlet conditions, and in fact remained subcooled from an equilibrium sense
for a number of points. Therefore, a method to establish a positive quality was necessary. The
literature search produced no model for true quality.
Hammouda (1996) presented a notional graph of the variation of true mass quality as a
function of length in IFB. Collier (1994) presents a similar graph on page 295. Hammouda' s
graph is not based on measurements, but instead from his interpretation of conditions based on
his observations of IFB. The slope of mass quality in IFB is positive at negative equilibrium
qualities. Near where equilibrium quality equals zero, the slope of mass quality with length
increases. At some low value of quality, mass and equilibrium qualities are equal, after which
equilibrium quality is greater. At an equilibrium quality of one, the mass quality is less than one
due to vapor superheating. This model encompasses the following three concepts in IFB: the
subcooled liquid experiences some sensible heating; vapor is present and accumulates while the
bulk flow is subcooled; due to vapor superheating, the flow will not be entirely vaporized when
the equilibrium quality equals one.
To complete the set of equations, the level of bulk vapor superheat, the amount of liquid
sensible heating, and the nature of the wall friction must be determined.
To quantify vapor superheat, several concepts were tried, including theory presented by
Burmeister (1993). He presents a theoretical derivation for the mixing cup temperature.
Following are the applicable energy equation and boundary conditions used:
aT 1 (vrq, )
8: r dr
subj ect to
dT (r = r z) qw
8 T ( = 0, z)
= 0 (3.20)
Following are the assumptions used in his development.
1. constant wall heat flux
2. circular duct
3. flow velocity and temperature profiles are fully developed
5. Pr is constant and 1
6. Prt = 1.0
7. Law of the wall applies, with sublayer, buffer, and core zones
8. u/UCL = 7/0)1 7, and radial temperature profile has the same form
The results of his analysis give a mixing cup temperature of the following form:
T,,, = 5 (T Tw) +T, (3.21)
The centerline temperature in this equation is the liquid temperature.
It seemed logical that the temperature profie in the vapor could be modeled as a turbulent
flow. Most of the assumptions listed above are well satisfied by these NASA data, and
arguments can be made for the remaining assumptions. Use of this model with the 5/6th
coefficient caused numerous energy balance errors, primarily in the 1961 data set. Various
coefficient values were tried between the theoretical 5/6th and the commonly-used '/. Energy
balance errors were minimized with the smallest coefficient of '/. Therefore, it was determined
to proceed with this value. This coefficient value is consistent with the analyses of Takenaka
(1989) in his IFB studies. Nijhawan et al. (1980) performed experiments in which they
measured vapor superheats in post-CHF flowing water. They observed significant superheating
of the vapor. Interestingly, their data strongly support the use in this effort of '/ for the vapor
Liquid energy flow
To complete the theory for a Hyve-equation model, an assumption must be made regarding
the energy state of the liquid. Theory regarding heat transfer to the liquid flow can be found in
six-equation models, also called two-fluid models. Hedayatpour et al. (1993), in their two-fluid
model of a vertical line cooling with liquid nitrogen, used theory for water droplet heating in
superheated steam from Lee and Ryley (1968). The Nusselt number is modeled as
Nu = 2 + 0. 74 Reo 5 Pr 0 33 (3.22)
where the Reynolds number is evaluated at droplet conditions, and the Prandtl number is
evaluated at film conditions. This model was used in a flow geometry identical to that used in
this dissertation liquid core flowing homogeneously inside an annulus of vapor.
Hammouda (1997) observed that the heat transfer coefficients for the wall-to-vapor and the
vapor-to-liquid can both be modeled as functions of Reynolds number to a power and Prandtl
number to a different power. With some assumptions, he concluded that the ratio of vapor-to-
interface and wall-to-vapor heat fluxes were controlled as follows:
q"- T~ Ts
c": Tw -T
He gave no experimental justification for this model except that he noted predictions from his
two-fluid model provide better prediction accuracy than other IAFB prediction methods he
The assumption used in this dissertation is that the liquid experiences no sensible heating.
It remains at its inlet temperature throughout the heated tube unless the local pressure drops to
the saturation pressure for the liquid temperature. From this point onward, the liquid temperature
assumes the saturation temperature at the local pressure.
Rationale for this assumption comes from the fact that vapor is definitely present during
IFB, even for subcooled flows. Therefore, the liquid certainly does not absorb and evenly
distribute 100% of the energy from the tube wall. That is, the fluid does not increase in
temperature to saturation before it starts to generate vapor. This observation easily extends into
the saturated condition in which it is logical to assume that a saturated liquid also does not
absorb 100% of and evenly distribute the energy input from the wall. The true nature of the
liquid heating almost certainly lies between the extremes of no sensible heating and
Using some assumptions, the exact theoretical time-dependent liquid temperature profie as
it flows through the core of the tube can be solved. The liquid core is modeled as an infinitely
long rod of constant radius R having a uniform initial temperature T1 and instantaneously
subj ected to a uniform temperature bath at temperature Tsat. It is assumed that the bath
temperature is the saturation temperature of the fluid at the local pressure. That is, any liquid
that rises above the saturation temperature evaporates and leaves the liquid core and does not
heat the remaining liquid. Only liquid that is at the saturation temperature or lower remains to
conduct heat from the liquid/vapor interface inwards. This model also assumes that the liquid is
at a uniform temperature across its radius at the time heat is applied (the test section inlet), that
liquid radial velocity gradients are unimportant to heat transfer, and properties are constant. That
is, heat transfer in the liquid can be modeled by conduction alone.
The mathematical model that captures the physics of this problem is
aT k 1d 8 T\
= r -(3.24)
at pC r dr dr
subject to the following boundary conditions:
T(R, t) = T, (3.25)
= 0 (3.26)
and the initial condition:
T(r,0)= T, (3.27)
The time-dependent solution of this problem is (Arpaci, 1966)
T~r~)= T +2(, -T,) *J f )(3.28)
n (A,, R) Jz(,, R)
where hnR are the characteristic roots of the Bessel function of the first kind of order zero. The
solution of interest from this model is the average temperature rise for a typical differential liquid
volume that passes through the heated test section. Following are the values that will be used for
* R = 2E-3 m
* k = 0.1 W/m-K
* p = 60 kg / m3
* C = 2E4 J / kg-K
* Ts = 28 K
* T1= 25 K
These values correspond to a liquid subcooling of 3 K, which is a typical value. Also, a typical
differential fluid volume residence time in the test section of 1/30th Of a Second will be used.
Figure 3-13 presents the results that strongly support the assumption to ignore sensible heating of
the liquid. This is the assumption that will be applied in the model.
The frictional losses are modeled with a Blasius-type relation for the friction factor and a
two-phase friction multiplier developed by Rogers (1968) at Los Alamos National Lab. His
model was developed for friction modeled as only the liquid component of the two-phase flow
flowing alone. Thus, the friction factor is
f =0.079Re 0 25 (3.29)
Rogers' model was developed specifically for two-phase internal flow hydrogen. Although his
model is largely theoretical with some data validation, it is applicable to the entire two-phase
hydrogen pressure range, and is presented in closed form as follows:
2x 1( 8-,c 0 8187 0.1324(12.759 P)+ 003966(1.759P'] P)3
f5 x 1 PE (330
where pressure is in atmospheres, and E is
E = 1.896x 2.646x2 + 1.695 x3 (3.31)
During the implementation of this theory, two observations directed the Einal form of the
algorithm. First, implementing the theory requires an iterative scheme with discretized quality
and void fractions. Each combination of quality and void fraction will result in errors in
predicted pressure drop and energy flow relative to measurement. Acceptable levels of error
must be defined, which results in a quality-void fraction pair domain of solutions from which a
Einal pair must be selected. Second, it was found that there are some points for which this model
will not simultaneously satisfy both momentum and energy conservation. This is due mostly to
the inaccuracies of the model, and probably to a lesser extent due to inaccuracies in experimental
measurements. For most points, momentum and energy conservation are satisfied with
negligible errors associated with the necessity of discretization.
It is for these two reasons that 'smart' iteration techniques failed. Several other methods of
Ending the correct quality-void fraction solution were implemented that relied on reducing the
error in energy and momentum by determining the correct direction to change each value.
However, these iteration methods were found to be inadequate due to the nature of the equations
in the problem and due to the fact that, in some cases, the solution of least error is greater than
the error limits for most other points.
Performing a 'dumb' progression of quality/void fraction pairs, while not conservative of
CPU time, was found adequate. Figure 3-14 presents the flow chart of the algorithm. Note that
the thermodynamic state of the vapor and liquid are known since liquid temperature, vapor
temperature (through the superheat equation), and local pressure are known.
The error limits place on calculated momentum and energy changes are 2% of
measurement. All quality-void fraction pairs that agreed with the measured pressure loss to
within 2% were saved for processing in the energy balance. This preliminary solution set was
then input to the energy balance. The solution domain is constrained by noting the contribution
of velocity to the total energy flow. It is significantly less than that of enthalpy even for the high
velocity flows. Therefore, the energy balance is a very weak function of void fraction and a very
strong function of quality. Thus, the quality range is always reduced to one or a few discretized
values, but with a range of void fractions that satisfy the momentum equation within the error
It is logical to use the liquid and vapor velocities to discriminate between the remaining
solutions. Various methods were tried. One method required the vapor velocity to be greater
than the liquid velocity at all points, but this did not work best for runs near the critical pressure.
A slip of less than one appears to satisfy these runs best. Another constraint that led to problems
for high pressure runs was to require the vapor velocity to increase monotonically up the tube. It
was finally determined to select the minimum vapor velocity from the set of solutions that
satisfied the energy balance within the specified error limit. This constraint eliminated extremely
high vapor velocities, some well over the sonic velocity, while giving reasonable results for high
To address points for which momentum and energy conservation can not simultaneously
be satisfied, it was determined to equally increase the accepted momentum and energy errors
until a solution was obtained for both. Note that increasing the acceptable range of errors on
momentum consistently decreases the calculated errors in energy balance, so this method
identified the lowest level of error for both quantities while giving preference to neither.
Table 3-1. Table of experimental conditions for the NASA data se .
set run G Pin q" dp Subcool T D, inner Wall thickness Material
kg/m2-s kPa kW/m2 kPa K cm cm
1 1.1146 327 759 1193 27 -0.1 1.288 0.025 Inconel X
2 2.1152 643 969 948 22 -3.1 1.288 0.025 Inconel X
3 3.1143 329 743 735 15 -0.1 1.288 0.025 Inconel X
4 4.1151 488 1023 768 13 -2.6 1.288 0.025 Inconel X
5 5.115 662 1045 752 12 -3.8 1.288 0.025 Inconel X
6 6.1142 630 733 719 20 -1.3 1.288 0.025 Inconel X
7 1.1246 873 1075 1324 45 -4.1 1.113 0.081 Inconel
8 2.1247 536 1103 1308 76 -2.7 1.113 0.081 Inconel
9 3.1248 895 889 1242 41 -3.3 1.113 0.081 Inconel
10 4.1251 531 868 817 20 -2.4 1.113 0.081 Inconel
11 1.542 1237 616 1357 211 -2.9 0.851 0.051 304 Stainless steel
12 2.541 1119 861 1324 86 -4.8 0.851 0.051 304 Stainless steel
13 3.539 892 984 703 32 -6.6 0.851 0.051 304 Stainless steel
14 4.538 906 982 425 12 -7.2 0.851 0.051 304 Stainless steel
15 5.2 1553 1251 1766 102 -3.8 0.851 0.051 304 Stainless steel
16 6.201 1286 1112 1733 113 -2.7 0.851 0.051 304 Stainless steel
17 7.54 1178 759 2093 250 -2.6 0.851 0.051 304 Stainless steel
18 8.203 1129 1221 1733 93 -3.3 0.851 0.051 304 Stainless steel
19 9.204 1121 812 1635 148 -0.2 0.851 0.051 304 Stainless steel
20 10.535 945 685 1798 201 -1 0.851 0.051 304 Stainless steel
21 11.536 932 746 2076 223 -0.9 0.851 0.051 304 Stainless steel
22 1.568 3444 1265 1128 274 -6.6 0.478 0.079 Inconel
23 2.577 1965 1141 1112 232 -4.9 0.478 0.079 Inconel
24 3.559 2466 1059 1112 272 -6.2 0.478 0.079 Inconel
25 4.558 2446 1072 981 272 -6.2 0.478 0.079 Inconel
26 5.562 3186 856 670 160 -6.2 0.478 0.079 Inconel
27 6.56 2383 823 654 201 -4.3 0.478 0.079 Inconel
28 7.561 2735 817 670 146 -5.3 0.478 0.079 Inconel
29 8.564 2669 594 294 104 -3.8 0.478 0.079 Inconel
30 9.565 3406 613 310 124 -4.6 0.478 0.079 Inconel
31 10.563 2165 561 294 109 -2.7 0.478 0.079 Inconel
32 1.1802 1617 310 376 117 0 0.795 0.079 Inconel
33 2.1803 1242 279 376 68 0 0.795 0.079 Inconel
34 3.1804 849 228 376 51 0 0.795 0.079 Inconel
35 4.1805 575 188 376 38 0 0.795 0.079 Inconel
36 5.1806 1653 359 621 79 -0.8 0.795 0.079 Inconel
37 6.1807 1123 311 637 87 0 0.795 0.079 Inconel
38 7.1808 804 259 637 71 0 0.795 0.079 Inconel
39 8.2001 1553 399 981 107 -2.4 0.795 0.079 Inconel
40 9.2002 1242 359 981 103 -1.6 0.795 0.079 Inconel
41 10.2 858 303 997 114 -0.2 0.795 0.079 Inconel
42 10.2 721 257 997 101 0 0.795 0.079 Inconel
43 12.22 1379 448 1144 132 -1.3 0.795 0.079 Inconel
44 13.201 1626 457 1357 136 -3.2 0.795 0.079 Inconel
45 14.201 1206 399 1373 141 -2.1 0.795 0.079 Inconel
46 15.201 849 339 1373 146 -0.7 0.795 0.079 Inconel
47 16.201 712 286 1373 126 0 0.795 0.079 Inconel
48 17.22 1297 490 1520 149 -2 0.795 0.079 Inconel
49 18.22 922 408 1520 153 -0.8 0.795 0.079 Inconel
50 19.22 621 335 1520 132 0 0.795 0.079 Inconel
51 20.201 1516 498 1651 165 -3.2 0.795 0.079 Inconel
Table 3-2. Comparison of Core et al. (1959) heat transfer coefficients with Hendricks et al.
0 q" Win I ~average h dteec
HO~fCO run yr-m vwr ka "= R %
nmcre et al. 1.3 ;u 4/4 921 1.4
MenarCKS et al. A b/ 3/6 1W 9.1U1
nmcre et al. 4.4 bb2 YYI ZU 2I
MenarlCKS et al. 4 / / 1U 3Ub-
MenarCKS et al. 4 / / 1U 31U-3
noemn~ et al. 6.3 MU1 6;tr 4.3
MenarCKS et al. Is b/b 3/6 14 1.U 2
n mcre et al. I. 494 1bb@ 4b .1/
nenarCKS et al. DU 1 ;4; lbZ 3b .1
n mcreet al. b. bb li dbbb
MenarCKS et al. DU 1 lb U 13 3954
MenarCKS et al. 42 12 I; W/3dU-
Table 3-3. Comparison of average heat transfer coefficients for similar runs in the Wright and
Walters (1959) data set and TN 765.
G q" Pin average h difference
Source run kg/(m2-sec) kW/m2 kPa kWI(m2-K) %
Wright & Walters (1959) 15 908 260 250 3.83
Hendricks et al. (1961) 34 849 376 228 3.08 24.3
Wright & Walters (1959) 23 522 385 167 3.11
Hendricks et al. (1961) 35 575 376 188 3.00 3.7
Wright & Walters (1959) 35 427 390 179 2.47
Hendricks et al. (1961) 35 575 376 188 3.00 -17.7
Table 3-4. Summary of test conditions for maj or forced convection internal tube flow boiling
Mass flux, [kg/m2-Sec] Heat flux, [kW/m ]1 System pressure, [kPa]
Source min max min max min max
Core et al. (1959) 322 10271 16 98121 193 1469
Lewis et al. (1962) 4 231 11 1261 207 510
Wright & Walters (1959) 410 11721 10 3901 138 276
Hendricks et al. (1961) 575 16531 376 16511 188 498
Hendricks et al. (1966) 327 34441 294 20931 594 1265
thermal cond. heat xfer coeff. wall thickness heat flux 95% length
set # W/m-K kW/m2-K 1E-4m MW/m2 ticknesse
5 15 2 5 19
10 45 2 5 1 15
15 15 1 5 1 12.8
20 15 4 5 16
25 15 2 2 1 13.8
30 15 2 8 17
35 15 2 5 0.3 8.3
40 15 2 5 2 8.5
Table 3-6. Distance into tube wall from start of heating at which tube metal temperatures are
reduced by at least 5% from the nominal level.
Inner DIT Thermal Heat xfer 95%
Runs Diameter Thickness ratio Material cond. Heat Flux Coeff. distance
[cm] [cm] [W/m2-K] [MW~m2] [kW/m-K] [cm]
1- 1.288 0.025 51.5 Inconel X 13 1 1 0.41
7-1C 1.113 0.081 13.7 Inconel 13 1 1 0.73
11-21 0.851 0.051 16.7 304 SS 20 1 1 0.71
22-31 0.478 0.079 6.1 Inconel 13 1 1 0.69
32-51 0.795 0. 079 10.1 Inconel 13 1 1 0.71
Result of parametric sensitivity study of end axial heat conduction.
Q Statie-ptrersse $p
DP States difeCNtrenia
pressure transtneerQ I r
RT eanlatance themsemter s
VP Vauu pm
traller, 1000-00L Dever
g~sase Illqu~li rt ll uply Jir
la Pree- home
Figure 3-1. NASA TN 765 experimental setup.
Vent line to
cotrlf valve -1
Vacuum pu mp exhaust
to roof stack
Entrance venturi '
Figure 3-2. NASA TN-3095 experimental setup.
p Presure tap
T Pulid tapp.
5'F Buae theicrmocupl
Figure 3-3. NASA TN-765 data test section.
Exhaust flow-metering orifice
enmitd for clarity
P Pressure tap
AP Differntial pressure tap
Figure 3-4. NASA TN-3095 test section
T~uomg 0.062 e.d. Irf 0.038 1. d.;
transition to 118 o. 1. pressure tubing
tube; typically 0.032a.d. by0.0195 Ld.J SIlver solder
stalic pressure port
(a) Pressure taps.
Ceramic cement an
tube for insulation
wire 130 gage I0.010~ diam)
Junction is flattened,
ontou red to tube
curvature and spot We
bld Welded thermocouple junctions.
f" Ceramic cement over
function as filler and bond
Glass tape as lining '
to protect Junchan/
from clamp -/
\-Thermocouple lunction flattened
pressed and contoured to tube
Figure 3-5. TN 3095 instrumentation.
I ~U I
Figure 3-7. Radial metal temperature profiles as a function of metal thermal conductivity.
0.0 .0 .1 15 00 0 00 5 00 0 0.3 .4 .4 .5
AXIAL POSITON. m
o.000o 0 DOS
0.035 D 040 0.045 0.05D
Figure 3-6. Nodal distribution and heat generation distribution used to model end effects at tube
15, INNUER WALL
45, INNI\ER WALL
15, OUTER WALL
45, OUTER WALL
0.2, INNIER WAL
0.8, INNIER WAL
0.2, OUTER WA
0.8, OUTER WA
A |-L POS T ON, THICKNIESSES
Figure 3-8. Radial metal temperature profiles as a function of metal thickness.
5 k=15 kW/m
|-L POS T ON, THICKNUESSES
Figure 3-9. Effect of specified parameters on tube end wall axial heat transfer.
I 'I I
r I ~B
I 3 t
a I L 1 3 1
161 D I
"I i I I
Figure 3-10. Difference in wall to liquid temperature for all data considered.
22 u=3444,q=1128 26 G=3186,q=670
G [kg m2-sec]
q [kW m2]
E E ATION
Figure 3-11. Wall to liquid hydrogen temperature differences for four runs that show at least one
very low difference. These small differences are associated with pre-CHF conditions.
Theoretical Radial Liquid Core Temperature Profile at Heated Test Section Exit, Based
on Typical Experimental Values
Core radius = 0.2 cm
k = 0.1 W /m-K
-p = 60 kg / m3
C = 2E4 J / kg-K
Tsat = 28 K
T, = 25 K
SLiquid residence time = 0.0333 seconds
0 0.05 0.1
Radial Position, cm
Figure 3-12. Theoretical liquid core temperature profile at the exit of the heated test section, to
support the assumption to ignore sensible heating of the liquid.
Figure 3-13. Flow diagram for momentum and energy analysis of data.
ANALYSIS AND VALIDATION OF MOMENTUM MODEL RESULTS
The NASA data set comprises 51 steady state runs in which there are 13 data points each.
The first point is at the heated test section inlet. For runs 1 31, the 13th point is 6.3 cm before
the heated test section exit. For runs 32 51, the 13th point is at the heated test section exit. The
runs fall naturally into five groups based on inner diameter. Table 4-1 lists the tube dimensions,
the run numbers associated with each tube, and a reference number that will be used for
convenience in later analyses.
It was determined through various means that the data set needed to be refined. Following
is a description of the approach to this process.
The points that are affected by inlet and end conditions, and any calculations that include
these affected points, should be excluded from analyses. For runs 1-31, point 1 at the test section
inlet falls into this category. Only results between point 2 at 6 cm and point 13 at 55 cm will be
considered. For runs 32-51, points 1, 2, and 13 at the inlet, 0. 1 cm, and at the test section exit
will be excluded. Only results between point 3 at 1.6 cm and point 12 at 29 cm will be
Validity of some data is questioned. The basis for questioning these points lies in apparent
discontinuities between adjacent values. Figure 4-1 presents several examples. Run 42 point 8
at 19 cm shows a rise in wall temperature of 40 50 K relative to adj acent wall temperatures.
This magnitude of temperature rise and fall over a 7 cm length, and the fact that the event is
exceptional in these data, begs an explanation. A similar effect is evident in run 32 at 27 cm. It
may be that a unique flow structure occurs for a short length in these runs. The computer model
is robust enough to accommodate many, but not all, of these changes and solve for the
momentum and energy balances within the specified limits.
For some points in which the wall temperature increases drastically from the previous
point, the model cannot satisfy the energy balance. This is because of the assumption that all of
the vapor at a point is at the calculated vapor temperature, which is a function of the wall
temperature. If there is a large increase in wall temperature, then the increase in mean vapor
enthalpy may require a larger energy addition than the energy added through heating from the
previous point, even with zero change in quality. To satisfy these points, the quality would have
to be reduced, which is assumed to not be possible in the model. This is why these points of high
increase in wall temperature are consistently associated with negative energy addition errors -
the measured added energy can not attain the increase is vapor energy.
Tube three exhibits a consistent decrease in wall temperature at the 34 cm location. Figure
4-2 presents the wall temperatures for all 11 runs on tube three. This is interpreted as a bias in
the measurement. Therefore, in making calculations using the wall temperature, these erroneous
experimental values have been replaced by a linear interpolation between adj acent points. The
only other wall temperature point that was deemed obviously out-of-family was run 42, at the 19
cm elevation. This point also was replaced by a linear interpolation between adj acent points.
While other points in the data set showed erratic trending, it was usually uncertain which points
should be modified. A common characteristic is for adjacent points to trend oppositely, e.g., one
low and the next high. Which point was biased was usually not determined. Therefore, no
modifications were made.
Also evident in figure 4-1 are the end effects in which heat is conducted axially within the
tube metal at the inlet and exit of the test section, as discussed in chapter 3. The steep gradient in
wall temperatures for most runs between points 2 and 3 at 0.14 cm and 1.5 cm at the inlet and
points 12 and 13 at 29 cm and 30.5 cm prove the end effect. What is of particular interest are the
several runs, 32 and 36 in figure 4-1, in which the inlet and exit temperatures are actually higher
than their adj acent measured temperatures inward from the ends. This indicates that there is an
end effect other than axial heat conduction influencing measured wall temperatures. This can be
explained by considering that the collars brazed onto the test section ends to apply a voltage will
not distribute the current absolutely evenly across the tube metal radius. The current flow will
distribute itself across the thickness of the metal over a finite distance, and will be concentrated
near the brazed collar at the ends. Therefore, the current density will be higher at the tube outer
wall where the collar is brazed and will therefore generate more heat towards the outer part of
the wall, where the thermocouple is attached.
The inlet and exit wall temperatures are lower than their adj acent wall temperatures for the
runs that have high wall temperatures, as runs 40, 42, 47, and 50 indicate, in spite of the
concentration of current near these end thermocouples. Comparing these trends with the rising
wall temperatures at the ends observed in runs 32 and 36 exemplify the relative impact of the
independent effects of axial heat conduction and current concentration. For runs with low wall
temperatures, the temperature rise due to current concentration is greater than the temperature
decrease due to axial heat conduction, with the net effect that the measured wall temperature
rises. The opposite net effect is evident in the high wall temperature runs. That is, axial heat
conduction has a greater effect on measured temperature than does current concentration.
The pressure data exhibited uneven trending, to varying degrees, in all runs. This
unevenness can present problems for a modeling algorithm using the pressure data to solve for
other flow conditions. Therefore, a smooth regression line was generated for each run to
represent the axial pressure profile. It was found that a third order least squares fit modeled all
runs very well, with correlation coefficients very near unity for most runs. Table (4-2) presents
these correlation coefficients, standard deviations, and normalized (to pressure drop across
length of tube) standard deviation.
In discussing the model results, momentum and energy balances are discussed in terms of
normalized values for each point. For example, the 'fractional pressure drop error' between two
points refers to the difference in the calculated pressure drop between the previous point and the
current point and the measured difference, divided by the measured difference. In the same way,
the 'fractional energy add error' refers to the difference in the calculated energy addition
between the previous point and the current point and the measured addition, divided by the
While run eight momentum and energy balances are satisfied, its void fraction is clearly
impossible. This is because, while the overall pressure drop of 76 kPa in run eight is not high,
the pressure gradient in the lower portion of the tube is extremely high, as shown in figure 4-3,
comparing runs seven and eight. These runs have similar system pressures and heat fluxes,
while their mass fluxes are 873 kg/m2-S in Tun SeVen and 536 kg/m2-S in run eight. For most
elevation intervals, there is a low and high void fraction solution for the associated pressure drop.
The pressure drop between the inlet and the 6.35 cm elevation in run eight does not allow for a
high void fraction solution given the energy addition and vapor temperature, leaving only the
low void fraction solution as an option. This is the only run in the data set in which this problem
arises. It is noteworthy that this run has always been unique and presented difficulties regardless
of the various modeling methods attempted. As will be shown in the results, run eight pressure
profile is significantly different from the predicted pressure profiles from this dissertation model
and the homogeneous equilibrium model, which happen to trend very closely with each other.
As a result of the obviously erroneous void profile, run eight will be excluded from further
analy si s.
Figure 4-4 shows that the fractional pressure drop and energy addition errors for run 14 are
numerous and relatively large. The algorithm cannot achieve such a low pressure drop for most
points given the vapor temperature and energy addition. Run 13 is very similar to 14 in mass
flow and system pressure, but with a heat flux of 703 kW/m2 VeTSes 425 kW/m2 for run 14. The
lower heat flux in run 14 will certainly cause lower pressure loss, but the overall pressure drop -
12 kPa verses 32 kPa is a little more than 1/3rd that of run 13, which seems to be a great
reduction given the relatively modest decrease in heat flux. The high system pressure of about
75% of the critical pressure should also mitigate the difference in pressure loss between the two
runs, since the difference is liquid and vapor density is not too great, and should thus be less
sensitive to the pressure loss associated with vaporization. Finally, the profile of the pressure
appears to have two inflection points, which, though possible, indicates a very complicated flow
structure. This reverse s-shape appears in other problem runs.
It is interesting to note that run 14 is the only run that is subcooled from an equilibrium
standpoint throughout the length of the heated test section. While the model appears to handle
other subcooled conditions adequately, the highly subcooled nature of run 14 may present a
special problem. As a result, run 14 will be excluded from the correlation process.
Tube 4 runs, 22-31, have by far the highest mass fluxes in the database. This is likely the
reason that four of these runs (22, 26, 29, and 30) exhibit low wall temperatures in the lower
portions of the test section. As discussed in chapter 3 and shown in figure 3-1 1, these runs are
associated with pre-CHF conditions and therefore will be excluded from consideration. Two
other runs, 28 and 31, produce poor energy balances. Thus, four runs (23, 24, 25, and 27) remain
for further consideration. As will be discussed later, these remaining four runs will be excluded
from the correlation process since the nature of their test conditions and resulting slips are
removed from the general body of data.
Three other runs (32, 36, 44) will also be excluded later, based on their high velocity slip
profiles that are inconsistent with all other slip profiles. This will be addressed later. While
there are occasional momentum and energy balance errors in other runs, it has been determined
that useful information can be obtained from them. Therefore, all other runs will be considered
further. Thus, after excluding pre-IFB runs (22, 26, 29, and 30), bad momentum and energy
balance runs (14, 28, and 3 1), and run eight with a bad void profile, there remains a total of 43
runs for further consideration.
It is determined from the above discussion that these runs that cannot be satisfactorily
processed by the model should be excluded from the analysis of results. This should not be
interpreted as a condemnation of the data from these runs, but rather an observation that the
model is not capable of resolving the data, and results using those data will be misleading.
The effect of vapor superheat on the energy balance was investigated by modifying the
coefficient to the parenthesized term in Burmeister' s (1993) model, presented in equation 3-21.
Figure 4-5 presents the results using run 39, in which the term 'C' in the legend refers to this
coefficient. Coefficient values of 5/6th (Burmeister' s theoretical result), 2/3rd, and V/2 (TOSulting in
the commonly used mean film temperature) were tried. Increasing the vapor superheat in
general improves the overall energy balance. The decision to model the vapor superheat as the
mean temperature of the wall and saturation temperature was based upon this analysis. As
discussed in chapter 3, using the mean film temperature is consistent with the experimental
findings of Nijhawan et al. (1980).
As figure 4-6 shows, the processing of data from these runs does not generate perfect
results. Momentum balance within 10% of measurement is achieved in all but 10 instances, and
to with 5% of measurement for all but 20 points. Energy balance is achieved within 10% of
measurement in all but 11 cases, and to within 5% of measurement for all but 44 cases. The
great maj ority of incremental pressure drops and energy additions for each run are well modeled.
Momentum and energy balances that fall outside the targeted 2% variance from measurement are
typically caused by steep changes in wall temperatures. The reason this causes problems lies in
the assumption that all of the vapor is at the mean temperature. If the wall temperature increases
markedly, then so too does the mean vapor enthalpy. The measured energy input to the flow is
less than the energy increase determined from the local pressure and mean vapor temperature.
The algorithm selects the energy solution closest to measurement from the quality/void fraction
domain generated by the momentum balance, but the energy balance error is still larger than the
targeted accuracy in these few cases. In these cases, the calculated quality does not change from
one point to the next.
Figure 4-7 presents the resulting void profiles from the model for the 43 runs. In general,
the void profiles rise steeply but smoothly. Where discontinuities occur, in general there are
steep increases in wall temperatures. In this figure, the four lowest void profiles correspond to
runs 23-25, and 27, all on tube 4. These are by far the highest remaining mass flux runs in the
culled data set.
Figure 4-8 presents the resulting velocity slip ratios from the model for the culled data set.
The trend of the slip profiles are in general smooth. The high mass flux runs 32, 36, and 44 on
tube 5 produce the highest slip ratios. The trends, while generally smooth for all runs, are
Validation of Model Results
Figure 4-7 shows that an extremely steep void fraction build-up occurs in IFB, and departs
markedly from the relatively shallow build-up predicted by models such as that resulting from
the Lockhart-Martinelli parameter. These results are consistent with findings of Per Ottosen
(1980) in which void profiles in IFB conditions were measured using y-ray scattering. Ottosen
published the first known results from the use of y-ray absorption to measure void fraction in low
velocity IFB nitrogen. Figure 4-9 presents results from three of his many runs. It is apparent
that void fractions versus equilibrium quality (he made no attempt to quantify true mass quality)
rise very steeply. He observed the transition from IAFB to DF at void fractions between 80-
90%. All of his experiments were at approximately constant wall temperature conditions.
Additionally, all his data represent much lower mass fluxes than these hydrogen data. Perhaps
most importantly, the mass velocities were at least an order of magnitude lower than those in
these hydrogen data.
While a fine quantitative comparison is not made here due to the differences in
experimental conditions, a qualitative comparison is reasonable. It is apparent that extremely
rapid void fraction build-up is a characteristic of IFB.
Rohsenow and coworkers (Dougall and Rohsenow, 1963; Laverty and Rohsenow, 1967;
Forslund and Rohsenow, 1969) used nitrogen in their studies of IFB. In their work, they
determined the actual mass quality. They observed that the transition from IAF to DF occurred
at a mass quality of about 10%. Combining this observation with Ottosen' s of the void fraction
at transition, it can be concluded that void fractions of 80 90% at a mass quality of 10% are
typical. These experimental observations agree well with the results of this model.
Range of Validity
To avoid the momentum and energy balances of the high mass flux runs on tube 4, the
range of validity of this model has been reduced in terms of mass flux only. A total of eight runs
have been excluded from further analysis due to the inability of this model to reproduce the
pressure drop and energy balances. Forty-three runs remain. The remaining data for which the
balances are acceptable have the following range: pressures from 180 kPa to the critical pressure,
mass fluxes from 300 kg/m2-S to 2500 kg/m2-S, and heat fluxes from about 370 kW/m2 to about
Table 4-1. List of tube numbers, dimensions, and runs executed with the tubes.
Tube ref # Inner diameter Length Run numbers
1 1.8 60.96 1-6
2 1.113 60.96 7-10
3 0.851 60.96 11-21
4 0.478 60.96 22-31
5 0.79 30.48 32-51
Table 4-2. Statistical analysis of pressure data show goodness of fit through R2, and relative
unevenness of data through normalized (by pressure drop across test section length)
standard deviation. Results are from least squares fit of third order.
oJ Norm a R2
set Pa a / P
1 248 9.05E-03 0.999
3 229 1.57E-02 0.998
5 232 1.90E-02 0.997
7 328 7.28E-03 1
9 2202 5.39E-02 0.981
11 229 1.09E-03 1
13 195 6.15E-03 1
15 232 2.28E-03 1
17 279 1.11E-03 1
19 288 1.95E-03 1
21 301 1.35E-03 1
23 408 1.76E-03 1
25 1697 6.24E-03 1
27 1907 9.48E-03 0.999
29 820 7.86E-03 1
31 301 2.76E-03 1
33 681 1.00E-02 0.999
35 466 1.21E-02 0.999
36 838 1.06E-02 0.999
37 1888 1.02E-02 0.999
39 1082 1.01E-02 0.999
40 875 7.31E-03 1
41 1050 9.19E-03 0.999
42 7040 7.32E-03 1
43 871 46.5E-03 1
46 2189 1.30E-02 0.999