Numerical Investigations of the Cooling Performance of an R410A Closed-Loop Spray Cooling System

: To investigate the spray cooling characteristics and the impact of spray parameters such as chamber pressure, spray height, and spray tilt angle on heat transfer e ﬃ ciency, a mathematical model based on the Eulerian ‒ Lagrangian frame was established for an R410A closed-loop spray cooling system. The results revealed that the spray pa tt ern is conical, with the center velocity significantly higher than the edge velocity. The temperature distribution of the cooling surface and liquid ﬁ lm height both exhibit a “W” shape, and the surface temperature is lower where the liquid ﬁ lm is thin. There is an optimal liquid ﬁ lm height of approximately 5 µm, at which the cooling surface temperature is the lowest. The surface temperature increases with an increase in the spray chamber pressure. Considering average cooling surface temperature, the optimal tilt angle is 40° with an average surface temperature of 330.1 K. When considering wall temperature and wall heat transfer coe ﬃ cient uniformity, however, the optimal tilt angle is 10°, leading to the average surface temperature of 332.6 K. When increasing the optimal spray height to 70 mm, the average surface temperature is 313.4 K.


Introduction
Currently, more and more compact size results in a high power density in electronic devices.Ten years ago, the heat flux density of electronic chips had reached up to 100 W/cm 2 [1,2], and it is expected to achieve over 1000 W/cm 2 [3][4][5] in the near future.Thermal management technologies like micro-channels, heat pipes, and spray cooling have great potential, among which spray cooling technology stands out due to its high heat dissipation capacity, low surface temperatures, system stability, and cost-effectiveness, attracting widespread attention from researchers [6].
Up to now, different cooling agents have been tested for spray cooling including distilled water, anhydrous ethanol, liquid ammonia, and various refrigerants like FC-77, R134a, R22, R32, R410A, etc. Water has a high latent heat of vaporization (2489.65 kJ/kg) [7], but its boiling point at atmospheric pressure is 100 °C, which is higher than the fault temperature of electronic chips such as CPU (75 °C).The surface tension of ethanol is small, which leads to a larger spreading area of droplets, so as to obtain higher heat transfer performance [8].However, anhydrous ethanol is flammable, yielding significant safety concerns in experiments.Liquid ammonia can achieve a high critical heat flux (750 W/cm 2 ) because of its low boiling point (−33 °C) and high latent heat of vaporization (1200 kJ/kg) [9], but its toxic and flammable characteristics hinder its implementation.Refrigerants like FC-77, R134a, R22, and R32 with boiling points of 97 °C, −26.1 °C, −40.8 °C, and −51.7 °C show better heat dissipation capacity, but FC-77 has the drawback of high wall superheat, R134a is susceptible to moisture absorption [10], R22 is destructive to the ozone layer (ODP:0.055)[11], and R32 is prone to flammability [12].In comparison, R410A offers advantages such as high latent heat of vaporization (256.7 kJ/kg), low boiling point (−51.4°C), non-toxicity, non-flammability, and ODP of 0 [13], which make it a promising candidate for heat dissipation in relation to electric chips.Some researchers have conducted experimental research on spray cooling with R410A, seeking its optimal heat transfer performance through the investigation of experimental parameters.Zhou et al. [14] constructed a closed spray cooling experimental system using R410A as the cooling medium.Together with the tetrahedral structure to enhance heat transfer, a critical heat flux (CHF) of 411.2 W/cm 2 was achieved.Wang et al. [15] increased the CHF to as high as 522.1 W/cm 2 with a working temperature < 55.6 °C, which is enhanced by 73.4% on surfaces with fins of 3 × 0.5 mm.Longo et al. [16] used R410A as the cooling medium and conducted a comparative analysis of the condensation process between microfin tubes and smooth tubes, and it was found that the heat transfer capacity of the microfin tube was 3.5-3.6times that of a smooth tube with a mass velocity of 200 kg/(m 2 •s).However, the parameters such as liquid film velocity and liquid film height cannot be quantitatively analyzed only through experiments.
Compared to experiments, numerical simulations can reveal the cooling mechanism of sprays with a low cost [17].Cai et al. [18] established a mathematical model of the spray cooling system from the single-phase region to nucleate boiling via water and ethanolwater mixed cooling media.They found that the liquid film in the nucleate boiling region is thinner than in the single-phase region, and the heat transfer coefficient decreases with increasing spray height and increases with increasing spray pressure.Tian et al. [19] used a one-dimensional homogeneous flow model to simulate refrigerant gas-liquid two-phase flow inside the nozzle, and a discrete phase model (DPM) to simulate spray atomization, droplet evaporation, and fragmentation.They found that in the high-concentration spray region near the nozzle (Z ≤ 50 mm), the droplet temperature exhibited a "W" distribution radially.Away from the nozzle region (Z > 70 mm), the "W" distribution of droplet temperature completely disappeared, and a relatively flat low-temperature zone appeared.However, R410A is rarely used as the refrigerant to analyze the heat transfer performance of spray cooling in the current numerical calculation of spray cooling.
Besides the limited experimental research on closed-loop spray cooling systems employing R410A [20], theoretical research is far behind, which necessitates the investigation of the relationship between the cooling preformation and spray characteristics.At present, there are few studies on the systematic quantitative description of liquid film velocity and height, and the cooling surface temperature distribution.In this sense, the numerical simulation for an R410A closed-loop spray cooling system was conducted, allowing for the calculation of phenomena that cannot be directly observed in spray cooling experiments, such as liquid film velocity and height, and the cooling surface temperature distribution as well.Through the investigations of the impact of different spray parameters on heat transfer performance, this study will provide guidance for future applications of spray cooling.

Physical Model
To simulate the experimental condition of the R410A closed spray cooling system [13], the spray chamber shown in Figure 1a is simplified into a cylindrical computation domain with a diameter of 70 mm and a height of 110 mm (Figure 1b).At the bottom of the spray chamber, there is a protrusion along the z-axis, with a height of 1 mm and a diameter of 15 mm, made of polyether ether ketone (PEEK).This protrusion serves to reduce liquid submersion and act as thermal insulation.A circular heating surface with a diameter of 15 mm is placed in the middle of the protrusion, represented in blue in Figure 1b.A solid conical nozzle is employed to investigate the influence of spray inclination angle θ between the nozzle and the z-axis on heat transfer performance, with a constant distance H (25 mm) from the nozzle outlet to the center of the cooling surface.

Mathematical Model
The spray cooling process is extremely complex, involving interactions between droplets with the circumambient air and liquid film on the cooling surface.In this study, a Eulerian-Lagrangian approach is employed to establish a mathematical model, with liquid R410A as the cooling medium.The model utilizes a discrete phase model for the liquid phase and a continuous phase model for the spray environment, which consists of gaseous R410A.

Continuous Phase Model
For the liquid film and gas in the computational domain, a continuous phase model is constructed, including equations for mass, energy, momentum conservations, as well as component transport.These equations are represented by Equations ( 1)-(4) [21] as follows: where  is the density of the continuous phase fluid,  is time,  is the velocity of the continuous phase fluid,  is thermodynamic temperature,  is pressure,  is gravitational acceleration,  is fluid thermal conductivity, μ is viscosity,  is the mass of the continuous phase gas due to the evaporation of droplets,  is the internal energy of the continuous phase fluid, and ∆ is the change in continuous momentum due to droplets.The Realizable k-ε turbulence model is selected for its capability to provide more accurate predictions of the spreading rate of plane and round jets compared to the Standard k-ε model, as shown in Equations ( 5) and ( 6) [22]: where  is the turbulent viscosity,  is the turbulent kinetic energy generated by mean velocity gradients,  is the turbulent kinetic energy generated by buoyancy forces, in numerical calculations,  ,  ,  are model constants, and their values are typically set to 1.9, 1.44, and 1.0 [23], respectively.

Discrete Phase Model
Cryogen droplets are subjected to various forces, such as drag force, gravity, and lift, during their motion within a fluid.In the discrete phase model, one can employ Newton's second law to establish the motion equation of droplets in a Lagrangian frame.Since the environmental temperature changes little and the lateral velocity gradient is small for droplets of micro size, the thermophoretic force, Brownian motion, and Saffman lift can be neglected.Only drag force is taken into account, as shown in Equations ( 7)-( 11): where  is the velocity vector of the droplet,  is the droplet density,   / is the drag force per unit mass, d is the droplet diameter,  is the drag coefficient for the droplet,  , is the drag coefficient for a spherical droplet,  is the relative Reynolds number,  is the shape correction factor for the droplet: when the droplet is disk-shaped,  = 1, and when the droplet is spherical,  = 0.

Coupling between the Continuous and Discrete Phase Models
When calculating the trajectories of discrete-phase particles, a two-way coupling approach is used to consider the interaction between the continuous phase and the discrete phase.The exchange of momentum, heat, and mass between the continuous phase and the discrete phase is described by Equations ( 12)-( 14): where Δ is the time step,  ˙ , is the initial mass flow rate of liquid droplets,  , is the initial mass of the droplet,  in is the mass of the droplet as it enters the continuous phase grid,  out is the mass of the droplet as it leaves the continuous phase grid,  is the specific heat of the particle,  in is the temperature of the liquid droplet upon entering the continuous phase grid,  out is the temperature of the liquid droplet upon leaving the continuous phase grid,  ref is the reference temperature for enthalpy (300 K), and  lat ref is the latent heat of vaporization under reference conditions.

Thermal Properties and Boundary Conditions
The thermal properties of the R410A coolant were obtained using the REFPROP software (9.1), as shown in Table 1: As the flowchart shows in Figure 2, firstly the continuous phase is solved and iterates to stability; then, the computation of the discrete phase model starts from the nozzle exit, to update the velocity, temperature, and position.The second-order upwind scheme and the second-order implicit discretization scheme are, respectively, used for the spatial and temporal discretization of the momentum, turbulent kinetic energy, and turbulent dissipation rate equations, and the Simplec algorithm is employed for the pressure-velocity coupling.According to the design of the R410A closed spray cooling experimental system, the bottom of the spray chamber serves as the pressure outlet.As mentioned in Section 2.1, the constant heat flux wall condition is implemented on the cooling surface, where the wall film model [25] is applied for droplets forming a liquid film.The other surfaces of the spray chamber are set as escape boundaries, meaning that droplets reaching these surfaces will leave the wall and are no longer considered in the computation.At the beginning of the simulation, the spray chamber is filled with R410A gas at an ambient temperature of 300 K.

Model Validation
This study employs a polyhedral meshing technique in the Fluent (2022R1) meshing.To ensure more accurate results, three boundary layers are set, as shown in Figure 3.The distribution of grid mass is shown in Figure 4.The minimum mass of grid cell is 0.29, and the maximum mass of grid cell is 0.99.The grid independent test is conducted by using different cell numbers: 5000, 35,000, 100,000, 200,000, 300,000, and 450,000 grids.As shown in Figure 5a, when the cell number exceeds 200,000, there is minimal change in the temperature and liquid film height on the cooling surface.Therefore, 200,000 cells are selected for subsequent numerical simulations.
Based on the experimental data of the R410A closed spray cooling system in our group, the following simulations were conducted with a spray height of 25 mm, spray inlet pressure of 2.05 MPa, chamber pressure of 0.54 MPa, spray mass flow rate of 6.

Cooling Characteristics
In this study, numerical simulations were initially conducted with the following experimental parameters to investigate the characteristics of closed spray cooling: spray height of 25 mm, nozzle outlet pressure of 1.8 MPa, spray chamber pressure of 0.8 MPa, heat flux density of 100 W/cm 2 , coolant mass flow rate of 4.9 g/s, coolant temperature of 300 K, and a spray angle of 33.4°.

Spray Morphology and Velocity
From Figure 6a, it can be observed that the overall spray morphology is conical, and its velocity distribution shows that the center velocity of the conical spray is significantly higher than its periphery.This is because when the spray enters the chamber from the nozzle in the form of fine droplets, it undergoes flash evaporation due to heat absorption, resulting in a lower number of droplets at the periphery of the conical spray, which leads to internal contraction of the spray.Due to its high velocity, there is almost no radial migration, resulting in a gradual decrease in spray velocity from the center to the periphery.At the bottom end of the spray, on the cooling surface, the overall spray velocity is lower than at the top of the cone, and there is a stagnation point at the center of the cooling surface.From Figure 6b, it is evident that there is a velocity stagnation point at the center of the cooling surface for the liquid film (highlighted in green at the center).Moving radially from the stagnation point to the outer edge of the cooling surface, the liquid film velocity gradually decreases.Analyzing Figures 7a and 8a in conjunction, it is observed that in regions with higher liquid film velocity, the film thickness is smaller, leading to lower temperatures on the cooling surface.Conversely, in areas with lower liquid film velocity, the film thickness is larger, resulting in higher temperatures on the cooling surface.This is attributed to the fact that in regions with faster liquid film velocity, there is a larger flow rate, causing the cooling medium to accumulate radially outward, leading to a smaller liquid film thickness in areas with higher velocity.Additionally, regions with higher liquid film velocity exhibit enhanced convective heat transfer, resulting in better heat exchange efficiency and lower temperatures on the cooling surface.

The Distribution of Liquid Film Height
Figure 7a,b provide the contour plot and radial distribution of liquid film height on the cooling surface (dashed line in Figure 7a).From these figures, the liquid film height exhibits a "W" shape, wherein at the center of the cooling surface, the presence of a velocity stagnation point leads to a relatively thicker liquid film (highlighted in green at the center).Moving radially from the stagnation point to the outer edge of the cooling surface, the liquid film thickness gradually increases.As we can see from Figure 7a, the liquid film formed by droplets on the cooling surface can be divided into the liquid film bulge region (green at the circle's center), the impingement liquid film flow region (blue), and the free liquid film flow region (green, yellow and red), which is consistent with the observation from reference [26].When droplets reach the wall, the velocity in the direction of vertical heat source becomes zero, and droplets accumulate in the liquid film bulge region.The liquid film formed by the droplet flows around the wall, forming the impingement liquid film flow region in the spray area.As the flow progresses, the liquid film outside the spray area decreases in velocity due to the disappearance of the pressure gradient, and the liquid film tends to become thicker.This part becomes the free liquid film flow region.

Temperature Distribution on the Cooling Surface
Figure 8a,b provide the temperature distribution contour and radial temperature distribution on the cooling surface (dashed line in Figure 8a).Similar to Figure 7b, the radial temperature distribution on the cooling surface also exhibits a "W" shape.According to the mechanism of spray cooling, the single-phase region is primarily associated with liquid film evaporation and forced convection heat transfer, while the two-phase region is mainly linked to nucleate boiling heat transfer.Through comparing Figures 6b, 7a, and 8a, it is evident that areas with higher liquid film velocity and thinner liquid film exhibit faster liquid film evaporation and lower surface temperatures.This suggests that in R410A closed-loop spray cooling, the dominant heat transfer mechanisms are surface nucleate boiling heat transfer and secondary nucleate boiling heat transfer.Therefore, to enhance the heat transfer performance of closed-loop spray cooling, it is essential to focus primarily on the two-phase region.

The Impact of Spray Pressure
This section investigates the influence of the pressure inside the spray chamber on spray cooling heat transfer.Keeping all other spray cooling parameters constant, the spray chamber pressures were set to 0.5, 0.6, 0.7, 0.8, 0.9, and 1.0 MPa, respectively.
From Figures 9 and 10, it can be observed that both the radial distribution of liquid film height and the surface temperature exhibit a "W" shape.Under the same pressure, the wall temperature is lower in the region where the liquid film is thinner.The differences in liquid film height distribution are not significant under different pressures.It is worth noting that compared with Figures 9 and 10, the temperature is higher in places where the liquid film is thin under different pressures.This is because with the increase of pressure in the spray chamber, the heat transfer performance deteriorates, the wall temperature increases, the evaporation is intense, and the liquid film becomes thin.The relationship between spray pressure and the average temperature of the cooling surface, as well as the average liquid film height, is shown in Figure 11.From Figure 11 it can be found that as the spray pressure increases, the average temperature of the heat exchange surface gradually rises.This is due to the fact that with the increase in spray pressure, the corresponding saturation temperature of the cooling medium increases, reducing the temperature difference between the wall surface temperature and the saturation temperature, which results in the decrease of heat transfer efficiency.At the same time, with the spray inlet pressure remaining constant and an increase in the spray pressure, the pressure difference decreases, leading to a deterioration in spray atomization and a decrease in heat transfer efficiency.This also results in an increase in the average wall temperature.

The Impact of Spray Tilt Angle
This section investigates the effect of spray tilt angle, which is the angle formed between the nozzle and the perpendicular direction of the cooling surface, on the heat transfer efficiency of spray cooling.The calculations ensure that the distance H (25 mm) from the nozzle outlet to the center of the heat exchange surface remains constant, while only changing the angle θ between the nozzle and the z-axis.The impact of spray tilt angles (represented by θ) in the range of 0° to 45° is studied, with intervals of 5°.
The distribution cloud charts for different spray tilt angles is shown in Figure 12.The figure reveals that as the spray tilt angle increases, the liquid film velocity increases towards one side gradually, causing liquid film accumulation on that side.In the areas where the liquid film accumulates, the wall temperature is higher, and the heat transfer coefficient is lower.Additionally, it is observed that in regions with high heat transfer coefficients and low wall temperatures, the liquid film height is not necessarily at its minimum.This further demonstrates the existence of an optimal value for liquid film height.The relationship between the radial distribution of film thickness and wall temperature with respect to the tilt angle is illustrated in Figure 13a,b.In these figures, blue represents the lowest film thickness and wall temperature, while red represents higher film thickness and wall temperature.From the figures, it is easy to observe that the radial distribution of film thickness and wall temperature changes consistently.When the tilt angle is 0°, both exhibit a symmetric "W" shape.The proportion of different colors reflects the local distribution of film thickness and wall temperature.As the tilt angle increases from 0° to 45°, the radial distribution of film thickness and wall temperature gradually increases on the left side of the cooling surface, while it decreases on the right side.This is because as the tilt angle changes from 0° to 45°, the liquid velocity is greater on the left side of the cooling surface, causing liquid droplets to accumulate on the left side.This, in turn, results in higher film thickness on the left side compared to the right side.As mentioned earlier in the distribution of wall temperatures, areas with higher film thickness have less effective heat transfer, resulting in higher wall temperatures, while areas with lower film thickness exhibit better heat transfer and lower wall temperatures.From Figure 14, which shows the relationship between the average temperature of the heat exchange surface and the tilt angle, it is evident that the cooling surface temperature reaches its lowest point as the tilt angle changes.At a tilt angle of 40°, the average temperature of the cooling surface is 330.1 K.It is apparent that changing the spray tilt angle can easily lead to an asymmetric and non-uniform distribution of cooling surface temperatures.In the cooling process of high-power chips, researchers aim to avoid local overheating.Generally, evaluating heat transfer performance requires examining parameters such as wall temperature and wall heat transfer coefficient.Therefore, this study introduces two evaluation criteria: wall temperature uniformity and wall heat transfer coefficient uniformity, as shown in Figure 15.It is clear that as the spray tilt angle increases, wall temperature uniformity remains relatively consistent.However, wall heat transfer coefficient uniformity remains stable between 0° and 10° but gradually deteriorates beyond 10°.Considering the results from Figure 14, the optimal spray tilt angle is 10°, with an average cooling surface temperature of 332.6 K.   16, it can be seen that as the spray height increases, the radial distribution of liquid film height gradually changes from a "W" shape to a more uniform distribution and eventually forms an "M" shape, with the liquid film thickness decreasing.Through a comparison between Figures 16 and 17, there is an optimal value for local liquid film height, which is around 5 µm.Excessively high or low local liquid film heights will lead to higher local wall temperatures.From Figure 18, it can be observed that the average temperature of the cooling surface decreases as the spray height increases, until it reaches a minimum value and then starts to increase.Therefore, there is an optimal spray height of 70 mm, at which the corresponding liquid film on the surface is the thinnest.Thinner liquid films are advantageous for evaporative heat transfer.When the spray height is less than 70 mm, the impact velocity of droplets on the heat exchange surface is higher, increasing the chances of droplets splashing and rebounding, which is not conducive to liquid film formation and thus reduces heat transfer efficiency.When the spray height is greater than 70 mm, on one hand, the interaction between droplets and air is enhanced, and the interactions between droplets themselves are reduced.This leads to a decrease in droplet momentum and attenuation.On the other hand, thicker liquid films are formed on the heat exchange surface, resulting in slower liquid film flow velocity, reducing the convective heat transfer in the liquid film, and worsening the heat transfer efficiency.The spray coverage gradually increases with the elevation of the spray height, expanding from 1.13 cm 2 to 28.27 cm 2 .Correspondingly, the spray area utilization rate at the optimal spray height ranges from 9.8% to 17.4%.

Conclusions
In this study, a mathematical model based on the Eulerian-Lagrangian frame was established for an R410A closed-loop spray cooling system.The conclusions are as follows: (1) The spray pattern is conical, with the center velocity significantly higher than that at the edge, and the liquid film velocity gradually decreases from a central stagnation point to the periphery of the cooling surface.(2) The temperature distribution of the cooling surface and liquid film height both exhibit a "W" shape.(3) The surface temperature gradually increases with an increase in spray chamber pressure, and considering the average cooling surface temperature, the optimal tilt angle is 40° with an average surface temperature of 330.1 K.When considering wall temperature and wall heat transfer coefficient uniformity, however, the optimal tilt angle is 10°, leading to the average surface temperature of 332.6 K. (4) There is an optimal liquid film height of approximately 5µm, at which the cooling surface temperature is the lowest.The optimal spray height is around 70 mm, where the average cooling surface temperature is 313.4K. Correspondingly, the spray area utilization rate at the optimal spray height ranges from 9.8% to 17.4%.
out Mass of a droplet as it leaves the continuous phase grid, kg  Specific heat of a particle, J•kg

Figure 1 .
Figure 1.Spray cooling experimental system, geometric model, and nozzle position: (a) Experimental system spray chamber; (b) Geometrical model.

Figure 2 .
Figure 2. Schematic diagram of the solution algorithm.

Figure 6 .
Figure 6.Velocity distribution contour in the spray chamber: (a) longitudinal droplet velocity distribution; (b) liquid film velocity distribution on the cooling surface.

Figure 7 .Figure 8 .
Figure 7. Liquid film distribution: (a) liquid film height distribution on the cooling surface; (b) Radial distribution of liquid film height.

Figure 9 .
Figure 9.The radial distribution of liquid film height at different spray chamber pressures.

Figure 10 .
Figure 10.The radial distribution of cooling surface temperature at different spray chamber pressures.

Figure 11 .
Figure 11.The relationship between wall temperature, film thickness, and spray chamber pressure.

Figure 13 .
Figure 13.The relationship between wall liquid film height, wall temperature radial distribution, and spray tilt angle: (a) Radial distribution of wall liquid film height; (b) Radial distribution of wall temperature.

Figure 14 .
Figure 14.The relationship between average cooling surface temperature and spray tilt angle.

Figure 15 .
Figure 15.The relationship between uniformity and spray tilt angle.4.2.3.The Impact of Spray Height This section investigates the effect of spray height on the heat transfer in spray cooling.Keeping other parameters constant, the spray height varies from 20 mm to 100 mm with steps of mm, and the results are shown in Figures 16-18.The radial distribution of liquid film height and wall temperature under different spray heights is shown in Figures 16 and 17, respectively.From Figure16, it can be seen that as the spray height increases, the radial distribution of liquid film height gradually changes from a "W" shape to a more uniform distribution and eventually forms an "M" shape, with the liquid film thickness decreasing.Through a comparison between Figures16 and 17, there is an optimal value for local liquid film height, which is around 5 µm.Excessively high or low local liquid film heights will lead to higher local wall temperatures.

Figure 16 .
Figure 16.The radial distribution of liquid film height at different spray heights.

Figure 17 .
Figure 17.The radial distribution of cooling surface temperature at different spray heights.

Figure 18 .
Figure 18.The relationship between wall temperature and spray height.

− 1 •K − 1
in Temperature of a liquid droplet upon entering the continuous phase grid, K  out Temperature of a liquid droplet upon leaving the continuous phase grid, K lat ref Latent heat of vaporization, J•kg −1 q Heat flux, W•cm −2 X Coordinate position, mm h Heat transfer coefficient, W•K −1 •m −2Greek letters The density of continuous phase fluid, kg•m −3  Fluid thermal conductivity, W•m −1 •K −1  Viscosity, kg•m −1 •s −1  Turbulent viscosity, kg•m −1 •s −1 dDroplet density, kg•m −3 θ Spray tilt angle, 0