Numerical Investigation on Spray Characteristics With Upstream Flow Pulsation of a Pintle Injector

The volume of fluid (VOF) model and the adaptive mesh refinement (AMR) method are used to study the spray characteristics of a gas–liquid pintle injector and the effects of mass flow pulsation of liquid on it. A pintle injector is a thrust adjusting device that changes the injection area using movable parts. Pressure pulsation in the supply pipeline is simulated by periodically changing the mass flow rate of the inlet. Spray characteristics with constant and pulsating upstream flows are compared with each other. The effect of frequency and amplitude of upstream liquid flow pulsation on the spray performance was studied. The results reveal that holding the mass flow rate of the upstream liquid flow constant, under the impact of gas flow, the liquid block, the liquid filament, and a large number of small droplets are peeled off from the liquid film. The film breakup position stays relatively fixed, and the spray has a conical shape. However, when the upstream liquid flow is pulsating, the film breakup position changes periodically, and the spray has a “Christmas tree”-shape. The pulsation frequency has little effect on the spray angle. But it strongly determines the droplet size and the spatial distribution of the spray. In addition, the pulsation amplitude can enhance the phenomenon of “Christmas tree.” With the increase in pulsation amplitude, the liquid film at the outlet of the pintle injector appears with a periodic phenomenon of “contraction–expansion.”


INTRODUCTION
Liquid rockets have the advantages of high specific impulse, high thrust, and repeatable ignition (Yang and Anderson, 1971;Ferraiuolo and Riccio, 2019;Khan and Qamar, 2020;Liang et al., 2021). In recent years, with the development of manned lunar exploration project and reusable rocket technology, large-scale variable-thrust liquid rockets have become a hot spot for investigation. The pintle injector is selected as the best choice for the injection system (Zhang et al., 2020). The most famous pintle injector engine is the Lunar Module Descent Engine (LMDE) used in the American "Apollo program," with a maximum thrust of 44.52 kN and a 10:1 thrust adjustment capability (Dressler, 2013). At present, several SpaceX engines also use pintle injectors, and its engine variable thrust performance is in a leading position in the world (Dek et al., 2020).
In liquid rocket engines, the atomization of the liquid propellant can promote evaporation, mixing, and combustion. Since the 1980s, researchers from the United States (Dressler and Bauer, 2000;Gamertsfelder et al., 2020), South Korea (Son et al., 2015;Son et al., 2016a;Radhakrishnan et al., 2018), China (Zhang et al., 2014;Fang and Shen, 2017;Chen et al., 2019;Wang et al., 2020), Japan (Sakaki et al., 2015;, and other countries have carried out a lot of investigations on the spray characteristics of pintle injectors, mainly focusing on the angle and the breakup process of the spray, together with size, velocity, mixing features, and other aspects of the droplets. The pintle injector originated from the Jet Propulsion Laboratory (JPL) of the California Institute of Technology. It was first used to investigate the reaction rate of the hypergolic propellant (Dressler and Bauer, 2000), and then the pintle injector was used in the liquid rocket engine by TRW. Wang et al. (2020) studied the breakup and mixing processes of the pintle injector unit. Two different spray patterns were observed. Chen et al. (2019) studied the spray boundary, droplet size distribution, and flow field characteristics under different momentum ratios and found that the spray characteristics were also affected by the local momentum ratio. Fang and Shen (2017) obtained the effects of the structure parameters on SMD, droplet size distribution uniformity index, and spray angle. Sakaki et al. (2015) and Sakaki et al. (2018) used the ethanol/liquid oxygen pintle injector to evaluate the performance of a rocket engine. Nardi et al. (2015) also designed a pintle injector engine and conducted experiments. Min et al. (Son et al., 2015) experimentally studied the effects of momentum ratio and Weber number on the spray angle. They found that in addition to the lower liquid injection pressure drop, the spray angle was proportional to the momentum ratio and Weber number under all experimental conditions. Subsequently, a combination of experiments and numerical simulations were used to analyze the velocity distribution, spray angle, and the droplet spatial distribution. The discrete phase model (DPM) was used to explore the effects of injection conditions and geometric parameters on spray characteristics, but the model has certain flaws and cannot capture a detailed breakup phenomenon (Zhang et al., 2020).
Combustion instabilities caused by the coupling of heat release and acoustic pressure in the combustor may generate intense pressure fluctuations and lead to excessive heat that transfers to combustor walls and injector plates Yang et al., 2021). Although tremendous human, material, and time resources have been invested by many countries since 1940s to determine ways to solve this problem, it has not yet been resolved. There are two significant characteristics of pressure fluctuations in the unstable combustion chamber: large amplitudes and obvious periodicity (Miller et al., 2007;Pomeroy and Anderson, 2016).
Propellant atomization is the basis of liquid rocket engine combustion. The length of liquid film breakup and droplet size distribution basically determine the combustion process, efficiency, and stability (Yang and Anderson, 1971). During the practical work of the liquid rocket engine, the liquid mass flow in the propellant supply pipeline has a certain form of pulsation, which directly affects the liquid film spray characteristics. The atomization subjected to pressure fluctuations may have a similar periodicity with which combustion instability occurs. It might be coupled with the combustion process under certain conditions, which may lead to pressure oscillations and unstable heat release in the combustion chamber (Yang et al., 2018).
The injection pressure oscillations accompany the mass flow oscillations that modulate the velocities of the liquid with time.
Meanwhile, the faster droplets chase the slower droplets. A large proportion of droplets are superposed together at a certain point. This phenomenon called the "klystron effect" also appeared in impinging jet injectors (Heister et al., 1997) and swirl injectors (Khil et al., 2012). Mulmule et al. (2010) studied the crushing process of the circular liquid film formed by the front collision of two jets under sound pressure pulsation and found that the diameter of the circular liquid film and the average droplet size decreased. Ahn et al. (2012) developed mechanical pulsators and obtained the periodically changing air core diameters inside the transparent vortex chamber and the nozzle, as well as the spray near the exit nozzle using a highspeed camera. Yang and Turan (2017) studied the effect of disturbance frequency and amplitude on jet atomization. It is found that the jet is insensitive to high-frequency disturbance. The amplitude of disturbance has an important effect on droplet number and droplet size.
In summary, the steady spray characteristics of the pintle injector have been studied sufficiently. Meanwhile, upstream flow pulsation has a great influence on spray. However, the unsteady spray characteristics of the pintle injector have not yet been analyzed fully. Therefore, it is necessary to deeply understand the spray characteristics of the pintle injector under the condition of mass flow pulsation from the upstream propellants. In this study, a numerical study on the unsteady atomization process of a pintle injector is carried out. Effects of pulsation in the supply system on the spray performances are deeply analyzed.

Numerical Methods
The liquid film breaks up near the pintle injector outlet, and the droplets are very dense. In this study, the commercial software FLUENT is used for numerical calculation. Therefore, this study adopts the realizable k-epsilon model and the VOF model, combined with the adaptive mesh refinement (AMR) method to capture the detailed gas-liquid interface. In this numerical study, the unsteady pressure-based solver is selected. The separation scheme, the pressure-implicit with splitting of operators (PISO), and the second-order spatial discretization are employed.

VOF Model
The VOF model is utilized in the article to capture the atomization process and the interaction of liquid with gas (Dinic and Sharma, 2019). The advantage of the VOF model is that there is no need to deal with physical phenomena such as interface fragmentation and aggregation because these topological changes are implicitly included in the VOF model. Therefore, the VOF model is very suitable for atomization numerical calculation of liquid film fragmentation and aggregation.
In the VOF model, the volume fraction of each phase is defined as α q , and the volume fraction of all phases always sum to unity in every control volume. In the current study, air and water are the primary and second phases, respectively. So, there are only three situations in the control volume: the cell filled with air α 1 1, the cell filled with water 0 < α 1 < 1, the cell filled with water and air, where α 1 is the liquid volume fraction in the control volume.
When α 1 is eligible 0 < α 1 < 1 in the control volume, the tracking of the gas-liquid interface is achieved by solving the continuity equation. For the water, the equation is written as follows: where ρ is the average density and u is the velocity in the control volume. The equation is solved by explicit time formulation.
The velocity u is shared among water and air in each control volume and is obtained from momentum equation, which is given by the following equation: where p is the pressure; g represents gravity, which is neglected in this article; and F is surface tension. ρ and μ are the average density and average dynamic viscosity in the control volume, respectively, which are listed as follows: where ρ l and ρ g are the density of water and air, respectively. μ 1 and μ g are the dynamic viscosity of water and air, respectively. The same as the momentum equation, the energy equation also shared among phases, which is given as follows: where k eff represents effective thermal conductivity. T and E are the temperature and specific total energy, respectively.

Adaptive Mesh Refinement Method
For flow problems involving complex shapes and free surfaces, sufficient grid resolution is required to ensure computational accuracy and stability. The significant advantage of the adaptive mesh refinement (AMR) method is to reduce calculations without loss of calculation accuracy (Xiang et al., 2020;Li and Soteriou, 2014). An octree is a tree-like data structure that describes a threedimensional space. The use of an octree adaptive mesh for discretization can make dynamic adaptive algorithms simple and flexible. The AMR can automatically refine or coarsen the mesh at the intersection of the two phases. When the liquid block leaves the grid position, the hanging grids are automatically removed, which greatly reduces the calculation time and storage space, and improves the calculation efficiency. Figure 1 is a schematic diagram of the encryption process, where the grid elements are represented by cube elements. In this article, the gradient of the liquid phase volume fraction is used to adapt the grid at the gas-liquid interface, using two levels of adaptive grids to capture more subtle fluctuations at the gas-liquid interface.

VOF-to-DPM
However, it is very challenging to sample droplets through the VOF model at a given computational domain. In order to facilitate the statistics of droplets, the VOF-to-DPM (Jie and Peng, 2021) is studied and analyzed to study droplets. Due to the transition mechanism of the VOF-to-DPM, droplets can also be sampled at any given location in the simulation. Moreover, droplet properties such as volume, equivalent diameter, velocity, and position can be determined with the DPM. Since there are very few studies on the atomization process under pulsation conditions using the VOF-to-DPM, which cannot verify the effectiveness of VOF-to-DPM in capturing pulsation atomization, the VOF model is still adopted to study the spray characteristics with mass flow pulsation of liquid, and the VOF-to-DPM is only used to count the droplet size. The near-spherical lumps are transferred from the VOF-to-DPM point masses based on two criteria on the lumps. The first criterion is based on the droplet size, and the droplets must be between the minimum and maximum volume-equivalent diameter. The specific formula is as follows: where V lump is the volume of the VOF liquid block (i.e., "lump") and R is the defined conversion radius. In this article, the VOF-to-DPM is mainly used to count the particle size of droplets. In order to convert droplets into particles as much as possible, the minimum and maximum conversion radius are set as 0 and 0.5 mm, respectively. The second criterion is based upon the droplet shape to support more coherent and smooth conversions. This criterion ensures that small irregular droplets or filaments are not transferred. The shape of the liquid lump is calculated in two ways: standard deviation of the normalized radius and average radius-surface orthogonality. The maximum asphericity S 1 is calculated from the standard deviation of the normalized radius. As shown in Figure 2, the distance between the facet center and the center of gravity of the lump is calculated for every facet of the lump surface. These distance values are weighted by the surface area of the boundary facet. The specific formula is as follows: where n represents the surface element, m is the number of facets, and r′ is the average radius of the lump. The other parameters are shown in Figure 2. Then, the standard deviation of these weighted distances is computed and normalized by the average radius of the lump. Only droplets with asphericity below the specified maximum value can be converted. The specific formula is as follows: The maximum asphericity S 2 is calculated based on average radius-surface orthogonality. Specifically, a vector is drawn from the center of gravity of each lump to the central point of the lump boundary, for every facet on the lump surface. Each vector is normalized and used to calculate a dot product with the unit normal of the facet, yielding a measure of the relative orthogonality that ranges from 0 to 1. In this article, these asphericity standards are both set to 0.5. The specific formula is as follows: If the liquid lump satisfies both size and shape criteria, the liquid lump is removed from the VOF model and converted into particle in the Lagrange system. All locally adaptive mesh refinements used to solve the gas-liquid interface are automatically restored so that the Lagrange particle packets can be placed in a single-large grid. If the volume of the liquid mass is much larger than that of the coarser mesh, the liquid mass is converted into as many Lagrange particle packets as possible. Kim et al. (2007) and Chadha et al. (2021) employed a coupled VOF and DPM approach to model the atomization of fuel within a gas turbine injector, and the numerical results fit well with the experiment and can be used to predict the experimental results.

Pintle Injector and Boundary Conditions
The structure of the pintle injector is shown in Figure 3, which is the same as that in the reference (Radhakrishnan et al., 2018). Specific dimensions are shown in Table 1. The liquid is injected from the inner channel. Then, the thin liquid film is formed and interacted with the high-speed gas in the outer channel. After that, the liquid film breaks up and tears down into small droplets.   Air and water are set to the primary and secondary phases in the simulation, respectively. The computational domain is shown in Figure 4. The minimum mesh cell size of the liquid out cross section is 40 μm. The pintle injector is axisymmetric. In order to reduce the calculation, a quarter of the spray area is selected as the computational domain. Meanwhile, the following mass flow rate of gas and liquid is set to a quarter of the true value. Air is treated as an ideal gas with a constant mass flow rate of 0.759g∕ s at its inlet. Water is considered an incompressible phase with a density of 998.2 kg/m 3 . The water inlet boundary conditions are compiled using UDF and set as follows: m m o [1 + m′sin(2πft)]. In the formula, m is the actual inlet mass flow rate of liquid, m 0 is the stable mass flow rate of liquid and set as 15g/s. m′ (the pulsation amplitude of mass flow rate of liquid) is defined as the ratio of mass flow rate of pulsation to the stable mass flow rate of liquid; f is the pulsation frequency of mass flow rate of liquid.
The outlet boundary is set as a pressure outlet, and the backflow volume fraction of water is set to 0, which means that the effect of backflow of water from outside the computational domain is not considered. The other boundaries are set as non-slip wall boundary conditions. The calculated Reynolds number and turbulence intensity are R e 3194 and I 0.0583, respectively, and the surface tension coefficient is 0.072, the hydraulic diameter of liquid inlet is 1.55, and the hydraulic diameter of gas inlet is 1.
In order to study the effect of pulsation frequency on spray, the pulsation amplitude is set to 0.15, and the specific frequency setting are shown in Table 2. In addition, the effect of pulsation amplitude on spray is studied. The frequency is set as 2,000 Hz, and the specific pulsation amplitude settings are shown in Table 2.

Physical Principles Validation
The pressure at the minimum liquid outlet cross section (as shown in Figure 3) is calculated by CFD. The local pressure is much higher than the "saturated vapor pressure (2.33KPA)" at room temperature. So, the "cavitation" phenomenon could not occur, and the VOF model could be used for numerical analysis.
At the same time, the velocity of water with the maximum mass flow rate of liquid is calculated. The velocity of water at the minimum liquid outlet cross section is 27.1 m/s, and the air velocity at the gas outlet of the pintle injector is 176.26 m/s. They are both less than the local sound speed, so the gas and liquid channels are believed to be unchoked, which is in line with the design principle of the pintle injector (Son et al., 2016b).

Numerical Model Validation
Since experimental research about the effects of upstream pulsation on spray performance of a pintle injector has not been published, in the current study, to verify the accuracy of the numerical model, two validation cases are carried out. The verification process is shown in Figure 5.
The model was verified by the pintle injector spray under steady flow conditions, namely, m l 29.5g/s and m g 3.036g/s (case 2 in the reference (Radhakrishnan et al., 2018)). Various convergence tests were carried out to evaluate the accuracy of the numerical results. Grid independence tests were carried out as well. Numerical model based on three different meshes (302,158 cells, 610,273 cells, and 926,574 cells) were compared with each other, as shown in Figure 6. Although the coarse grid can simulate the appearance of the spray field, it cannot capture the details of the liquid film breakup, which cannot meet the research needs. The intermediate and fine grids can accurately capture liquid filament peeling and a large number of droplets, improve the fineness of the atomization process, and obtain a more realistic droplet distribution. The spray angle is measured and compared with the results in the references (Radhakrishnan et al., 2018). The error for all three angles are less than 10%, which shows that the model can simulate the spray angle effectively. The droplet particle size is counted, and the SMD is calculated with an error of less than 15%. Considering the actual computing resources and ensuring the accuracy of the calculation, a mesh with 610,273 grids is finally selected. Comparing the results of the steady spray with the experimental and numerical results of the reference (Radhakrishnan et al., 2018), as shown in Figure 6A, the analysis of the overall spray shape and the spray angle is basically consistent, which can verify the accuracy of the gas-liquid interaction model under steady flow condition.
Then, the upstream pulsation model is verified by using the pulsation boundary conditions with a cylindrical jet, which is the same as that in reference (Srinivasan et al., 2011). The numerical result is shown in Figure 7B, and the overall performance of the jet surface fluctuation pattern is in good agreement with results as shown in Figure 7A from reference (Srinivasan et al., 2011). It shows that the boundary condition of liquid inlet compiled by UDF can effectively simulate the flow pulsation of liquid. So it could be concluded that both the VOF method and the pulsation model in the current article are validated.

Spray Characteristics With Steady Flow
Case 1 (f 0Hz) is selected to study the spray characteristics without upstream pulsation. The liquid film breakup process is shown in Figure 8. The numerical results effectively capture the liquid film fluctuation process and the whole breakup process. Due to the action of the aerodynamic shear, the liquid film falls off the large liquid block and filament. Small droplets gather around the large liquid block, and the liquid blocks further breakup to a large number of liquid droplets under the action of aerodynamic force again. The overall spray pattern is a stable cone-shaped structure. Figure 9 shows the streamline and velocity contour of case 1(f = 0). There are recirculation zones in the central region of the spray field. The outside droplets enter the central recirculation zone. After the gas and liquid collide, the velocity rapidly decays from 184 m/s at the gas outlet of the pintle injector to about 70 m/s. During this gas-liquid colliding process, there is a momentum exchange between the gas and the liquid. The impact of the gas on the liquid film makes the liquid film fluctuate. As the liquid film develops downstream, the thickness of the liquid film rapidly decreases, which promotes the breakup of the liquid film, as shown in Figure 10.

Spray Characteristics With Unsteady Flow
Case 3 (f = 1,000 Hz) is selected to reveal the spray characteristics under unsteady flow conditions, and the liquid film breakage pattern is shown in Figure 11. Different from the liquid film breakup with steady upstream, upstream pulsation promotes the liquid film breakup. There are two main reasons of liquid film breakage. First, the gas force and surface tension work together. The liquid film produces some irregular holes, and small droplets are generated near the holes (Chu et al., 2020). The size and the number of holes are further developed. Finally, the liquid film breaks up; meanwhile, the upstream fluctuation increases the instability of the liquid film and promotes the breakup, resulting in liquid filaments and liquid blocks earlier. Under the action of aerodynamic force, the liquid filaments and liquid blocks broke up further and produced small droplets. Statistics of the spray angle at different moments show that the spray angle does not change much in one cycle. The spatial distribution of the spray stays relatively stable. Compared with the spray under steady conditions, upstream pulsation can enlarge the spatial distribution of spray.  Frontiers in Aerospace Engineering | www.frontiersin.org May 2022 | Volume 1 | Article 876191 8 Figure 12 shows the evolution of the spray in one cycle for case 3 (f = 1,000 Hz). Three obvious periodic droplet groups appear. The droplet group presents a shape of circular ring in space. With the downstream development of spray, the ring gradually becomes larger, so the number of droplets per unit space volume decreases and becomes sparse. The entire spray shows like a" Christmas tree." The main reason is that the upstream flow pulsation causes periodic change of liquid film velocity at injector exit. The faster droplets can catch up slower droplets, so the droplet concentration region appears. Moreover, because of the radial velocity, the droplet continues to radially diffuse, and finally present the shape of the "Christmas tree." This phenomenon is the specific embodiment of "klystron" (Kang et al., 2016a;Kang et al., 2016b) effect on the pintle injector.
Auxiliary lines are drawn on the spray cycle images, as shown in Figure 12. It can be seen that the length of the complete liquid film above the red line remains nearly fixed, while the broken position of the liquid film below the red line changes periodically. This is obviously different from the basically fixed position of liquid film breakup under steady upstream. The main reason is that the liquid film just coming out of the injector is still relatively thick and cannot be broken immediately. As the liquid film develops downstream, the thickness becomes thin, and when the broken boundary conditions of the liquid film are reached, the fragmentation occurs. At the same time, due to upstream pulsation, the liquid film velocity also changes periodically. When the liquid film velocity is small, the momentum of the liquid film is small and the liquid film is easy to be broken by the gas impact, so the breakup position of the liquid film is relatively close to the upper. When the liquid film speed is high, the momentum of liquid film is big and the liquid film moves rapidly downward, so the breakup position is relatively later. When the liquid film velocity changes periodically, the breakup position of the liquid film also changes periodically. In the spray region between the red line and yellow line, it is obvious that the droplet space quantity density changes periodically, which is also caused by upstream pulsation. The periodic change of space quantity density of droplet might lead to the periodic change of combustion heat release, which may be not conducive to the combustion stability of the liquid rocket engine.

Spray Pattern
When the liquid rocket engine is working, the violent vibration of the rocket causes the vibration of the propellant pipe, which can  disturb the stability of mass flow rate of upstream propellant, affecting the spray characters (Yang et al., 2018). This part analyzes the effects of the pulsation frequency of upstream liquid propellant on the spray characteristics. The pulsation amplitude m′ is set to 0.15, and the pulsation frequency f is set to 500, 1,000, 2,000, and 3,000 Hz in order. Figure 13 shows the corresponding spray patterns under different frequencies. It can be seen that there are periodic accumulation of droplets, indicating that the atomization of the pintle injector can respond to a wide range of pulsation frequency of upstream liquid propellant. With the increase of the pulsation frequency, the axial distance of the droplet groups decreases, and the number of droplet groups appearing increases. This indicates that the frequency of liquid block and droplets falling off the liquid film is increasing, and the falling off frequency is consistent with the pulsation frequency of upstream liquid propellant.

Cross-Sectional Mass Flow Rate
The mass flow rate of the liquid through the surface is monitored at different cross-sectional positions as shown in Figure 4, and the result is shown in Figure 14. From the aspect of waveform of the flow curve of the monitoring surface, under steady flow conditions, the mass flow rate at A, B, and C surfaces is almost constant, and the mass flow rate at D and E surfaces fluctuates evenly up and down the average value. It can be concluded that there is no significant difference in the velocity of droplets passing through the same location, the droplet is not locally gathered, and the droplet distribution in the spray field is relatively uniform. When there is a flow pulsation, the flow curve on A surface is sinusoidal, the flow curve on B and C surfaces is slightly deformed, and the flow curve on D and E surfaces turns into a steep peak waveform. The fluctuation value above the sinusoidal curve is much larger than that below that, and the value According to the Rayleigh criterion (Rayleigh, 1878), when the phase of combustion chamber pressure and the combustion exothermic oscillation are consistent, the instable combustion may be excited. The droplet accumulation in the spray field may cause spatial oscillation of the combustion exotherm, so the mass flow pulsation of upstream liquid propellant at the special pulsation frequency may stimulate the combustion instability. When the liquid rocket is unstable at low and intermediate frequencies, the pressure curve at the position of the injector of the combustion chamber is sawtooth, and there is a correlation between the steep peak wave flow curve and the sawtooth pressure curve. Therefore, mass flow pulsation may be one of the reasons for the unstable combustion of liquid rockets at low and intermediate frequencies.
From the aspect of the fluctuation frequency of the flow curve of the monitoring surface, the frequency response of the mass flow rate of the passing monitoring surface to the upstream pulsation is consistent with that. However, there is a phase difference (or a time lag τ) between the response and the upstream pulsation, as shown in Figure 14D. It can be seen from the figure that the pulsation frequency has great effects on the phase difference. Take the phase difference of the flow curves of the A and C surfaces as an example. When the frequency is 500 and 1,000 Hz, the phase difference is small; when the frequency is 2,000 Hz, the phase difference is about 45°; and when frequency is 3,000 Hz, the phase difference increases to around 160°. In conclusion, the magnitude of the phase difference is positively correlated with the pulsation frequency.
The pulsation frequency of upstream liquid not only affects the phase difference but also affects the peak value of the flow curve. Taking the C surface flow curve as an example, as the pulsation frequency increases, the peak value of the flow curve becomes significantly larger. Compared with the frequency of 1,000 Hz, the peak value under the conditions of 2,000 and 3,000 Hz is increased by 27.8 and 73%, respectively. It shows that the higher the pulsation frequency, the more obvious the droplet groups aggregation phenomenon, which can also be seen from Figure 13. Most of the flow curves on the D and E surfaces also conform to this trend. However, the peak value of the E surface in Figure 14E is smaller than the upstream flow peak. The main reason is that the liquid film produces droplet groups earlier under high-frequency conditions. The droplet groups change Frontiers in Aerospace Engineering | www.frontiersin.org May 2022 | Volume 1 | Article 876191 12 sparsely before reaching the E surface so the peak flow of the mass rate on the E surface is reduced.

Droplet Size
The droplet size distribution is analyzed, and the quantitative effect of pulsation frequency on droplet size is shown in Figure 15. Macroscopically, upstream pulsation can significantly reduce the droplet size. As the frequency of upstream pulsation increases, SMD, average droplet size, and median droplet size all increase. The main reason is that with the increase of the frequency, the same spray area can produce more droplet groups. After the same gas impact, small droplets with high pulsation frequency are more likely to gather and produce large droplets. What is more special is that the SMD value with a pulsation frequency of 3,000 Hz is greater than all other cases, probably because a few large droplets appear, so the SMD value increases sharply. Figure 16 shows the effects of pulsation frequency on the probability distribution of the droplet size. When the upstream pulsates, the proportion of large droplets decreases while the proportion of small droplets increases, and the mean value of the probability distribution curve is moving to the right. It also can be concluded that under the action of upstream pulsation, the average droplet size decreases. Meanwhile, the droplet size increases with the increase of the pulsation frequency, but it does not exceed the droplet size under steady state flow condition. The main reason is that upstream pulsation enhances the instability of the liquid film. After the gas impact, liquid film breakup earlier, and the time of droplets subjected to gas-liquid shear becomes longer, and the droplets become smaller. The numerical results show that when the upstream pulsates, the number of droplets in the local area of the atomization field increases significantly and the droplet size decreases.

Droplet Spatial Position
The spatial position of the droplets is displayed in the coordinate system and projected. The effects of pulsation frequency on the droplet space distribution are shown in Figure 17. In the spray region, the number of "bulging" structures is 0, 1, 2, 4, and 7 in order. The number of "bulging" structures is linearly related to the pulsation frequency, indicating that the number of "bulging" structures response to the upstream pulsation frequency is linear. It can be clearly seen that the central region of the spray zone has smaller droplets, and the outer droplets are larger, indicating that the internal recirculation region plays a role in improving the quality of atomization.

Spray Pattern
This part analyzes the effects of the pulsation amplitude of upstream liquid propellant on the spray characteristics. The pulsation frequency f is set to 2,000 Hz, and the pulsation amplitude m′ is set to 0.01, 0.05, 0.15, 0.2, 0.3, and 0.5 in order. Figure 18 shows the corresponding spray patterns. It can be seen that when the pulsation amplitude is very small (m′ 0.01), there is no droplet group in the spray area. There is no significant difference between the spray field with upstream pulsation and that under steady upstream. It can be considered that the pintle injector works under "steady" condition at this situation. Only when the pulsation amplitude is higher than a critical value, the spray field can change significantly. It can be seen from Figures 18B-E that with the increase in pulsation amplitude, the phenomenon of periodic droplet groups enhances gradually, and the liquid film deformation near the injector exit is more severe, indicating that the increase of pulsation amplitude can make "klystron" effect more significant. The number of droplet groups remains constant, indicating that the frequency of liquid block and droplet falling off from the liquid film is not affected by the pulsation amplitude.  In Figure 18, the auxiliary line with a fixed angle is made, and it can be seen that the overall spray region does not change significantly, and the spray angle is almost constant. However, the spray shrink toward the central region with the pulsation amplitude increases, resulting in the increase of droplet space density in the central region. The main reason is that the liquid  velocity is relatively smaller when the liquid flow is less than the average flow in the cycle. After being impacted by gas, the liquid film is quickly blown down to the center region of the spray field, resulting in an increase in droplets in the center region, and then the overheating of the pintle injector tip may occur in the actual work of the rocket engine. In previous experiments (Vasques and Haidn, 2017), the problem of erosion of the pintle injector tip had occurred.

Liquid Film Breakup Process
The dominant feature in the early stages of "Christmas tree"shaped spray is the periodic "contraction and expansion" of the liquid film near the pintle injector exit. This phenomenon becomes more pronounced as the amplitude increases. Case 10 (m′ 0.5, f 2000Hz) with more obvious phenomenon was selected for in-depth analysis. A complete sequence of the formation of "contraction and expansion" phenomenon is explained in Figure 19. The period from 0 to 4/8T is a "contraction" process, and the period from 4/8 to 1T is an "expansion" process. When the mass flow rate of liquid is in a decreasing stage, the liquid velocity also decreases, and the liquid film would be squeezed by the gas flow more, resulting in the liquid film to shrink inward. When the liquid flow is in the rising stage, the liquid film speed increases, the liquid film expands outward, and the broken position of the liquid film also moves downward. The formation of the next "contraction and expansion" process of the liquid film begins to repeat. This section elaborates on the pulsating circulation of the liquid film located near the needle bolt injector outlet, as it has an importanan effect on the downstream spray pattern.

Cross-Sectional Mass Flow Rate
The effect of pulsation amplitude on the mass flow rate of liquid at the monitoring surface is shown in Figure 20. Under the condition of small amplitude (m′ 0.01), the mass flow curve of liquid has no obvious fluctuation period, and the flow curve is similar to that under steady condition. As the amplitude (m′ 0.05 − 0.5) increases, the pulsation cycle of the flow curve becomes apparent, yet each change trend of the flow curve at different monitoring surfaces is different. The peak value of the flow curve of A, B, and C surfaces increased with the increase in pulsation amplitude, while the peak value of D and E surfaces increases first and then decreases with the increase in pulsation amplitude, and the detailed analysis of each flow curve of the mass flow rate at different monitoring surfaces is presented in the following.
A surface calculates the liquid mass flow at the liquid inlet, so the flow curve is sinusoidal with different amplitudes. The mass flow rate of liquid at the injector liquid outlet was calculated through B surface. When the amplitude was between 0.01 and 0.2, the flow curve of B surface was very similar to that at A surface. When the amplitude is 0.3 the waveform of the flow curve at B surface changes partially. When the amplitude is 0.5, the waveform at B surface changes obviously, and the peak value increases 66.7% compared with that of the amplitude of 0.3. It shows that the fluctuation of the liquid film at the injector outlet becomes more obvious with the increase in amplitude. Compared with the flow curve at B surface, the flow curve at C surface is deformed more strongly with the same amplitude.
In the amplitude range of 0.01-0.5, the peak value of flow curves of D and E are different from those of A, B, and C surfaces. The peak flow at D surface increases between amplitudes of 0.01 and 0.2, and decreases between amplitudes of 0.3 and 0.5. The flow curve of E surface also conforms to this trend, but the amplitude corresponding to the peak transition is different from that of D surface. The flow peak decreased, indicating that the degree of droplet aggregation decreased. In the decreasing range of the peak value, the higher the amplitude is, the more obvious the liquid film swings at the injector outlet, and the more difficult it is for droplet groups to maintain the form of aggregation in the downstream. E surface is located downstream of D surface, so the droplet group distribution at E surface is more dispersed, and the amplitude corresponding to the drop point of peak flow at E surface is smaller than that in D surface.
The changing trend of the flow curve at D and E surfaces of case 10 (m = 0.5) conforms to the previous analysis, but the value of the mass flow rate is smaller than that of other amplitudes. To explain this phenomenon, a parameter is defined: capturing rate of the mass flow rate η average mass flow rate at C (D or E) surface average mass flow rate at A surface, where the average mass rate at C (D or E) is obtained by calculating the average mass flow rate at C, D, and E surfaces with different amplitudes in Figure 20. The results are shown in Figure 21. It can be seen from Figure 21 that η at E surface of case 10 (m′ 0.5) is small and much liquid is lost. The main reason for mass flow loss of liquid is that the mesh size of the downstream part is larger than that of the upper part, and the droplet size of the downstream part is smaller, so the mesh at the downstream part cannot capture too small droplets. The following part analyzes the particle size, and it is found that the droplet size of case 10 is the smallest.
The reason why case 10 (m′ 0.5) tends to appear as smaller droplets is mainly due to the large amplitude of liquid mass flow. When it is near the trough of the mass flow curve, the mass flow rate of liquid is very small, and it directly breaks up into small FIGURE 21 | Capture rate of the mass flow rate of the liquid at C, D, and E surfaces with different amplitudes (the pulsation frequency is set as 2,000 Hz).
Frontiers in Aerospace Engineering | www.frontiersin.org May 2022 | Volume 1 | Article 876191 droplets due to the impact of gas flow. Very small droplets cannot be captured effectively, so part of the liquid mass is lost.
The η at C surface is lower than that of other working conditions, and it can be inferred that the position of small droplets produced when the amplitude is larger is closer to the upstream position.

Droplet Size
The quantitative effects of pulsation amplitude on particle size are shown in Figure 22. Macroscopically, pulsation amplitude has great effects on droplet size. With the increase in pulsation amplitude, SMD, mean particle size, and median particle size all decrease. The reason is that liquid film instability is enhanced, and liquid film breakup is more likely to occur because of mass flow pulsation. Moreover, the droplet size decreases with the increase in pulsation amplitude. When the mass flow rate is greater than the average mass flow rate, the liquid film velocity is bigger and the gas-liquid interaction is stronger, which is conducive to the liquid film breakup. When the mass flow rate is lower than the average mass flow rate, the liquid film velocity is smaller, and it is easier to be impacted by the external airflow and broken into droplets. Therefore, the presence of upstream pulsation reduces the droplet size, and the larger the pulsation amplitude is, the more obvious the droplet size decrease is. Figure 23 shows the effects of pulsation amplitude on the probability distribution of droplet size. With the increase in pulsation amplitude, the proportion of large droplets decreases and that of small droplets increases. It shows that under the action of upstream pulsation, the average droplet size decreases. On the other hand, with the increase in pulsation amplitude, the distribution curve becomes thin and high, indicating that the droplet particle size distribution becomes concentrated and uniform.
In conclusion, the increase in pulsation amplitude is beneficial to the decrease in droplet size, and the droplet size distribution is relatively concentrated. In terms of droplet combustion heat release, this may be conducive to the uniform heat release of droplet combustion and reduces the pressure difference at different positions in the combustion chamber.

CONCLUSION
This investigation uses numerical methods to study the spray characteristics of the pintle injector with and without pulsation in the upstream liquid propellant, and the main conclusions are as follows: (1) Using the volume of fluid (VOF) model and the adaptive mesh refinement (AMR) method to effectively capture the liquid film fluctuation process and the entire atomization process of the pintle injector under steady conditions, the overall spray shape is a cone structure, and the existence of the inner central recirculation region improves the droplet mixing effect. (2) Under unsteady flow conditions, the droplet morphology is more diversified, the breakup position of the liquid film changes periodically, the spray angle does not change much in one cycle, the overall spray region is relatively stable, the central spray area expands outward, and the overall spray shape appears as a "Christmas tree." (3) With the increase in the pulsation frequency, the axial distance of the droplet groups decreases, and the number of droplet groups appearing increases, and the falling frequency of large liquid blocks is consistent with the pulsation frequency. The response of the mass flow rate of the monitoring face to the upstream flow pulsation is FIGURE 22 | SMD, mean, and median diameter at different pulsation amplitudes (the pulsation frequency is set as 2,000 Hz).
FIGURE 23 | Probability density function (PDF) distribution of droplet diameter at different pulsation amplitudes (the pulsation frequency is set as 2,000 Hz).