3. Materials and Methods
In this work, we performed an analysis of credible abnormal situations and highly unlikely accidents, ensuring that the UCN source’s deuterium chamber performs safely even with conservative approaches. We conducted analytical calculations using Python-based iterative schemes to gain initial insights into the container’s behavior under various accidental conditions. These conservative calculations, grounded in fundamental engineering principles and existing material property data, addressed factors such as pressure build-up, thermal loads, and associated stresses. Although these analyses do not incorporate complex validated simulation tools or the latest nuclear data sets at this stage, they lay the groundwork for more comprehensive studies. Subsequent investigations will leverage advanced computational modeling, enhanced nuclear data inputs, and robust validation efforts to further refine and optimize the container design, ultimately contributing to safer and more efficient UCN source operations. The properties of deuterium employed in our calculations are present in the
Table 1. Notably, the high value of the specific heat of deuterium together with its high mass provides us with relatively high time of response to accidental scenarios. In cases without available deuterium chemical properties’ data, hydrogen properties were adopted.
In this work we consider all the possible scenarios, including system’s cooling failure, loss of isolation vacuum, which with potential fracture of deuterium chamber is leading to the maximum hypothetical accident - mixture of air and deuterium causing deflagration or detonation. In all the scenarios, first we track the temperature evolution of the deuterium, which we find from the enthalpy balance. These scenarios will be discussed in following paragraphs, highlighting the temperature and pressure evolutions, and resulting mechanical stresses experienced by the structure of the UCN source.
Strong assumption of zero temperature gradients within the deuterium moderator and container walls was adopted, hence, all reported results represent average values. Also, we assume the container is sealed and rigid, so the volume
V is constant. Because the system can reach relatively high densities at cryogenic temperatures, applying the ideal gas law (
) may lead to non-negligible inaccuracies. In particular, strong intermolecular forces and significant molecular volumes invalidate the assumptions of an ideal gas under these conditions [
25]. Therefore, to account for the real-gas behaviors, we use the Redlich–Kwong equation of state (EoS), which introduces empirical corrections for molecular volume and intermolecular attractions:
where
a and
b account for intermolecular forces and nonzero molecular volumes [
27]. Compared with widely known Van der Waals EoS, Redlich–Kwong typically improves accuracy at elevated temperatures, while retaining simplicity for computational purposes. Redlich-Kwong EoS requires iterative or numerical methods to solve for
P at the final temperature
and known moles
n (or molar volume
), especially since these expressions are nonlinear. We implemented these calculations in Python, using standard root-finding routines.
With the internal pressure determined, the principal mechanical stresses in the container walls are evaluated using the thin-walled approximation - Mariotte model. The deuterium premoderator chamber consists of:
Cylindrical walls (Deuterium–Vacuum and Deuterium–Helium interfaces),
Spherical forward cap (Deuterium–Vacuum interface),
Annular plate (Deuterium–Helium coolant interface).
The theory behind how stresses are computed for each geometry using appropriate forms of the Mariotte equations is presented in
Appendix A. Failure is assumed to occur if
exceeds the allowable stress
, where
is the 0.2 % yield stress of Aluminum Alloy 5056 and
is the safety factor. Due to all the aforementioned assumptions, a conservative safety factor of
is adopted (ASME BPVC Section VIII, Div. 1 (UG-23, UG-25)), which also accounts for uncertainty in material properties at cryogenic and elevated temperatures (20–600 K), and lack of fatigue or stress concentration analysis in this preliminary phase. The goal is to design the container that will plastically deform, but not catastrophically fail under extreme impulsive loading.
Finally, although this analysis focuses on the deuterium chamber, the helium coolant envelope acts as a secondary containment structure. Its walls serve as an additional safety barrier, enhancing the mechanical resilience of the overall cryogenic module in the event of partial failure.
3.1. Power Control and Cooling Failure
The uncontrolled power increase of the WWR-K reactor for more than 20 % (7.2 MWth) is leading to automatic control rods insertion, shutting down the reactor, according to WWR-K reactor’s safety guides. If the helium refrigerator fails, the helium coolant will no longer be available to recondense deuterium vapor, leading to a gradual rise in pressure within the moderator tank.
In cryogenic systems, three primary heat transfer mechanisms are typically considered: (1) thermal radiation from surrounding structures, (2) thermal conduction along mechanical supports and piping (“highways” and suspensions), and (3) energy deposition from neutron and gamma radiation originating from the reactor core. In the present analysis, conductive heat transfer via supports and piping is neglected, as it would require a comprehensive thermal model of the entire facility and is expected to contribute significantly less than the other two mechanisms. The radiative heat load was explicitly calculated and amounts to 24 W for the current configuration.
The contribution from reactor-induced radiation was estimated for the highest possible power using the Monte Carlo code MCNP 6.2 [
29]. For the baseline configuration, where both the deuterium and helium moderator shells are made from Aluminum Alloy 5056, the following energy deposition rates were obtained:
Eventually, the total heat load sums up to 171.5 W, which needs to be removed by a cooling system.
An alternative design was also evaluated, in which the walls of the helium moderator were assumed to be constructed from Zircaloy-4, while the deuterium moderator retained aluminum alloy walls. This change resulted in a reduction in heat generation within both the deuterium and helium volumes. However, the heat deposition in the Zircaloy-4 walls increased compared to aluminum, partially offsetting the gains in thermal performance. For this reason, the evaluation of the Zircaloy-based configuration has been deferred to future studies.
In the present work, deuterium is initially stored as a liquid at about 20 K, which is well below its normal boiling point. The overall heat rate that warms the liquid deuterium is specified as an external input. We adopt a piecewise approach to account for phase change and subsequent vapor heating:
The total heat
needed to bring the entire mass
from
K to some final temperature
is thus obtained by summing these contributions:
where each term is activated only when its respective temperature or phase-change regime is reached. In our simulations, a user-defined heat rate
is specified, and we integrate or step through time to determine how
accumulates and in which phase the deuterium resides at each stage. This analysis assumes a worst-case scenario in which pressure cannot be released due to simultaneous failure of both the pressure relief valve and deuterium return channel.
3.2. Loss of Isolation Vacuum
Loss of the isolation vacuum is among the most critical hazards in cryogenic systems, as any influx of ambient gas substantially increases heat transfer to cold surfaces. The severity of such an event depends on the rate of temperature rise within the cryogenic container: faster warming accelerates pressure buildup and raises the risk of mechanical failure. Initially the vacuum vessel has a pressure equal to Pa, reaching the atmospheric pressure at the end of the accident.
Two principal routes can lead to a total loss of insulating vacuum: (i) failure of the vacuum pump or the protective circuitry, allowing air to enter the vacuum container, or (ii) rupture of the premoderator chamber filling the vacuum vessel with deuterium. In the latter case, leaking D2 would expand into a larger vessel volume, quickly reducing the likelihood of overpressurization. Moreover, the absence of oxygen in the vacuum eliminates deflagration risks. Notably, reaching the vacuum space requires a simultaneous breach of both the deuterium chamber and its thin helium coolant envelope—multiple failures that further decrease the probability of such an incident. Even if a small leak occurs, the large vacuum volume and higher wall thickness lead to lower pressures and lesser stresses than in other scenarios.
Because the more consequential hazard is air ingress, this study focuses on the potential for ambient air to leak into the vacuum vessel. We analyze both the cylindrical sections and spherical cap of the moderator assembly, where radiation, conduction, and natural convection act with geometry-specific parameters. The annular plate surface is omitted from these calculations because it only faces helium coolant channels at low temperature, which would further reduce the deuterium temperature and thus yield a conservative overestimate of heat transfer when excluded, however stress calculations that the annular plate experiences are included in this work.
3.2.1. Cylindrical Geometry
We assume the deuterium container is coaxial with the vacuum vessel, forming a cylindrical gap whose pressure P increases over time from to . At each time step, the net heat flux into the container wall, , is computed as follows:
Radiative Heat Transfer between concentric cylinders of emissivities
and
:
where
is the Stefan–Boltzmann constant,
F is a view factor close to unity for nested surfaces, and
is the inner cylinder’s external area.
Conductive Heat Transfer through gas:
where
L is the cylinder’s axial length,
and
are the vacuum vessel’s inner radius and the container’s outer radius, respectively, and
depends on the instantaneous pressure in the gap.
-
Convective Heat Transfer, which becomes significant as the gap pressure rises: The convective heat transfer rate between the cylindrical surfaces is determined using the Rayleigh number (
) and Nusselt number (
). The Rayleigh number is calculated as:
where
g is the gravitational acceleration,
is the thermal expansion coefficient,
and
are the temperatures of the interacting surfaces,
L is the characteristic length,
is the kinematic viscosity, and
is the thermal diffusivity. For natural convection, the Nusselt number is estimated based on the Rayleigh number:
The convective heat transfer coefficient is then given by:
The convective heat transfer rate is finally computed as:
where
is the interacting surface area.
3.2.2. Spherical Cap Geometry
At the vessel’s forward region, the deuterium container features a spherical cap rather than a cylindrical surface. The gap between this cap and a larger enclosing dome is relatively small (on the order of 15 mm), so a parallel-plate approximation is employed to characterize radiative, conductive, and convective fluxes. Denoting the cap’s temperature by
and the surrounding structure’s temperature by
, we have:
Here,
is the emissivity,
is the Stefan–Boltzmann constant,
is the spherical cap’s external surface area,
is the effective thermal conductivity of the gas, and
is the small radial distance of the gap. As in the cylindrical case, the total heat flux into the spherical cap is:
3.2.3. Container–Deuterium Coupling
By summing the heat flux contributions from the cylindrical and spherical sections, the overall thermal load on the premoderator container is determined. Helium cooling of deuterium and deuterium container’s wall was neglected to make the calculations more conservative, hence radiation heat rate coming from the WWR-K reactor was added (). The integrated model then updates the deuterium temperature and phase fraction over each time step, allowing an estimation of how quickly the system warms and how rapidly internal pressures rise.
At each time step the container wall temperature
is first advanced with the net heat:
with
and
as the wall’s mass and specific heat, respectively. A separate conduction path through the wall thickness governs heat flow into the deuterium. Depending on phase-change behavior and fluid properties,
and the vapor fraction are updated each step. The updated wall then conducts heat through the metal-helium–metal stack into the deuterium volume. We treat that conductive flux as:
where
is the total wetted area (cylindrical + spherical cap) and
are the individual layer thicknesses. The energy that actually reaches the deuterium over a time step
is:
where
stands for any direct radiation from the reactor that heats the deuterium.
Because the sign of
in Eq. (
17) can reverse, energy can flow
from deuterium
back to the wall whenever
, ensuring the two bodies equilibrate instead of diverging. This explicit energy-balance update is repeated for every
until the prescribed simulation time
is reached. These transient temperature and pressure profiles are subsequently used with the mechanical stress calculations (see
Appendix A) to verify if the Aluminum 5056 chamber remains within safe operational limits.
Ultimately, this unified approach provides a conservative estimate of heat ingress during a vacuum-loss event, whether in a cylindrical segment or at the spherical cap. It thus captures the critical pathways for radiative, conductive, and convective heat transfer under realistic leak conditions, enabling a thorough assessment of the deuterium system’s thermal and structural margins.
3.3. Liquid Deuterium Deflagration/Detonation Scenario
It is essential to evaluate whether a detonation—characterized by supersonic flame propagation and destructive blast waves—could feasibly occur in the proposed system. For cryogenic flammable gas mixtures such as deuterium–air, detonation requires conditions that are rarely achievable in practice.
The initiation energy required for detonation is orders of magnitude higher than for deflagration. Moreover, geometric constraints play a key role: for a self-sustaining detonation wave to form and propagate, the smallest dimension of the vessel or piping must exceed a critical threshold, typically greater than 12 times the detonation cell size
. Experimental data [
26] indicate that for stoichiometric deuterium–air mixtures,
is small compared to other fuels but still requires channel diameters larger than those present in our system. For the most conservative stoichiometric case, the detonation critical tube diameter equals to 30 cm. In the AlSUN configuration, both the deuterium chamber and associated cryogenic piping are designed to be narrower than this critical limit, and no flame-accelerating obstructions (e.g., baffles or sudden contractions) are present.
Even the possibility of a deflagration-to-detonation transition (DDT) is highly unlikely. While DDT has been observed in large unconfined clouds under specific conditions—such as in petrochemical accidents involving multi-ton vapor clouds—these conditions include high turbulence, delayed ignition, and favorable geometry. None of these are present in the AlSUN moderator system. Additionally, any air that could leak into the cryogenic system would immediately freeze on contact with cold surfaces, preventing homogeneous mixing and flame propagation. Consequently, the possibility of detonation is dismissed on both physical and geometrical grounds, however it was assessed in this work for the completeness of analysis.
3.3.1. Deuterium Deflagration Analysis
First, we start with a more plausible deflagration scenario, where combustion occurs subsonically in a partially confined space. The primary hazard associated with deflagration is the rapid pressure rise resulting from exothermic chemical reactions in a closed or semi-closed system. This scenario considers how a known mass of cold liquid deuterium (), initially at approximately 20 K, may undergo deflagration upon contact with ambient air. Under normal operation, both deuterium and air components remain frozen, and gas-phase combustion is only possible during cryostat warm-up. To define bounding conditions for safety assessment, the container is modeled as adiabatic and isochoric, assuming thermally insulating walls and neglecting mechanical deformation over the timescale of the event. Adiabatic Isochoric Complete Combustion (AICC) pressures, , are therefore adopted as conservative estimates for the maximum overpressure that could arise.
The deflagration accident could originate from one of two principal initiating events. The first is a simultaneous rupture of the vacuum vessel and deuterium chamber shell, which would allow significant air ingress into the moderator volume. Although unlikely—given the robust shielding of the reactor thermal column and structural isolation of these components—this worst-case condition must be evaluated to set conservative design limits. The second, more credible pathway is a rupture of the deuterium feed or return lines while the fill-valve is open. In such a case, the ingress of air is limited to the volume of the piping network, and it is negligible compared to the 239.7 kg required for a stoichiometric mixture, which has volume concentrations of fuel and air such that no excess fuel or oxygen remains at the completion of the chemical reaction.
In the considered design, the deuterium container holds a prescribed quantity of
(14 kg), while during a hypothetical accidental scenario air enters at a specified mass and temperature, simplified to 23.2% of
and 76.8% of
by mass, ignoring minor species like Ar,
, etc. The combustion of deuterium is represented by a single-step reaction:
where partial or complete consumption of deuterium depends on whether the oxygen is below stoichiometric (
), exactly stoichiometric (
), or in excess (
). Unreacted deuterium or oxygen, therefore, remains if one reactant is limiting, and the nitrogen remains inert.
Because the container volume is fixed, and no heat is assumed to leave or enter, the system evolves under adiabatic, constant-volume conditions. The final temperature is found by balancing enthalpies under the assumption that the gaseous species (, , , and ) behave ideally above 298 K. Liquid deuterium below its boiling point near 24 K is heated using piecewise constants for heat capacity and an appropriate latent heat term at vaporization. All chemical dissociation at high temperature and wall heat losses are neglected, producing a conservative upper-bound temperature estimate.
To formalize this, we write a typical adiabatic energy balance:
where
and
are the moles and enthalpies of each species
i, respectively, and
is the net enthalpy of reaction for the fraction of deuterium that reacts. Below 298 K, we account for liquid-vapor transitions of deuterium and any sensible heating to 298 K; above 298 K, we use Shomate polynomial [
20] expressions to calculate temperature-dependent enthalpies. The Shomate polynomial for the molar enthalpy
of each gas-phase species typically takes the form:
where
is often given in units of J/mol (or kJ/mol), t = T(K)/1000 and
A through
F are fitted coefficients valid over a particular temperature range (e.g., 298–6000 K) [
20]. For each species
i:
and the Shomate form provides a closed-form expression for
. The standard enthalpy of formation at 298 K,
, is included in
for non-elemental species (e.g.,
).
In practice, the Python-based iterative solver code first brings each reactant from its actual initial temperature to a standard reference (298 K), applies the reaction enthalpy at 298 K, and then raises the products from 298 K to the final temperature
. Because the net system enthalpy must remain constant in an adiabatic enclosure, the solver iteratively adjusts
until the overall enthalpy balance is satisfied. The final state represents a conservative temperature
that may exceed real-world conditions if partial dissociation, incomplete combustion, or heat losses occur. The defined
allows to compute the
employing the Ideal Gas Law:
where
is the initial pressure of the gas mixture before combustion,
is the initial temperature of the gas mixture before combustion,
is the total number of moles of product gasses after combustion,
is the total number of moles of reactant gasses before combustion. As a bounding estimate, this provides critical insight into maximum pressures and, by extension, into stresses on container walls. Thus, the model serves as a straightforward yet sufficiently comprehensive framework for preliminary safety assessments of a liquid-deuterium moderator vessel subjected to a sudden influx of air and subsequent deflagration.
3.3.2. Deflagration to Detonation Transition Analysis
While deflagration is the more likely combustion regime under accidental air ingress, a conservative safety analysis must also consider the possibility of a transition to detonation. In particular, deflagration-to-detonation transition (DDT) is a well-known phenomenon in confined geometries containing reactive mixtures. Under certain conditions—especially when turbulence, confinement, and geometric obstacles are present—a slow flame can accelerate, resulting in a shock-coupled reaction front and abrupt pressure escalation. Although the likelihood of DDT in the current system is considered low due to cryogenic conditions, limited air ingress volume, and minimal obstruction in flow paths, the structural impact of such a regime must be assessed. Detonation leads to sharply elevated pressures and impulse loads that far exceed those of steady deflagration. As such, it defines the upper bound of mechanical demand on the vessel walls and surrounding components.
The detonation of a deuterium–air mixture is approximated as an impulsive internal load applied to the vessel walls. In such scenarios, the stand-off distance between the explosive front and the structure—as well as the degree of confinement—play crucial roles in determining the intensity and spatial distribution of the pressure loads. In general, gaseous explosions result in a combination of step-loading (rapid but sustained pressure rise) and impulse loading (transient pressure pulse), and are typically modeled under the assumption that the gas completely fills the available vessel volume.
Detonations often arise from uncontrolled ignition events where deflagration accelerates and transitions into detonation. This process can lead to the development of a high-speed shock wave strongly coupled to the structural response of the confining vessel. As the detonation propagates through the reactive medium, the proximity of the wavefront to structural walls leads to transient loading with high spatial and temporal gradients in pressure.
Experimental observations [
21] show that a detonation wave consists of a leading shock front closely followed by a thin reaction zone in which the chemical transformation of fuel and oxidizer occurs. This zone produces combustion products at extremely high temperatures (2000–3000 K), resulting in a steep rise in pressure and density. For stoichiometric hydrogen–air mixtures—a well-characterized surrogate for deuterium–air systems—the detonation travels at the Chapman–Jouguet (CJ) velocity, approximately
m/s. Due to limited experimental data on deuterium detonations, hydrogen parameters are adopted here as a conservative approximation.
The reaction zone thickness in gaseous detonations is typically below 10 mm for fuel–air mixtures and up to 100 mm for fuel–oxygen mixtures. However, hydrodynamic instabilities and turbulence increase the effective width to 10–100 times larger, introducing fluctuations in pressure and flow fields. Nonetheless, for structural calculations, these transient fluctuations are often neglected in favor of time-averaged pressure values. The ideal peak pressure immediately behind the detonation front, denoted as , is approximately 15.6 bar for hydrogen–air systems.
From a structural mechanics perspective, the detonation imposes a propagating, spatially nonuniform load. There is no pressure acting ahead of the detonation front due to its supersonic nature. At any fixed location, the pressure rises abruptly upon the arrival of the detonation front and then decays through an expansion wave. This expansion zone typically reaches halfway between the leading front and the initiation point. In the trailing region behind the expansion wave, the gas becomes stationary, and the pressure decreases to around .
Several mechanisms can amplify the pressure beyond :
Shock reflections from end caps or curved walls can produce secondary shock waves, raising the local pressure to 2–3 .
Geometric focusing—especially in annular or curved geometries—may intensify the pressure by constructive interference of wavefronts.
Strong confinement prevents rapid venting of combustion products, allowing pressure to accumulate.
Transient overshoot during DDT, due to rapidly accelerating turbulent flames, can generate brief spikes above the steady-state CJ value.
As supported by experimental and modeling studies (e.g., [
21]), the combined effect of these phenomena can produce peak pressures as high as 5 times
in confined vessels. Therefore, in this analysis, we adopt a conservative design pressure of
to ensure a robust safety margin in the subsequent structural stress evaluation. However, the factor is much lower and equals 2–3
for typical pipe runs or vented enclosures, and since we have an air ingress, the detonation can’t be considered fully confined. The structural system is idealized as one-dimensional with a simplified multi-layered geometry, consisting of concentric spherical or annular walls. The detonation wave is assumed to propagate in the direction of the reactor’s active zone, which represents the most conservative trajectory in terms of structural loading and confinement effects.