Numerical simulation of tubes-in-tube heat exchanger in a mixed refrigerant Joule–Thomson cryocooler

Mixed refrigerant Joule-Thomson (MRJT) cryocoolers can produce cryogenic temperatures with high efficiency and low operating pressures. As compared to the high system pressures of around 150–200 bar with nitrogen, the operational pressures with non-azeotropic mixtures (e.g., nitrogen-hydrocarbons) come down to 10–25 bar. With mixtures, the heat transfer in the recuperative heat exchanger takes place in the two-phase region. The simultaneous boiling and condensation of the cold and hot gas streams lead to higher heat transfer coefficients as compared to single phase heat exchange. The two-phase heat transfer in the recuperative heat exchanger drastically affects the performance of a MRJT cryocooler. In this work, a previously reported numerical model for a simple tube-in-tube heat exchanger is extended to a multi tubes-in-tube heat exchanger with a transient formulation. Additionally, the J-T expansion process is also considered to simulate the cooling process of the heat exchanger from ambient temperature conditions. A tubes-in-tube heat exchanger offers more heat transfer area per unit volume resulting in a compact design. Also, the division of flow in multiple tubes reduces the pressure drop in the heat exchanger. Simulations with different mixtures of nitrogen-hydrocarbons are carried out and the numerical results are compared with the experimental data.


Introduction
Mixed refrigerant Joule-Thomson (MRJT) cryocoolers can produce cryogenic temperatures with high efficiency and low operating pressures. With MRJT cryocoolers, temperatures in the range of 80-230 K can be achieved by varying the percentage of the different components in a mixture. This temperature range can serve applications like cooling of infrared sensors, cryosurgery probes, gas chillers, etc. As with any J-T cryocooler, precooling of the high pressure stream in MRJT cryocoolers is not only essential for obtaining cryogenic temperatures but also for limiting the maximum pressure after compression. This precooling is achieved with the recuperative heat exchanger. The operational pressures with non-azeotropic mixtures (e.g., nitrogen-hydrocarbons) come down to 10-25 bar as compared to system pressures of 150-200 bar with nitrogen alone. Also, the heat transfer in the recuperative heat exchanger with mixtures occurs under two-phase conditions due to simultaneous boiling and condensation of the cold and hot fluid streams. Due to this reason, the heat transfer coefficients are high as compared to single phase heat transfer. As a result, the size of the heat exchanger reduces. Further reduction in size is possible with a multi tubes-in-tube heat exchanger which offers more heat transfer area per  [2], [3]) reported work on optimization of the mixtures and performance of the MR J-T systems. Alexeev et al. [4] simulated tubes-in-tube heat exchanger with different mixtures. However, the numerical predictions were not compared with experimental data except for the pressure drop on the shell side. Ardhapurkar et al. [5] predicted the hot side temperatures with simple energy balance equations from the measured values of temperatures on the cold side for a multi tubes-in-tube heat exchanger. Thermal losses and pressure drop in the heat exchanger were not considered. In general, numerical simulation of tubes-in-tube recuperative heat exchanger for MRJT cryocooler has received less attention in the literature. This may be due to the lack of generic heat transfer correlations for condensation and boiling in the cryogenic range for gas mixtures.
The objective of this work is to develop a numerical model for the simulation of a tubes-intube heat exchanger for a MRJT cryocooler. A numerical model is useful for optimizing the geometrical and operating parameters of the cryocooler. Previously, the authors have reported a numerical model for a simple tube-in-tube heat exchanger for a MRJT cryocooler [6]. Modified correlations were used for the flow boiling and condensation of fluid streams. Additionally, the J-T expansion process is also considered to simulate the cooling process of the heat exchanger from ambient temperature conditions. In this work, the developed numerical model is extended to tubes-in-tube heat exchanger and the numerical results are compared with experimental data for three different mixtures of nitrogen-hydrocarbons.

Numerical model
The fluid streams and solid tubes of the heat exchanger are divided into a series of control volumes (CVs) along its length. The governing equations of mass, momentum and energy are then written in a discrete form over these CVs. The CV arrangement for the inner tubes and the fluid inside shown in Figure 1 and Figure 2 shows the same for the external fluid. The numerical model is also briefly described in this section.

Assumptions
The following assumptions are made in model derivation for the tubes-in-tube heat exchanger: i) heat transfer and fluid flow only along the heat exchanger length i.e., 1D; ii) flow divides uniformly in all the inner tubes; iii) pressure drop is calculated with the homogeneous model; iv) inner and outer tubes are assumed to be adiabatic at ends; v) emissivity of the outer tube surface is constant and receives radiation from ambient temperature; vi) axial conduction in the fluid is negligible; vii) body forces and axial stresses in the fluid are negligible; viii) effect of helical coil is negligible. 2

Governing equations
The mass conservation equation over a fluid CV is given by: The momentum conservation gives: The shear stress is calculated with the homogeneous flow model. The two-phase mixture velocity at a given cross-section is calculated as The energy equation is written in terms of enthalpy due to two-phase flow conditions. The energy balance is given by: where, t=time, ρ= density (kg/m 3 ),ṁ=flow rate (kg/s), H=enthalpy (J/kg), T w =wall temperature (K),T = mean temperature of CV (K), p=pressure (N/m 2 ), V=velocity (m/s), h=heat transfer coefficient (W/m 2 K), A=cross-sectional area (m 2 ), τ w = wall shear stress (N/m 2 ), l p = wetted-perimeter (m) and x g is the gas mass fraction. Subscripts g and l indicate gas and liquid states respectively, while w denotes tube wall. An overbar denotes average over a CV. A general energy equation for the solid tubes is given below: Q conv represents the heat transfer per unit length due to convection from the surfaces of the solid elements.Q rad is the heat transfer per unit length due to radiation considered only for the outer tube surface. k and Cp denote the thermal conductivity (W/mK) and the specific heat (J/kgK) of the tube material respectively.

Boundary conditions
The inlet temperature (T), pressure (p) and mass flow rate (ṁ f ) are known for both hot and cold fluid streams. The inner and outer solid tubes are assumed to be adiabatic at ends. Thus, The initial temperature in the solid tubes and fluid streams is initialized to ambient temperature. T a,e is the temperature of the fluid after isenthalpic expansion, to the pressure p c at the inlet of the external annulus. T a,e goes on reducing from its initial value (ambient temperature) to cryogenic temperature T c at steady state. The subscripts h and c denote the inlet variables for the hot fluid in the inner tubes and the cold fluid in the external annulus respectively.

Heat transfer correlations
The validity of the existing correlations for condensation and boiling over the cryogenic range is not well established. Recently, Ardhapurkar et al. [7] assessed the existing flow boiling heat transfer correlations with the experimental data reported by Nellis et al. [8] for mixtures (nitrogen-hydrocarbons). The Granryd correlation [9] was modified and recommended for flow boiling of mixtures at cryogenic temperatures. The heat transfer coefficient (h m ) for flow boiling calculated according to [7] is given as: where, h lo is the liquid only heat transfer coefficient calculated from the modified Dittus-Boelter equation with properties of mixture as given below.
Here, G= mass flux (kg/m 2 s), d= characteristic length, P r l = μ l Cp l k l is the liquid Prandtl number with μ=viscosity(N-s/m 2 ) and Cp=specific heat (J/kg K). F p , the parameter for flow boiling of pure refrigerants, is given by: X tt is the Martinelli parameter for turbulent-liquid and turbulent-vapour flow calculated as: The parameter A G in the above equation is: where, C lg is the enhancement factor to account for the gas and liquid interface effects. For evaporation of refrigerants, Granryd [9] recommended C lg = 2. Ardhapurkar et al. [7], in the modified Granryd correlation, proposed C lg = 1.4 for G > 500 kg/m 2 s. Cp w is the apparent local specific heat for a non-azeotropic mixture and is defined as The correlation of Cavallini and Zecchin [10], recommended by [11], is employed for calculating the condensation heat transfer coefficients (h cond ). This is given as where Re eq =G eq d/μ l is the equivalent Reynolds number for two-phase flow and the equivalent mass flux G eq = G((1− x g )+x g (ρ l /ρ g ) 0.5 ). The above condensation correlation is corrected with the Silver [12] and Bell and Ghaly [13] method to account for the non-isothermal condensation process of mixtures as in [11].

Experimental set-up and heat exchanger configuration
The experimental set-up comprising of the heat exchanger amongst other devices and instrumentation is depicted in Figure 3. The details of the experimental set-up and measurement of mass flow rates, temperatures and composition can be found in the work reported by Ardhapurkar et al. [5]. Temperature sensors were installed on the cold side to measure the temperatures along the length of the heat exchanger. The tubes-in-tube heat exchanger is helically coiled with an overall length of 6.75 m. The heat exchanger dimensions are given in Table 1. Usually, the high pressure gas mixture enters the inner tubes at a pressures of 10-20 bar with a temperature of around 300 K. The low pressure gas mixture flows through the external annulus at lower pressure of around 4-6 bar in the opposite direction. Its inlet temperature is dependent on the lowest temperature that can be attained after expanding a given mixture composition from high pressure to low pressure. This temperature can be around 100-150 K. During its travel through the heat exchanger, the high pressure stream condenses inside the inner tubes while the low pressure stream evaporates in the external annulus. This forms a counter-flow heat exchanger with two-phase heat transfer.

Results and discussion
The aforementioned numerical model is employed for simulating the tubes-in-tube heat exchanger described in Table 1. Ardharpurkar et al. [5] carried out experiments with this heat exchanger for three different mixtures of nitrogen-hydrocarbons. They measured temperatures at eight locations including inlet and outlet on the low pressure (cold) side. Also, inlet and outlet temperatures on the high pressure (hot) side were measured. The intermediate temperatures on the hot side were obtained with energy balance. With the current model, the numerically obtained temperatures are compared with those reported by [5]. The operating conditions and  Table 2. The variation of thermophysical properties of the gas mixture, with temperature and pressure along the heat exchanger, is also taken into account. AspenOne software [14] is used for calculating these properties. The no load temperatures obtained during the experiments are also given in Table 2. With Mix#2, a lowest temperature of 98 K is achieved due to the higher percentage of nitrogen. Obviously, the lowest temperature is a function of the mixture composition.  Figure 4 shows the comparison between the measured and predicted temperature profiles of the hot and cold fluid streams for Mix#1 [15]. On the hot side, an outlet temperature of 132 K is predicted by the numerical model. The corresponding gas fraction is zero which is corresponds to the liquid state of the fluid with a bubble point of 136.7 K. This conforms with the experimental observation of liquid at the hot side outlet. The cold stream enters the heat exchanger in two-phase condition at a temperature of 127.89 K. Towards the hot end of the heat  exchanger the single phase heat transfer zone is observed. This is depicted by the linear and almost parallel temperature profiles. In the middle portion of the heat exchanger, the numerical predictions are on a higher side for both the hot and cold streams. The numerical values match well at the cold end of the heat exchanger. The change in slope of the temperature profiles, due to phase change, on both the hot and cold sides is also captured by the numerical model. The temperature profiles of the hot and cold fluids for Mix#2 are shown in Figure 5. The hot fluid enters as single phase gas at ambient temperature and leaves the heat exchanger at a temperature of 108.39 K. This corresponds to a two-phase state as its bubble point is around 104.48 K. The cold side fluid enters in a two-phase state at a temperature of 98 K which is well 6

ICECICMC
IOP Publishing IOP Conf. Series: Materials Science and Engineering 171 (2017) 012070 doi:10.1088/1757-899X/171/1/012070 above its bubble point temperature of 86.9 K. It leaves the heat exchanger at a temperature of 299.138 K in a single phase gas state as its dew point is around 248.72 K. In the single phase region near the hot end of the heat exchanger, the numerical and experimental values agree well with each other. In the middle portion of the heat exchanger, which is the transition region,  the rise in the temperature on the hot and cold sides is also well predicted by the numerical model. However, the numerically obtained temperatures are on the lower side as compared to the experimental values. The numerical and experimental profiles match again at the cold end of the heat exchanger. The temperature profile comparisons for Mix#3 are shown in Figure  6. In this case, the high pressure fluid entering the heat exchanger at ambient temperature leaves at a temperature of 117.275 K in a two-phase state. The low pressure fluid entering the heat exchanger in a two-phase state at a temperature of 110.9 K exits in a gas state at a temperature of 301.606 K. The numerical predictions and experimental observations match well at both the hot and cold ends of the heat exchanger for Mix#3. The rise in temperature of both the hot and cold fluid streams close to the single phase region is also well reflected in the numerical predictions. Also, the change of slope at the cold end is also replicated by the numerical model. Table 3 shows the outlet temperature on both the hot and cold sides. In all the three cases, the maximum relative differences between the numerical and experimental values of outlet temperatures are below 2.1%. Considering the lack of general correlations for condensation and boiling in the cryogenic range and the relative simplicity of the developed model, the qualitative trends and overall heat transfer in the heat exchanger are predicted reasonably well by the numerical model.

Conclusions
Numerical simulation of a tubes-in-tube heat exchanger in a MRJT cryocooler is presented in this paper. One dimensional transient model is employed for the two-phase heat transfer in the heat exchanger. Simulations are carried out for three different mixtures of nitrogen-hydrocarbons and the results at steady state are compared with the experimental data. Modified Granryd correlation is used for evaluating the heat transfer coefficients for flow boiling of the low pressure fluid. The condensation heat transfer coefficients are estimated with the correlation of Cavallini and Zecchin with the Silver-Bell-Ghaly correction. Thermophysical properties are evaluated at the local conditions of temperature and pressure along the heat exchanger length. The J-T expansion process is also taken into account to simulate the cooling in the heat exchanger from ambient temperature condition. The numerically predicted temperature profiles agree reasonably well with the experimental data and their qualitative trends are also reproduced very well. The maximum relative differences between the experimental and numerical values of outlet temperatures are less than 2.1%.