Numerical Investigation of the E ff ect of Sudden Expansion Ratio of Solid Fuel Ramjet Combustor with Swirling Turbulent Reacting Flow

In this paper, the effect of sudden expansion ratio of solid fuel ramjet (SFRJ) combustor is numerically investigated with swirl flow. A computational fluid dynamics (CFD) code is written in FORTRAN to simulate the combustion and flow patterns in the combustion chamber. The connected-pipe facility is used to perform the experiment with swirl, and high-density Polyethylene (HDPE) is used as the solid fuel. The investigation is performed with different sudden expansion ratios, in which the port and inlet diameters are independently varied. The results indicated that the self-sustained combustion of the SFRJ occurs around the reattachment point at first, and then the heat released in reattachment point is used to achieve the self-sustained combustion in the redevelopment zone. The average regression rate is proportional to the sudden expansion ratio for the cases with a fixed port diameter, which is mainly dominated by the enhancement of heat transfer in backward-facing step. However, the average regression rate is inversely proportional to the sudden expansion ratio for the cases with fixed inlet diameter, which is influenced by the heat transfer mechanism of developed turbulent flow in the redevelopment zone.


Introduction
Nowadays, solid fuel ramjet (SFRJ) has been widely used and researched in the military and aerospace fields, due to the advantages of simple structure, high reliability, and easy long-term storage.Thus, how to improve combustion performance is a hotspot of current researches.As can be seen from previous researches, the combustion performance of SFRJ could be highly determined by the regression rate, and it is easily influenced by the geometry of SFRJ combustor.
The regression rate and flow field characteristics of SFRJ and scramjets has been extensively investigated with experimental and theoretical studies.An experimental investigation on SFRJ is performed by Schulte [1], with different altitudes, flight speeds, and mass fluxes of in-coming air, to determine the fuel regression behavior.Polyethylene (PE) was used as the standard fuel.The results show that decreasing the port diameter, chamber pressure, and inlet temperature will increase the regression rate.Gany et al. [2] investigated the combustion that is affected by the geometry of Polymethylmethacrylate (PMMA) fueled SFRJ.The research founded that the local regression rate could be sensitive to the local convective heat flux and it was significantly affected by the port diameter.Ferreira et al. conducted experimental research on Polyethylene (PE) [3], at which the port diameter is varied from 20 mm to 45 mm.They found that changing port diameter has a significant effect on the regression rate.Reference [4] numerically investigated the combustion performance of the  The experimental procedure is described in detail in the paper [22].The repeatability of the connected pipe facility and the accuracy of the experimental result are shown and discussed in paper [22].

RANS Equations
The Reynolds-averaged Navier-Stokes equations for two-dimensional (2D) axisymmetric compressible flow can be written, as follows [23]: In which the conservative vector Q, flux vectors E, F and EV, FV.The axisymmetric source terms of H and HV and the source term S generated by the pyrolysis of fuel could be expressed in the following form:  The experimental procedure is described in detail in the paper [22].The repeatability of the connected pipe facility and the accuracy of the experimental result are shown and discussed in paper [22].

RANS Equations
The Reynolds-averaged Navier-Stokes equations for two-dimensional (2D) axisymmetric compressible flow can be written, as follows [23]: In which the conservative vector Q, flux vectors E, F and EV, FV.The axisymmetric source terms of H and HV and the source term S generated by the pyrolysis of fuel could be expressed in the following form:  The experimental procedure is described in detail in the paper [22].The repeatability of the connected pipe facility and the accuracy of the experimental result are shown and discussed in paper [22].

RANS Equations
The Reynolds-averaged Navier-Stokes equations for two-dimensional (2D) axisymmetric compressible flow can be written, as follows [23]: In which the conservative vector Q, flux vectors E, F and EV, FV.The axisymmetric source terms of H and HV and the source term S generated by the pyrolysis of fuel could be expressed in the following form: The experimental procedure is described in detail in the paper [22].The repeatability of the connected pipe facility and the accuracy of the experimental result are shown and discussed in paper [22].

RANS Equations
The Reynolds-averaged Navier-Stokes equations for two-dimensional (2D) axisymmetric compressible flow can be written, as follows [23]: In which the conservative vector Q, flux vectors E, F and E V , F V .The axisymmetric source terms of H and H V and the source term S generated by the pyrolysis of fuel could be expressed in the following form: In Equations ( 2) and (3), u, v, w, E, p, ρ, and T represent axial, radial and tangential velocity, energy, pressure, density, and temperature, respectively.λ is thermal conductivity, τ is viscous stress, ρ i , c i , w i , D i , and h i (i = 1, 2, . . ., 6) represent the density, mass fraction, mass source that is generated by the gas-phase chemical reaction, diffusion coefficient, and enthalpy of component i of unit mass.
ρD i is calculated by Equation (4).Where S ct (0.8) and S c (0.5) represent turbulent and laminar Schmidt number.
Equation ( 5) is used to calculate the laminar viscosity [24]: in which X i/j and µ li/j are mole fractions and laminar viscosity, respectively.M i/j represent the molecular weight.In addition, the subscripts i and j stand for components i and j.
The turbulent viscosity is calculated by Equation ( 6): where a 1 is assumed to be constant a 1 = 0.31 and ω and k represent the turbulent dissipation rate and turbulent kinetic energy.Ω represents the absolute value of vorticity.f 2 is calculated by Equation (7), in which y represents the distance to the wall surface.
Equation ( 3) is employed to calculate the thermal conductivity λ, where the turbulent and laminar Prandtl numbers are constant, P rt = 0.9 and P rl = 0.72.
Equation ( 9) is applied to define the pressure of gas mixture: The Newton iteration method can be used to evaluated the temperature T, as shown in Equation (10): In Equations ( 9) and ( 10), the universal gas constant is represented by R u .h 298 i is enthalpy at temperature 298 K. C pi used in Equation (10) and defined by Equation (11), the coefficients α ki (k = 1, 2, 3, 4, 5) are obtained from the chemical kinetics package [25].

Governing Equations of Solid Domain
The two-dimensional/axisymmetric Fourier's equation of solid domain is given by [26]: where ρ s is density of solid fuel, c s is heat capacity of solid domain, and h s , is the enthalpy of formation of unit mass.The energy change that is caused by the reaction in solid domain is represented by ρ s .r∂h s /∂y [27].

Numerical Method
The Reynolds-averaged Navier-Stokes equations have been solved by density-based, cell-centered, and finite-volume method with multi-block and structured grids.In this research, the simulation requires the high accuracy of predicting the heat transfer and flow conditions of the shear layer led by the sudden expansion combustor, thus, based on the flow conditions in this model, the k-ω SST turbulence model [28] is employed as Reference [29] described.The convective flux through cell face is computed by AUSMPW+ (Advection Upstream Splitting Method by Pressure Based weight Function) scheme [30].The third-order monotone upstream centered scheme for conservation laws (MUSCL) is employed for the convection term.Moreover, Van Albada limiting function is adopted to eliminate non-physical oscillations in the regions of large gradients.For temporal discretization, the lower upper symmetric Gauss-Seidel implicit method [31] is employed.The second-order central differences is used to evaluate the fluxes of the viscous terms [32].More information regarding this model was described in the paper [22].
In this simulation, at the first 0.5 s of the global time, the ignition gas with total temperature 2500 K, mass flow rate 0.3 kg/s, and the components of which are N 2 , H 2 O, and CO 2 flows into the combustor together with incoming air and lasts for 0.5 s.As the ignition time promoting, when the combustion became steady until the regression rate is almost constant, then the simulation is stopped.

Chemical Reaction Model
In our previous research, two models (eddy-dissipation and finite rate model) were used in the simulation to check the differences in simulating the combustion process of swirl intake SFRJ [33].The results indicate that there are no significant differences between the two models.The finite rate model is used due to its excellent behavior in the over fast chemistry simulation, accurate simulation of chemical kinetics, and the saving of computational resources in order get more insight into the combustion characteristics [20].While in the chamber of SFRJ, the combustion is very complex, and it is difficult to come up with an accurate model for high-density Polyethylene (HDPE) combustion.Thus, to simplify the simulation of the chemical reaction, the assumption is made that C 2 H 4 is the only pyrolysis product of polyethylene, and N 2 , CO, H 2 O, CO 2 , and O 2 are regarded as the main products in the chamber of SFRJ [34].As Reference [35] presents, Table 1 shoes the reaction model of the gas phase.In addition, the regression rate of HDPE is calculated by the Arrhenius empirical equation [36], as shown by Equation (13). .
In Equation ( 13), the value of A py and E apy could be found in Reference [37], and the temperature of solid fuel surface could be represented by T ω .

Reaction
A/(cm

Boundary Conditions
At the inlet boundary, the swirl velocity could be defined by using profile, in which the profiles of axial and swirl velocities are obtained from the experiment purposed by Dellenback [38].The axisymmetric boundary condition is used at the axis, a no-slip adiabatic boundary condition is applied for the wall boundaries.At the pressure outlet, all axial gradients of flow quantities are set to zero, and the energy balance equation (Equation ( 14)) is applied to calculate the wall temperature T w for the interface of the solid and fluid domains.

Case Description
A number of unsteady simulations are performed for different port diameters with a fixed inlet diameter (cases 1-5), and different inlet diameters with fixed port diameter (cases 3 and 6-9).As listed in Table 2, the air inlet temperatures (T in ) are 540 K and the inlet mass flow rates ( .m air ) are 0.3 kg/s for each case.The port diameters (d p ) range from 70 mm to 90 mm, at which the inlet diameters (d in ) are fixed at 40 mm.The inlet diameters range from 30 mm to 50 mm, at which the port diameter is fixed at 80 mm. and d p /d in define the port-to-inlet diameters.Meanwhile, the diameters of throat (d t ) for each case are 28.5 mm, to avoid the effect of throat diameter on combustion.In addition, the swirl intensity (s) for each case is fixed at 0.6 and the length of HDPE (L) is 500 mm.

Computational Models
Three different cell numbers (185,633, 103,321, and 67,482) are carried out for case 3 in order to check the independency of the simulation results, and the y + is shown in Figure 4. Figure 5 shows the regression rate.The axial and tangential velocities of three cases are described in Figures 6 and 7, respectively.As shown in Figures 6 and 7, the value of velocity was performed on x-axis, and the y-axis was normalized as r/R, in which r represents the radial position, and R means the radius of the SFRJ chamber.As could be seen from Figures 6 and 7, the tangential and axial velocities show little difference between the cases with maximum and medium number of cells; however, for the case with minimum cell number, the result has an obvious gap when compared with other two cases.In addition, as for the regression rate profiles (Figure 5), the case with the minimum cell number, the regression rate has significant difference with the other two cases and the experimental results (see Table 3).Therefore, when considering the high simulating accuracy and saving the resources of calculation, the mesh number around 10 5 is used in this present study.
is fixed at 80 mm. and dp/din define the port-to-inlet diameters.Meanwhile, the diameters of throat (dt) for each case are 28.5 mm, to avoid the effect of throat diameter on combustion.In addition, the swirl intensity (s) for each case is fixed at 0.6 and the length of HDPE (L) is 500 mm.

Computational Models
Three different cell numbers (185633, 103321, and 67482) are carried out for case 3 in order to check the independency of the simulation results, and the y + is shown in Figure 4. Figure 5 shows the regression rate.The axial and tangential velocities of three cases are described in Figures 6 and 7, respectively.As shown in Figures 6 and 7, the value of velocity was performed on x-axis, and the y-axis was normalized as r/R, in which r represents the radial position, and R means the radius of the SFRJ chamber.As could be seen from Figures 6 and 7, the tangential and axial velocities show little difference between the cases with maximum and medium number of cells; however, for the case with minimum cell number, the result has an obvious gap when compared with other two cases.In addition, as for the regression rate profiles (Figure 5), the case with the minimum cell number, the regression rate has significant difference with the other two cases and the experimental results (see Table 3).Therefore, when considering the high simulating accuracy and saving the resources of calculation, the mesh number around 10 5 is used in this present study.In this research, the mesh number of all the cases is 124602 the minimum and 148873 the maximum.For the simulations of SFRJ, to ensure the accuracy of boundary temperature and heat-flux, the mesh is clustered near the solid fuel surface and the wall boundary (as shown in Figure 8) to resolve the turbulent kinetic energy, temperature, and other flow variables with large gradients.In order to meet the requirement of SST turbulence model, the y + has to be kept at the value of 1-5.In this research, the mesh number of all the cases is 124,602 the minimum and 148,873 the maximum.For the simulations of SFRJ, to ensure the accuracy of boundary temperature and heat-flux, the mesh is clustered near the solid fuel surface and the wall boundary (as shown in Figure 8) to resolve the turbulent kinetic energy, temperature, and other flow variables with large gradients.In order to meet the requirement of SST turbulence model, the y + has to be kept at the value of 1-5.In this research, the mesh number of all the cases is 124602 the minimum and 148873 the maximum.For the simulations of SFRJ, to ensure the accuracy of boundary temperature and heat-flux, the mesh is clustered near the solid fuel surface and the wall boundary (as shown in Figure 8) to resolve the turbulent kinetic energy, temperature, and other flow variables with large gradients.In order to meet the requirement of SST turbulence model, the y + has to be kept at the value of 1-5.

Model validation
In order to confirm that the in-house code that was developed in this research has sufficient accuracy in predicting the combustion and flow in SFRJ, the swirl model, chemical reaction model, and heat diffusion of the solid domain were validated by using benchmark cases.The swirl model was validated by using Nejad's experimental result [39]; the swirl intensity for the benchmark case is 0.5 and the swirl intensity is defined as: where R is the inlet radius, w and u represent the tangential and axial velocities, while r is the radius positions of the flow area.Figure 9 shows the predicted results of axial and tangential velocities, where y-axis was normalized as y/H, x-axis was normalized as W/U0 and U/U0, in which y represents the radial position and H means the diameter of the tube that used for calculation; moreover, W, U, and U0 represent the tangential, axial, and maximum axial velocity, respectively.It could be seen from Figure 9 that the agreement between the simulation and experimental result could meet the requirements of the study that are presented in this paper.

Model Validation
In order to confirm that the in-house code that was developed in this research has sufficient accuracy in predicting the combustion and flow in SFRJ, the swirl model, chemical reaction model, and heat diffusion of the solid domain were validated by using benchmark cases.The swirl model was validated by using Nejad's experimental result [39]; the swirl intensity for the benchmark case is 0.5 and the swirl intensity is defined as: where R is the inlet radius, w and u represent the tangential and axial velocities, while r is the radius positions of the flow area.Figure 9 shows the predicted results of axial and tangential velocities, where y-axis was normalized as y/H, x-axis was normalized as W/U 0 and U/U 0 , in which y represents the radial position and H means the diameter of the tube that used for calculation; moreover, W, U, and U 0 represent the tangential, axial, and maximum axial velocity, respectively.It could be seen from Figure 9 that the agreement between the simulation and experimental result could meet the requirements of the study that are presented in this paper.A famous benchmark case [40] is employed to validate the chemical reaction model.A diameter of 15 mm spherical projectile under the operating condition of the pressure of 42662 Pa, temperature of 250 K, and the supersonic speed (1685 m/s) of the stoichiometric mixture (H2/Air) flow through the computational domain.An eight element reactions and seven species chemistry model is used for the simulation of the benchmark case.Moreover, Figures 10 and 11 show the results.As in these figures described, the simulation the results that were obtained from the in-house code reproduce the experimental data very well.It becomes evident that the accuracy of the combustion model could meet the requirements of this research.A famous benchmark case [40] is employed to validate the chemical reaction model.A diameter of 15 mm spherical projectile under the operating condition of the pressure of 42662 Pa, temperature of 250 K, and the supersonic speed (1685 m/s) of the stoichiometric mixture (H2/Air) flow through the computational domain.An eight element reactions and seven species chemistry model is used for the simulation of the benchmark case.Moreover, Figures 10 and 11 show the results.As in these figures described, the simulation the results that were obtained from the in-house code reproduce the experimental data very well.It becomes evident that the accuracy of the combustion model could meet the requirements of this research.A semi-infinite plate case is applied to validate the heat diffusion of the solid domain.A rectangular plate with the size of 2 mm * 15 mm is utilized in this validation.In this simulation, the heat transfer was calculated by using Equation (10), and the plate has a density of 7840 kg/m 3 and specific heat of 465 J/(kg•K).The initial temperature is set to 300K, the thermal conductivity is 49.8 W/(m•K).In addition, the boundary condition at x = 0 has two different conditions: temperature Tw = 900 K for case a and heat flux qw = 70 kw for case b.To be compared with numerical solutions, the analytical solution that is given by Equations ( 16) and ( 17) is conducted for the two cases.
where erf () and α = λ/qc represent the error function and thermal diffusivity of the material, respectively.
The results of numerical and analytical solutions are presented in Figure 12.It could be seen that excellent agreements are achieved.A semi-infinite plate case is applied to validate the heat diffusion of the solid domain.A rectangular plate with the size of 2 mm * 15 mm is utilized in this validation.In this simulation, the heat transfer was calculated by using Equation (10), and the plate has a density of 7840 kg/m 3 and specific heat of 465 J/(kg•K).The initial temperature is set to 300K, the thermal conductivity is 49.8 W/(m•K).In addition, the boundary condition at x = 0 has two different conditions: temperature T w = 900 K for case a and heat flux q w = 70 kw for case b.To be compared with numerical solutions, the analytical solution that is given by Equations ( 16) and ( 17) is conducted for the two cases.
where erf () and α = λ/qc represent the error function and thermal diffusivity of the material, respectively.The results of numerical and analytical solutions are presented in Figure 12.It could be seen that excellent agreements are achieved.A semi-infinite plate case is applied to validate the heat diffusion of the solid domain.A rectangular plate with the size of 2 mm * 15 mm is utilized in this validation.In this simulation, the heat transfer was calculated by using Equation (10), and the plate has a density of 7840 kg/m 3 and specific heat of 465 J/(kg•K).The initial temperature is set to 300K, the thermal conductivity is 49.8 W/(m•K).In addition, the boundary condition at x = 0 has two different conditions: temperature Tw = 900 K for case a and heat flux qw = 70 kw for case b.To be compared with numerical solutions, the analytical solution that is given by Equations ( 16) and ( 17) is conducted for the two cases.
where erf () and α = λ/qc represent the error function and thermal diffusivity of the material, respectively.
The results of numerical and analytical solutions are presented in Figure 12.It could be seen that excellent agreements are achieved.In order to validate the numerical accuracy of the in-house code, the experiment on the case of d p 70 mm d in 40 mm was carried out and compared with the simulation results, and the discussions are as follows.Figure 13 described the internal surface of HDPE after combustion, and Figure 14 shows the regression rate contour of internal surface measurement by 3D scanner.According to Figures 13  and 14, the reattachment point is located at the place with the maximum regression rate.In addition, it could be seen from Figure 13 that the carbon particle that was attached to the HDPE surface showed the obvious characteristic of swirl flow.shows the regression rate contour of internal surface measurement by 3D scanner.According to Figure 13 and Figure 14, the reattachment point is located at the place with the maximum regression rate.In addition, it could be seen from Figure 13 that the carbon particle that was attached to the HDPE surface showed the obvious characteristic of swirl flow.Figure 15 shows the comparison of numerical and experimental regression rate.According to Figure 15, the simulation result and experiment result have the same trend, and the regression rate reaches the maximum value at the reattachment point.However, when compared with the result of experiment, the reattachment point of the simulation moves towards the downstream of the flow field; this may due to the higher velocities loss and more complicated flow conditions in experiment when compared with the simulation conditions.In the paper, this code could safely predict the regression rate for a small combustion time about 2s, due to neglecting the effect of moving boundary.All of the results that are presented in this paper are based on working time within 1.5s only.Thus, the moving mesh could be safely neglected.Figure 13 and Figure 14, the reattachment point is located at the place with the maximum regression rate.In addition, it could be seen from Figure 13 that the carbon particle that was attached to the HDPE surface showed the obvious characteristic of swirl flow.Figure 15 shows the comparison of numerical and experimental regression rate.According to Figure 15, the simulation result and experiment result have the same trend, and the regression rate reaches the maximum value at the reattachment point.However, when compared with the result of experiment, the reattachment point of the simulation moves towards the downstream of the flow field; this may due to the higher velocities loss and more complicated flow conditions in experiment when compared with the simulation conditions.In the paper, this code could safely predict the regression rate for a small combustion time about 2s, due to neglecting the effect of moving boundary.All of the results that are presented in this paper are based on working time within 1.5s only.Thus, the moving mesh could be safely neglected.Figure 15 shows the comparison of numerical and experimental regression rate.According to Figure 15, the simulation result and experiment result have the same trend, and the regression rate reaches the maximum value at the reattachment point.However, when compared with the result of experiment, the reattachment point of the simulation moves towards the downstream of the flow field; this may due to the higher velocities loss and more complicated flow conditions in experiment when compared with the simulation conditions.In the paper, this code could safely predict the regression rate for a small combustion time about 2 s, due to neglecting the effect of moving boundary.All of the results that are presented in this paper are based on working time within 1.5 s only.Thus, the moving mesh could be safely neglected.Table 4 shows the average regression rate and characteristic velocities that were obtained from numerical and experimental studies.According to the deviation calculated in Table 4, it can be concluded that the accuracy of this code for engineering investigations is acceptable, and the in-house code that is used in this research has enough precision in predicting the combustion and flow characteristics in SFRJ.Table 4 shows the average regression rate and characteristic velocities that were obtained from numerical and experimental studies.According to the deviation calculated in Table 4, it can be concluded that the accuracy of this code for engineering investigations is acceptable, and the in-house code that is used in this research has enough precision in predicting the combustion and flow characteristics in SFRJ.

Flow Field Characteristics
In the swirl flow condition, the combustion characteristics of SFRJ is mainly determined by the flow field characteristics.In this section, case1 is used for analyzing the flow field characteristic; the combustion and regression rate become steady at the time of 1.352 s, thus, all of the results that are presented in following figures are the state in 1.352 s.According to Figure 16, the temperature reaches 2800 K around the recirculation zone and the reattachment point.In the redevelopment zone, the reaction of gaseous fuel and incoming air occurs, where the temperature reaches 2400 K.It could be seen from Figures 16 and 17 that the interface of the region with high CO 2 ratio and the region with high O 2 ratio has the maximum temperature.As described in Figure 17b, the region with higher mass fraction of CO coincides with the high-temperature region.Thus, it is obvious that the chemical reaction mainly took place in the region with more CO 2 and CO.The interface between these two regions is named as the flame surface, which was also experimentally identified in Reference [34].There are two combustion mechanisms (diffusion process and chemical process) in the combustion chamber.It is believed that the diffusion-controlled reaction takes place in the shear layer and the redevelopment zone with higher velocity (see Figure 18a-c).Meanwhile, the particles of propellant gas have longer residence time in recirculation zone due to their relatively low velocity.Thus, the chemical process controls the reaction of the recirculation zone.The axial and tangential velocities at different radial position are introduced in Figure 18a,b.The maximum value of velocity was found at the place of the shear layer in the recirculation zone.At the reattachment point, the place of maximum velocity approaches toward the solid surface.The place of maximum velocity approaches toward to the core flow area with the increasing of X/L.The axial and tangential velocities at different radial position are introduced in Figure 18a,b.The maximum value of velocity was found at the place of the shear layer in the recirculation zone.At the reattachment point, the place of maximum velocity approaches toward the solid surface.The place of maximum velocity approaches toward to the core flow area with the increasing of X/L.Moreover, at the core flow field in the chamber, the axial velocity increases with an increasing of X/R, and the tangential velocity decreases with an increasing of X/R, owing to the attenuation of the swirl intensity along the axial direction.

Theoretical Analysis of SFRJ Performance
Figure 19 shows the radial profiles of temperature.In recirculation zone, the flame surface is at the shear layer, thus, the temperature distribution in the recirculation zone is more uniform than that of reattachment point and redevelopment zone.At reattachment point, the flame surface approaches

Theoretical Analysis of SFRJ Performance
Figure 19 shows the radial profiles of temperature.In recirculation zone, the flame surface is at the shear layer, thus, the temperature distribution in the recirculation zone is more uniform than that of reattachment point and redevelopment zone.At reattachment point, the flame surface approaches to the fuel surface, and the distance of maximum value of temperature to the solid surface is about 3 mm.At the end of the combustor, the flame surface is gradually approaching to the core flow area with the increasing of X/L. Figure 20 describes the regression rate of HDPE at different times.As shown in Figure 20, the first 0.5 s of the global time in the simulation, the regression rate rapidly grows, due to the ignition gas with a total temperature of 2500 K and mass flow rate 0.3 kg/s.After 0.5 s, the sustained combustion is obtained and the regression rate of reattachment point reduces to 0.4 mm/s.When the combustion time lasts 1 s, the regression rate of SFRJ is almost constant, and it is believed that combustion of SFRJ reaches the steady state.
Figure 22b, the value of turbulent viscosity (μt) is decreased with the distance approaching the fuel surface.Therefore, as could be seen from the profiles, the μt reaches to the maximum value around reattachment point.While in redevelopment zone, the μt decreases with increasing of X/L.
Through the analysis above.In the reattachment point, the flame surface shifts towards the fuel surface and μt reaches the maximum value, which leads to the maximum value of regression rate and heat flux.However, the temperature and μt near the fuel surface gradually decrease downstream of the reattachment point, which result in the decreasing of regression rate and heat flux.The analysis above shows that the regression rate could be sensitive to the variation of turbulent viscosity and temperature distribution.It becomes evident that, in Figure 21, the regression rate and heat flux introduce a similar trend, both are rapidly increased in the recirculation zone and decrease gradually downstream of flow field.This phenomenon indicates that the heat flux of solid fuel surface and regression rate are proportional to each other.
There is a strong relationship between the heat flux and the thermal conductivity of the solid domain, and the thermal conductivity is proportional to the turbulent viscosity [36].As Figure 22a introduced, the temperature gradient reaches the maximum value at the reattachment point.In addition, in the redevelopment zone, the temperature gradient gradually decreases.According to Figure 22b, the value of turbulent viscosity (µ t ) is decreased with the distance approaching the fuel surface.Therefore, as could be seen from the profiles, the µ t reaches to the maximum value around reattachment point.While in redevelopment zone, the µ t decreases with increasing of X/L.   Figure 23 shows the heat transfer coefficient along the HDPE surface, which could be calculated by the equation below. ( ) Through the analysis above.In the reattachment point, the flame surface shifts towards the fuel surface and µ t reaches the maximum value, which leads to the maximum value of regression rate and heat flux.However, the temperature and µ t near the fuel surface gradually decrease downstream of the reattachment point, which result in the decreasing of regression rate and heat flux.The analysis above shows that the regression rate could be sensitive to the variation of turbulent viscosity and temperature distribution.
Figure 23 shows the heat transfer coefficient along the HDPE surface, which could be calculated by the equation below.
Energies 2019, 12, x FOR PEER REVIEW 19 of 30 temperature of the inlet and Prandtl numbers, respectively.γfluid is specific heat ratio of propellant and  represents the Mach number at inlet. Figure 23, when the value of heat transfer coefficient is negative, the solid domain is in exothermic status, conversely, the solid domain is in endothermic status.During the ignition period (0-0.5 s), the heat transfer coefficient is positive, then a short exothermic status on the surface of the fuel occurs at the time of 0.502 s, and the heat transfer coefficient reaches the maximum value of exothermic stage.At the time of 0.6 s, around the region of reattachment point, the heat transfer coefficient is negative.However, for other regions of the fuel surface, the heat transfer coefficient is positive.This indicates that the combustion is achieved around the reattachment point, and the heat released in reattachment point will be applied to achieve combustion in the redevelopment zone.As time progresses, the heat transfer coefficient in the redevelopment zone has gradually changed from positive to negative, which implies that the self-sustained combustion is successfully achieved.

Effects of Sudden Expansion Ratio on Regression Rate
Figure 24 introduces the local regression rate of different port-to-inlet diameters.In this figure, the inlet diameter ranges from 30 mm to 50 mm and the port diameter fixed at 80 mm.Further, the dp/din (sudden expansion ratio) ranges from 1.6-2.67.As the analysis above, the regression rate rapidly grows in the recirculation zone and it gradually decreases downstream of the flow field.However, for different inlet diameters, the regression rates in the reattachment point and redevelopment zone for each case are increased with increasing of dp/din.Figure 25 indicates the heat flux for different cases, the heat flux introduced the same trend as the regression rate, in the area of reattachment point and redevelopment zone, and the heat flux increases with increasing of dp/din.
The distributions of temperature and temperature around the fuel surface are shown in Figure 26a,b, respectively.As could be seen from Figure 26a, the flame surface is getting closer to the fuel surface with the decrease of dp/din in the recirculation zone.The value of maximum temperature reaches 2800 K.In the reattachment point and redevelopment zone, the distance between the flame surface and fuel surface does not change with dp/din.In the core flow area, the temperature of the downstream flow field increased with an increasing of dp/din.Figure 26b displays the distributions of temperature near the fuel surface.According to this figure, in recirculation zone, the temperature was irregular for different cases.While in reattachment point, the temperature is proportional to In Equation ( 18), λ fluid is the thermal conductivity of the fluid domain, ∂T/∂n represents the temperature gradient, and T w and T aw stand for the temperature of the fuel surface and recovery temperature, respectively.T aw can be calculated by Equation ( 19) and T ∞ and P rl represent the temperature of the inlet and Prandtl numbers, respectively.γ fluid is specific heat ratio of propellant and Ma ∞ represents the Mach number at inlet.
Figure 23, when the value of heat transfer coefficient is negative, the solid domain is in exothermic status, conversely, the solid domain is in endothermic status.During the ignition period (0-0.5 s), the heat transfer coefficient is positive, then a short exothermic status on the surface of the fuel occurs at the time of 0.502 s, and the heat transfer coefficient reaches the maximum value of exothermic stage.At the time of 0.6 s, around the region of reattachment point, the heat transfer coefficient is negative.However, for other regions of the fuel surface, the heat transfer coefficient is positive.This indicates that the combustion is achieved around the reattachment point, and the heat released in reattachment point will be applied to achieve combustion in the redevelopment zone.As time progresses, the heat transfer coefficient in the redevelopment zone has gradually changed from positive to negative, which implies that the self-sustained combustion is successfully achieved.

Effects of Sudden Expansion Ratio on Regression Rate
Figure 24 introduces the local regression rate of different port-to-inlet diameters.In this figure, the inlet diameter ranges from 30 mm to 50 mm and the port diameter fixed at 80 mm.Further, the d p /d in (sudden expansion ratio) ranges from 1.6-2.67.As the analysis above, the regression rate rapidly grows in the recirculation zone and it gradually decreases downstream of the flow field.However, for different inlet diameters, the regression rates in the reattachment point and redevelopment zone for each case are increased with increasing of d p /d in .Figure 25  dp/din.Moreover, the same trend was found for the temperature in the redevelopment zone, but the differences among the five cases are quite small.dp/din.Moreover, the same trend was found for the temperature in the redevelopment zone, but the differences among the five cases are quite small.The distributions of temperature and temperature around the fuel surface are shown in Figure 26a,b, respectively.As could be seen from Figure 26a, the flame surface is getting closer to the fuel surface with the decrease of d p /d in in the recirculation zone.The value of maximum temperature reaches 2800 K.In the reattachment point and redevelopment zone, the distance between the flame surface and fuel surface does not change with d p /d in .In the core flow area, the temperature of the downstream flow field increased with an increasing of d p /d in .Figure 26b displays the distributions of temperature near the fuel surface.According to this figure, in recirculation zone, the temperature was irregular for different cases.While in reattachment point, the temperature is proportional to d p /d in .Moreover, the same trend was found for the temperature in the redevelopment zone, but the differences among the five cases are quite small.Figure 27 shows the distributions of turbulent viscosity; it can be seen that, in recirculation zone, the turbulent viscosity was irregular for different cases.However, in the reattachment point, with the increase of dp/din, the turbulent viscosity significantly increased.In the redevelopment zone, the turbulent viscosity is lower than that of the recirculation zone and reattachment point overall and it increases with the decrease of dp/din in the redevelopment zone.
In summary, for the cases with different inlet diameters, the turbulent viscosity could be significantly affected by the change of sudden expansion ratio.The variation of turbulent viscosity will lead to the changes of surface thermal conductivity and temperature, and will ultimately affect the regression rate.Figure 27 shows the distributions of turbulent viscosity; it can be seen that, in recirculation zone, the turbulent viscosity was irregular for different cases.However, in the reattachment point, with the increase of d p /d in , the turbulent viscosity significantly increased.In the redevelopment zone, the turbulent viscosity is lower than that of the recirculation zone and reattachment point overall and it increases with the decrease of d p /d in in the redevelopment zone.For the regression rate affected by port diameter with fixed din, as described in Figure 28, the port diameter ranges from 70 mm to 90 mm and the inlet diameter is fixed at 40 mm.The dp/din ranges from 1.75-2.25 in the recirculation zone, the regression rate increased rapidly, while unlike the regression rate with different din, at the place of reattachment point, the regression rate of different cases is almost the same.Moreover, in the redevelopment zone, the local regression rate is inversely proportional to dp/din.As paper [17] described, the increasing of the regression rate may be due to the enhancement of mass flux that is caused by the decreasing of port diameter.According to the analyses above, similar to the cases with different din, the heat flux shows the same trend with the regression rate.(See Figure 29)  In summary, for the cases with different inlet diameters, the turbulent viscosity could be significantly affected by the change of sudden expansion ratio.The variation of turbulent viscosity will lead to the changes of surface thermal conductivity and temperature, and will ultimately affect the regression rate.
For the regression rate affected by port diameter with fixed d in , as described in Figure 28, the port diameter ranges from 70 mm to 90 mm and the inlet diameter is fixed at 40 mm.The d p /d in ranges from 1.75-2.25 in the recirculation zone, the regression rate increased rapidly, while unlike the regression rate with different d in , at the place of reattachment point, the regression rate of different cases is almost the same.Moreover, in the redevelopment zone, the local regression rate is inversely proportional to d p /d in .As paper [17] described, the increasing of the regression rate may be due to the enhancement of mass flux that is caused by the decreasing of port diameter.According to the analyses above, similar to the cases with different d in , the heat flux shows the same trend with the regression rate.(See Figure 29).For the regression rate affected by port diameter with fixed din, as described in Figure 28, the port diameter ranges from 70 mm to 90 mm and the inlet diameter is fixed at 40 mm.The dp/din ranges from 1.75-2.25 in the recirculation zone, the regression rate increased rapidly, while unlike the regression rate with different din, at the place of reattachment point, the regression rate of different cases is almost the same.Moreover, in the redevelopment zone, the local regression rate is inversely proportional to dp/din.As paper [17] described, the increasing of the regression rate may be due to the enhancement of mass flux that is caused by the decreasing of port diameter.According to the analyses above, similar to the cases with different din, the heat flux shows the same trend with the regression rate.(See Figure 29)   As Figure 30a and b introduced, the distance between the fuel surface and flame surface does not change with dp/din, due to the centrifugal force that is generated by swirl flow.In the recirculation zone and reattachment point, the temperature reaches the maximum value of 2800 K, in the redevelopment zone, the temperature reduces to 2400 K.As analyzed above, the flame surface moves far from the fuel surface.However, in the core flow area, the temperature does not change with dp/din.As Figure 30b described, in the recirculation zone and in reattachment point, the temperatures around the surface are irregular with dp/din, while, for the temperature of redevelopment zone, it is inversely proportional to dp/din.Moreover, as Figure 30b introduced, in the reattachment point, the temperature gradient is higher than that of other place, and it became steady at the end of the combustor.
(a) As Figure 30a,b introduced, the distance between the fuel surface and flame surface does not change with d p /d in , due to the centrifugal force that is generated by swirl flow.In the recirculation zone and reattachment point, the temperature reaches the maximum value of 2800 K, in the redevelopment zone, the temperature reduces to 2400 K.As analyzed above, the flame surface moves far from the fuel surface.However, in the core flow area, the temperature does not change with d p /d in .As Figure 30b described, in the recirculation zone and in reattachment point, the temperatures around the surface are irregular with d p /d in , while, for the temperature of redevelopment zone, it is inversely proportional to d p /d in .Moreover, as Figure 30b introduced, in the reattachment point, the temperature gradient is higher than that of other place, and it became steady at the end of the combustor.As Figure 30a and b introduced, the distance between the fuel surface and flame surface does not change with dp/din, due to the centrifugal force that is generated by swirl flow.In the recirculation zone and reattachment point, the temperature reaches the maximum value of 2800 K, in the redevelopment zone, the temperature reduces to 2400 K.As analyzed above, the flame surface moves far from the fuel surface.However, in the core flow area, the temperature does not change with dp/din.As Figure 30b described, in the recirculation zone and in reattachment point, the temperatures around the surface are irregular with dp/din, while, for the temperature of redevelopment zone, it is inversely proportional to dp/din.Moreover, as Figure 30b introduced, in the reattachment point, the temperature gradient is higher than that of other place, and it became steady at the end of the combustor.In summary, for the cases with different port diameters, the heat flux of the fuel surface and the regression rate perform irregular changes in the area of reattachment point and recirculation zone with the variation of dp/din.However, the regression rate and heat flux reproduce inversely proportional to the sudden expansion ratio in the redevelopment region.
The regression rate that was obtained from simulation for different cases are shown in Figure 32 and Figure 33.The air-fuel equivalence ratio is shown in Figures 34 and 35, which is calculated by Equation (20).Where the AFR and AFRstoich represent the Air-fuel ratio and Air-fuel ratio at stoichiometry, respectively, in which AFR was calculated by Equation ( 21    Figure 31 shows the distributions of turbulent viscosity, in recirculation zone, the μt of different cases are almost the same, and at the position of reattachment point, the μt are irregular with dp/din.However, the μt is inversely proportional to dp/din in the redevelopment zone.In summary, for the cases with different port diameters, the heat flux of the fuel surface and the regression rate perform irregular changes in the area of reattachment point and recirculation zone with the variation of dp/din.However, the regression rate and heat flux reproduce inversely proportional to the sudden expansion ratio in the redevelopment region.
The regression rate that was obtained from simulation for different cases are shown in Figure 32 and Figure 33.The air-fuel equivalence ratio is shown in Figures 34 and 35, which is calculated by Equation (20).Where the AFR and AFRstoich represent the Air-fuel ratio and Air-fuel ratio at stoichiometry, respectively, in which AFR was calculated by Equation ( 21  In summary, for the cases with different port diameters, the heat flux of the fuel surface and the regression rate perform irregular changes in the area of reattachment point and recirculation zone with the variation of d p /d in .However, the regression rate and heat flux reproduce inversely proportional to the sudden expansion ratio in the redevelopment region. The regression rate that was obtained from simulation for different cases are shown in Figures 32  and 33.The air-fuel equivalence ratio is shown in Figures 34 and 35, which is calculated by Equation (20).Where the AFR and AFR stoich represent the Air-fuel ratio and Air-fuel ratio at stoichiometry, respectively, in which AFR was calculated by Equation ( 21) [41].
AFR stoich is calculated by a thermochemical code, Chemical Equilibrium with Applications (CEA), as proposed by NASA, based on the result from Chemical Equilibrium with Applications, the AFR stoich = 14.8.
Figure 32 introduces the maximum and average regression rates for different d p /d in , in which d p is fixed as constant.Moreover, an approximately linear relationship is established based on the distributions of average regression rate and maximum regression rate.According to the profiles, a similar tendency was found between maximum regression rate and the average regression rate, but the slope of the fitting curve of maximum regression rate is higher than that of average regression rate.
From the analyses above, we can conclude that, the regression rate is proportional to the heat flux of fuel surface.In the reattachment point, around the fuel surface, the heat flux reaches the maximum value, and it increased with decreasing of d in , which indicates that, for the cases with a fixed d p and inlet mass flow rate, the overall heat transfer behavior of the fuel surface is mostly dominated by the heat transfer around the reattachment point.
AFRstoich is calculated by a thermochemical code, Chemical Equilibrium with Applications (CEA), as proposed by NASA, based on the result from Chemical Equilibrium with Applications, the AFRstoich = 14.8.
Figure 32 introduces the maximum and average regression rates for different dp/din, in which dp is fixed as constant.Moreover, an approximately linear relationship is established based on the distributions of average regression rate and maximum regression rate.According to the profiles, a similar tendency was found between maximum regression rate and the average regression rate, but the slope of the fitting curve of maximum regression rate is higher than that of average regression rate.
From the analyses above, we can conclude that, the regression rate is proportional to the heat flux of fuel surface.In the reattachment point, around the fuel surface, the heat flux reaches the maximum value, and it increased with decreasing of din, which indicates that, for the cases with a fixed dp and inlet mass flow rate, the overall heat transfer behavior of the fuel surface is mostly dominated by the heat transfer around the reattachment point.As Figure 33 describes, an approximately linear relationship is established based on the distributions of average regression rate and maximum regression rate.Unlike the cases with fixed dp, the fitting curve of maximum regression rate is almost smooth and the average regression rate is a monotony decrease function of port-to-inlet ratio.
Based on the above analysis, an irregular behavior of the regression rate and heat flux of solid fuel surface is performed in the recirculation zone and reattachment point.When compared with the obvious regularity of downstream regression rate and heat flux distribution and overall average regression rate shown in Figure 33, it could be concluded that, in the redevelopment zone, the heat transfer enhanced by the backward-facing step could be negligible.Consequently, it is believed that the average regression rate is influenced by the heat transfer mechanism of developed turbulent flow in the redevelopment zone.As Figure 33 describes, an approximately linear relationship is established based on the distributions of average regression rate and maximum regression rate.Unlike the cases with fixed d p , the fitting curve of maximum regression rate is almost smooth and the average regression rate is a monotony decrease function of port-to-inlet ratio.
Based on the above analysis, an irregular behavior of the regression rate and heat flux of solid fuel surface is performed in the recirculation zone and reattachment point.When compared with the obvious regularity of downstream regression rate and heat flux distribution and overall average regression rate shown in Figure 33, it could be concluded that, in the redevelopment zone, the heat transfer enhanced by the backward-facing step could be negligible.Consequently, it is believed that the average regression rate is influenced by the heat transfer mechanism of developed turbulent flow in the redevelopment zone.As shown in Figure 34, the κ is increased with an increasing of dp/din, and in the case of 1.6 and 1.8 sudden expansion ratio, the value of κ is less than one, which indicated that the chemical reactions took place in a fuel-rich state.However, the cases where the value of κ is exceeded one indicated that the reaction occurs in an oxide-rich state.
Figure 35 shows the air-fuel equivalence ratio for different cases with din fixed as constant.It could be seen in Figure 35, that, when compared with the cases of fixed din, the trend is flat and it does not change with dp/din.This may lead to the irregular behavior that took place in the reattachment point and the recirculation zone.Moreover, it is believed that the chemical reactions are all taking place in an oxide-rich state, as judged by the value of κ.As shown in Figure 34, the κ is increased with an increasing of d p /d in , and in the case of 1.6 and 1.8 sudden expansion ratio, the value of κ is less than one, which indicated that the chemical reactions took place in a fuel-rich state.However, the cases where the value of κ is exceeded one indicated that the reaction occurs in an oxide-rich state.
Figure 35 shows the air-fuel equivalence ratio for different cases with d in fixed as constant.It could be seen in Figure 35, that, when compared with the cases of fixed d in , the trend is flat and it does not change with d p /d in .This may lead to the irregular behavior that took place in the reattachment point and the recirculation zone.Moreover, it is believed that the chemical reactions are all taking place in an oxide-rich state, as judged by the value of κ.As shown in Figure 34, the κ is increased with an increasing of dp/din, and in the case of 1.6 and 1.8 sudden expansion ratio, the value of κ is less than one, which indicated that the chemical reactions took place in a fuel-rich state.However, the cases where the value of κ is exceeded one indicated that the reaction occurs in an oxide-rich state.
Figure 35 shows the air-fuel equivalence ratio for different cases with din fixed as constant.It could be seen in Figure 35, that, when compared with the cases of fixed din, the trend is flat and it does not change with dp/din.This may lead to the irregular behavior that took place in the reattachment point and the recirculation zone.Moreover, it is believed that the chemical reactions are all taking place in an oxide-rich state, as judged by the value of κ.

Conclusions
The present study investigated the combustion characteristics effect by the sudden expansion ratio of solid-fuel ramjet combustor by the simulation approach.The simulations with unsteady, reacting, and turbulent flow with swirl in SFRJ were conducted via a CFD code that was written in FORTRAN.First, the experimental measurements of SFRJ is used to validate the accuracy of solid fuel decomposition, the chemical reaction kinetics has been validated by a benchmark case, the sudden expansion combustor case is applied to validate the swirl flow simulation, and heat diffusion in the solid domain is validated by a semi-infinite plate case.
Secondly, the effect of sudden expansion ratio on the flow field behavior and combustion phenomena of SFRJ has been investigated by a series of simulations.The relationship between the sudden expansion ratio and average regression rate have been proposed in this investigation.The main conclusions are summarized below:

•
The physical reasons for the regression rate that are affected by the sudden expansion ratio were obtained.The variation of turbulent viscosity due to the sudden expansion ratio changes could significantly affect the solid fuel regression rate and the fuel surface heat transfer.

•
It becomes evident that the combustion process is closely related to the heat transfer process of solid fuel surface.Based on the analysis of the heat transfer coefficient, the self-sustained combustion occurs around the reattachment point at first, and then gradually spread to the redevelopment zone.The heat released in the reattachment point will be used to achieve the self-sustained combustion in the redevelopment zone.

•
The linear relationship between the sudden expansion ratio with fixed port diameter and average/maximum regression rate was obtained.The result indicates that the average regression rate and maximum regression rate are more sensitive to the change of the port-to-inlet diameter.Additionally, the overall heat transfer behavior of the fuel surface was mostly dominated by the maximum heat transfer around the backward-facing step.

•
Based on the analysis of linear relationship between regression rate and sudden expansion ratio with fixed inlet diameter.The average regression rate is mainly affected by heat transfer mechanism in a fully developed turbulent flow in rthe edevelopment zone, and it decreased with the increasing of the sudden expansion ratio.

Conclusions
The present study investigated the combustion characteristics effect by the sudden expansion ratio of solid-fuel ramjet combustor by the simulation approach.The simulations with unsteady, reacting, and turbulent flow with swirl in SFRJ were conducted via a CFD code that was written in FORTRAN.First, the experimental measurements of SFRJ is used to validate the accuracy of solid fuel decomposition, the chemical reaction kinetics has been validated by a benchmark case, the sudden expansion combustor case is applied to validate the swirl flow simulation, and heat diffusion in the solid domain is validated by a semi-infinite plate case.
Secondly, the effect of sudden expansion ratio on the flow field behavior and combustion phenomena of SFRJ has been investigated by a series of simulations.The relationship between the sudden expansion ratio and average regression rate have been proposed in this investigation.The main conclusions are summarized below:

•
The physical reasons for the regression rate that are affected by the sudden expansion ratio were obtained.The variation of turbulent viscosity due to the sudden expansion ratio changes could significantly affect the solid fuel regression rate and the fuel surface heat transfer.

•
It becomes evident that the combustion process is closely related to the heat transfer process of solid fuel surface.Based on the analysis of the heat transfer coefficient, the self-sustained combustion occurs around the reattachment point at first, and then gradually spread to the redevelopment zone.The heat released in the reattachment point will be used to achieve the self-sustained combustion in the redevelopment zone.

•
The linear relationship between the sudden expansion ratio with fixed port diameter and average/maximum regression rate was obtained.The result indicates that the average regression rate and maximum regression rate are more sensitive to the change of the port-to-inlet diameter.Additionally, the overall heat transfer behavior of the fuel surface was mostly dominated by the maximum heat transfer around the backward-facing step.

•
Based on the analysis of linear relationship between regression rate and sudden expansion ratio with fixed inlet diameter.The average regression rate is mainly affected by heat transfer mechanism in a fully developed turbulent flow in rthe edevelopment zone, and it decreased with the increasing of the sudden expansion ratio.

Figure 4 .
Figure 4. y + (y-axis) along the fuel surface (x-axis) of different mesh number.

Figure 4 .
Figure 4. y + (y-axis) along the fuel surface (x-axis) of different mesh number.

Figure 5 .
Figure 5. Axial distribution of local regression rate of different mesh number.

Figure 6 .
Figure 6.Axial velocity in different axial locations for different mesh number.

Figure 5 .
Figure 5. Axial distribution of local regression rate of different mesh number.

Figure 5 .
Figure 5. Axial distribution of local regression rate of different mesh number.

Figure 6 .
Figure 6.Axial velocity in different axial locations for different mesh number.

Figure 6 . 30 Figure 7 .
Figure 6.Axial velocity in different axial locations for different mesh number.Energies 2019, 12, x FOR PEER REVIEW 9 of 30

Figure 7 .
Figure 7. Tangential velocity in different axial locations for different mesh number.

Figure 7 .
Figure 7. Tangential velocity in different axial locations for different mesh number.

Figure 9 .
Figure 9. Normalized axial velocity (a) and tangential velocity (b) at different axial position.

Figure 10 .
Figure 10.The temperatures and pressures of simulation and experiment.

Figure 9 . 30 Figure 9 .
Figure 9. Normalized axial velocity (a) and tangential velocity (b) at different axial position.A famous benchmark case[40] is employed to validate the chemical reaction model.A diameter of 15 mm spherical projectile under the operating condition of the pressure of 42,662 Pa, temperature of 250 K, and the supersonic speed (1685 m/s) of the stoichiometric mixture (H 2 /Air) flow through the computational domain.An eight element reactions and seven species chemistry model is used for the simulation of the benchmark case.Moreover, Figures10 and 11show the results.As in these figures described, the simulation the results that were obtained from the in-house code reproduce the experimental data very well.It becomes evident that the accuracy of the combustion model could meet the requirements of this research.

Figure 10 .
Figure 10.The temperatures and pressures of simulation and experiment.Figure 10.The temperatures and pressures of simulation and experiment.

Figure 10 .
Figure 10.The temperatures and pressures of simulation and experiment.Figure 10.The temperatures and pressures of simulation and experiment.

30 Figure 11 .
Figure 11.The composition distributions of simulation and experiment. 0

Figure 11 .
Figure 11.The composition distributions of simulation and experiment.

30 Figure 11 .
Figure 11.The composition distributions of simulation and experiment.

Figure 12 .
Figure 12.Temperature profiles case (a) T w = 900 K and case (b) heat flux q w = 70 kw.

Figure 13 .
Figure 13.Internal surface of solid fuel after combustion.

Figure 14 .
Figure 14.Regression rate contour of HDPE solid fuel grain.

Figure 13 .
Figure 13.Internal surface of solid fuel after combustion.

Figure 13 .
Figure 13.Internal surface of solid fuel after combustion.

Figure 14 .
Figure 14.Regression rate contour of HDPE solid fuel grain.

Figure 14 .
Figure 14.Regression rate contour of HDPE solid fuel grain.

Energies 2019 , 30 Figure 15 .
Figure 15.The regression rate of simulation and experiment.

Figure 15 .
Figure 15.The regression rate of simulation and experiment.

Figure 16 .Figure 17 .
Figure 16.The contour of temperature and streamline.

Figure 17 .
Figure 17.The contours of (a) mass fraction of CO2, (b) mass fraction of CO, and (c) mass fraction of O 2 .

Figure 18 .
Figure 18.Radial profiles of axial (a) and tangential velocity (b) at different axial positions along the combustor.

Figure 18 .
Figure 18.Radial profiles of axial (a) and tangential velocity (b) at different axial positions along the combustor.

Figure 19 .
Figure 19.Temperature at different axial positions in the combustor.Figure 19.Temperature at different axial positions in the combustor.

Figure 19 .
Figure 19.Temperature at different axial positions in the combustor.Figure 19.Temperature at different axial positions in the combustor.Energies 2019, 12, x FOR PEER REVIEW 17 of 30

Figure 20 .
Figure 20.The local regression rate of different time.Figure 20.The local regression rate of different time.

Figure 20 .
Figure 20.The local regression rate of different time.Figure 20.The local regression rate of different time.

Figure 20 .
The local regression rate of different time.

Figure 21 .
Figure 21.The profiles of local heat flux and regression rate along the solid surface.

Figure 21 .Figure 22 .
Figure 21.The profiles of local heat flux and regression rate along the solid surface.Energies 2019, 12, x FOR PEER REVIEW 18 of 30

Figure 23 .
Figure 23.The heat transfer coefficient along high-density Polyethylene (HDPE) surface of different time for swirl flow.

Figure 23 .
Figure 23.The heat transfer coefficient along high-density Polyethylene (HDPE) surface of different time for swirl flow.
Figure24introduces the local regression rate of different port-to-inlet diameters.In this figure, the inlet diameter ranges from 30 mm to 50 mm and the port diameter fixed at 80 mm.Further, the d p /d in (sudden expansion ratio) ranges from 1.6-2.67.As the analysis above, the regression rate rapidly grows in the recirculation zone and it gradually decreases downstream of the flow field.However, for different inlet diameters, the regression rates in the reattachment point and redevelopment zone for each case are increased with increasing of d p /d in .Figure25indicates the heat flux for different cases, the heat flux introduced the same trend as the regression rate, in the area of reattachment point and redevelopment zone, and the heat flux increases with increasing of d p /d in .Energies 2019, 12, x FOR PEER REVIEW 20 of 30

Figure 24 .
Figure 24.Local regression rate obtained from simulation for different cases.

Figure 25 .
Figure 25.The local heat flux for different cases.

Figure 24 .
Figure 24.Local regression rate obtained from simulation for different cases.

Figure 24 .
Figure 24.Local regression rate obtained from simulation for different cases.

Figure 25 .
Figure 25.The local heat flux for different cases.

Figure 25 .
Figure 25.The local heat flux for different cases.

Figure 26 .
Figure 26.Temperature (a) and temperature around fuel surface (b) at different axial positions of different cases.

Figure 26 .
Figure 26.Temperature (a) and temperature around fuel surface (b) at different axial positions of different cases.

Figure 28 .
Figure 28.Local regression rate obtained from simulation for different cases.

Figure 27 .
Figure 27.Turbulent viscosity (µ t ) at different axial positions of different cases.

Figure 28 .
Figure 28.Local regression rate obtained from simulation for different cases.Figure 28.Local regression rate obtained from simulation for different cases.

Figure 28 .
Figure 28.Local regression rate obtained from simulation for different cases.Figure 28.Local regression rate obtained from simulation for different cases.

Figure 29 .
Figure 29.The local heat flux along the solid surface.

Figure 29 .
Figure 29.The local heat flux along the solid surface.

Energies 2019 , 30 Figure 29 .
Figure 29.The local heat flux along the solid surface.

Figure 31
Figure31shows the distributions of turbulent viscosity, in recirculation zone, the μt of different cases are almost the same, and at the position of reattachment point, the μt are irregular with dp/din.However, the μt is inversely proportional to dp/din in the redevelopment zone.

Figure 30 .
Figure 30.Temperature (a) and temperature around HDPE grain surface (b) at different axial positions of different cases.

Figure 31
Figure31shows the distributions of turbulent viscosity, in recirculation zone, the µ t of different cases are almost the same, and at the position of reattachment point, the µ t are irregular with d p /d in.However, the µ t is inversely proportional to d p /d in in the redevelopment zone.

Figure 30 .
Figure 30.Temperature (a) and temperature around HDPE grain surface (b) at different axial positions of different cases.

Figure 31 .
Figure 31.Radial profiles of the turbulent viscosity (μt) at different axial positions in , solid fuel ramjet (SFRJ) combustor of different cases.

Figure 31 .
Figure 31.Radial profiles of the turbulent viscosity (µ t ) at different axial positions in, solid fuel ramjet (SFRJ) combustor of different cases.

Figure 32 .
Figure 32.Regression rate for different sudden expansion ratio with dp fixed as constant.

Figure 32 .
Figure 32.Regression rate for different sudden expansion ratio with d p fixed as constant.

Figure 33 .
Figure 33.The regression rate for different sudden expansion ratio with din fixed as constant.

Figure 34 .
Figure 34.The regression rate for different sudden expansion ratio with dp fixed as constant.

Figure 33 .
Figure 33.The regression rate for different sudden expansion ratio with d in fixed as constant.

Energies 2019 , 30 Figure 33 .
Figure 33.The regression rate for different sudden expansion ratio with din fixed as constant.

Figure 34 .
Figure 34.The regression rate for different sudden expansion ratio with dp fixed as constant.Figure 34.The regression rate for different sudden expansion ratio with d p fixed as constant.

Figure 34 .
Figure 34.The regression rate for different sudden expansion ratio with dp fixed as constant.Figure 34.The regression rate for different sudden expansion ratio with d p fixed as constant.

Figure 35 .
Figure 35.The regression rate for different sudden expansion ratio with din fixed as constant.

Figure 35 .
Figure 35.The regression rate for different sudden expansion ratio with d in fixed as constant.

Table 3 .
Deviation analysis of average regression rate of different mesh number.

Table 3 .
Deviation analysis of average regression rate of different mesh number.

Table 3 .
Deviation analysis of average regression rate of different mesh number.

Table 4 .
The characteristic velocities and average regression rate that obtained from numerical and experimental result.
Author Contributions: Formal analysis, Investigation, Methodology, Writing-original draft, W.L.; Supervision, X.C.; Validation, W.C.; Writing-review & editing, O.M.This research was funded by NATIONAL NATURAL SCIENCE FOUNDATION OF CHINA, grant number 51306092; Jiangsu Postdoctoral Research Foundation, grant number AD41872.The authors declare no conflict of interest. Funding: