NASA/TM–2011–215988

Issues of Long-Term Cryogenic Propellant

Storage in Microgravity

C. B. Muratov

Department of Mathematical Sciences, New Jersey Institute of Technology, Newark, NJ 07102

V. V. Osipov, V. N. Smelyanskiy

Applied Physics Group, Intelligent Systems Division

NASA Ames Research Center Moffett Field, CA 94035-1000

October 2011
NASA STI Program . . . in Profile

Since its founding, NASA has been dedicated to the

advancement of aeronautics and space science. The

NASA scientific and technical information (STI)

program plays a key part in helping NASA maintain

this important role.

The NASA STI Program operates under the

auspices of the Agency Chief Information Officer. It

collects, organizes, provides for archiving, and

disseminates NASA’s STI. The NASA STI Program

provides access to the NASA Aeronautics and Space

Database and its public interface, the NASA

Technical Report Server, thus providing one of the

largest collection of aeronautical and space science

STI in the world. Results are published in both

non-NASA channels and by NASA in the NASA

STI Report Series, which includes the following

report types:

• TECHNICAL PUBLICATION. Reports of

completed research or a major significant phase of

research that present the results of NASA

programs and include extensive data or

theoretical analysis. Includes compilations of

significant scientific and technical data and

information deemed to be of continuing reference

value. NASA counterpart of peer-reviewed formal

professional papers, but having less stringent

limitations on manuscript length and extent of

graphic presentations.

• TECHNICAL MEMORANDUM.

Scientific and technical findings that are

preliminary or of specialized interest, e.g., quick

release reports, working papers, and

bibliographies that contain minimal annotation.

Does not contain extensive analysis.

• CONTRACTOR REPORT. Scientific and

technical findings by NASA-sponsored contractors

and grantees.

• CONFERENCE PUBLICATION.

Collected papers from scientific and technical

conferences, symposia, seminars, or other

meetings sponsored or co-sponsored by NASA.

• SPECIAL PUBLICATION. Scientific, technical,

or historical information from NASA programs,

projects, and missions, often concerned with

subjects having substantial public interest.

• TECHNICAL TRANSLATION. Englishlanguage translations of foreign scientific and

technical material pertinent to NASA’s mission.

Specialized services also include creating custom

thesauri, building customized databases, and

organizing and publishing research results.

For more information about the NASA STI

Program, see the following:

• Access the NASA STI program home page at

http://www.sti.nasa.gov

• E-mail your question via the Internet to

help@sti.nasa.gov

• Fax your question to the NASA STI Help Desk at

443-757-5803

• Phone the NASA STI Help Desk at 443-757-5802

• Write to:

NASA STI Help Desk

NASA Center for AeroSpace Information

7115 Standard Drive

Hanover, MD 21076–1320
NASA/TM–2011–215988

Issues of Long-Term Cryogenic Propellant

Storage in Microgravity

C. B. Muratov

Department of Mathematical Sciences, New Jersey Institute of Technology, Newark, NJ 07102

V. V. Osipov, V. N. Smelyanskiy

Applied Physics Group, Intelligent Systems Division

NASA Ames Research Center Moffett Field, CA 94035-1000

National Aeronautics and

Space Administration

NASA Ames Research Center Moffett Field, CA 94035-1000

October 2011
The use of trademarks or names of manufacturers in this report is for accurate reporting and does not constitute an

offical endorsement, either expressed or implied, of such products or manufacturers by the National Aeronautics and Space

Administration.

Available from:

NASA Center for AeroSpace Information

7115 Standard Drive

Hanover, MD 21076-1320

443-757-5802
Abstract

Modern multi-layer insulation (MLI) allows to sharply reduce the heat leak into cryogenic propellant

storage tanks through the tank surface and, as a consequence, significantly extend the storage

duration. In this situation the MLI penetrations, such as support struts, feed lines, etc., become one

of the most significant challenges of the tank’s heat management. This problem is especially acute

for liquid hydrogen (LH2) storage, since currently no efficient cryocoolers exist that operate at very

low LH2 temperatures (∼20K). Even small heat leaks under microgravity conditions and over the

period of many months give rise to a complex slowly-developing, large-scale spatiotemporal physical

phenomena in a multi-phase liquid-vapor mixture. These phenomena are not well-understood nor

can be easily controlled. They can be of a potentially hazardous nature for long-term on-orbital

cryogenic storage, propellant loading, tank chilldown, engine restart, and other in-space cryogenic

fluid management operations. To support the engineering design solutions that would mitigate these

effects a detailed physics-based analysis of heat transfer, vapor bubble formation, growth, motion,

coalescence and collapse is required in the presence of stirring jets of different configurations and

passive cooling devices such as MLI, thermodynamic vent system, and vapor-cooled shield. To

develop physics-based models and correlations reliable for microgravity conditions and long-time

scales there is a need for new fundamental data to be collected from on-orbit cryogenic storage

experiments. Our report discusses some of these physical phenomena and the design requirements

and future studies necessary for their mitigation. Special attention is payed to the phenomena

occurring near MLI penetrations.

Contents

1 Introduction 2

2 Background 4

3 Physics of long-term cryogenic storage 6

3.1 Basic thermodynamics . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 8

3.1.1 Time to saturation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 8

3.1.2 Time to complete evaporation . . . . . . . . . . . . . . . . . . . . . . . . . . . 9

3.1.3 Self-pressurization and liquid heating due to mechanical work . . . . . . . . . 9

3.1.4 Vapor-cooled shield and broad area cooling . . . . . . . . . . . . . . . . . . . 10

3.1.5 Penetration cooling . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 10

3.1.6 Helium dissolution hazard . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 15

3.2 Pressure control . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 16

3.2.1 TVS hazard induced by ullage motion . . . . . . . . . . . . . . . . . . . . . . 16

3.2.2 Boiling near hot spots . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 16

3.2.3 Formation of a multi-phase liquid-vapor foam . . . . . . . . . . . . . . . . . . 18

3.3 Dynamics of bubble growth and collapse . . . . . . . . . . . . . . . . . . . . . . . . . 19

3.3.1 Onset of bubble nucleations at hot spots . . . . . . . . . . . . . . . . . . . . . 19

3.3.2 Explosive nucleate boiling hazard . . . . . . . . . . . . . . . . . . . . . . . . . 21

3.3.3 Bubble growth over a hot spot . . . . . . . . . . . . . . . . . . . . . . . . . . 21

3.3.4 Bubble collapse and accumulation in the subcooled liquid . . . . . . . . . . . 22

1
3.3.5 Bubble collapse by pre-start pressurization . . . . . . . . . . . . . . . . . . . 23

4 Mathematical models of bubble dynamics 24

4.1 Conductance-based models . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 24

4.2 Bubble growth as function of bubble location . . . . . . . . . . . . . . . . . . . . . . 25

4.3 Interaction of bubbles . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 26

4.4 Ostwald ripening . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 28

4.5 Motion of bubbles suspended in the liquid . . . . . . . . . . . . . . . . . . . . . . . . 29

4.6 Bubble creation and departure from the wall . . . . . . . . . . . . . . . . . . . . . . 30

5 Computational analysis 30

6 Summary and conclusions 34

A Analysis of heat exchange between the cooling tube and the strut 42

B Analysis of heat spreading in the tank wall 43

C Effect of heat loss into the liquid on the heat conduction in the tank wall 44

D Spreading of heat from the wall into the liquid 44

E Bubble collapse in the subcooled liquid via conduction through thermal boundary layer 45

1 Introduction

This paper aims at discussing some of the basic physics issues associated with long-term storage

of cryogenic liquids in zero gravity or microgravity environments. By “long-term” we mean, for

example, the durations of the currently envisioned extended storage periods in the low earth orbit

(LEO), which range from months to years. For NASA’s present and future space exploration

missions, understanding the behavior of cryogenic liquids over long periods of storage is of crucial

importance, because of the fundamental role played by cryogenic propellants, primarily liquid

hydrogen (LH2) and liquid oxygen (LOX), in rocket propulsion, specifically for long-range missions

[1]. The very feasibility of using liquid propellant engines based on LH2 and LOX in the future

long-range missions depends on the success of storing these propellants under microgravity or zero

gravity for extended periods of time. One of the currently considered exploration strategies calls

for the development of propellant storage and transfer facilities in LEO [2]. These “fuel depots”

will need to be able to spend significant amounts of time (at least on the order of several months) in

LEO without any substantial propellant losses due to boil-off [2–5]. With the current passive heat

insulation technologies, it is theoretically possible to reduce the cryogen boil-off rate to below 3%

per month [6]. Even so, this issue becomes a challenge when the required storage duration exceeds

6 months, and yet a greater one for manned missions to Mars [7].

Cryogenic fluid management (CFM) in microgravity provides a number of fundamental physical

challenges, many of which were previously discussed in the literature [6, 8–13]. This is especially

relevant to storage of LH2 because of its low critical temperature. One of the main features of

2
microgravity environments is that due to much-reduced levels of g-forces and their generally timevarying character, the vapor bubbles that form as a result of boil-off at the tank walls mainly near

MLI penetrations (hot-spots) will not rise quickly towards the ullage space, as they do under normal

gravity. Instead, they may slowly grow to very large sizes (tens of centimeters), or they may detach

from the wall, migrate toward the stagnation areas of the stirring flow and accumulate there, forming

regions of saturated liquid and complicated foam-like vapor-liquid structures whose properties may

be not easy to control. These processes are governed by the complex heat transfer mechanisms in

the near-wall region; capillary and g-forces; complex dynamics of nucleate boiling; bubble growth,

detachment and collapse; chemical traces that can accumulate in the liquid with time and affect

its properties. In long-term storage missions foam or bubble colonies can grow at the expense of

the single ullage space. They may not be easily removed by tank pressurization because the heat

released from vapor condensation may raise the temperature of the liquid surrounding the bubbles

to the saturation temperature at the higher pressure. Capillary forces may be sufficiently strong

to prevent the detachment of the foam from the tank walls by any realistic stirring flow that keeps

the ullage intact. Basic challenges, therefore, include control of the tank pressure, temperature,

ullage space size and location, boil-off venting, and work of liquid acquisition devices (LAD) that

can be clogged by the foam. Similarly, since the role of buoyancy-driven convection, which is

the main mechanism of heat transfer on earth, is greatly reduced in microgravity, heat transfer

mechanisms will be significantly altered. Vapor and fluid motion, in turn, will be dominated by the

capillary forces, heat transfer-mediated bubble dynamics, bubble coalescence, Ostwald ripening and

the induced thermocapillary convection. The resulting bubble patterns and near-wall dynamics,

especially around the MLI penetrations can substantially depend on the type of the wall material,

chemical traces, vibrations and other external factors.

In view of these complications, the basic technical issues that need to be dealt with in today’s

design of successful cryogenic storage and transfer devices for long-term operation in microgravity:

heat transfer management, pressure control, design of tank stirring, mass gauging, liquid acquisition, and fluid transfer are much more challenging [11], compared to the Apollo era short duration

missions, in which a low level of gravity was propulsively maintained [14–16].

A successful treatment of the pressing technical issues of cryogen management in microgravity

is impossible without a thorough mechanistic understanding of the underlying physical processes of

nucleate boiling. Surprisingly, detailed physical understanding of nucleate boiling phenomena is still

lacking today (see e.g. [17]). This may be due to the fact that boiling is a strongly non-equilibrium

phenomenon in which an interplay between stochastic nucleation events at the micro-scale and

complicated deterministic nonlinear dynamics at macro-scale takes place (for reviews, see [18–21]).

At the same time, in view of the fact that microgravity presents quite a different environment

compared to the usual environment on earth, one should exercise caution in applying the engineering correlations developed under earth gravity conditions to the design of cryogenic systems

to be operated in space [21, 22]. To address the above technology gaps, it is necessary to collect

fundamental data on liquid-vapor structure and dynamics during long-term storage in microgravity

from carefully designed in-space long-duration experiments. Experimental work should be done in

combination with a detailed physics analysis, mechanistic modeling, first principles computational

and multi-scale approaches.

Here we perform some basic physical estimates in order to evaluate the relative importance of

different physical processes during long-term cryogenic storage. We concentrate our efforts on LH2,

since it is the cryogen of primary importance to rocket propulsion and is also the most difficult

3
low gravity and the uncertainty of the position of the liquid-vapor interface if the propellants are not settled. 

Cryogenic Fluid Management (CFM) technology development is addressing these issues through ground testing and 

analytical model development, while having crosscutting applications and benefits to virtually all missions requiring 

in-space operations with cryogens. Liquid hydrogen (LH2) is the most challenging of the three propellants but has 

the larger technology database since it has been used as the test fluid for many CFM experiments since the 1960’s.1 2

There is less CFM test experience and data with liquid methane (LCH4) and liquid oxygen (LO2) and these 

propellants are the primary focus of the current development activity. 

Thermal Control

- Insulation (launch environments

and in-space)

- Vapor cooled shields

- Low conductivity/ cooled support

structure

Liquid Acquisition

- Capillary retention devices for 

OMS/RCS

Pressurization

- Cold helium

Pressure Control

- Zero-g venting 

(thermodynamic vent and 

heat exchanger)

Lightweight Cryogenic Tank

- Metallic (Al-Li)

- Composite overwrap (pressure 

fed system)

Propellant Gauging

- Settled propellant

- Inventory (Bookkeeping)

- Pressure-volume-temperature (PVT)

- High accuracy zero-g techniques

Vent or to vapor 

cooled shields

E F

STRUT

Liquid Propellant

To Main 

Engine/ RCS 

Manifold

Compartment Divider &

Fluid Communication

Window (screen)

Spherical or 

cylindrical tanks

Pink - flight or extensive ground 

demonstration

Blue - key components ground 

tested

Orange - technology 

development required

Color Code for Text

Thermal Control

- Insulation (launch environments

and in-space)

- Vapor cooled shields

- Low conductivity/ cooled support

structure

Liquid Acquisition

- Capillary retention devices for 

OMS/RCS

Pressurization

- Cold helium

Pressure Control

- Zero-g venting 

(thermodynamic vent and 

heat exchanger)

Lightweight Cryogenic Tank

- Metallic (Al-Li)

- Composite overwrap (pressure 

fed system)

Propellant Gauging

- Settled propellant

- Inventory (Bookkeeping)

- Pressure-volume-temperature (PVT)

- High accuracy zero-g techniques

Vent or to vapor 

cooled shields

E F

STRUT

Liquid Propellant

To Main 

Engine/ RCS 

Manifold

Compartment Divider &

Fluid Communication

Window (screen)

Spherical or 

cylindrical tanks

Pink - flight or extensive ground 

demonstration

Blue - key components ground 

tested

Orange - technology 

development required

Color Code for Text

Pink - flight or extensive ground 

demonstration

Blue - key components ground 

tested

Orange - technology 

development required

Color Code for Text

Figure 1. Conceptual exploration vehicle cryogenic propellant tank elements. 

The primary CFM technology elements are: passive and active thermal control, pressure control, liquid acquisition, 

mass gauging and propellant feed line conditioning. Other CFM technology areas affecting in-space propulsion 

systems include leak-free couplings and disconnects, light weight composite tanks and support structure, leak 

detectors and component and subsystem integration. The conceptual tank schematic in Figure 1 illustrates how these 

elements are combined into a cryogenic propellant storage and delivery system. 

The CFM development goal of allowing fluid management functions without settling propellants enables the 

exploration architecture and provides major propulsion system benefits by simplifying vehicle operations, reducing 

system mass, and expanding operational and architectural options. Performing these functions with settled 

propellants provides a potential backup mode, reducing the overall CFM risk. 

NASA/TM—2007-214810 2

Figure 1. The schematics of a cryogenic storage tank (from [6]).

in terms of CFM due to its low boiling point. Let us emphasize that we aim at obtaining only

relatively rough estimates that take into account the long-term nature of storage. Thus, our main

tool will be dimensional analysis, with minimal reference to more advanced mathematical tools.

Once the main physical processes acting on the considered long timescales are identified, relevant

space experiments can be designed and more precise calculations may be made using advanced

mathematical and high-fidelity computational tools. In short, our main goal is to identify these

processes and the issues, such as safety hazards and design optimization parameters, which arise

specifically during extended periods in microgravity.

2 Background

We start with some basic considerations relevant to large-scale cryogenic storage in microgravity.

A conceptual representation of a cryogenic fuel tank [6] is shown in Fig. 1. To fix ideas, let us

consider one of the proposed designs for the LH2 tank of the Earth Departure Stage (EDS) for

moon missions [23], in which the tank has the shape of a rounded cylinder with height H0 = 12 m

and radius R0 = 2.5 m, containing 15 tons of LH2. The scale of the tank is similar to that used in

4
the S-IVB stage of the Saturn V rocket in the 1960’s and 70’s, and a brief comparison is, therefore,

appropriate.

Indeed, NASA’s most comprehensive experience with large-scale cryogenic storage tanks in

orbital conditions goes back to the Apollo moon missions1[15, 16, 25] (see also reviews of other

CFM experimental activities in [9, 14]). The third stage of the Saturn V rocket was propelled by

LH2 and LOX, containing 19,800 kg and 88,800 kg, respectively (here and below the numbers are

from the Apollo 17 mission [26, 27]). About a quarter of the propellants (5,000 kg of LH2 and

25,200 kg of LOX) was utilized for orbital insertion, and the remaining amount was used for the

translunar injection burn. Once in orbit, the third stage spent about 3 hours in LEO. The tank

insulation employed (polyurethane foam attached to the interior side of the tank wall) brought

the LH2 boil-off amount down to about 1,000 kg, still an acceptable margin of under 10% of LH2

available for the second burn. The flow of continuously vented hydrogen vapor (GH2) was used to

provide enough thrust (on the order of 10−5g − 10−4g) to ensure that the propellants were settled

at the bottoms of the tanks at all times during the orbital coast phase. The presence of small

upward g-force ensured that the boil-off bubbles rose in the thin convective layer along the tank

walls without entering the bulk liquid. Immediately prior to the second burn, the LH2 tank was

rapidly pressurized by the stored heated helium gas (GHe), raising the tank pressure from about

1.5 to 2 atm. This short-time pressurization, followed by firing of the ullage motors must have

condensed the smaller vapor bubbles and made the larger bubbles move towards the ullage. The

resulting LH2 liquid at the bottom of the tank should have, therefore, contained little or no bubbles,

allowing to safely fire the engine. After that, the thrust of the engine would maintain the gas-free

liquid at the LH2 intake, with screens adding an extra protection.

While hugely successful in bringing man to the surface of the moon and back, this approach

may not be applied to the missions currently under consideration. They key reason why the Apollo

CFM approach worked for the lunar program was that the required in-orbit storage time for the

cryogenic propellants was short enough, so it was possible to tolerate a large boil-off rate and, as a

consequence, avoid microgravity conditions altogether during the time in orbit. Thus, the Apollo

approach carefully avoided dealing with long-term CFM issues associated with microgravity. Any

kind of large modern long-range mission, however, would require storing cryogens for extended

periods of time. To achieve this, one would need to drastically reduce the amount of boil-off and

work in micro-g or zero-g environments. New approaches are, therefore, needed to answer these

emerging challenges.

Modern multi-layer insulation (MLI) allows to dramatically reduce the boil-off rate compared

to the Apollo missions. Let us assume that the tank is wrapped in an MLI blanket of 50 layers.

Using the Lockheed correlation [28], we can estimate the heat flux through the MLI, given the outer

environment temperature T0 ' 240K, to be q0 ' 7.4 × 10−2 W/m2[29]. Taking for simplicity the

tank area to be S0 = 2πR0H0 ' 200m2, we obtain a lower bound of Q0 = q0S0 ' 15 W for the total

heat flow into the tank, with the corresponding boil-off rate of at least 90 kg/month, or 0.6% of

the propellant mass per month. Of course, these numbers must underestimate the actual heat flow,

since they do not take into account heat leaks through various MLI penetrations by struts, feed

lines, etc., as well as imperfections in the MLI itself. Let us also note that the heat flux depends

1We note that liquid helium (LHe) has flown for extended periods of time on a number of scientific missions [13].

Let us point out, however, that in these missions LHe is chilled down to the superfluid state. This makes the case of

LHe storage very different from all other cryogens, since, in contrast to all other cryogens, superfluid LHe has infinite

heat conductivity [24], which prevents it from thermal stratification and nucleate boiling.

5
very sensitively on the outer environment temperature T0. In the extreme case of T0 ' 400K, we

find q0 ' 0.8 W/m2, an order of magnitude higher than the one computed previously, resulting in

the heat flow of Q0 ' 160 W and an unacceptably high boil-off rate of 6% per month.

Even with very efficient MLI insulation, the loss of propellant becomes prohibitive for extended

missions. Therefore, active boil-off reduction techniques are necessary to improve retention of

usable propellants. One idea that has been developed over recent years is to employ zero-boil-off

technology (ZBOT) [3, 30–32]. While ZBOT approach was demonstrated to be successful in the

case of cryogens with higher boiling points, e.g. LOX, no cryocoolers enabling ZBOT yet exist that

could operate at LH2 temperatures [33]. Another idea to further decrease the heat inflow into the

LH2 tank is to use the concept of broad area cooling (BAC), whereby the tank is surrounded by a

network of tubes carrying a coolant fluid [6, 12, 31, 33, 34]. Circulating the fluid through the tubes

with the subsequent heat removal by cryocoolers operating at higher temperature may then allow

to significantly reduce heat penetration. Use of vapor-cooled shield (VCS) [6, 7, 13, 35] would be

particularly efficient, since thermalizing GH2 with the outer layers of the MLI could increase the

storage time up to a factor of 6 (see Secs. 3.1.4 and 3.1.5 for more details). We note, however,

that strong localized heat leaks through MLI penetrations (the main subject of the present paper)

provide one of the greatest CFM challenges for long-term cryogenic storage.

3 Physics of long-term cryogenic storage

Let us now perform some basic estimates for the sample 15-ton LH2 tank whose dimensions were

introduced in Sec. 2. Since the precise parameters of the MLI performance admit significant

variation, we will take an overall heat leak per unit area with an ample margin: q0 = 0.4 W/m2,

giving a total heat influx of 80W through the MLI (see also [6,29]). In addition, we will assume that

another 40W of heat enters the tank via penetrations in the form of localized heat sources, giving

the total incoming heat flow of Q0 = 120W. We note that MLI penetrations, such as support struts,

feed lines, etc., may provide the greatest challenge in the tank’s heat management. For example,

taking the characteristic parameters of the orbiter support strut from the Space Shuttle external

tank, which is a tubular structure of radius R ' 20 cm, thickness d ' 5 mm and length L ' 1

m, with thermal conductivity κ ' 7 W/(m·K) of Inconel 718 alloy at T ' 100K [36], we find that

the conductive heat leak into the tank can be estimated as Qstrut ' 2πκRdT0/L = 11W. Note

that this formula may significantly underestimate Qstrut, since it does not take into account the

additional heat entering the strut through its own thermal insulation. Similarly, for a titanium

strut with d ' 1 cm and κ ' 15 W/(m·K) [37] we find Qstrut ' 45W, and for Al 2219 strut with

d ' 1.5 cm and κ ' 70 W/(m·K) [38] we find a prohibitively high Qstrut ' 320W. In view of

the preceding considerations, however, the conductive heat leaks must not exceed several Watts

per penetration in order for the tank to remain within the acceptable thermal budget. There is,

therefore, a significant trade-off between the structural and thermal properties of the materials

used, requiring strong materials with low thermal conductivity and a possible need for external

penetration cooling [6].

To proceed, we need to specify the operating parameters for LH2 in the tank. We will assume

that the tank is initially at pressure p0 = 1.6 atm, corresponding to the saturation temperature

Ts0 = 22K of LH2 (parahydrogen), and that LH2 is subcooled to TL0 = 20.3K, corresponding

to saturation temperature at 1 atm. Here and everywhere below the definitions and values of

the parameters of hydrogen used are are listed in Table 1 (for pressure p0 at saturation [39]).

6
Table 1. Physical parameters used in the estimates.

Parameter Value Meaning

DHe 5 × 10−9 m2/s Diffusivity of helium in LH2

H0 12 m Tank height

P 5 - 15 W Local heat leak power

Q0 120 W Total heat leak

R0 2.5 m Tank radius

Rg0 1.8 m Ullage bubble radius

RH2 4,124 J/(kg·K) Gas constant of hydrogen

RHe 2,077 J/(kg·K) Gas constant of helium

S0 200 m2 Tank surface area

T0 240K Exterior environment temperature

TL0 20.3K Subcooling temperature

Ts0 22K Saturation temperature of LH2

V0 240 m3 Tank volume

Vg0 24 m3 Ullage volume

cL 10,820 J/(kg·K) Specific heat of LH2 at constant pressure

cw 10 J/(kg·K) Specific heat of aluminum

g0 9.81 m/s2 Earth’s acceleration of gravity

g 0 − 10−6g0 Microgravity acceleration

h 1 cm Tank wall thickness

p0 1.6 atm Operating pressure

q0 0.4 W/m2 Heat flux through the MLI

qL 4.35 × 105 J/(kg·K) Latent heat of LH2 vaporization

βL 0.0192 K−1 Thermal expansion coefficient of LH2

κL 0.101 W/(m·K) Heat conductance of LH2

κv 0.019 W/(m·K) Heat conductance of GH2

κw 20 – 200 W/(m·K) Heat conductance of aluminum

μL 1.16 × 10−5 Pa·s Viscosity of LH2

ρL 68.7 kg/m3 Density of LH2

ρv 2.07 kg/m3 Density of GH2

ρw 2,700 kg/m3 Density of aluminum

σL 1.65 × 10−3 N/m Surface tension of LH2

Initially, a 10% by volume ullage space with volume Vg0 = 24 m3is present, pressurized by cold

GHe. The required mass of GHe to produce the excess pressure of 0.6 atm is equal to MHe =

(p0 − patm)Vg0/RHeTL0 ' 35 kg. We note that at large ullage volumes (as LH2 is lost due to

boil-off or transfer from the tank), large amounts of cold helium gas are required to pressurize the

tank. For example, when the ullage occupies 50% of the tank volume, one would need to supply

175 kg of GHe, respectively. In practice, even greater amounts of helium may be required due to

dissolution of GHe in LH2 on long storage timescales (see Sec. 3.1.6). The liquid is subcooled in

order to avoid the presence of vapor bubbles in the bulk LH2, which are a potential hazard for

engine restart, etc.

7
Let us briefly comment here on the nature of microgravity environment experienced by the

cryogenic tank in LEO. For simplicity, we will consider the tank whose axis is oriented along its

velocity vector in a circular orbit. The contributions to g-forces can be separated into three parts.

The first part has to do with the spatial non-uniformity of the earth’s gravitational field and the

centrifugal force. The resulting effective g-force will be a linear function of the distance to earth,

vanishing at the tank axial mid-plane. The apparent gravity g will point away from earth on the

side of the tank farthermost from earth, and towards the earth at the side closest to earth. The

maximum value of |g| can be shown to be given by |g| ' 3g0R0/Rorbit ' 10−6g0, where g0 is the

acceleration of gravity on earth. This gives the Bond number Bo = ρL|g|R2

0

/σL ' 2.5, indicating

that capillary forces balance the apparent gravity forces in the case of the single ullage space.

In this situation the equilibrium ullage shape should be approximately an oblate spheroid, whose

center is located in the plane passing through the tank’s axis and is normal to the direction towards

Earth. Note that more generally the apparent g-forces will be time-dependent on the time scale of

the orbital period. From the balance of inertia and capillary forces, we find that the characteristic

timescale for the ullage to settle is of the order of tullage ∼

p

ρLR3

0

/σL ' 15 mins.

The second contribution to the g-forces is the g-jitter, caused by the motion of the astronauts.

These are impulses of acceleration with amplitude roughly of order 10−4g0 and duration 1 s in

random directions (see e.g. the discussion in [40]). It provides the random vibrational background

to the motion of fluid. Note that this contribution would be absent in an unmanned vehicle, or

would be the only significant g-force on space missions beyond LEO. Also, both the g-jitter and the

space-dependent apparent gravity may assist bubble coalescence. Finally, the third contribution to

g-forces is due to the Coriolis force, which arises due to motion of the liquid inside the tank.

3.1 Basic thermodynamics

We begin by evaluating the heat budget of the tank and related issues.

3.1.1 Time to saturation

First, let us calculate the time needed for LH2 to come from the subcooled condition to the saturation temperature under an assumption of perfect mixing (e.g. by an active mixer inside the

tank) and in the absence of any boiling and active cooling. This time is given by equating the total

amount of heat that entered the tank to the increase in the LH2 heat content:

tsat,mixed =

cLρL(V0 − Vg0)(Ts0 − TL0)

Q0

' 26 days. (1)

On the other hand, in the absence of mixing, boiling, and any g-forces the heat will only penetrate

from the tank wall to the depth equal to the thermodiffusion length L =

p

κLt/(cLρL) of LH2 in

time t. Then the same balance leads to

tsat,thermodiffusive =

cLρLκL(Ts0 − TL0)

2

q

2

0

' 16 days. (2)

Note that, replacing Ts0 in Eq. (2) with 23K, one can see that in order to achieve a superheat of

∼1K at the tank wall (at which nucleate boiling normally occurs in LH2 under normal conditions on

earth [41–43]), one would need to wait t ' 40 days. Finally, for a tank in LEO, taking into account

8
convective transport in microgravity with g = 10−6

g0, we obtain for the tank wall temperature Tw

an estimate Tw−TL0 ' 0.22K, which is based on the Nusselt number NuR0 = q0R0/(κL(Tw−TL0)) '

46 and the Rayleigh number for these parameters Ra = gβL(Tw −TL0)cLρ

2

LR3

0

/(μLκL) ' 2.8×107,

and we used the correlation NuR0 = 0.15Ra1/3[44, Eq. (9.31)]. Thus, the timescale on which the

bulk of LH2 heats to the boiling point under the considered heat loads in the absence of any other

sources and sinks of heat is about 1 month. Note that on such a long timescale heat conductance

alone is sufficient to carry the heat into the tank interior.

3.1.2 Time to complete evaporation

On the other hand, assuming that the tank is maintained at constant pressure p0, we can find the

LH2 storage time by equation the total heat that entered the tank to the heat needed to vaporize

all the LH2:

tstorage =

(qL + cL(Ts0 − TL0))ρL(V0 − Vg0)

Q0

' 22 months. (3)

In particular, about 700 kg or ∼5% by volume, of LH2 will be lost to boil-off in one month. This

also means that cold GHe needs to be supplied at a rate of ∼16 kg/month to enable maintaining

LH2 at subcooled conditions.

3.1.3 Self-pressurization and liquid heating due to mechanical work

While it will take on the order of one month in the absence of any other heat sinks for the bulk of

LH2 to be heated to the saturation temperature, local boiling may begin soon due to the localized

heat sources through penetrations. With the power from those sources equal to Q0 − q0S0 = 40W

and assuming that the boil-off bubbles remain attached to the hot spots on the tank walls, we can

estimate the the rate of boil-off of GH2 as

M ̇ GH2,boil-off =

Q0 − q0S0

qL + cL(Ts0 − TL0)

' 8 kg/day. (4)

This results in an increase of the tank pressure at the rate of  ̇p = M ̇ GH2,boil-offRH2Ts0/Vg0 ' 0.3

atm/day.

Let us also point out that GH2 forming as a result of boil-off performs mechanical work against

the liquid in order to expand. The mechanical work done on LH2 will be converted into heat

through viscous dissipation in LH2, resulting in a further temperature increase in the liquid bulk.

If the boil-off rate M ̇ GH2,boil-off is given by Eq. (4), then the work done by the vapor on the liquid

is

W =

p0M ̇ GH2,boil-off

ρv

' 7W. (5)

In other words, about 20% of the heat that goes into boil-off ends up heating the bulk liquid. To

this, one should also add the amount of work supplied to the liquid through stirring.

9
3.1.4 Vapor-cooled shield and broad area cooling

The use of VCS can significantly reduce the amount of boil-off by utilizing sensible heat of the

boil-off vapor [6, 7, 13, 35]. If the GH2 is simply vented overboard, the amount of heat removed

from the tank by the vapor will be Qboil-off ' qLm ̇ GH2, where  ̇mGH2 is the mass flow of GH2

through the vent. If, on the other hand, the vapor is allowed to thermalize with the outer surface

of the MLI at temperature T0 ' 240K before being vented, the amount of heat removed equals

QV CS ' (qL +cL(T0 −TL0))  ̇mGH2 (taking into account that the specific heat of GH2 does not vary

significantly with temperature in the considered interval [39]). The ratio of the two is given by

QV CS

Qboil-off

' 1 +

cL(T0 − TL0)

qL

' 6.5. (6)

This means that the use of VCS may reduce boil-off rate for a fixed heat load by over a factor of 6

by intercepting the heat entering into the tank. Further increase in the efficiency may be achieved

by passing the warm GH2 through para-ortho converters [7].

Broad area cooling (BAC) is another promising concept that takes advantage of actively cooled

GHe running through tubes within the MLI to capture the heat entering from the tank environment

and allows to significantly reduce the heat flow per unit area through the MLI [33] (see also [34]

for a related concept of active co-storage). One of the main difficulties in applying the BAC

technology is efficient thermal bonding of the BAC tubing to the tank structures [6, 31]. Recently,

both the “tube-to-tank” and the “tube-to-shield” concepts, whereby the BAC tubing is bonded

either directly to the tank wall or to an intermediate layer of the MLI, respectively, have been

successfully demonstrated [12, 31, 33, 45]. Let us point out that these concepts may also be used as

part of a VCS in conjunction with TVS [35]. Note, however, that in the absence of boiling inside

the tank the VCS/TVS system may not be able to intercept the incoming heat, since the heat

entering the tank may raise the LH2 temperature locally near the tank walls without producing a

pressure or bulk LH2 temperature rise. In this case, the TVS will not operate, and the propellant

heating near the tank walls may lead to the potentially dangerous explosive nucleate boiling hazard

(see Sec. 3.3.2 for more detail). We also note that in order for the tube-to-tank BAC design to

work in an LH2 tank, the cold GHe must be circulated at the temperature of about 20K to avoid

boiling of LH2 at the tank walls. Since no efficient cryocoolers currently exist operating at these

temperatures, the tube-to-tank concept is not currently applicable to LH2 storage.

3.1.5 Penetration cooling

One of the most challenging problems in the tank heat management is to control the heat leaks

from various MLI penetrations, such as the tank structural supports (struts) and propellant feed

lines [6]. As was discussed in Sec. 3.1.4, using either VCS or BAC for passive or active cooling

of the heat shield, respectively, one could all but eliminate the heat leaks through the MLI into

the tank. In this case, the penetrations will provide most of the heat load to the propellant. Note

that in contrast to the MLI, where radiative heat transfer dominates, it is not possible to adapt

the active tube-to-shield BAC concept operating at intermediate temperatures (see Sec. 3.1.4) to

intercept most of the heat coming through penetrations because of the dominant role of conductive

heat transfer there. A natural idea (an extension of the VCS concept) is, therefore, to sacrifice

some of the stored hydrogen to passively cool the penetrations at their points of contact with the

tank walls and, in particular, to suppress possible boiling in those areas. The hydrogen is most

10
conveniently supplied by TVS, providing a regulatory feedback between the tank heat load and the

amount of LH2 used for cooling.

In view of the above discussion, to emphasize the physical phenomena under consideration,

in this section we will assume that all heat entering into the tank comes from the penetrations.

To make our considerations more quantitative, let us assume for simplicity that the tank contains

Nstrut = 8 identical penetration structures (struts). Consistently with our previous assumptions in

Sec. 3, each strut is assumed to carry heat in the amount of Qstrut = Q0/Nstrut = 15W conductively

into the tank. Under our assumptions about the storage times for the tank, the maximum available

LH2 budget (not including the effect of VCS or BAC) that can be used to cool the penetrations is

Jstrut =

Qstrut

qL

= 3.3 × 10−5kg/s. (7)

This mass flow of cold GH2 and/or LH2 can be used, e.g., to locally cool the MLI penetrations in

the shape of thick pipes. We consider a helical coil in perfect thermal contact with the strut surface,

which is formed by winding a thin hydrogen-carrying tube around the strut (see Fig. 2). We now

analyse efficiency of such a method of local cooling of the penetrations. In the following, we assume

that the cooling tube has characteristic length L = 10 m. We also assume that TVS operates in

the regime in which LH2 at p0 = 1.6 atm and TL0 = 20.3K passes through a Joule-Thompson valve

with pressure p1 = 0.2 atm at the other end. By considering an isenthalpic process, we then find

that hydrogen will be exiting the valve in the form of a two-phase mixture at T1 ' 16K containing

the mass fraction η ' 0.9 of LH2 [39]. We will consider the situation in which hydrogen from TVS

is fully vaporized and is then thermalized with the LH2 propellant before being supplied into the

cooling tube (see Fig. 2). In order to estimate the minimum tube radius Rtube, consider that GH2

is supplied in the form of gas at T = TL0. Using the Darcy-Weisbach equation [46, 47]

∆p =

f ρu2L

4Rtube

(8)

for the pressure drop ∆p across the tube, where f is the Darcy friction factor, ρ is the fluid density,

u is the average fluid velocity, together with Haaland’s approximation to the solution of Colebrook’s

equation for f [48]:

f =

 

1.8 log10 "

drough

7.4Rtube 1.11

+

6.9

Re #!−2

, (9)

where drough is the tube surface roughness and Re is the flow Reynolds number, we can estimate

the minimum inner tube radius Rtube needed to supply the necessary LH2 flow Jstrut given by Eq.

(7). Taking drough = 1.5 μm [44] and ∆p = p1 = 0.2 atm, we find that the minimum tube radius

is Rtube ' 1.1 mm. Note that this estimate is quite insensitive to various assumptions on the

parameters, such as the tube roughness or the pressure drop.

Let us now estimate how effective these tubes will be in cooling the penetrations. Assume that

the temperature of GH2 flowing inside the coil tube that winds around the hot strut changes from

Tv = T1 to Tv = T2. By thermodynamic considerations, we must have

Qstrut = cp(T2 − T1)Jstrut. (10)

11
!"#$%&"#$

!"#$%&''(#$)%*+,-%

!"#$%&''(#$)%*+,-%

./0%1#*"%!.2034%

5-*6(%*6$7%16((%

5.8%

/'*%

9:6&";-$*%

!20<34%

!=004%

Cold 

gas flow 

from TVS

Exit flow 

a) b) Figure 2. Cooling a penetration by GH2 from TVS. (a) The tank schematics and (b) the close-up

of the cooled penetration.

It then follows that the GH2 supplied at T1 = TL0 = 20.3K will have the capacity to remove

Qstrut = 15W from the strut, if it is heated to T2 ' 65K. Let us now estimate length of the tube

which can remove this heat flow. For steady turbulent flow, the heat transfer coefficient htube can

be obtained from the Dittus-Boelter correlation [44]

htube =

κv

2Rtube

Nutube, Nutube ' 0.023 × Re4/5

tube Pr1/3tube, (11)

where Nutube, Retube and Prtube are the Nusselt, Reynolds and the Prandtl number associated

with the gas flow in the tube, respectively. Note that the formula in Eq. (11) can be equivalently

rewritten as

htube = 0.023cp



Jstrut

πR2

tube

4/5 

μv

2Rtube 1/5

Pr−2/3, (12)

i.e, at fixed Jstrut the heat transfer coefficient is proportional to R

−9/5

tube and, therefore, grows rapidly

with decrease of Rtube. Considering a tube of radius Rtube = 2 mm, we find that for the mass flow

given by Eq. (7) we have Retube ' 104 and Prtube ' 0.8, giving the Nusselt number Nutube ' 32

and the heat transfer coefficient htube ' 130 W/(m2·K) (using the parameters for GH2 at T = TL0

and p = p1 [49]). This corresponds to the power Ploop ' 4π

2RtubeRstruthtube(Ts0 − TL0) = 1.8W

taken away by the first loop of the tube, provided the strut at the base is at saturation temperature

Ts0. Therefore, several loops are needed to take away the heat from the strut. Let us note that the

pressure drop ∆p = p1 ensures the flow Jstrut in a tube of length L ' 200 m (see Eq. (8)).

The analysis of heat exchange between the strut maintained at T = T0 at the warm end and in

thermal contact with cold GH2-carrying tubes is rather involved and is presented in Appendix A.

The results depend significantly on many factors, including the dimensions and the material of the

penetration. Similarly, the ability of the tube to take away the heat from the penetration depends

in a non-trivial way on these factors, as we demonstrate below. For concreteness, in the following

we will assume (nominal regime) that a thermally insulated tubular strut of radius Rstrut = 10 cm,

12
a) b)

c) d)

0 20 40 60 80 100

20.5

21.0

21.5

22.0

22.5

23.0

23.5

x, mm

T, K

0 20 40 60 80 100

20

22

24

26

28

30

32

x, mm

T, K

0 10 20 30 40 50 60

20

25

30

35

x, mm

T, K

0 20 40 60 80 100

20

25

30

35

x, mm

T, K

Figure 3. Temperature distributions in the strut (blue line) and the tube (purple line) as a function

of the distance along the strut to the tank wall. (a) Nominal regime; (b) The strut radius is

increased to Rstrut = 15 cm; (c) The strut is shortened to Lstrut = 30 cm; (d) The GH2 flow rate

is decreased to Jstrut = 2 × 10−5 kg/s.

thickness dstrut = 5 mm, length Lstrut = 50 cm is made of pure titanium and take the value of

κstrut ' 10 W/(m·K) for thermal conductivity, corresponding to the low temperature end of this

material parameter [37]. Under these assumptions, the heat leak through the strut in the absence

of vapor cooling is equal to Qstrut = 15W, consistent with the earlier discussion. We will also take

the pitch of the helix dtube = 2 cm which is smaller than the strut radius and bigger than the tube

diameter in order to ensure high surface area of thermal contact and absence of thermal shorts

through the tubes.

The solution of the governing equations for the temperature distribution in both the strut and

the tube as a function of the distance from the tank along the strut under an assumption that the

tube captures all the heat entering into the strut at the warm end is presented in Fig. 3(a). One can

see that a temperature boundary layer develops near the tank wall, in which the strut temperature

fails to follow cold GH2 temperature on the length scale lstrut ' 2.5 cm (see also the discussion in

Appendix A). The maximum temperature at the strut base reaches 21.6K, which is still below the

saturation temperature Ts0 = 22K of LH2. Therefore, this rather crude analysis demonstrates that

boiling at the attachment point of the strut may be suppressed by TVS-produced vapor cooling

under suitable conditions.

Let us now see how sensitive this conclusion is to variations of various parameters of the considered cooling system. For example, consider the case in which the strut radius is increased by

50%. The solution of the governing equations is then presented in Fig. 3(b). One can see that in

13
LH2/GH2 flow from 

JT valve at T=16K 

ls

h

GH2 flow at T=20K 

Ltube

!"#$%

!&'(%)&**%

a) b) 

Figure 4. Cooling a penetration by LH2/GH2 mixture from TVS. (a) The tank schematics and (b)

the details of the tube-to-tank scheme.

this case the strut base temperature increases to T = 26.5K, corresponding to a 4.5K superheat.

Similarly, decreasing the strut length by 40% results in a strut base temperature T = 31K, a 9K

superheat. Finally, a decrease of the GH2 supply rate by 40% results in the strut base temperature

T = 29K, a 7K superheat. In all cases, the obtained superheat is big enough to initiate boiling.

One should, therefore, carefully consider the heat conduction problem associated with penetrations

in order to asses the feasibility of the proposed cooling strategy. One important factor to keep

in mind is that, according to the model predictions, the efficiency of the proposed vapor cooling

system decreases with an increase in the effective heat conductance coefficient of the strut. The

latter must also include the effect of the additional heat conduction pathways introduced by the

highly conducting thermal bonding material and the tubes themselves. It is not clear at this point

whether this approach may always guarantee boiling suppression near the MLI penetration points.

A more detailed and accurate analysis of the system is required.

Before concluding this section, we briefly comment on the possibility to use LH2 from the JouleThompson valve directly rather than GH2 from the TVS heat exchanger to cool the penetrations

(see Fig. 4). This can be accomplished, e.g., by the tube-to-tank configuration, carrying the

two-phase LH2/GH2 mixture from the Joule-Thompson valve to the penetrations. Let us point

out in the first place that the added cooling power from the latent heat of LH2 is in fact only

a small fraction of the total cooling power of the GH2 stream. Indeed, for the same reason that

the use of VCS may potentially increase the heat leak capacity of the LH2 tank by a factor of 6

(see Sec. 3.1.4), the fraction of the latent heat of LH2 to the sensible heat of GH2, when heated

from TL0 = 20.3K to T0 = 240K is estimated to be qLJstrut/(cp(T0 − TL0)) ' 0.2. Therefore, the

advantage of using LH2 rather than GH2 to cool the penetrations is rather questionable.

At the same time, in order to supply LH2 from the TVS output, the supply tubes must be

sufficiently isolated to prevent boil-up and thermalization of the hydrogen flowing through the

tubes with LH2 inside the tank. If not, LH2 will be transformed into GH2 due by the heat entering

the tubes that carry the liquid from the Joule-Thompson valve towards the tubes wound around

the strut (see Fig. 4). Indeed, the amount of heat entering a tube segment of length L may be

estimated as κwLh(TL0 − T1)/ls, where κw is the tank wall thermal conductivity, h is the wall

14
thickness, and ls is the length of temperature equilibration due to convection. In the presence of

natural convection alone in microgravity the value of ls may be estimated to be ls ' 1 m (see

Sec. 3.3.1 below for more detail), so for a highly conductive tank wall material with κw = 200

W/(m·K) and h = 1 cm we get that LH2 in the supply tube will boil completely on the length

L ' ηqLJstrutls/(ηκwh(TL0 − T1)) ' 1.7 m. Therefore, bringing LH2 from TVS to a penetration

without thermal insulation of the tubes on the way is problematic. On the other hand, suppose

that by using suitable thermal insulation the heat entering the supply tube on the way to the

penetration to be cooled is reduced sufficiently, so that LH2 at T1 = 16K is present at the first

point of contact between the cooling tube and the penetration (see Fig. 2(b)). Then, since LH2

in the tube is significantly cooler than the surrounding LH2 in the tank and the tank walls, the

tube may draw a substantial amount of heat from the tank wall rather than the penetration. It

is possible to adapt the analysis of Appendix C (leading to Eq. (18) below) to show that a 4K

subcooling of a strut base of radius Rstrut = 10 cm would result in a power P ' 17 W drawn from

the wall, for the same tank wall parameters as above, in addition to Pstrut = 15W coming from the

strut. This contradicts conservation of heat at the strut base. Therefore, the tube would not be

able to cool the penetration, indicating that the tube temperature cannot be maintained at T = T1.

This, in turn, implies the onset of film boiling in the cooling tube, which would greatly reduce its

heat transfer coefficient to the values characteristic of the use of GH2 vapor. Thus, the advantage

of using LH2 is once again greatly diminished.

3.1.6 Helium dissolution hazard

We point out that at TL0 = 20.3K the solubility limit of GHe in LH2 is ∼0.5% by weight [50].

Therefore, all LH2 in the tank is capable to absorb up to 75 kg of GHe, which is about 2 times

more than the mass MHe ' 35 kg of GHe in the ullage at the beginning. Let us estimate the

time in which GHe may dissolve in LH2. Taking the diffusion coefficient of the dissolved helium

DHe = 5 × 10−9 m2/s (assumed to be of the same order as the available value for neon [51], see

also [52]), assuming that the ullage has the shape of a spherical bubble of radius Rg0 = 1.8 m and

estimating the diffusive flux at the ullage boundary to be DHeρHe,sat/lHe, where ρHe,sat = 0.005ρL

is the helium saturation density and lHe =

√

DHet is the helium diffusion length, respectively, in

LH2, we obtain that in the absence of active mixing the dissolved mass of GHe in time t is

MHe,dissolved ' 4πR2

g0

lHeρHe,sat, MHe,dissolved ' 1.6 kg, t = 1 month. (13)

Let us note that in the presence of active mixing the rate of GHe dissolution may be significantly higher. For example, if LH2 is circulated with an average axial velocity ubulk = 1 mm/s

across the tank, then the Peclet number Pe = 2Rg0ubulk/DHe ' 7 × 105for the ullage bubble.

Therefore, the dimensionless Sherwood number Sh = 0.65Pe1/2 ' 550 ( [53, Eq. (3.52)], assuming spherical ullage in a background flow with velocity ubulk), and the dissolution rate becomes

M ̇ He = 2πRg0DHeρHe,satSh ' 30 kg/month. Thus, the entire mass of GHe may dissolve in LH2

in only one month, leading to the collapse of the ullage pressure. The dissolution rate may further

increase due to ullage bubble motion. Therefore, on the long-term storage timescales one needs to

evaluate the potential ullage collapse hazard due to GHe dissolution in LH2. One also needs to

take into consideration the possible effect of dissolved non-condensible helium gas on the boiling

characteristics [21, 54, 55].

15
3.2 Pressure control

When the sufficient level of superheat is reached at the tank walls and/or enough heat enters

through penetrations, nucleate boiling will start. As vapor is created due to boil-off, the tank

pressure rises, requiring venting in the absence of heat removal. As was shown in Sec. 3.1, with

the considered heat leak the excess pressure p0 − patm will increase by a factor of two in just 2 days

in the absence of venting.

In microgravity, venting becomes and issue, since the location of the ullage space in the tank

is not known and, therefore, it is not possible to guarantee that vapor, not liquid, is vented in

the process. Note that venting LH2 caused the vehicle to tumble out of control during the AC4 test flight in 1964 [8, 9]. To circumvent this problem, the thermodynamic vent system (TVS)

technology is proposed, in which the two-phase mixture of liquid and vapor may be expanded in

a Joule-Thompson device and then passed through a heat exchanger to take away heat from the

warmer fluid in the tank, see Fig. 5 for schematics [6, 10, 31, 56]. We note that ideally the work of

TVS would result in the LH2 loss rate equal to that of a tank with the same heat budget venting the

vapor directly overboard. In practice, increased losses are inevitable due to finite TVS efficiency.

The basic physical principle of TVS operation is to control tank pressure by keeping the bulk

LH2 at subcooled conditions. Therefore, successful operation of TVS relies crucially on its ability

to efficiently move heat from the regions of nucleate boiling to the TVS heat exchanger inside the

tank by utilizing fluid mixers, such as an axial jet or spray bars. We note that, on one hand, LH2

circulation must be sufficiently strong in order for the heat through penetrations to be removed

from vapor bubbles. On the other hand, if the circulation is too strong, then the ullage bubble may

be distorted or fragmented, resulting in its possible capture by the TVS intake, see Fig. 6.

3.2.1 TVS hazard induced by ullage motion

In microgravity the ullage containing the saturated vapor and helium can drift to the TVS intake

and become captured by it. This dangerous effect is unexplored. It may cause the TVS to cease

normal operation. Bubbles containing helium may come out of the TVS, decreasing the helium

mass in the main ullage. These effects may cause generation of many new bubbles with GHe and

the formation of the multi-phase liquid-gas foam instead of the regular ullage. It maybe then be

impossible to annihilate this foam by pre-start pressurization, since the collapsing bubbles will

contain non-condensible helium gas.

To estimate the critical LH2 circulation velocity ubulk, above which the ullage motion hazard

becomes significant, we compute the flow Weber number We = 2ρLu

2

bulkRg0/σL. The flow will

begin to affect the shape of the ullage bubble when the Weber number We & 1; for example, for

ubulk = 0.5 cm/s we will have We ' 4, indicating the onset of ullage motion [19].

3.2.2 Boiling near hot spots

Let us now estimate the rate of heat removal from a bubble growing near a hot spot at the tank

wall (Fig. 5(b) and Fig. 7). A detailed treatment of this question is difficult without specifics of

the LH2 flow pattern, which depend on the tank design (see [32, 57] for recent numerical studies).

Therefore, for simplicity we will assume that LH2 flow is generated by the axial jet which creates

a counterflow pattern (see Fig. 5), in which LH2 moves with average velocity ubulk in the region

r < R0 − L, where r is the distance to the tank axis, and with average velocity ucounter in the

16
Cold LH2 with temperature TL=20.3K 

and pressure p=1.6 atm

Interface temperature Ts=22K 

Hot spot with 

power P = 5W 

!"#$%&''()$

External heat flow q0 = 0.4 W/m2

*(&+,-&+$.(/0)$$

!"#

ullage:

GH2+GHe

cold liquid 

jet

mixer

LAD

warm 

fluid

1)/0$

TVS

ubulk

L

2R0

uloc

Hot spot

$"#

!"#$$

%&''()$

Wall temperature Tw > 22K 

Wall temperature 

Tw = 22K 

Figure 5. (a) A Tank/TVS configuration with an axial jet mixer. (b) Schematics of bubble growth

on the tank wall. In (a), blue arrows indicate LH2 circulation; red arrows indicate heat influx

through the MLI; on the left, a vapor bubble is attached to a hot spot (red) at the point of a

localized heat leak through a penetration.

opposite direction in a layer of thickness L next to the tank wall. Since at the wall the flow velocity

is zero, in a laminar flow the magnitude of the LH2 velocity next to a bubble of radius R may

be estimated as uloc ' 4Rucounter/L, where we used linear interpolation from u = 0 at r = R0 to

u = 2ucounter at r = R0 − L for the axial liquid velocity u. From mass conservation we find that

uloc '

4R(R0 − L)

2

(2R0 − L)L2

ubulk. (14)

Taking L ' 1 m and the value of ubulk = 0.5 cm/s obtained in the preceding paragraph, one can

see that a localized heat source with power P = 5W may produce a steady bubble of radius R ' 13

cm. This follows from the balance of heat, taking into consideration that the Nusselt number

NuR = P/(2πRκL(Ts0 − TL0)) ' 35 for the bubble can be correlated as NuR = 0.65Pe1/2

R with the

Peclet number PeR = 2ulocRcLρL/κL ' 3000 in this case (see [53, Eq. (3.52)]), assuming small

contact angle and, hence, an almost spherical bubble, and ignoring the presence of the wall for the

flow for simplicity). Similarly, a bulk velocity ubulk = 1 mm/s would result in a steady bubble of

radius R = 20 cm. We note that even with this smaller value of the circulation speed ubulk one

would need to supply LH2 through the axial jet at the rate GLH2 = π(R0 − L)

2ρLubulk ' 0.5 kg/s,

which, for an axial jet with an orifice of 10 cm diameter would require the fluid velocity of about

1 m/s. The corresponding circulation time for the whole tank would be 8 hours in this case. At

the same time, for a given leak power the system will not be able to control the growth of bubbles

attached to hot spots whose radii are smaller than the ones estimated above.

17
Ullage: 

vapor+ 

helium 

!"#$%&

'('')*&

Hot side 

of HEX. 

Ullage: 

vapor+ 

helium 

+$,&

-#$,& +$,&

-#$,&

Figure 6. An illustration of the TVS hazard induced by ullage motion (spray bar TVS schematics

adapted from [56]).

3.2.3 Formation of a multi-phase liquid-vapor foam

The above estimates assumed that the bubble receiving power P from the heat leak remains attached to the hot spot at all times. In microgravity, bubbles will have a much reduced tendency to

detach from the wall and rise towards the ullage [21, 58–60]. Using the correlation of Fritz [17, 61],

we find that the bubble radius at departure in microgravity can be estimated at Rd ' 30 cm,

assuming, e.g., a contact angle γ ' 20◦(small contact angles are characteristic of LH2, a strongly

wetting liquid [61, 62]). Note that the Weber number for the bubble WeR = 2ρLu

2

locR/σL ' 0.005

is small in this case, and so the bubble may not be blown off easily from the hot spot by the flow.

If the bubble departure radius Rd exceeds the steady bubble radius obtained above, the bubble will

stop growing and assume its steady-state radius, and the heat will be transferred from the hot spot

through the bubble to the liquid, as desired.

The departure radius is proportional to the contact angle [17, 61] and becomes smaller than

the steady-state bubble radius estimated in the preceding paragraph for contact angles γ . 10◦.

Thus, for these smaller contact angles bubbles will be leaving the area of the hot spot and entering

into the bulk liquid. Other mechanisms of bubble departure, such as bubble coalescence, g-jitter,

and the effect of the liquid flow [17, 63, 64], which are found to be important under reduced gravity

conditions [58–60], may further reduce the bubble departure radius. As a result, bubbles may be

injected into the liquid and start to move with the flow, reaching, in particular, the flow stagnation

regions. Continuously arriving and collapsing bubbles will raise the temperature in the stagnation

regions above the subcooling temperature TL0 = 20.3K. As a result, the bubbles will stop collapsing

and will accumulate in those regions. When bubbles arrive there, a complex multi-phase liquidvapor foam-like mixture may form with temperature at saturation, Ts0 = 22K. The bubble colonies

18
liquid

warm

heat 

cold

cold

warm solid

hot

liquid

vapor

Figure 7. The schematics of the tank wall in contact with a heat leak. A vapor bubble grows by

partially absorbing heat from the hot wall.

forming in this way may further thermally insulate the tank walls from the subcooled liquid and

result in the formation of more bubbles through nucleate boiling, making TVS cooling ineffective.

3.3 Dynamics of bubble growth and collapse

We now discuss the dynamics of vapor bubbles inside LH2 in more detail.

3.3.1 Onset of bubble nucleations at hot spots

Let us begin by estimating the time needed for a single bubble to appear near a localized heat leak

with power P = 5W. We assume that the heat enters the tank through a penetration connected

to the exterior side of an aluminum LH2 tank. We note that in contrast to the Apollo design, in

which thermal insulation was on the interior tank surface [16], the MLI has to be installed on the

tank’s exterior surface, since it needs to operate in vacuum. Hence, in the absence of any special

inner surface coating, the stored LH2 will be in contact with a highly conductive metal surface.

At TL0 = 20.3K, heat conductance κw of the tank wall lies in the range of 20 – 200 W/(m·K),

depending on the composition of the aluminum alloy used [13, 65]. Also note that bubble inception

superheat will be higher than the one needed to maintain nucleate boiling and will vary depending

on the tank wall finish [21, 66–68].

Suppose first that LH2 is at the subcooled temperature TL0 when the heat leak is applied.

Because of the much higher heat conductance of aluminum, this heat will first spread into the

tank wall. If lw =

p

κwt/(cwρw) is the thermodiffusion length of aluminum, then one can roughly

estimate the temperature increase in the tank wall by equating the amount of heat P t entering

aluminum in time t to the heat content πl2

whcwρw(T − TL0) of a cylindrical section of the tank wall

with radius lw and thickness h. The timescale of temperature spreading may also be estimated by

equating lw to r0:

tw =

cwρwr

2

0

κw

' 0.3 s ÷ 3 s, (15)

where r0 ' 5 cm is the radius of the penetration. A more precise analysis (see Appendix B) gives

19
an extra logarithmic factor in the expression for the maximum temperature in the hot spot:

T ' TL0 +

P

4πhκw

ln(at/tw), (16)

where a ' 6.1. Taking h = 1 cm, we find that for t = 3tw the value of T varies in the range 21K

– 26K, depending on the heat conductance κw. Therefore, three distinct scenarios are possible,

depending on the values of κw and P in the absence of boiling.

The first scenario is realized, if the heat leak power P is sufficiently high and the heat conductance κw is sufficiently low. For the considered value of P and κw = 20 W/(m·K), which is

at the low end of the range of κw, by Eq. (16) we get T ' 24K already at t = tw = 3 s. This

corresponds to superheat of 2K, already above the nominal 1K superheat for the nucleate boiling

onset in LH2 [41–43]. In this situation, nucleate boiling will start immediately, with all the heat

going into a single vapor bubble.

The second scenario is realized when the heat conductance κw is sufficiently high. Then in

the presence of losses through convection the heat from the leak spreads to distances up to (see

Appendix C)

ls =

s

κwhR0

κLNuR0

, (17)

where NuR0is the Nusselt number associated with convection. The maximum temperature in the

hot spot can be estimated as (see Appendix C)

max T ' TL0 +

P

2πκwh

ln 

bls

r0



, (18)

where b ' 1.85. For example, with κw = 200W/(m·K), corresponding to the high end of the range

of κw, and NuR0 ' 46 estimated in Sec. 3.1 for free convection in microgravity, we obtain ls ' 1

m. About the same Nusselt number is also obtained for forced convection with average velocity

ubulk = 1 mm/s considered in Sec. 3.1, using the correlation [44, Eq. (7.23)]

NuR0 = 0.332 ×



ρLubulkR0

μL

1/2 

cLμL

κL

1/3

' 44. (19)

Then, according to Eq. (18), we get max T ' 21.8K, so the temperature in a hot spot remains

below the saturation temperature, and the heat is removed convectively from the leak, as desired.

Note that in zero gravity and absence of any boiling or mixing the value of ls in Eq. (18) should

be replaced with R0, and the value of TL0 should be replaced with the spatially averaged tank

temperature Tw0. The latter will be increasing on the timescale tsat,thermodiffusion (see Sec. 3.1),

eventually leading to nucleate boiling. On the other hand, an addition of an active mixer next to

a hot spot will further reduce the value of ls and, therefore, further suppress boiling. This may be

a better strategy for mitigating the effect of heat leaks through MLI penetrations (compare with

the strategies discussed in Sec. 3.1.5).

The third scenario is realized when the wall heat conductance is low, and the heat leak power

is also sufficiently low, provided that convective heat transfer is negligible. The latter takes place,

e.g., in zero gravity and in the absence of active mixing. We consider this scenario in more detail

in Sec. 3.3.2.

20
3.3.2 Explosive nucleate boiling hazard

If the power P is not sufficient to initiate nucleate boiling quickly in time tw, then, according to

Eq. (16) it may take an exponentially long time for the wall temperature to reach the required

superheat, resulting in a long delay in the onset of nucleate boiling. Consider, for example, the

case of κw = 20 W/(m·K) and P = 1W. Then, by Eq. (16) a temperature T ' 24K corresponding

to a 2K superheat will only be reached at t = 1.4 hours. In this time the heat will spread to

lw ' 2 m along the tank wall and lL '

p

κLt/(cLρL) ' 2.6 cm into LH2. When nucleation

occurs, the heat stored in this superheated layer of LH2 will be used to convert a mass mH2 into

vapor, with mH2 ' πl2

wlLcLρL(T − Ts0)/qL = 1 kg for the considered parameters. The resulting

vapor volume V ' 0.5 m3 will then be violently released into the tank in a short time. This is an

example of the phenomenon of explosive nucleate boiling [67–69] that was ubiquitously observed in

the microgravity boiling experiments on board the Space Shuttle Columbia [40]. We note that the

accompanying pressure spike may present a potential hazard for the operation of the storage tank.

Moreover, many vapor bubbles may be injected into LH2 as a consequence of explosive boiling,

contributing to the formation of liquid-vapor foam.

Let us note that in the absence of mixing the size of the hot spot is limited by the length scale

(see Appendix D)

L =

κw

κL

h. (20)

Assuming that κw ' 20 W/(m·K), as before, we find that L ' 2 m. Correspondingly, the maximum

temperature in the hot spot is limited by the expression given by Eq. (18) with ls replaced by L

(see Appendix D). However, in this case the steady state superheat will extend into the liquid also

to the length L, creating a much larger mass of superheated LH2, whose explosive boiling may lead

to catastrophic consequences.

3.3.3 Bubble growth over a hot spot

Once a vapor bubble is nucleated, it will grow by drawing the heat from surrounding superheated

liquid which, in turn, receives heat from hot spots on the tank wall. We note that for small heat

fluxes considered here the dominant heat transfer mechanism will be transient heat conduction (for

a recent discussion of different growth mechanisms, see [70]).

Assuming that all the power from a localized heat leak is used to convert liquid into vapor in

a single bubble, the bubble radius as a function of temperature, or, equivalently, the growth time

for a bubble of a given radius, are given by

R =



3P t

4πρvqL

1/3

, tgrowth =

4πρvqLR3

3P

. (21)

With P = 5W, we then find

R ' 20 cm, tgrowth = 1.7 hour, R ' 2 cm, tgrowth = 6 seconds. (22)

Note the difference in the dependence of the bubble radius R on t with the classical t

1/2 dependence

obtained in [71, 72] (see also [62, 64]). This is due to the fact that the bubbles under consideration

grow near a hot spot on a thin strongly conductive tank surface. Therefore, heat flux into the

21
bubble is mediated by the high heat conductance through the tank wall (for a recent discussion of

the importance of heat conductance in the heater during nucleate boiling, see [73–75]; see also Sec.

4.2).

Note that far from the hot spots this formula remains valid, if one sets P = q0A, where A is

the wall area per single bubble. The latter is valid in the case of high heat conductance of the tank

wall, when on average each bubbles will be able to intercept all the heat entering through the area

A of the tank wall. For example, for A = L

2

, with L = 10 cm, we find that the bubbles will reach

the radius of R = 2 cm in t ' 2 hours. When neighboring bubbles at the tank wall grow large

enough to come in contact with each other, more complicated dynamics involving coalescence and

detachment from the tank wall will occur.

3.3.4 Bubble collapse and accumulation in the subcooled liquid

As a result of several possible bubble departure mechanisms, vapor bubbles may detach from the

tank wall and enter into the bulk liquid. Several scenarios are possible here, depending on the

size of the departing bubbles, the level of microgravity, and the liquid flow rate. Larger bubbles

may rise toward the area of zero gravity (the tank mid-plane for a tank in a circular LEO) under

the action of buoyancy forces. The rise time may be estimated by balancing the buoyancy force

4

3

πR3ρLg with viscous drag 4πμLRu, where u ' R0/trise is the bubble velocity relative to the liquid

(recall that R0 is the tank radius) [53], to obtain

trise =

3μLR0

ρLgR2

' 5 min, R = 2 cm, (23)

which is decreasing with the increase of bubble radius R. These bubbles may then be swept towards

the ullage bubble by the flow generated by the mixer and coalesce with it (see Fig. 5). The timescale

of this process is given by

tflow =

H0

ubulk

' 3 hours, ubulk = 1 mm/s. (24)

Note that by Eq. (22) it takes much longer for a bubble to reach the ullage than to grow to the

considered size. Hence bubbles may accumulate in the liquid and be carried along by the flow and

deposited in its stagnation areas.

Let us now estimate the bubble lifetime, assuming that upon departure it enters the region of

the subcooled LH2. Assuming first that the heat escapes the bubble via steady conduction and

equating the conductive heat flow 4πR2κL(Ts0 − TL0)/R to the heat release rate 4πR2(dR/dt)ρvqL

due to condensation, and then solving the obtained differential equation, we obtain the time for

the bubble of radius R to collapse into the subcooled liquid:

tcollapse1 =

qLρvR2

2κL(Ts0 − TL0)

. (25)

Note that this equation is asymptotically exact in the limit of vanishing subcooling, but underestimates tcollapse1 for larger subcoolings due to the condensation blocking effect [76]. Indeed, assuming

that the condensation rate is dominated by the conduction through the thermal boundary layer of

width l =

p

κLt/(cLρL) and, hence, one should replace the factor 1/R by 1/l in the expression for

22
the heat flux (a more precise analysis, also giving the coefficient of proportionality is presented in

Appendix E)

tcollapse2 =

πρ2

v

q

2

LR2

4cLρLκL(Ts0 − TL0)

2

. (26)

The blocking mechanism is not important when tcollapse1 . tcollapse2, which is the case in the

examples that follow.

From Eq. (25), one can see that a bubble of radius R = 20 cm will have a lifetime tcollapse1 ' 30

hours (tgrowth ' 1.7 hour). For smaller bubbles with, say, R = 2 cm, we have tcollapse1 ' 20 mins

(tgrowth ' 6 sec). Note that these estimates agree very well with the results of direct numerical

simulations of the full system of hydrodynamic equations describing the bubble collapse which are

presented in Sec. 5. Once again, in the considered case the bubbles are produced faster than they

are collapsing in the subcooled liquid, i.e. tgrowth  tcollapse. Note, however, that as the bubbles

collapse, they release heat into the liquid, so that the subsequent bubbles are exposed to a smaller

degree of subcooling. Therefore, as more bubbles arrive, their lifetime may steadily increase and

the liquid around them reach saturation, resulting in the formation of stable bubble colonies. These

colonies then contribute to the formation of the multiphase liquid-vapor foam.

3.3.5 Bubble collapse by pre-start pressurization

The presence of vapor bubbles in LH2 is highly undesirable for the engine restart, since those

bubbles may enter into engine feed lines and result in cavitation of the turbopumps [1]. A way to

reduce or eliminate the vapor bubbles prior to engine restart is to pressurize the tank with GHe

and, at the same time, apply a small settling thrust from ullage engines. Pressurization changes

the LH2 saturation temperature relative to the bulk liquid temperature, making vapor bubbles to

condense.

Several issues arise in the course of pre-start pressurization which may result in an incomplete

vapor bubble collapse, making the procedure inefficient. First, bubble condensation releases heat

into the bulk liquid, increasing its temperature and potentially bringing it to the new saturation

temperature and stopping further condensation. This is particularly relevant to the vapor bubble

colonies. Consider, for example, the case in which the tank originally at pressure p0 = 1.6 atm

is pressurized to a new pressure p = 2 atm. The corresponding new saturation temperature

is Ts = 23K. Now, suppose a bubble colony, which remains at the old saturation temperature

Ts0 = 22K has vapor volume fraction f. Then the total amount of heat this foam-like multiphase

fluid can absorb is

Q = (1 − f)cLρLV (Ts − Ts0), (27)

where V is the colony volume. This amount of heat, in turn, can condense only the volume

Vcond = Q/(qLρv) of the vapor. Comparing the value of Vcond with the total vapor volume fV , we

can see that by purely thermodynamic considerations all the bubbles will not be able to condense,

if

f ≥



1 +

ρvqL

cLρL(Ts − Ts0)

−1

' 0.45. (28)

23
In other words, it will not be thermodynamically possible to condense all the vapor, if the volume

fraction of vapor exceeds a critical value given by Eq. (28).

On the other hand, even if the volume fraction f of vapor is below the critical value, vapor

condensation may not occur during the time interval of pressurization due to the general slowness

of the condensation process and, in particular, due to the condensation blocking phenomenon for

larger bubbles. One can once again use Eqs. (25) and (26) in these two regimes to estimate the

collapse time for bubbles of different size, provided that Ts0 − TL0 is replaced with Ts − Ts0. Using

these formulas, we now find that

tcollapse1 = 30 mins, R = 2 cm, tcollapse1 = 50 hours, R = 20 cm. (29)

It is clear that larger bubbles may not be eliminated by pressurization in a reasonable time. However, by a settling acceleration g ∼ 10−4

g0 it is possible to move larger bubbles toward the ullage.

Using Eq. (23) with R0 replaced by H0, we find that bubbles of radii R ≥ 5 mm will be able to

move towards the ullage in time trise = 5 mins. At the same time, for smaller bubbles the collapse

time is bounded above by tcollapse1 ' 2 mins, assuming the best case scenario given by Eq. (25).

This indicates that there are rather tight constraints for achieving the desired result from pre-start

pressurization. Also note that during pre-start pressurization bubbles will be continuously generated at the hot spots. Provided their departure radius is below 5 mm, these bubbles will not be

eliminated at the moment of engine start. In addition, larger bubbles forming inside the engine’s

start box, if any, will be trapped by the capillary screens of the liquid acquisition device (LAD)

and will not be able to rise to the ullage. Finally, helium-filled bubbles forming as a result of a

possible ullage capture by the TVS intake cannot be eliminated by pre-start pressurization.

4 Mathematical models of bubble dynamics

Here we present a discussion of more detailed mathematical models that can be used for the analysis

of heat transfer mechanisms and bubble growth and interaction in microgravity under heat loads

relevant for long-term cryogenic storage.

4.1 Conductance-based models

Consider a situation in which bubbles form due to boil-off at the tank wall in the vicinity of a

heat leak in a quiescent saturated liquid (Fig 7). In the absence of vapor bubbles the heat from

the localized heat leak will enter into the tank wall and rapidly spread along the wall. When the

heat spreads sufficiently far from the source, it will then start entering into the liquid inside the

tank. Mathematically, this can be described by the following equations (for simplicity, we neglect

the temperature dependence of all the material parameters and assume that the liquid occupies a

semi-infinite space next to the flat tank wall):

cwρw

∂T

∂t = κw



∂

2T

∂x2

+

∂

2T

∂y2



+

κL

h

∂T

∂z +

1

h

q(x, y, t), z = 0, (30)

cLρL

∂T

∂t = κL



∂

2T

∂x2

+

∂

2T

∂y2

+

∂

2T

∂z2



, z > 0, (31)

where the first equation is the heat conductance equation in the wall, averaged over the wall

thickness h (see [77] for a similar treatment), assuming the wall coincides with the plane z = 0, and

24
the second equation is the equation of heat conductance in the liquid. The two equations above

are coupled via the boundary condition, namely, via continuity of temperature T at z = 0, and

the derivative ∂T /∂z is evaluated as z → 0

+. The quantity q appearing in Eq. (31) represents the

heat flux from a localized heat source. In writing the above equations we neglected the possible

convective terms (including thermocapillary forces [54, 55, 58, 59, 78]), which are expected to be

quite small under the considered heat loads. For example, in the case of a bubble of radius R = 2

cm and heat penetration length l ' 2R, the associated Rayleigh number is

Ra =

gβL(Ts0 − TL0)cLρ

2

L

l

3

μLκL

' 900, (32)

which is not sufficient for the onset of microconvection [79].

4.2 Bubble growth as function of bubble location

Vapor bubbles play the role of heat sinks for the heat conductance problem in Eqs. (30) and (31),

and will, therefore, screen the heat coming from the wall into the bulk liquid (for a recent review of

single bubble heat transfer mechanisms, see [70]). To incorporate the bubbles, we need to introduce

the boundary condition T = Ts(p), where Ts is the saturation temperature at the ambient pressure

p, at the boundary of each bubble. As a consequence, heat will flow in or out of each bubble, making

them grow or shrink, respectively. Under the considered conditions of very slow dynamics (on the

time scale of hours or greater), the growth of a vapor bubble will be thermally-limited [71, 72],

with transient conduction dominating other possible heat transfer mechanisms [70]. For a spherical

bubble (assuming zero contact angle for simplicity) this leads to

4πρvqLR

2

dR

dt = Q, Q = κL

Z

∂B

∂T

∂ν dA, (33)

where ∂B denotes the boundary of the bubble, ∂/∂ν is the derivative in the direction of the outward

normal to ∂B. Moreover, since the underlying conductance-limited bubble growth dynamics will

be slow, it is reasonable to adopt a quasi steady-state approximation in solving Eqs. (30) and

(31). This approximation will be valid when the time scale of thermal diffusion exceeds that of

conductance-mediated growth. Estimating these timescales:

tthermodiffusion =

cLρLR2

κL

, tconductance =

qLρvR2

κL(T − Ts)

, (34)

one can see that the condition tconductance & tthermodiffusion is roughly satisfied for superheats under

1K in LH2.

Consider now an isolated bubble with the center r0 which is at least distance d > R away from

the wall. It is then possible to use the method of images [80] to express the solution of the obtained

elliptic boundary value problem in the form of an infinite series:

T(r) = Ts(p) + 1

2πκL

X∞

n=1

Z

wall

qn(r

0

)

|r − rn(r

0)|

d

2

r

0

. (35)

25
Here

q1(r

0

) =  ̃q(r

0

), r1 = r

0

, q ̃ = −κL

∂T

∂z









z=0

, (36)

qn+1(r

0

) = −

Rqn(r

0

)

|rn(r

0) − r0|

, rn+1(r

0

) = r0 +

R2(rn(r

0

) − r0)

|rn(r

0) − r0|2

, n odd, (37)

qn+1(r

0

) = qn(r

0

), rn+1(r

0

) = Prn(r

0

), n even, (38)

where P is the mirror reflection with respect to the xy-plane. In view of the fact that qn+2/qn ≤

R/d < 1, the series in Eq. (35) converges exponentially fast. Furthermore, by Gauss theorem, the

total heat flow Q into the bubble is given explicitly by

Q = −2

X∞

n=1

Z

wall

q ̃2n(r) d

2

r. (39)

This expression may then be substituted into Eq. (33), which gives an equation for the bubble

radius in terms of the flux  ̃q from the wall to the liquid. This flux, in turn, can be found by solving

the integral equation obtained by plugging the wall temperature, obtained from Eq. (35), into Eq.

(30) with the time derivative set to zero.

To get a sense of the phenomena described by the equations obtained above, let us consider a

situation, in which a large bubble is suspended in the liquid at a distance d from the wall, such

that R  d  L, where L is defined in (20), directly over a point heat source q(x, y) = P δ(x)δ(y).

Then it is possible to retain only the first term in the series in Eq. (35). This allows one to estimate

Q as follows

Q ∼

R

d

qd ̃

2 ∼

P Rd

L2

, (40)

where we took into account that  ̃q ∼ P/L2 and that only the region of size of order d will give a

significant contribution to the integral in (39). According to Eqs. (33) and (20), this implies that

the bubble will grow on the time scale t ∼ ρvqLR2κ

2

wh

2/(P dκ2

L

). Taking R ∼ d ∼ 1 cm and P ∼ 1

W, we can estimate t ∼ 10 hours for κw = 20 W/(m·K). On the other hand, if d → R, the entire flux

P will be absorbed by the bubble, giving a much faster time scale of growth t = 4πρvqLR3/(3P)

(i.e. Eq. (21)), and for R = 1 cm we now have t ∼ 4 s. Thus, the timescales associated with the

bubble dynamics depend sensitively on the bubble location.

4.3 Interaction of bubbles

From the arguments leading to Eq. (40), it is clear that many bubbles may need to be created to

take away the heat entering the tank from a localized heat leak. In principle, it is possible to use the

approach of Sec. 4.2 to incorporate the effect of multiple bubbles on heat transfer. However, with

the number of bubbles rapidly increasing, the problem of calculating the series representation of the

temperature distribution quickly becomes intractable. Instead, a homogenization approach which

reduces the bubble-filled liquid to an effective medium is advantageous (for a general reference, see

e.g. [81]). Note that a related class of problems was recently treated by homogenization techniques

in [82].

26
a) b) cold heat heat

warm

cold

warm

Figure 8. Multiple bubbles attached to the wall (a) or inside the liquid (b).

There are two main situations that need to be considered (see Fig. 8). In the first, it is assumed

that only bubbles attached to the wall are present (Fig. 8(a)). In this case each bubble will work

as a sink absorbing a heat flow QR = CRκwh(T − Ts)/ ln(l/R), where l  R is the characteristic

distance between bubbles and CR is a constant that needs to be obtained from the solution of

the homogenization cell problem. Introducing the distribution function f2(R, r, t) which gives the

number of bubbles of radius between R and R + dR in the wall area element d

2

r around r at time

t, we then find that the averaged heat flow into the bubble per unit wall area is

Q ̄(r, t) = K2(r, t)(T(r) − Ts), K2(r, t) = κwh

Z ∞

0

CR

ln(l/R)

f2(R, r, t)dR, (41)

where K2 is the effective (homogenized) heat transfer coefficient in two dimensions. Substituting this formula into Eq. (30), we now obtain the respective homogenized equation for the wall

temperature

cwρw

∂T

∂t = κw



∂

2T

∂x2

+

∂

2T

∂y2



−

K2(r, t)

h

(T − Ts) + 1

h

q(x, y, t). (42)

If all the bubbles have the same radius R, and the spatial density of bubbles is constant equal to

c2 = 1/l2, then the coefficient K2 becomes

K2 =

2c2κwhCR

| ln(c2R2)|

. (43)

From this equation, one can obtain the screening length λ2 ∼

p

hκw/K2 ∼ l ln1/2(l/R). This will

be the length scale of spreading of heat from a heat leak in the presence of bubbles attached to the

wall.

In the second scenario one needs to consider bubbles suspended in the liquid. The homogenization problem in this situation (with fixed bubbles) was solved in Ref. [83]. The resulting

homogenized version of Eq. (31) becomes

cLρL

∂T

∂t = κL



∂

2T

∂x2

+

∂

2T

∂y2

+

∂

2T

∂z2



− K3(r, t)(T − Ts), (44)

27
where, introducing the distribution function f3(R, r, t) of bubbles of radius between R and R + dR

in the volume element d

3

r, we get

K3(r, t) = 4πκL

Z ∞

0

Rf3(R, r, t)dR. (45)

Once again, if c3 = 1/l3is the constant density of bubbles inside the liquid, we find

K3 = 4πκLRc3, (46)

and the respective screening length is λ3 ∼

p

κL/K3 ∼ l

p

l/R.

One can further use the expressions obtained above to find the effective heat transfer coefficient

K2 appearing in Eq. (42). Assuming that T varies on a much longer spatial scale than λ3, we can

reduce Eq. (44) to a one-dimensional boundary value problem (also taking advantage of the quasi

steady-state approximation). As a result, the solution for T (homogenized) in the liquid may be

written as

T(x, y, z) ' Ts + (T(x, y, 0) − Ts)e

−z/λ3

, λ3 =



4π

Z ∞

0

Rf3(R)dR−1/2, (47)

where we assumed for simplicity that the bubble distribution is constant in space and time. Plugging

this expression into Eq. (30), we obtain Eq. (42) with

K2 =

κL

λ3

, (48)

giving the screening length in the plane λ2 ∼

p

hλ3κw/κL =

√

Lλ3. The problem needs to be

studied further when both mechanisms of heat transfer act simultaneously.

4.4 Ostwald ripening

Let us now discuss a possibility of Ostwald ripening in the suspension of bubbles inside the liquid or

a collection of bubbles attached to the wall [66]. In microgravity, an intriguing new mechanism may

result in Lifshitz-Slyozov-type dynamics of vapor bubbles, even in the absence of mass transport

between bubbles. Instead, the mechanism can be mediated via temperature, as we demonstrate

below. Due to surface tension, the vapor pressure inside a bubble of radius R is equal to p =

p0 + 2σL/R [84]. At the same time, in order to achieve evaporation-condensation equilibrium, the

vapor pressure must be equal to the saturation pressure at temperature on the bubble wall. This

leads to a peculiar feedback mechanism: The temperature on the surface of a smaller bubble with

radius R1 is higher than that on the surfaces of a larger bubble with radius R2:

T1 > T2 if R1 < R2, T1,2 = Ts +

2TsσL

qLρvR1,2

, (49)

where we used the Clapeyron-Clausius relation [84]. As a consequence, heat will flow from the

smaller bubble to the larger one, resulting in condensation in the smaller bubble and evaporation

in the larger bubble under suitable conditions. It is possible to show that the resulting problem

is mathematically equivalent to the one describing Ostwald ripening during solute precipitation in

a supersaturated solution and can, therefore, be described within the Lifshitz-Slyozov theory [85].

28
Estimating the time scale of ripening by balancing the amount of heat 4πR2ρvqLdR needed to

increase the bubble radius by dR with the inflow of heat 4πR2κLTsσL/(qLρvR2)dt in time dt, we

obtain

tripening1 '

ρ

2

v

q

2

LR3

κLσLTs

, (50)

and the usual Lifshitz-Slyozov R ∼ t

1/3 growth law. For LH2, we find that this process is quite

slow in the case of bubbles suspended in the liquid, with tripening1 ∼2.5 days for bubbles of average

radius R ∼ 1 mm. At the same time, if the bubbles are attached to the highly conducting tank wall,

then the heat conductance between the bubbles is enhanced, with the heat flow into the bubble of

order 2πRhκwTsσL/(qLρvR2)dt, resulting in an estimate

tripening2 '

ρ

2

v

q

2

LR4

hκwσLTs

, (51)

and an unusual dependence R ∼ t

1/4

. Once again, using h = 1 cm and κw = 200 W/(m·K) for

highly conductive tank wall material, we obtain t ' 20 days for R = 2 cm. Therefore, the latter

ripening mechanism may play an important role in selecting the bubble sizes.

4.5 Motion of bubbles suspended in the liquid

Let us now discuss the motion of bubbles under the action of microgravity. Let g(r, t) denote the

apparent acceleration of gravity in the tank (which includes forces of inertia). Note that both the

direction and the magnitude of the vector g may vary in space and time. Since the magnitude of

g is very small, the movement of bubbles will be dominated by viscous drag. An isolated bubble

away from the walls will, therefore, move with velocity [53, Eq. (3.15)]

v = u −

ρLR2

3μL

g, (52)

where u is the background liquid velocity, and we neglected the vapor density compared to ρL.

Note that when |g| ∼ 10−5 m/s2 and R ∼ 1 cm, the bubble would move with speed |v| ∼ 2 mm/s

in LH2, comparable to the characteristic fluid circulation velocities considered in Sec. 3.2.2. In the

case of many hydrodynamically non-interacting bubbles, one can use this equation for v, together

with Eq. (33), to write the continuity equation for the bubble density f3(R, r, t):

∂f3

∂t + ∇ · u −

ρLR2

3μL

g



f3



+

∂

∂R 

Q

4πρvqLR2

f3



= D∇2f3. (53)

The heat flow Q, in turn, depends on R and the temperature distribution in the liquid, obtained,

e.g. by solving the homogenized Eq. (44). In writing Eq. (53) we also introduced an effective

diffusion constant D depending on R which possibly arises due to turbulence or vibrations (the

g-jitter is not included in g). Finally, note that the bubbles, in turn, will exert a body force onto

the liquid and may, therefore, affect the flow.

29
4.6 Bubble creation and departure from the wall

Equation (53) needs to be supplemented with boundary conditions at the tank wall. Let us note

that because of microgravity, the bubbles that form at the surface by heterogeneous nucleations

will not easily detach from the wall. One source of such detachments are vibrations. At the same

time, another important mechanism was recently identified in the experiments conducted in the

low gravity environment of the NASA’s KC-135 aircraft [58–60]. In this mechanism, two bubbles

have to be sufficiently close to each other, so that at a certain moment their coalescence occurs.

The energy released in this process may be sufficient to break the resulting large bubble free from

the wall [60]. Taking these two effects into account, we can write the following boundary condition

for f3 in Eq. (53):

−D

∂f3

∂z









wall

= k1(R)f2(R, r, t)

+

Z ∞

0

Z ∞

0

k2(R, R0, R00)f2(R

0

, r, t)f2(R

00

, r, t)dR0dR00, (54)

where k1 and k2 are the corresponding single and double bubble detachment rates, and f2 is the

distribution function of the bubbles attached to the tank wall. The latter, in turn, is governed by

its own balance equation which contains growth terms and sink terms due to bubble departure:

∂f2

∂t +

∂

∂R 

Q

4πρvqLR2

f2



= Jδ(R − R0)



1 − A0

Z ∞

0

f2(R, r, t)dR

−k1(R)f2(R, r, t) −

Z ∞

0

Z ∞

0

k2(R, R0, R00)f2(R

0

, r, t)f2(R

00

, r, t)dR0dR00, (55)

where J is the nucleation rate per unit area, which depends on the local temperature T(r, t), A0

is the area per nucleation site, and R0 is the radius of newly nucleated bubbles. The first term in

the right-hand side of Eq. (55) states that all empty nucleation sites will be activated with rate J.

5 Computational analysis

We now illustrate how the basic physics analysis discussed in the preceding sections can be supplemented by high-fidelity computational studies of the underlying fluid dynamics equations in order

to obtain accurate information about the bubble dynamics. For definiteness, we consider the problem of vapor bubble collapse in subcooled LH2 upon bubble departure from a localized hot spot.

To simplify matters, we assume that upon departure a spherical bubble is carried along with the

flow of the subcooled liquid, where it shrinks uniformly by vapor condensation.

Let us begin with the governing equations for the liquid phase in the vicinity of the bubble.

The motion of the liquid phase can be described by incompressible Navier-Stokes equation, written

under the assumption of spherical symmetry [79]:

∂u

∂r +

2u

r

= 0, ρL



∂u

∂t + u

∂u

∂r 

= −

∂p

∂r + μL



∂

2u

∂r2

+

2

r

∂u

∂r −

2u

r

2



, r > R, (56)

where r is the radial coordinate, R = R(t) is the bubble radius, u = u(r, t) is the radial liquid

velocity for r > R, and p = p(r, t) is the pressure inside the liquid. The liquid velocity u vanishes

30
at infinity. Similarly, the thermal diffusion equation in the liquid has the form [79]

cLρL



∂T

∂t + u

∂T

∂r 

= κL



∂

2T

∂r2

+

2

r

∂T

∂r 

, r > R, (57)

where T = T(r, t) is the liquid temperature for r > R, and we neglected the heat generation due

to viscous dissipation. From the conservation of mass, the corresponding boundary condition for

the velocity u at r = R are

j(t) = ρL(u(R + 0, t) − R ̇(t)), (58)

where j is the condensation mass flux at the bubble interface, and here and below the dot denotes

time derivative.

Since the considered bubble dynamics is slow compared to the timescale of sound propagation

in the vapor phase, with very good accuracy the pressure pv inside the bubble is independent of

space (a similar situation is realized in models of gaseous combustion [86]) and can be obtained

from the condition of mechanical equilibrium at the interface:

pv(t) = p(R, t) − 2μL

∂u(R + 0, t)

∂r +

2σL

R(t)

− j

2

(t)(ρ

−1

v

(R, t) − ρ

−1

L

), (59)

where the second term accounts for viscous forces [53], the third term is due to a pressure jump

due to surface tension [53], and the fourth term is the recoil force due to evaporation/condensation

flux [75, 87]. With this assumption in mind, let us now turn to the governing equations for the gas

phase. The conservation of vapor mass is described by [79]

∂ρv

∂t +

∂(ρvu)

∂r +

2ρvu

r

= 0, r < R, (60)

where ρv = ρv(r, t) is the vapor mass density and u = u(r, t) is the radial vapor velocity for r < R.

Furthermore, neglecting the contributions of the gas velocity and treating the vapor as an ideal

gas, we can write the energy conservation in the form [79]

∂(cvρvT)

∂t +

∂(cpρvuT)

∂r +

2cpρvuT

r

= κv



∂

2T

∂r2

+

2

r

∂T

∂r 

, r < R, (61)

where T = T(r, t) is the vapor temperature for r < R and cv is the specific heat of vapor at constant

volume. Note that the ideal gas equation of state used in Eq. (61):

pv(t) = ρv(r, t)RvT(r, t), r < R, (62)

relates the vapor density ρv and temperature T at every point inside the bubble. The continuity

of mass and energy at the interface read

j(t) = ρv(R, t)(u(R − 0, t) − R ̇), (63)

and

qLj(t) = κv

∂T(R − 0, t)

∂r − κL

∂T(R + 0, t)

∂r , T(R − 0, t) = T(R + 0, t), (64)

31
respectively. Finally, to close the obtained system of equations, we relate the condensation flux

j with the conditions at the interface via the Hertz-Knudsen relation [88, 89], using a simplified

equation for saturated vapor pressure ps as a function of the temperature at the liquid-vapor

interface [90]:

j(t) = pv(t) − ps(T(R, t))

p

2πRvT(R, t)

, ps(T) = pc



T

Tc

λ

, (65)

where Rv = RH2, pc = 1.315 × 106 Pa and Tc = 33.2K are the critical pressure and temperature,

respectively, and λ = 5 for hydrogen.

The continuity equation in both the gas and the liquid phase allows one to solve explicitly for

the fluid velocity u at every point in space. For the liquid phase, which is incompressible, the

solution is well-known:

u(r, t) = a(t)

r

2

, r > R, (66)

where a(t) is some function of time only. Then, integrating Eq. (56) over r, with the help of Eq.

(66) we obtain

p(R, t) = p0 + ρL



a ̇

R

−

a

2

2R4



. (67)

Also, substituting Eq. (66) into Eq. (57), we have

cLρL



∂T

∂t +

a

r

2

∂T

∂r 

= κL



∂

2T

∂r2

+

2

r

∂T

∂r 

, r > R. (68)

In the gas phase, we can use the assumption that the vapor pressure pv is space-independent

to obtain a similar equation for u. Indeed, from Eqs. (60) and (62), we find that

∂u

∂r =

1

T

∂T

∂t +

u

T

∂T

∂r −

p ̇v

pv

−

2u

r

, r < R. (69)

Substituting this expression into Eq. (61) and using Eq. (62) yields

cpρv



∂T

∂t + u

∂T

∂r 

= κv



∂

2T

∂r2

+

2

r

∂T

∂r 

+  ̇pv. (70)

Now, let us use Eqs. (69) and (62) in Eq. (60), multiply the obtained equation by r

2 and integrate

over r. After some algebra, we obtain a formula for u:

u(r, t) = 1

3γpv(t)



3(γ − 1)κv

∂T(r, t)

∂r − rp ̇v(t)



, r < R, (71)

where γ = cp/cv. Combining Eq. (71) with (70) then yields

cppv

RvT

∂T

∂t +

κv

T



∂T

∂r 2

−

rp ̇v

3(γ − 1)T

∂T

∂r = κv



∂

2T

∂r2

+

2

r

∂T

∂r 

+  ̇pv. (72)

32
a) b)

0 200 400 600 800 1000 1200

0.000

0.005

0.010

0.015

0.020

t, s

R, m

0 200 400 600 800 1000 1200

0

1

2

3

4

t, s

j, 10

-5 kgHm

2

sL

Figure 9. The dynamics of bubble collapse in the subcooled liquid obtained by solving the governing

fluid dynamics equations: (a) the bubble radius R and (b) the condensation flux j, as functions of

time.

We now substitute the obtained solutions for u into the boundary conditions at r = R. First of

all, from Eqs. (58) and (66), we find that

R ̇ =

a

R2

−

j

ρL

. (73)

On the other hand, from Eqs. (62), (63) and (71) we obtain

p ̇v =

3

R



j



γpv

ρL

− (γ − 1)cpT(R, t)



−

γpva

R2

+ (γ − 1)κv

∂T(R, t)

∂r 

, (74)

where j is related to other variables via Eq. (65). Now, using Eq. (59) together with Eqs. (66)

and (67), we also have

a ̇ =

a

2

2R3

+

R

ρL



pv − p0 −

4μLa

R3

−

2σL

R

+ j

2



(γ − 1)cpT(R, t)

γpv

−

1

ρL

 . (75)

This equation closes the system of equations describing T(r, t), R(t), a(t) and pv(t). Thus, in order

to obtain the full solution to the bubble dynamics, one needs to solve Eqs. (64), (65), (68), (72) –

(75). For the problem of bubble collapse, the following initial conditions may be used:

R(0) = R0, pv(0) = p0, a(0) = 0, (76)

T(r, 0) = Ts0, r < R0, T(r, 0) = TL0 +

(Ts0 − TL0)R0

r

, r > R0. (77)

The results of the numerical solution of the above equations under a simplifying assumption

that pv(t) = ps(T(R(t), t)), justified in the considered situation [76], and with R0 = 2 cm are

presented in Fig. 9. As expected, the bubble shrinks to zero radius in finite time due to vapor

condensation on the liquid-vapor interface. The heat released from the condensation is carried

away into the subcooled liquid, and the temperature inside the bubble stays essentially constant,

see Fig. 10. One can also see that the bubble collapse time is in very good agreement with the

estimate obtained in Sec. 3.3.4.

33
!"#$"%&'(%")*+,

-.

--

/.01

.

/20.

.

/202

.

32./4 32./4 32./4

5"6789':)*#;<,

'=. '=>..< '=?..< '=1..<

*@,) '=>..< '=?..< '=1..<

2)8#

2)8#

*&,)

Figure 10. Temperature (a) and velocity (b) distributions in and around the collapsing bubble at

times t = 0, 300, 600, 900 s.

6 Summary and conclusions

In conclusion, we have presented an overview of the physics-related issues associated with largescale, long-term storage of LH2 in microgravity. A special difficulty associated with LH2 storage, as

opposed to other cryogenic liquids, such as liquid oxygen or liquid methane, is its very low boiling

point (∼20K). Currently, no efficient cryocoolers exist that can operate at LH2 temperatures. This

presently puts the promising ideas of ZBOT technologies out of reach in the case of LH2. For the

same reason, other active cooling technologies, such as BAC, are also limited to the “tube-to-shield”

concept, which can only partially address the problem of LH2 heat load reduction. Importantly, in

the case of LH2 storage heat leaks through MLI penetration may become dominant challenges for

enabling long-term storage in microgravity.

Below we first list the main physics-related issues of LH2 storage and then discuss potential

challenges identified by our analysis.

In microgravity, the location, structure, and shape of the ullage space is generally unknown.

The dominant forces determining the ullage bubble shape are capillary forces. The location of the

ullage bubble is determined by a combination of the g-forces and the forces exerted on the ullage

bubble by the liquid stirring flow. The ullage bubble may drift inside the tank due to the generally

time-varying character of the g-forces.

Since the ullage bubble position is unknown and may be variable, boil-off vapor cannot be vented

directly overboard. A combination of TVS and a fluid mixer may be used to control the boil-off

pressure by maintaining the liquid in a subcooled state. This requires an addition of non-condensible

cold pressurant gas to avoid ullage bubble collapse. Further heat load reduction strategies, both

34
active (BAC applied to an intermediate layer of the MLI) and passive (VCS utilizing the boil-off

vapor from TVS to cool the MLI and penetrations) may also be necessary.

The stirring flow velocity of the liquid inside the tank is subject to competing constraints.

Stirring has to be sufficiently strong to ensure efficient heat transfer from the tank walls through

the liquid into the TVS heat exchanger. Stirring has to be sufficiently weak to avoid deforming or

breaking the ullage bubble. Excessive stirring may result in the injection of unwanted energy into

the liquid.

Ullage motion driven by a strong stirring flow may result in a TVS hazard, whereby the ullage

bubble may be captured by the TVS intake and block normal operation of TVS. As a consequence,

TVS may fail to control the tank pressure, or bubbles filled with non-condensible pressurant gas

may come out of the TVS recirculation outlet.

This implementation faces a number of challenges related to fundamental physics.

1. On the basis of correlations derived under earth gravity, it is conceivable that both natural convection in microgravity and stirring may keep the tank wall temperature below the

saturation point away from the possible temperature hot-spots due to localized heat leaks

through MLI penetrations. Efficient heat removal from the penetrations may require highly

conductive tank wall material of sufficient thickness and vapor cooling. On-orbit experiments

and further physics-based modeling is necessary to establish new reliable correlations which

are needed to design a proper strategy for tank temperature control.

2. Unpredictable changes of the flow conditions may result in overheating of the liquid near the

tank wall. For high wall conductivities this may lead to the explosive boiling hazard, whereby

a large mass of liquid becomes superheated in a rather large area around a temperature hotspot on the tank wall or near flow stagnation points after a long delay time. The ensuing

fast nucleate boiling may create a pressure jump and a strong liquid flow affecting the ullage

bubble.

3. Single large vapor bubbles attached to the tank walls or many small bubbles injected into

the liquid may form at the points of strong localized heat leaks. The choice of the scenario

depends sensitively on the mechanism of bubble departure from the tank wall during nucleate

boiling near hot-spots. Key factors determining this mechanism are the liquid contact angle,

the tank surface properties, microgravity magnitude and direction, g-jitter, Ostwald ripening

and bubble merging.

4. The bubbles inside the liquid may drift with the stirring flow and accumulate in its stagnation

areas, resulting in the formation of complex multi-phase liquid-vapor foam-like structures.

The foam may further thermally insulate parts of the tank wall from the subcooled liquid,

increasing the boiling rate in those regions. The associated pressure increase cannot be well

controlled by TVS.

5. Boil-off bubbles may not be easy to remove by an overpressurization, since the collapsing

bubbles release heat into the surrounding liquid and raise its temperature to saturation. In

the presence of TVS hazard by the ullage motion, overpressurization may be ineffective in

the presence of non-condensible pressurant gas inside the vapor bubbles.

6. The multi-phase liquid-vapor foam may grow at the expense of shrinking of the ullage bubble,

when the vapor mass is transferred from the ullage into the foam. This phenomenon is

35
complicated by the pressurant gas dissolution hazard, which is further amplified by strong

stirring. The dissolved non-condensible pressurant gas may enter into the vapor bubbles of

the foam. As a result, these bubbles may be difficult or impossible to eliminate by applying

overpressurization.

References

1. G. P. Sutton and O. Biblarz. Rocket Propulsion Elements. John Wiley and Sons, Inc., Hoboken,

NJ, 2001.

2. D. J. Chato. Cryogenic fluid transfer for exploration. Cryogenics, 48:206–209, 2008.

3. D. Plachta and P. Kittel. An updated zero boil-off cryogenic propellant storage analysis applied

to upper stages or depots in an LEO environment. Technical Report TM—2003-211691, NASA,

2003.

4. B. F. Kutter, F. Zegler, G. O’Neil, and B. Pitchford. A practical, affordable cryogenic propellant

depot based on ULA’s flight experience. In AIAA SPACE 2008 Conference and Exposition,

2008.

5. J. A. Goff, B. F. Kutter, F. Zegler, D. Bienhoff, F. Chandler, and J. Marchetta. Realistic nearterm propellant depots: Implementation of a critical spacefaring capability. In AIAA SPACE

2009 Conference and Exposition, 2009.

6. S. M. Motil, M. L. Meyer, and S. P. Tucker. Cryogenic fluid management technologies for

advanced green propulsion systems. Technical Report TM-2007-214810, NASA, 2007.

7. L. J. Salerno and P. Kittel. Cryogenics and the human exploration of Mars. Cryogenics,

39:381–388, 1999.

8. R. F. Lacovic, F. C. Yeh, S. V. Szabo, R. J. Brun, A. J. Stofan, and J. A. Berns. Management

of cryogenic propellants in a full-scale orbiting space vehicle. Technical Report TN D-4571,

NASA, 1968.

9. D. Glover. NASA cryogenic fluid management space experiment efforts 1960-1990. Technical

Report TM-103752, NASA, 1991.

10. F. T. Dodge. The New “Dynamic Behavior of Liquids in Moving Containers”. Southwest

Research Institute, San Antonio, TX, 2000.

11. D. J. Chato. Low gravity issues of cryogenic fluid management technologies enabling exploration, 2009.

12. M. P. Doherty, J. D. Gaby, L. J. Salerno, and S. G. Sutherlin. Cryogenic fluid management

technology for moon and mars missions. Technical Report TM-2010-216070, NASA, 2010.

13. M. Donabedian, editor. Spacecraft thermal control handbook. Volume 2, Cryogenics. The

Aerospace Press, El Segundo, CA, 2003.

36
14. D. J. Chato. The role of flight experiments in the development of cryogenic fluid management

technologies. Cryogenics, 46:82–88, 2006.

15. NASA Technical Report HSM-R421-67. Evaluation of AS-203 Low Gravity Orbital Experiment,

1967.

16. NASA-TM-X-881. Apollo Systems Description. Volume II. Saturn Launch Vehicles, 1964.

17. V. K. Dhir. Mechanistic prediction of nucleate boiling heat transfer–achievable or a hopeless

task? J. Heat Transfer, 128:1–12, 2006.

18. E. I. Nesis. Boiling of liquids. Sov. Phys. – Uspekhi, 8:883–907, 1965.

19. C. E. Brennen. Fundamentals of Multiphase Flows. Cambridge University Press, 2005.

20. N. I. Kolev. Multiphase flow dynamics, volume 2. Springer-Verlag, Berlin, 2nd edition, 2005.

21. V. K. Dhir. Boiling under microgravity conditions. In Proceedings of 12th International Heat

Transfer Conference, Grenoble, France, 2002.

22. J. Kim. Review of reduced gravity boiling heat transfer: US research. J. Jpn. Soc. Microgravity

Appl., 20:264–271, 2003.

23. NASA’s Exploration System Architecture Study (Final Report). NASA-TM-2005-214062, 2005.

24. E. M. Lifshits and L. P. Pitaevskii. Course of Theoretical Physics, volume 8. Pergamon Press,

London, 1980.

25. NASA-TM-X-68725. Apollo Program Summary Report, 1975.

26. NASA-TM-X-69534. Saturn V Launch Vehicle Flight Evaluation Report – AS-512. Apollo 17

Mission, 1973.

27. R. W. Orloff. Apollo by the numbers: A statistical reference. Technical report, NASA SP2000-4029, 2000.

28. C. W. Keller, G. R. Cunnington, and A. P. Glassford. Thermal performance of multilayer

insulations. Technical Report CR-134477, NASA, 1974.

29. J. R. Feller. Some simple calculations to assess the effectiveness of passive and active boil-off

reduction techniques to EDS. NASA-ARC Cryogenics Group, 2009.

30. L.J. Hastings, D.W. Plachta, L. Salerno, and P. Kittel. An overview of NASA efforts on zero

boiloff storage of cryogenic propellants. Cryogenics, 41:833–839, 2002.

31. M. E. Moran. Cryogenic fluid storage technology development: Recent and planned efforts at

NASA. Technical Report TM—2009-215514, NASA, 2009.

32. C. Panzarella, D. Plachta, and M. Kassemi. Pressure control of large cryogenic tanks in microgravity. Cryogenics, 44:475–483, 2004.

37
33. J. Feller, L. Salerno, A. Kashani, J. Maddocks, B. Helvensteijn, G. Nellis, and Y. Gianchandani.

Distributed cooling techniques for cryogenic boil-off reduction systems. In The Proceedings of

the 15th International Cryocooler Conference, 2008.

34. E. R. Canavan, R. F. Boyle, and S. Mustafi. Active costorage of cryogenic propellants for

exploration. In Mohamed S. El-Genk, editor, Space Technology and Applications International

Forum–STAIF 2008, volume 969, pages 733–740. AIP, 2008.

35. J. F. LeBar and E. C. Cady. The advanced cryogenic evolved stage (ACES) – a low-cost,

low-risk approach to space exploration launch. In AIAA 2006-7454, 2006.

36. K. R. Hanby. Handbook on Materials for Superconducting Machinery. Battelle Columbus Labs,

Ohio Metals and Ceramics Information Center, Columbus, OH, 1974.

37. M. P. Malkov. A reference on physics and engineering fundamentals of cryogenics (in Russian).

Energatomizdat, Moscow, 1985.

38. A. Woodcraft. Predicting the thermal conductivity of aluminium alloys in the cryogenic to

room temperature range. Cryogenics, 45:421–432, 2005.

39. R.D. McCarty, J. Hord, and H.M. Roder. Selected properties of hydrogen (engineering design

data), volume 168 of Natl. Bur. Stand. Monographs. National Bureau of Standards, 1981.

40. M. M. Hasan, C. S. Lin, R. H. Knoll, and M. D. Bentz. Tank pressure control experiment:

Thermal phenomena in microgravity. Technical Report TP-3564, NASA, 1996.

41. L. Bewilogua, R. Knoner, and G. Wolf. Heat transfer in boiling hydrogen neon nitrogen and

argon. Cryogenics, 6:36–39, 1966.

42. Y. A. Kirichenko, K. V. Rusanov, and E. G. Tyurina. Heat-transfer in subcooled liquid cryogens. Cryogenics, 23:209–211, 1983.

43. L. Bewilogua, R. Kn ̈oner, and H. Vinzelberg. Heat transfer in cryogenic liquids under pressure.

Cryogenics, 15:121–125, 1975.

44. F. P. Incopera and D. P. DeWitt. Fundamentals of Heat and Mass Transfer. John Wiley and

Sons, Inc., New York, 2007.

45. J. R. Feller, D. W. Plachta, G. Mills, and C. McLean. Demonstration of a cryogenic boiloff reduction system employing an actively cooled thermal radiation shield. Technical report,

NASA, 2010.

46. F. W. Schmidt, R. E. Henderson, and C. H. Wolgemuth. Introduction to Thermal Sciences.

John Wiley and Sons, Inc., 1993.

47. V. L. Streeter, E. B. Wylie, and K. W. Bedford. Fluid Mechanics. McGraw-Hill, 1998.

48. S. E. Haaland. Simple and explicit formulas for the friction factor in turbulent flow. J. Fluids

Eng., 102:89–90, 1983.

38
49. P.J. Linstrom and W.G. Mallard, editors. Thermophysical Properties of Fluid Systems. NIST

Chemistry WebBook, NIST Standard Reference Database Number 69. National Institute of

Standards and Technology, Gaithersburg MD, 2010.

50. W. B. Streett, R. E. Sonntag, and G. J. Van Wylen. Liquid—vapor equilibrium in the system

normal hydrogen—helium. The Journal of Chemical Physics, 40:1390–1395, 1964.

51. G. Cini-Castagnoli, A. Giardini-Guidoni, and F. P. Ricci. Diffusion of neon, ht, and deuterium

in liquid hydrogen. Phys. Rev., 123:404–406, 1961.

52. W. P. A. Hass, G. Seidel, and N. J. Poulis. Nuclear spin relaxation and molecular diffusion in

liquid hydrogen. Physica, 26:834 – 852, 1960.

53. R. Clift, J. R. Grace, and M. E. Weber. Bubbles, Drops, and Particles. Academic Press, New

York, 1978.

54. R. Marek and J. Straub. The origin of thermocapillary convection in subcooled nucleate pool

boiling. Int. J. Heat Mass Transfer, 44:619 – 632, 2001.

55. M. Barthes, C. Reynard, R. Santini, and L. Tadrist. Non-condensable gas influence on the

marangoni convection during a single vapour bubble growth in a subcooled liquid. Europhys.

Lett., 77:14001, 2007.

56. R. H. Flachbart, L. J. Hastings, A. Hedayat, S. L. Nelson, and S. Tucker. Testing the effect

of helium pressurant on thermodynamic vent system performance with liquid hydrogen. In

J. G. Weisend II, editor, CP985, Advances in Cryogenic Engineering: Transactions of the

Cryogenic Engineering Conference — CEC, volume 53, 2008.

57. S. H. Ho and M. M. Rahman. Three-dimensional analysis for liquid hydrogen in a cryogenic

storage tank with heat pipe-pump system. Cryogenics, 48:31–41, 2008.

58. J. Kim, J. F. Benton, and D. Wisniewski. Pool boiling heat transfer on small heaters: effect of

gravity and subcooling. Int. J. Heat Mass Transfer, 45:3919–3932, 2002.

59. C. D. Henry and J. Kim. A study of the effects of heater size, subcooling, and gravity level on

pool boiling heat transfer. Int. J. Heat Fluid Flow, 25:262–273, 2004.

60. UCLA Heat Transfer Lab. Nucleate boiling in the low gravity environment of the KC-135

aircraft. http://boiling.seas.ucla.edu/multimedia gallery.htm, 2009.

61. Y. A. Kirichenko, L. A. Slobozhanin, and N. S. Shcherbakova. Analysis of quasi-static conditions

of boiling onset and bubble departure. Cryogenics, 23:110–112, 1983.

62. W. B. Bald. Bubble-growth constants for liquid-hydrogen and liquid-helium. Cryogenics,

16:709–712, 1976.

63. A. Mukherjee and V. K. Dhir. Study of lateral merger of vapor bubbles during nucleate pool

boiling. Journal of Heat Transfer, 126:1023–1039, 2004.

39
64. S. Siedel, S. Cioulachtjian, and J. Bonjour. Experimental analysis of bubble growth, departure

and interactions during pool boiling on artificial nucleation sites. Exp. Therm. Fluid. Sci.,

32:1504–1511, 2008.

65. G. M. Zlobintsev, V. V. Kozinets, B. A. Merisov, E. M. Ogneva, V. I. Sokolov, and A. V.

Blinova. Thermal conductivity of structural aluminum alloys at low temperatures. Metal Sci.

and Heat Treat., 17:436–437, 1975.

66. J. Straub. Boiling heat transfer and bubble dynamics in microgravity. In J. P. Hartnett, editor,

Advances in Heat Transfer, volume 35. Academic Press, San Diego, CA, 2001.

67. M.-C. Duluc, B. Stutz, and M. Lallemand. Boiling incipience in liquid nitrogen induced by a

step change in heat flux. Int. J. Heat Mass Transfer, 51:1738 – 1750, 2008.

68. M. E. Bland, C. A. Bailey, and G. Davey. Boiling from metal-surfaces immersed in liquidnitrogen and liquid-hydrogen. Cryogenics, 13:651–657, 1973.

69. J. Li and G. P. Petersen. Microscale heterogeneous boiling on smooth surfaces—from bubble

nucleation to bubble dynamics. Int. J. Heat Mass Transfer, 48:4316–4332, 2005.

70. J. Kim. Review of nucleate pool boiling bubble heat transfer mechanisms. Int. J. Multiphase

Flow, 35:1067–1076, 2009.

71. M. S. Plesset and S. A. Zwick. The growth of vapor bubbles in superheated liquids. J. Appl.

Phys., 25:493–500, 1954.

72. L. E. Scriven. On the dynamics of phase growth. Chem. Eng. Sci., 10:1–13, 1959.

73. V. S. Nikolayev, D. A. Beysens, G. L. Lagier, and J. Hegseth. Growth of a dry spot under

a vapor bubble at high heat flux and high pressure. International Journal of Heat and Mass

Transfer, 44:3499 – 3511, 2001.

74. V. Nikolayev, D. Beysens, Y. Garrabos, C. Lecoutre, and D. Chatain. Bubble spreading during

the boiling crisis: modelling and experimenting in microgravity. Microgravity Sci. Technol.,

18:34–37, 2006.

75. V. S. Nikolayev. Dynamics of the triple contact line on the non-isothermal heater. Phys. Fluids,

22:082105, 2010.

76. V. V. Osipov and C. B. Muratov. Dynamic condensation blocking in cryogenic refueling. Appl.

Phys. Lett., 93:224105, 2008.

77. M. Pˇribyl, C. B. Muratov, and S. Y. Shvartsman. Discrete models of autocrine cell communication in epithelial layers. Biophys. J., 84:3624–3635, 2003.

78. Y. Abe and A. Iwasaki. Microgravity experiments on dual vapor bubbles of self-wetting fluid.

AIP Conf. Proc., 608:189–196, 2002.

79. L. D. Landau and E. M. Lifshits. Course of Theoretical Physics, volume 6. Pergamon Press,

London, 1987.

40
80. L. D. Landau and E. M. Lifshits. Course of Theoretical Physics, volume 8. Pergamon Press,

London, 1984.

81. V. A. Marchenko and E. Ya. Khruslov. Homogenization of partial differential equations, volume 46 of Progress in Mathematical Physics. Birkh ̈auser, Boston, MA, 2006.

82. C. B. Muratov and S. Y. Shvartsman. Boundary homogenization for periodic arrays of absorbers. Multiscale Model. Simul., 7:44–61, 2008.

83. G. C. Papanicolaou and S. R. S. Varadhan. Diffusion in regions with many small holes. In

Stochastic differential systems (Proc. IFIP-WG 7/1 Working Conf., Vilnius, 1978), volume 25

of Lecture Notes in Control and Information Sci., pages 190–206. Springer, Berlin, 1980.

84. L. D. Landau and E. M. Lifshits. Course of Theoretical Physics, volume 5. Pergamon Press,

London, 1980.

85. E. M. Lifshits and L. P. Pitaevskii. Physical kinetics. Pergamon Press, Oxford, 1981.

86. F. A. Williams. Combustion Theory. The Benjamin/Cummings Publishing Company, Inc.,

Menlo Park, CA, 1985.

87. Andrea Prosperetti and M. S. Plesset. The stability of an evaporating liquid surface. Physics

of Fluids, 27:1590–1602, 1984.

88. Y.-P. Pao. Application of kinetic theory to problem of evaporation and condensation. Phys.

Fluids, 14:306–312, 1971.

89. L. D. Koffman, M. S. Plesset, and L. Lees. Theory of condensation and evaporation. Phys.

Fluids, 27:876–880, 1984.

90. J. A. Clark. Universal equations for saturation vapor pressure. In IAA/ASME/SAE/ASEE

Joint Propulsion Conference and Exhibit, AIAA 2004-4088, 2004.

91. L. C. Evans. Partial Differential Equations, volume 19 of Graduate Studies in Mathematics.

American Mathematical Society, Providence, RI, 1998.

92. M. Abramowitz and I. Stegun, editors. Handbook of mathematical functions. National Bureau

of Standards, 1964.

93. C. B. Muratov. Droplet phases in non-local Ginzburg-Landau models with Coulomb repulsion

in two dimensions. Comm. Math. Phys., 299:45–87, 2010.

41
Appendix A

Analysis of heat exchange between the cooling tube and the strut

Here we derive a set of approximate differential equations that describe the distribution of

temperature in a tubular strut of material with thermal conductivity κstrut, radius Rstrut, thickness

dstrut, and length Lstrut, with a helical coil made of a GH2-carrying tube wound around (see Fig. 2).

The pitch of the tube helix is denoted by dtube and is assumed to be much smaller than Rstrut. Let

x > 0 be the coordinate along the strut and let x = 0 correspond to the point of the strut attachment

to the tank wall. Consider the heat balance in the strut on the interval [x0 −

1

2

dtube, x0 +

1

2

dtube].

Assuming that the strut temperature T(x) varies appreciably only in the x-direction, to a first

approximation the steady heat balance reads

2πRstrutdstrutκstrut

dT

dx









x=x0+ 1

2

dtube

x=x0− 1

2

dtube

= 4π

2RstrutRtubehtube(T(x0) − Tv(x0)), (A78)

where Tv(x) is the vapor temperature at position x along the strut in the tube, htube is the heat

transfer coefficient obtained in Eq. (12), and we assumed that the tube is in perfect thermal

contact with the strut. Using Eq. (11) to express htube in terms of the tube Nusselt number and

approximating the jump in the derivative of the temperature above by dtubed

2T /dx2

evaluated at

x = x0, after some algebra we obtain the following differential equation for T(x):

l

2

strut

d

2T

dx2

= T − Tv, lstrut =

s

κstrutdstrutdtube

πκvNutube

(A79)

with boundary conditions

dT

dx









x=0

= 0, T







x=Lstrut

= T0, (A80)

where we took into account that all the heat from the strut is captured by the tube and thus

does not enter into the tank at x = 0, and that the strut is thermalized with the environment at

x = Lstrut. Note that the parameter lstrut has the dimension of length and characterizes the length

scale at which the strut temperature would equilibrate to the tube temperature, if the latter were

kept constant. For the parameters of Sec. 3.1.5, we find lstrut ' 2.5 cm.

On the other hand, the heat balance for GH2 in the interval [x0 −

1

2

dtube, x0 +

1

2

dtube] reads

cpJstrutTv(x)







x=x0+ 1

2

dtube

x=x0− 1

2

dtube

= 4π

2RstrutRtubehtube(T(x0) − Tv(x0)). (A81)

Once again, approximating the jump in Tv(x) by dtubedTv/dx at x = x0, we can write this as a

differential equation

cpJstrutdtube

dTv

dx = 4π

2RstrutRtubehtube(T(x0) − Tv(x0)). (A82)

42
Then using Eq. (A78) and integrating the obtained equation, we get

Tv = TL0 + lv

dT

dx , lv =

2πRstrutdstrutκstrut

cpJstrut

, (A83)

where we took into account that GH2 enters with Tv = TL0 at x = 0. Here the parameter lv

also has a dimension of length and characterizes the length scale on which the vapor temperature

approaches that of the strut, if the latter follows that of the vapor (when cooling is efficient). For

the parameters of Sec. 3.1.5, we find lv ' 10 cm.

Finally, substituting Eq. (A83) into (A79), we obtain a single linear equation

l

2

strut

d

2T

dx2

+ lv

dT

dx − T + TL0 = 0. (A84)

Together with the boundary conditions in Eq. (A80) this simple boundary value problem admits

an exact solution, which is too cumbersome to be written down as an explicit formula. We omit

the final answer, which is plotted in Fig. 3 for various parameters.

Appendix B

Analysis of heat spreading in the tank wall

Consider a flat infinite layer of aluminum of thickness h, to which heat power P is injected

uniformly into a disk of radius r0 & h at one of the surfaces. Averaging the heat conductance

equation written in cylindrical coordinates over the film thickness, we obtain an equation

cwρw

∂T

∂t =

κw

r

∂

∂r 

r

∂T

∂r 

+

P

πr2

0h

H(r0 − r), (B85)

where H(x) is the Heaviside function. The solution T(r, t) of this equation with initial condition

T(r, 0) = TL0 can be written in terms of the two-dimensional heat kernel (see e.g. [91]). From that

solution, we find that the maximum temperature, which is equal to the value of T(0, t) is

Tmax(t) = TL0 +

P

2πr2

0hκw

Z t

0

dt0

t − t

0

Z r0

0

rdr exp(−cwρwr

2

/(4κw(t − t

0

)))

= TL0 +

P

4πcwρwκwhr2

0



4tκw



1 − e

−

cwρwr

2

0

4κwt



+ cwρwr

2

0Γ



0,

cwρwr

2

0

4κwt



= TL0 +

P

4πhκw

ln 

4κwt

e

γ−1cwρwr

2

0



+ O(t

−1

), (B86)

where Γ(a, z) is the incomplete Gamma-function and γ ' 0.5772 is the Euler constant [92]. Approximating the solution in the tank wall by Eq. (B86) and evaluating the coefficient inside the

logarithm numerically, we obtain the estimate in Eq. (16).

Appendix C

43
Effect of heat loss into the liquid on the heat conduction in the

tank wall

Consider an infinite layer of aluminum of thickness h injected with a heat leak P uniformly

through a disk of radius r0 & h, as in Appendix B, in contact with a layer of liquid of thickness

R0 undergoing free convection. If NuR0is the Nusselt number associated with convection and

κLNuR0/R0 is the corresponding heat transfer coefficient, then the equation of heat conductance

in the layer, averaged over the layer thickness (see e.g. [77, Appendix] for a similar treatment) and

written in polar coordinates is

cwρw

∂T

∂t =

κw

r

∂

∂r 

r

∂T

∂r 

−

κLNuR0(T − TL0)

hR0

+

P

πr2

0h

H(r0 − r). (C87)

In contrast to Eq. (B85), this equation admits steady state solutions as t → ∞, due to the screening

effect of the heat loss, characterized by the screening length given by Eq. (17) (similar equations

arise in the studies of Debye-H ̈uckel theory of electrolytes, in which ls is the Debye radius [84]).

Indeed, the stationary solution of Eq. (C87) is [93, Eq. (3.39)]:

T(r) = TL0 +

P

πr2

0κwh

×

(

l

2

s − r0lsK1(r0/ls)I0(r/ls), r ≤ r0,

r0lsI1(r0/ls)K0(r/ls), r ≥ r0,

(C88)

where In(x) and Kn(x) is the modified Bessel functions of order n of the first and second kind,

respectively [92], and so the solution decays exponentially at lengths of order ls. The maximum of

T(r) is achieved at r = 0, and we have

max T(r) = TL0 +

P

2πκwh

"

ln 

2 exp 

1

2 − γ



ls

r0

!

+ O



r

2

0

l

2

s

ln ls

r0

#

, (C89)

where γ ' 0.5772 is the Euler constant [92]. Evaluating the numerical factor inside the logarithm

and dropping higher order terms, we obtain Eq. (18).

Appendix D

Spreading of heat from the wall into the liquid

Consider a flat infinite layer of aluminum at z = 0 covered by a semi-infinite layer of LH2 at

z > 0. The steady heat equation for aluminum averaged over the layer thickness reads

κw



∂

2T

∂x2

+

∂

2T

∂y2



+

κL

h

∂T

∂z









z=+0

+

P

πr2

0h

H



r0 −

p

x

2 + y2



= 0, (D90)

where the last term is as before and the boundary term comes from the continuity of heat flux

from LH2 into aluminum. The corresponding heat conductance problem for LH2 is given by the

following boundary value problem:

∂

2T

∂x2

+

∂

2T

∂y2

+

∂

2T

∂z2

= 0, z > 0, (D91)

44
and continuity of T is assumed down to z = 0. Applying Fourier Transform in the x and y variables

to Eq. (D91):

Tˆ

q(z) = Z +∞

−∞

Z +∞

−∞

e

iq1x+iq2y

(T(x, y, z) − TL0) dx dy, (D92)

where q = (q1, q2), and solving the obtained equation, we get

Tˆ

q(z) = Tˆq(0)e

−qz, q =

q

q

2

1 + q

2

2

. (D93)

In particular, the Fourier Transform of the Dirichlet-to-Neumann map associated with (D91) is

∂Tˆ

q

∂z











z=+0

= −qTˆ

q(0). (D94)

Plugging this expression into the Fourier-transformed Eq. (D90) and using the expression of [93, Eq.

(B6)] for the transform of the characteristic function of a disk, after some algebra we obtain

Tˆ

q(0) = 2P J1(qr0)

r0κLq

2(Lq + 1), (D95)

where Jn(x) is the Bessel function of the first kind of order n, and we used Eq. (20). Inverting

the Fourier transform in Eq. (D95) and performing the integrations, after some algebra we find

max T(x, y, z) = T(0, 0, 0), with

max T = TL0 +

P

2κLr0



H1

r0

L



− Y1

r0

L



−

2L

πr0



, (D96)

where H1(z) and Y1(z) are the Struve function and the modified Bessel function of the second kind,

respectively [92]. Finally, expanding the expression in Eq. (D96) in a series in r0, we find that

max T = TL0 +

P

2πκwh

"

ln 

2 exp( 1

2 − γ)L

r0

!

+ O



r

2

0

L2

ln L

r0

#

, (D97)

which coincides with Eq. (C89) with ls replaced with L, to the leading order.

Appendix E

Bubble collapse in the subcooled liquid via conduction through

thermal boundary layer

Consider a bubble of radius R(t), whose interface is at saturation temperature Ts0 suspended in

a subcooled liquid with temperature TL0. In the thermal boundary layer approximation the heat

conduction near the interface may be approximated by a one-dimensional initial boundary value

problem for the liquid temperature T(r, t)

cLρL

∂T

∂t = κL

∂

2T

∂r2

, T(R, t) = Ts0, T(∞, t) = TL0, T(r, 0) = TL0. (E98)

45
The solution of this problem is

T(r, t) = TL0 + (Ts0 − TL0) erf





s

cLρL(r − R)

2

4κLt



 , (E99)

where erf(z) is the error function [92]. From this formula, we find that the heat flux through the

interface is

q(t) = −κL

∂T(R, t)

∂r =

r

cLρLκL(Ts0 − TL0)

2

πt

. (E100)

The substituting this expression to the heat balance equation 4πR2dRρvqL = −4πR2q(t)dt, whose

solution is given by Eq. (26).

46

REPORT DOCUMENTATION PAGE Form Approved

OMB No. 0704–0188

The public reporting burden for this collection of information is estimated to average 1 hour per response, including the time for reviewing instructions, searching existing data sources,

gathering and maintaining the data needed, and completing and reviewing the collection of information. Send comments regarding this burden estimate or any other aspect of this collection

of information, including suggestions for reducing this burden, to Department of Defense, Washington Headquarters Services, Directorate for Information Operations and Reports

(0704-0188), 1215 Jefferson Davis Highway, Suite 1204, Arlington, VA 22202-4302. Respondents should be aware that notwithstanding any other provision of law, no person shall be

subject to any penalty for failing to comply with a collection of information if it does not display a currently valid OMB control number.

PLEASE DO NOT RETURN YOUR FORM TO THE ABOVE ADDRESS.

Standard Form 298 (Rev. 8/98)

Prescribed by ANSI Std. Z39.18

1. REPORT DATE (DD-MM-YYYY)

01-10-2011

2. REPORT TYPE

Technical Memorandum

3. DATES COVERED (From - To)

4. TITLE AND SUBTITLE

Issues of Long-Term Cryogenic Propellant Storage in Microgravity

5a. CONTRACT NUMBER

5b. GRANT NUMBER

5c. PROGRAM ELEMENT NUMBER

5d. PROJECT NUMBER

5e. TASK NUMBER

5f. WORK UNIT NUMBER

6. AUTHOR(S)

V. B. Muratov, V. V. Osipov and V. N. Smelyanskiy

7. PERFORMING ORGANIZATION NAME(S) AND ADDRESS(ES)

NASA NASA Ames Research Center Moffett Field, CA 94035-1000

8. PERFORMING ORGANIZATION

REPORT NUMBER

L–

9. SPONSORING/MONITORING AGENCY NAME(S) AND ADDRESS(ES)

National Aeronautics and Space Administration

Washington, DC 20546-0001

10. SPONSOR/MONITOR’S ACRONYM(S)

NASA

11. SPONSOR/MONITOR’S REPORT

NUMBER(S)

NASA/TM–2011–215988

12. DISTRIBUTION/AVAILABILITY STATEMENT

Unclassified-Unlimited

Subject Category

Availability: NASA CASI (443) 757-5802

13. SUPPLEMENTARY NOTES

An electronic version can be found at http://ntrs.nasa.gov.

14. ABSTRACT

Modern multi-layer insulation (MLI) allows to sharply reduce the heat leak into cryogenic propellant storage tanks through

the tank surface and, as a consequence, significantly extend the storage duration. In this situation the MLI penetrations, such

as support struts, feed lines, etc., become one of the most significant challenges of the tank’s heat management. This problem

is especially acute for liquid hydrogen (LH2) storage, since currently no efficient cryocoolers exist that operate at very low

LH2 temperatures (∼20K). Even small heat leaks under microgravity conditions and over the period of many months give

rise to a complex slowly-developing, large-scale spatiotemporal physical phenomena in a multi-phase liquid-vapor mixture.

These phenomena are not well-understood nor can be easily controlled. They can be of a potentially hazardous nature for

long-term on-orbital cryogenic storage, propellant loading, tank chilldown, engine restart, and other in-space cryogenic fluid

management operations. To support the engineering design solutions that would mitigate these effects a detailed

physics-based analysis of heat transfer, vapor bubble formation, growth, motion, coalescence and collapse is required in the

presence of stirring jets of different configurations and passive cooling devices such as MLI, thermodynamic vent system, and

vapor-cooled shield. To develop physics-based models and correlations reliable for microgravity conditions and long-time

scales there is a need for new fundamental data to be collected from on-orbit cryogenic storage experiments. Our report

discusses some of these physical phenomena and the design requirements and future studies necessary for their mitigation.

Special attention is payed to the phenomena occurring near MLI penetrations.

15. SUBJECT TERMS

16. SECURITY CLASSIFICATION OF:

a. REPORT

U

b. ABSTRACT

U

c. THIS PAGE

U

17. LIMITATION OF

ABSTRACT

UU

18. NUMBER

OF

PAGES

46

19a. NAME OF RESPONSIBLE PERSON

STI Help Desk (email: help@sti.nasa.gov)

19b. TELEPHONE NUMBER (Include area code)

(443) 757-5802