2D and 3D numerical modelling of internal flow of Pressure-swirl atomizer

This paper compares 2D axisymmetric and 3D numerical models used to predict the internal flow of a pressure-swirl atomizer using a commercial software Ansys Fluent 18.1. The computed results are compared with experimental data in terms of spray cone angle (SCA), discharge coefficient (CD), internal aircore dimensions and swirl velocity profile. The swirl velocity was experimentally studied using a Laser Doppler Anemometry in a scaled transparent model of the atomizer. The internal air-core was visualized at high temporal and spatial resolution by a high-speed camera with backlit illumination. The internal flow was numerically treated as transient two-phase flow. The gas-liquid interface was captured with Volume of Fluid scheme. The numerical solver used both laminar and turbulent approach. Turbulence was modelled using kε, k-ω, Reynolds Stress model (RSM) and coarse Large Eddy Simulation (LES). The laminar solver was capable to predict all the parameters with an error less than 5% compared with the experimental results in both 2D and 3D simulation. However, it overpredicted the velocity of the discharged liquid sheet. The LES model performed similarly to the laminar solver, but the liquid sheet velocity was 10% lower. The twoequation models k-ε and k-ω overpredicted the turbulence viscosity and the internal air-core was not predicted.


Introduction
The Pressure swirl atomizers (PS) have an irreplaceable role in many industrial applications including combustion, spray cooling, spray drying etc. In a typical PS atomizer, the pumped liquid is fed via tangential ports into a swirl chamber where it gains high angular velocity and creates a low-pressure zone along a centre line of the swirl chamber. Subsequently, the air is pulled inside the low-pressure zone, so an air-core is formed. The swirling liquid is discharged from the exit orifice in a form of a conical liquid sheet at a certain spray cone angle (SCA) and consequently disintegrates due to aerodynamical forces into filaments and ligaments. Velocity and thickness of the liquid sheet, as well as the SCA, affect the size of resulting droplets.
Despite the simple geometry, the internal flow behaviour is complex, mainly due to the dominant swirling velocity component and the induced internal aircore which blocks a portion of the exit orifice. The internal vortex behaves as a Rankine vortex as the swirl velocity has its maximum located at the air-core surface [1]. Secondary flow effect as Görtler vortices could be also presented in a boundary layer inside the swirl chamber [2,3].
In open literature, many correlations are available (for a general overview see [4]), which allows us to estimate the atomizer performance such as droplet sizes, SCA, liquid sheet thickness, discharge coefficient, etc. However, the most of those correlations were made by fitting the experimental results and they were developed for limited geometry variations and operating conditions.
To better understand the link between the atomizer performance and design, the internal flow must be examined. Some authors studied the internal flow analytically. Simple non-viscous treatment, reviewed in [5,6], can be used for a basic insight into the flow behaviour, but it is not accurate enough to predict the discharge parameters. The better agreement can be achieved when the viscous flow is assumed [7,8], but some aspects like a temporal stability or secondary flow effects are unable to be resolved.
The number of papers related to the numerical simulation of the PS atomizers had arisen in recent years due to an increase in a computational performance. One of the first numerical studies of the PS atomizer was conducted in 1997 by Yule and Chinn [9]. They used a 2D axisymmetric model with the laminar solver and reported a deviation from the experimental data to be less than 3%. Similar numerical setup was used by Amini [7] and showed better agreement with the experimental data than the analytical viscous models. He compared those models at Reynolds number inside the inlet ports in a range of Re p = 11,000-122,000. The difference between 2D and 3D computational models was examined by Sumer et al [10]. They used the laminar solver and found that the aircore diameter was about 5% smaller in the case of the 3D model. However, the frequency of waves on the air-core surface was predicted by both models closely and with good agreement to the experimental data. A comparison of Large Eddy Simulation (LES) and the laminar models was performed by Madsen et al. [11]. They used a scaled atomizer and operated it in range of Re p = 12,000-41,000. At those operating regimes, the laminar model had a slightly better agreement to the experimental data than LES. The authors also examined simple turbulence models represented by the RNG and realizable k-ε models. However, these models were unable to predict the internal air-core. Galbiati et al. [12] compared the LES simulation with RNG k-ε and RSM models. They found that the variation among used models was insignificant when compared to deviation in results from the empirical correlations. They also noted that the flow field was consistent for the LES and k-ε model, while RSM had some discrepancies. Baharanchi et al. [13] examined several schemes to capture the liquid-air interface using a 2D simulation with the RNG k-ε turbulence model. A georeconstruct scheme was found to be the most suitable for capturing the air-core. They also discussed if it is necessary to model a surface tension and found that the surface tension had effect only if the Weber number is smaller than 204. From a practical point of view, it had a negligible effect on the developed flow.
In this paper, 2D and 3D numerical models of the PS atomizer are compared. The numerical solver used both laminar and turbulent approach. The numerical results are also compared with experimental data from [1].

Experimental setup
The experiments were performed on a cold test bench at the Brno University of Technology.

The Atomizer
A small-sized Simplex atomizer used in a combustion chamber of a turbojet aircraft engine was investigated. Due to small proportions, see the dimension in Fig. 1, the optical measurement inside the atomizer would be difficult. To overcome this issue the transparent scaled model was manufactured. The scaled model was 10 times larger than the original atomizer and it was made from cast PMMA. To maintain the same flow behaviour, dimensionless numbers such as Reynolds number, Swirl number, Weber number and Froude number must be kept the same.
Reynolds number related to the inlet ports is defined: where d p is the hydraulic diameter of the inlet port, w p is the mean velocity inside the port and is the liquid kinematic viscosity.
The swirl number S 0 is the ratio of swirl momentum to the axial momentum, and for the Simplex atomizers it can be defined as [14]: where r o = d o /2, d o and r s are in Fig. 1 and A i is the total area of the inlet ports. The S 0 is identical for both the original and the scaled atomizer. Froude number, Fr, compares the effect of gravity with the energy of the bulk flow and can be calculated as: where r oa is the radius of air-core in the exit orifice, Q is the volume flow rate. To minimalize the effect of gravity, it is necessary to keep Fr >> 1. Weber number, We, is defined as the ratio of inertia and the surface tension force. The surface tension force usually has negligible effect inside the swirl chamber. According to [13] and [14], matching the Weber number can be safely ignored. Similarly, no effect of surface tension was observed in our initial CFD simulations.
Both the original and scaled atomizers were operated at room temperature of 23 °C using JET A1 with physical properties: σ = 0.029 kg/s 2 , μ l = 0.0016 kg/(m·s) and ρ l = 795 kg/m 3 .
The original atomizer was operated at the inlet pressure Δp i = 1 MPa, the inlet mass flow rate was ̇ = 7.3 kg/h. The mean velocity inside the inlet ports was w i = 5.16 m/s, which resulted in Re p = 1021. The scaled atomizer was operated at the same Re p thus the inlet pressure had to be reduced to 10 kPa. The Fr for original atomizer was 293 but it decreased to 9.3 for the scaled atomizer.

Experimental setup
The high-speed camera Photron SA-Z was used to document the spatial and temporal behaviour of the aircore inside the scaled transparent model of the atomizer. The camera was set to 4,000 frames per second with a frame resolution of 1024×1024 px. The camera was used a backlit illumination with background LED panel. The typical result is shown in Fig. 2.  Fig. 2. The results were corrected according to Zhang [15] due to the different refractive index of the liquid and atomizer body. For detailed description of the experimental setup see [1].

Numerical setups
The transient CFD simulations were made using commercial software Ansys Fluent 18.1. A 2D axisymmetric model with a swirl, shown in Fig. 3, and a three-dimensional 3D model with periodic boundary conditions, as shown in Fig. 4, were used. As the atomizer had three inlet ports, a 120° section was modelled. Pressure-velocity coupling was done using PISO scheme. The liquid-air interaction was captured by a Volume of Fluid model with a geo-reconstruct scheme or Compressive scheme was used in the case of LES. The inlet boundary was set to the velocity inlet with w i = 5.16 m/s. In the case of the 2D model, it was set to conserve the mass flow rate in the radial direction and conserve the angular momentum in the tangential direction. The pressure outlet was set to the outer boundaries. The no-slip condition was used on the internal wall boundaries. The air was treated as both with constant density and as an ideal gas with energy equation.
Variable time stepping was used with a Courant number 0.5. A typical time step size was approximately 2×10 -8 s. As the flow was treated as transient, after reaching a quasi-static solution, time averaging was applied.
Several turbulence models based on a Reynoldsaveraged Navier-Stokes (RANS) equation, LES and laminar solver were used and compared to achieve results comparable with the experiment. Detailed mathematical description of the used models can be found in [16].
Simple two-equation models represented by k-ε and kω were chosen for their good accuracy for industrial applications. Those models determine a turbulent length scale and a time scale by solving two separate transport equations. The k-ε model is based on a transport equation for kinetic energy k and dissipation rate ε. In this paper, the RNG k-ε model for swirl dominant flows was used as it designed for swirling flows. The k-ω SST model with low-Re correction was chosen as the second two-equation model. This model combines the standard k-ω model for near wall treatment and the k-ε model in the free stream flow.
Reynolds Stress model (RSM) is the most advanced RANS model for the swirl dominant flows as it solves all the transport equations for the Reynolds stresses. This model was used with a low-Re and shear flow correction.
Large Eddy Simulation (LES) was represented by Wall-Adapting Local Eddy Viscosity model (WALE). However, it was performed on the same mesh as the RANS models, thus it is referred here as a Coarse LES. There were set no perturbations in the velocity inlet.

Mesh independence test
The all quad/hex structural meshes were created in Ansys Meshing for both 2D and 3D models. The 3D mesh used match control on the periodic boundaries. An average skewness of 2D mesh was 0.058 and an average aspect ratio of 1.18. The mesh independency test was made for the 2D model. The inlet pressure was found to be the most sensitive parameter to the number of cells. As the internal air-core blocks part of the exit orifice, it is necessary to capture the air-core boundary with high accuracy. Thus, the mesh needed to be fine near the interface. The effect of mesh resolution through the radius of the exit orifice is shown in Fig. 5. The error in the prediction of the inlet pressure is rapidly increasing if the number of cells is lower than 20. For the 3D simulation, a compromise between calculation speed and accuracy was chosen as the mesh with a resolution of 25 cells through the radius of the exit orifice was used. This 3D mesh had a total of 459,591 cells while the final 2D mesh contained 30 cells through the radius of the exit orifice with a total of 46,645 cells.

Results and discussion
In this chapter, the CFD simulations are compared with the experimental data in terms of the discharge parameters such as the discharge coefficient C D , SCA and atomizer efficiency , the air-core dimensions and the velocity profiles inside the swirl chamber. In the inlet port, the laminar flow can be assumed as Re p = 1021. However, it may change to transient or even turbulent flow inside the exit orifice as the velocity increases. Nevertheless, Chinn [14] foretold the laminar flow even for very high Re p as the swirl dominant flow tends to laminarise the flow itself and the atomizer is too small to develop a fully turbulent flow.
Initial CFD simulations were performed for both, the original and scaled atomizer. The dimensionless air-core diameters, dimensionless velocity profiles, SCA and C D were identical for both atomizers, therefore in this paper, the results are shown here in a dimensionless form and they are based on the CFD simulations of the originally sized atomizer.
The difference between incompressible and compressible gas phase was also investigated but no evident effect was found. This is in agreement with [17], where Mach number smaller than 0.3 ensures that the air phase can be safely handled as incompressible. Considering the maximal air velocity in our simulations, Mach number was much smaller than 0.1; thus, in this paper, the simulations are done with the constant air density.

Discharge parameters
An accurate prediction of the air-core dimensions is the crucial aspect as it significantly affects the discharge parameters. The discharge coefficient C D is an important discharge parameter defined as a ratio of the inlet mass flow to the theoretical mass flow rate: where A o is the cross-sectional area of the exit orifice. For the PS atomizers, the C D is relatively low as a portion of the exit orifice cross section is occupied by the air-core. The efficiency of conversion of the inlet pressure energy into kinetic energy at the atomizer exit is called the atomizer efficiency, and it is defined as: where v o is the velocity of the liquid sheet at the atomizer discharge. It was not possible to measure v o directly, thus the liquid sheet velocity had to be estimated from the high-speed records and it is measured approximately 1 mm downstream from the atomizer exit orifice. PIVlab software was used to estimate the velocity of the liquid sheet as it is shown in Fig. 6. The v o was obtained from the simulations in the same downstream distance as in the experiment. The swirl velocity component was neglected in both cases, only the axial and radial velocity components were considered.  In the experiment, the internal air-core had an almost constant diameter along the swirl chamber and extends in diameter in the exit orifice.
The two-equation turbulence models k-ε and k-ω were unable to predict the air-core inside the swirl chamber. The typical flow pattern from those models is shown in Fig. 7. As the air-core is not presented, all the discharge parameters are predicted with enormous error -see table 1. These models probably overpredicted the turbulent eddy viscosity which results in a decrease in the swirl velocity and consequently the air-core decayed. We exclude those models from further evaluations.
The RSM model was able to capture the air-core shape correctly. However, it underestimates its dimension. The smaller air-core leads to the higher C D as the smaller aircore block a smaller portion of the exit orifice. In the case of 2D RSM simulation, the air-core diameter was underestimated by 20% and by 10% in the case of the 3D model.
Both 2D and 3D laminar models and LES model were able to predict the inlet pressure, C D and the air-core dimensions closely. However, the 2D laminar model significantly overpredicted the velocity of the discharged liquid sheet v o and the atomizer efficiency. The LES predicted the atomizer efficiency most closely but it was still overpredicted about 20%. The reason may be that the CFD model neglected the surface roughness of the real atomizer.

Velocity profiles
The swirl velocity profiles are shown in Fig. 8. The twoequation models are excluded. The swirl velocity reaches its maximum at the air-core surface, which is typical for a Rankine vortex. In the experiment, it was six times greater than the mean velocity inside the inlet port. All the numerical models were able to predict the position of velocity maximum. Nevertheless, the 2D models were unable to capture the magnitude of the velocity peak accurately as the laminar 2D model overpredicted the swirl velocity maximum of about 15% while the 2D RSM underpredicted it about 15%. This correlates well with the dimension of the air-core from Table 1 where the 2D laminar model showed the largest dimensions of the internal air-core, while the 2D RSM performed opposite. Only small variations in the swirl velocity profiles were found at different axial distances from the top of the swirl chamber as it was shown in [1].
The radial and axial velocity profiles are shown in Fig  8 and Fig 9 respectively. The radial velocity is very low compared to the other velocity components. The 3D RMS and 3D laminar models predicted a positive local maximum inside the air-core while other models not. The 2D RSM, on the other hand, showed a negative local maximum in the position of r/r s = 0.2.
The axial velocity was almost zero up to the position of r/r s = 0.2 where it starts rapidly grow and reaches a maximum of 4×w i at r/r s = 0.05. Similarly, as in the case of swirl velocity, the 2D laminar model overpredicted the axial velocity magnitude while 2D RSM underpredicted it compared to 3D models.

Conclusions
The atomizer discharge parameters and the velocity profiles were predicted using various CFD models and the results were compared with the experimental data.
The surface tension and gas compressibility were neglected as they have no practical impact on the results.
The simple two-equation models k-ε and k-ω were unable to predict the air-core and the results are worthless. The 2D laminar model was able to closely predict C D and SCA but overestimated swirl velocity, liquid sheet velocity v o and η n .
The best results were obtained using the 3D laminar and 3D LES model, which give the same C D and SCA, but the LES model had a better prediction of η n in comparison with the experiment.
Further work will be related to the validation of timedependent variables as air-core instabilities and surface waves and also full 3D simulation will be included as the periodic boundary conditions assume that the centre of the air-core is fixed along the centreline of the swirl chamber.