Three-dimensional fluid-structural interaction and thermal analysis of a large diameter horizontal cryogenic transfer line

This paper focuses on the thermal stress analysis and the thermal bowing phenomenon in a large 12-inch diameter cryogenic transfer line during chill-down with liquid oxygen. To investigate the chill-down characteristics and corresponding thermal stress generation, a threedimensional fluid-structural numerical model has been adopted in this work. The numerical simulation is performed on a 2 m long transfer line having an internal diameter of 30.48 cm (12 inch). Three dimensional modelling has been performed to numerate the channel wall temperature distribution along with the flow, at the top and the bottom portions. Based on the temperature data from CFD calculations, transient structural analysis and deformation of the transfer line are studied using a finite element model. Results show that the bowing in the transfer line occurs when there is a strong circumferential variation in temperature. Furthermore, this study provides a numerical simulation that predicts the chill-down time and the time when a stratified flow occurs in the line. This numerical method has been found to be a useful tool for the evaluation, design, and development of a large diameter cryogenic transfer line for transportation of large quantity of cryogenic fluid.


Introduction
Transfer lines, that are used for supplying cryogenic liquids to different locations from a source, are typical component of cryogenic systems used for thermal and fluid management in space flight, superconducting magnet cooling, cryogenic liquefaction plant etc. [1]. With increase in demand of cryogenic fluids like LH 2 , LOX and LNG etc. in large volumetric flow rate for applications like on ground cryogenic rocket propellant loading operation, natural gas transportation etc., transfer of these liquids through large diameter transfer line at sub-ambient temperature becomes a requirement [2]. However, the transient operation of bringing the transfer line channel wall to cryogenic fluid temperature (known as chill down operation) before vapor free liquid can flow through a large diameter transfer line involves non-uniform cooling of feed line systems with severe circumferential temperature gradient [3]. Excessive boiling of the incoming liquid is expected during this period as initially the channel wall will be at high temperature compared to the liquid. The purpose of the present work is to characterize such phenomenon during large diameter feedline chill-down, from thermal and structural perspectives.
In general, cryogenic transfer line chill-down characterizes flow boiling heat transfer stages such as forced convection film, transition and nucleate boiling with various dispersed as well as separated twophase flow patterns [1]. A number of experimental as well as numerical studies have been reported on this subject. For instance, Hartwig et al. [1], Darr et al. [4], and Jin et al. [5] have performed extensive investigation to characterize the chill-down phenomenon, development of accurate correlations for heat transfer coefficient and boiling regime transition criteria, and also methods for chill-down heat transfer enhancement. Other major contributions are [6][7][8]. Further, [9][10] are examples for recent numerical IOP Publishing doi:10.1088/1757-899X/1240/1/012012 2 studies on chill-down. These literatures have provided several insights into the complex flow phenomenon occurring during chill-down. For a detailed description of previous experimental investigations, predictive models, flow phenomenon, various boiling stages, corresponding flow patterns, reader is referred to the recent review type articles [1], [4].
Several investigations have been reported for cool down of small diameter transfer lines; however, research on a large diameter feed line of the order of 12-inch is comparatively limited. Bronson et al. [2] investigated pipe bowing and pressure oscillations during LH 2 chill-down. The author suggested maintaining a minimum mass flow rate to avoid the occurrence of stratified flow and thus the consequent bowing phenomenon. One of the initial studies related to very large diameter feed line chill-down was carried out by Schwartz et. al [3] in which, chill-down methods, time estimation, effects of cold shocking, and thermal stress analysis were reported for 18-inch diameter pipe line using both LH 2 and LOX. The method of precooling of the pipe line was used to avoid severe stress that can generate by direct cold shocking with liquid. Other notable studies are related to chill-down problems are Liebenberg et al. [11] and Tummescheit et al. [12].
Bowing in Cryogenic pipelines is mainly due to the induced transient temperature gradients over the cross-section of the pipe. Flieder et al [13] analyzed the nature and magnitude of the bowing in the transfer line and also the thermal stresses generated due to the temperature gradient. Due to the effects of gravitational force under horizontal pipe configuration, occurrence of stratified flow is expected much predominantly in large diameter transfer lines with severe circumferential temperature gradient. Such flow pattern can generate undesirable thermal stresses which can cause bowing, contraction, and yielding of the piping materials. It may also cause detrimental displacement of the system in relation to supports, guides, and restraints. Therefore, a detailed investigation of chill-down phenomenon in large diameter transfer line is required and is the objective of the present work.

Flow field analysis
A non-homogeneous single fluid approach is used for simulating the two-phase flow during chill-down with phase velocity deference accounted by a slip model. The various governing equations and constitutive relations of the model are described below.

Governing equations:
The continuity equation for the mixture is, where, → is the mass-average velocity, is the mixture density, = ∑ =1 and is the volume fraction of phase .
The equation of volume fraction for the secondary phase: The momentum equation: Where is the number of phases, → is a body force, and is the viscosity of the mixture: , is the drift velocity for secondary phase: The energy equation: In the above equation, is the effective conductivity (∑ ( + )) and is the turbulent thermal conductivity, according to the sst k-w model. consists of any other volumetric heat sources, while the first term on the right-hand side of equation is energy transfer due to conduction.
The above equation is for compressible phase and ℎ is the sensible enthalpy. For an incompressible phase = ℎ IOP Publishing doi:10.1088/1757-899X/1240/1/012012 4 Shear stress transport (SST) k-w model is used as a turbulence model. This model gives a better prediction of flow separation and behaves well in adverse pressure gradients. Near the wall, Wilcox model and in the free stream k-epsilon model get activated. This way it utilized the appropriate model throughout the flow field. In addition, near to the stagnation points, the shear stress helps the k-omega model just to avoid the rise of huge turbulent kinetic energy.

Boundary and initial conditions:
In the present work LOX is used as working fluid. At inlet, the mass flow rate of 20 kg/s and 0 kg/s are specified for liquid oxygen and oxygen gas respectively. At the outlet of the pipeline, the pressure-outlet boundary condition is specified. To avoid excessive generation of thermal stresses in the pipe wall, it is assumed that the precooling is performed by passing gaseous nitrogen and the transfer line is brought to 150K from ambient condition. Therefore initially (t = 0 s) the complete domain is considered to be at 150 K. Standard earth gravity is imposed on whole geometry during the entire period. Further the pipe material is given a temperature-varying property of SS304 and the simulation is performed with the varying time step (10 -6 -10 -3 s).

Discretization procedure:
In the present work, the control volume method is used to discretized the governing equations. The Pressure-Implicit with Splitting of Operators (PISO) algorithm is used to solve the Navier-Stokes equations which uses pressure-velocity calculation procedure. Along with the pressure and velocity corrections, the PISO algorithm performs neighbour correction and skewness correction, as well. This algorithm uses more CPU time per iteration but it can significantly decrease the number of iterations required for transient problems. Moreover, second order upwind scheme is used which uses two points for the computation. In addition, the PRESTO tie is used in pressure discretization, because compared to standard pressure interpolation scheme, PRESTO discretization gives more accurate results since interpolation errors and pressure gradient assumptions on boundaries are avoided. In addition, it calculates pressures on the face instead of interpolating the pressure on the faces using the cell centre values.

Grid independence test:
To get the more accurate results using CFD software, it is important to carried out the grid independence test. The structured hexahedral mesh elements are used with sizes varying in circumferential as well as in axial directions. Throughout the entire domain, the density of mesh is kept high to get the acceptable interface, and to get the highly accurate boundary layer effects. A gridindependent test is carried out for a time period of 60 s, for different grid sizes by calculating the temperature at a point where the maximum fluctuation in cooling rate is observed, (see Fig.3). The optimal number of nodes is found to be 6 million from grid independent study conducted using models having number of nodes from 0.2 million to 7.5 million nodes. It is found that refinement after 6 million showed no further change in the temperature profile. The simulations are performed using an intel server with 20 core processor and 128 GB RAM and takes one week to finish one simulation.

Results and discussion
The cool-down phenomenon is explained for a 30.48 cm (12 inch) inner diameter pipe with the help of temperature profile in circumferential as well as in the axial directions. In fig.4, the contours of volume fraction of liquid oxygen are taken at different flow times and it can be observed that with time the volume fraction of LOX increases. Fig.5 shows the temperature profiles at various locations in the pipe along the axial direction (shown in fig. 1). These points are located in the pipe domain at the topmost and bottom-most part. It is observed that the points located in the top portion of the pipe that is T1, T2, T3 and T4 takes almost 1200 s to reach 90 K; however, the points situated at the bottom part of the pipe that is B1, B2, B3 and B4 are having comparatively higher cooling rates and it hardly takes 100 s to reach 90 K. It is mainly due to the gravitational force induced stratified flow that causes the direct contact of LOX with the bottom wall surface.  fig. 1 for location). It is observed that as we proceed from top to bottom, the cooling rate gradually increases. Moreover, it can be seen that the points near to the bottom like I3 and I4 follows the same curve as of B1 and this can be better understood with the help of velocity streamlines (this is explained later). However, points located at the same potential height like I1 and I6 or I2 and I5 are having same rate of cooling. It is now clear from fig.5 and fig.6 that there is huge variation in the cooling rate between the top and bottom points of the pipe. For a better understanding of the same, the temperature difference of the top and bottom points is plotted in fig.7, for all the different axial location discussed earlier (see fig. 1). From fig. 7, the maximum circumferential temperature gradient noted is about 55 K. Further, in fig. 7, at time t = 0 s and t = 1200 s, the temperature difference is zero, which is the initial (150 K) and final conditions (90 K). Fig.7 also tells the total cool-down time which happened when all the curve touches the time axis (x-axis).   Fig. 8 shows the velocity streamlines at plane 1 for different flow times and these streamlines provides more insights in to the chill-down flow phenomenon. At time t = 1 min and 3 min it is observed that even though the points I3 and I4 are located above the point B1, these three points are having the same cooling rate since the fluid velocity close to these locations are same. Furthermore, at times, t=5 min and t=10 min, the velocity plot shows an interface which can be reasonably explained based on velocity difference between liquid and vapor oxygen phases. It can also be observed that this interface is moving from the bottom towards the top because of an increase in LOX volume fraction. There is no presence of interface in streamlines at 15 min and 20 min, however there is a red point in the cross section which is moving upward from bottom to the centre of cross-section. This red point indicates that the points located near this area are having almost the same velocities and it is maximum. At time t = 20 min, as the pipe is completely filled with the liquid oxygen, it can be seen that the flow is like a fully developed flow having maximum velocity at the centre and close to zero velocity near the wall.

Fluid structural interaction
For structural analysis, the temperature data from Fluent results is exported to Ansys mechanical transient structural module and the analysis is performed for complete chill-down time. The fine meshing of the pipe is done again here and SS304 material is assigned to the pipe with temperature varying material properties [16]. In actual industrial practice, the pipe is allowed to float in the vacuum chamber and there is no complete restriction to the pipe deformation during chill-down. Therefore, while doing the structural analysis the pipe boundary condition is provided in such a way that it can deform freely in the axial direction while its motion is restricted in the radial direction. The average equivalent Von-Mises stress and maximum deformation are calculated at every time step and the condition of the pipe at time t = 120 s when the maximum deformation happens is shown in in fig. 9. It is observed from fig.9 that because of the temperature gradient, the pipe is contracting in the axial direction as well as in the circumferential direction. Nevertheless, the bowing in the pipe is mainly due to the temperature gradient in the circumferential direction.  Fig.9 Bowing and axial contraction of the pipe due to circumferential temperature gradient. (Total deformation in mm is shown) In fig.10, the average equivalent thermal stress and maximum deformation are plotted with respect to time for the entire period. It is found that the average thermal stress rises dramatically to about 50 MPa in 200 s followed by a decrease in its value to 27.5 MPa at time t = 1200 s. In addition, the curve for maximum deformation reaches at its peak at 3.73 mm in just 120 s and afterwards it drops to around 3.35 mm.
Considering the size of the pipeline considered for present investigation, the obtained deformation is not severe. However, it needs to be mentioned here that the present investigation assumed that the axial deformations of the pipe are fully accounted using flexible joints such as bellows. If such provisions are not provided sufficiently, the stress generations may be huge and can even end up in transfer line failure.

Conclusions
Present work investigated the chill-down of a 30.48 cm internal diameter 2 m long cryogenic transfer line from 150 K initial temperature. It is assumed that precooling to 150 K has been performed to avoid excessive thermal stress in the feed line that may arise if direct cold shocking with liquid is allowed. Major outcome of the investigation is listed below. 1. Large diameter cryogenic transfer line chill-down shows majority of the chill-down operation in stratified flow regime with circumferential temperature gradient of the order of 55 K. 2. When LOX is flowing at a mass flow rate of 20 kg/s, chill-down took 20 min to finish. 3. There is a huge difference in cooling rate between the top and the bottom points of the pipe, because of stratified flow. Further, there is variation in the cooling rate in the axial direction as well. Hence, contraction of the pipe in both these directions are observed. 4. The bowing of the pipe is due to the circumferential temperature gradient. When the pipe is chilled from 150K to 90K, the maximum average thermal stress noted is about 50 MPa and the maximum deformation over the entire time period is 3.73 mm.