Maximum thrust nozzle based on height constraint

1. Introduction

A hypersonic vehicle refers to an aircraft that cruises in or across the atmosphere at a speed of more than Ma5. Scramjet, one of the main engine forms of the hypersonic vehicle, has a simple structure, does not need to carry oxidant, and can take off and land horizontally[1,2]. Scramjet has four parts: inlet, isolator, combustor, and nozzle. The nozzle is the main component that provides thrust and lift. Given that a hypersonic vehicle has a wide flight envelope, the nozzle performs well in a large range of nozzle pressure ratio (NPR), and the area ratio of nozzle outlet and inlet usually reaches 10. Researchers [3] showed that a symmetric nozzle cannot satisfy such application requirements; thus, the nozzle with an asymmetric geometry structure is developed. In addition, the nozzle’s geometry should be constrained because of the highly integrated characteristics and design requirements of a hypersonic vehicle [4]. The height of the nozzle outlet affects the NPR at the design point, the aircraft’s windward drag, and the nozzle’s net thrust. Therefore, the nozzle height should be constrained directly.

Rao [5] derived the maximum thrust theory with the method of calculus of variations to improve the nozzle aerodynamic performance. This method is widely used in the nozzle design of rocket engines. Lu [6] proposed a design method for an asymmetric nozzle based on the streamline tracing method and maximum thrust theory. Results showed that compared with the nozzle with straight walls, the stream-traced nozzle obtains a 2.7% increase in thrust and a 69.5% increase in lift. Lv [7] proposed an asymmetric nozzle design method based on the method of characteristics (MOC) under geometric constraints and studied the influence of flap length on nozzle performance. Their proposed method can increase the axial thrust coefficient, lift, and pitching moment of the nozzle by 5.5%, 1098.2%, and 20.3%, respectively. Yu [8] proposed an inverse design method of the nozzle based on maximum thrust theory and MOC, thereby improving the thrust coefficient, lift, and pitching moment by 31.8%, 201%, and 56.6%, respectively. Subsequently, Yu [9] developed a nozzle design method based on the inverse design method, which can constrain the absolute size and position of the inlet and outlet. Chen [10] proposed a design method for the 2D thrust-optimized single expansion ramp nozzle(SERN), which increased the axial thrust by 6.38% and the lift by 180%, further revealing the relationship between aerodynamic performance and geometric design conditions.

The nozzle designed using MOC usually has good aerodynamic performance. However, its length is extremely large, which is not conducive to the integrated designs of the aircraft and the engine. Researchers proposed the theory of minimum length nozzle and truncated the nozzle to shorten the nozzle length further. Argrow [11] developed minimal length nozzle theory. However, the ideal expansion nozzle designed by this method is still extremely long. Thus, the geometric constraints of the thrust nozzle cannot be satisfied. Shyne [12] proposed an improved method, enabling the nozzle to satisfy the size constraints and reduce the weight. The results showed that the truncation requirements should make a trade-off between the weight and performance loss of the nozzle. Hoffman [13] proposed a design method for the compressed truncated nozzle, which can effectively shorten the nozzle length. The performance difference between the compressed truncated nozzle and the Rao nozzle is only 0.04%–0.34%.

In the present study, a design method for 2D SERN is proposed under height constraint. The accuracy and rationality of CFD simulation are confirmed in Section 0. The design method for SERN is introduced in Section 3. The design parameters are studied parametrically, and the influences of their mechanisms on nozzle geometry and performance are studied in Section 4. In Section 5, the performance of the nozzle designed using the proposed design method and that designed using the traditional truncated design method are compared under the same design and off-design points to confirm the superiority of the proposed design method.

2. Numerical Approaches

The commercial software ANSYS Fluent is used to obtain the flow fields and the aerodynamic performance of the nozzles. The governing equations are briefly described as follows:

$$\begin{array}{c}\frac{\partial Q}{\partial t}+\frac{\partial E}{\partial x}+\frac{\partial F}{\partial y}=0\\ Q=\left[\begin{array}{l}\rho \\ \rho u\\ \rho v\\ \rho e\end{array}\right],E=\left[\begin{array}{l}\rho u\\ \rho u^2+p-{\tau }_{xx}\\ \rho uv-{\tau }_{xy}\\ u(\rho e+p)-u{\tau }_{xx}-v{\tau }_{xy}\end{array}\right],F=\left[\begin{array}{l}\rho v\\ \rho uv-{\tau }_{yx}\\ \rho v^2+p-{\tau }_{yy}\\ v(\rho e+p)-u{\tau }_{xy}-v{\tau }_{yy}\end{array}\right]\\ {\tau }_{xx}=-\frac{2}{3}\mu (\nabla ·\overrightarrow{V})+2\mu \frac{\partial u}{\partial x},{\tau }_{yy}=-\frac{2}{3}\mu (\nabla ·\overrightarrow{V})+2\mu \frac{\partial v}{\partial y}\\ {\tau }_{xy}={\tau }_{yx}=\mu (\frac{\partial u}{\partial y}+\frac{\partial v}{\partial x})\\ e=\frac{u^2+v^2}{2}+\frac{p}{\rho }\end{array}$$

(1)

where ρ, u, v, p, τ, and e represent the density, x-direction velocity component, y-direction velocity component, pressure, shearing stress, and total energy per unit mass, respectively. The inviscid flux is computed using the Roe flux-difference splitting scheme on the control surfaces. The second-order upwind scheme is used to discretize the governing equations, and the method based on the least squares cell is used to compute the gradient. During the calculation, the mass flow through the inlet and the force on the ramp and flap are used to monitor the convergence. The calculation can be considered convergent when the above calculation results remain constant with the increase in iteration steps. In the numerical simulation, the boundary conditions (i.e., pressure inlet, pressure outlet, and pressure far field) are adopted.

2.1. Validation of numerical method

The SERN experimental model of NASA is used to verify the reliability and accuracy of the numerical simulation. The nozzle named GR3 in reference [14] is selected. The area ratio of the outlet to the inlet is 1.993, and the NPR is set to 8.571. The geometric model of the nozzle is shown in Figure 1.

Figure 2 shows the computational grid, and Figure 3 shows the comparison of the numerical results of the normalized pressure distribution on the ramp and flap between CFD and the experiment. Two turbulence models, k-ω SST and realizable k-ε, are used in CFD simulation. The numerical results show that when x/h_t>8, the realizable findings of the k-ε model are closer to the experimental results than those of the k-ω SST model. Thus, all numerical simulations in the following work use the realizable k-ε turbulence model.

2.2. Grid resolution independence study

A mesh independence study is conducted to eliminate the influence of the mesh size on the results. The computational model is the nozzle designed using the proposed method in Section 4.1, as shown in Figure 4. The computational grids of the nozzle internal flow field with different scales are used, namely, coarse (75000 cells), medium (150000 cells), and fine (300000 cells). Figure 5 shows the normalized pressure distribution on the ramp with different grid scales under the same working conditions. In general, the pressure distributions obtained by the three grid sizes are the same. However, the maximum relative error between the coarse grid and the fine grid is 0.25% at x/h_t≈1.8 through the amplification at x/h_t≈1.8 and 5.1. Moreover, the maximum relative error between the medium grid and the fine grid is 0.18% at x/h_t≈5.1. A medium-scale computational grid is used in the following numerical simulation work to consider the numerical accuracy and cost comprehensively.

3. Nozzle Design Methodology

3.1. Maximum Thrust Nozzle Theory Under Height Constraint

The proposed design method draws on Rao’s theory of maximum thrust nozzle and derives the optimal design conditions of the thrust nozzle under height constraints. According to the characteristic line diagram of the maximum nozzle under the height constraint given in Figure 6, the mathematical expressions of the two constraints along the control surface DE are as follows:

$$H_{DE}=y_D+{\displaystyle {\int }_D^E\tan \phi dx}$$

(2)

$${\dot{m}}_{DE}={\displaystyle {\int }_D^E\rho V\frac{\sin (\phi -\theta )}{\cos \phi }dx}$$

(3)

where ρ, θ, and V refer to density, velocity angle, and velocity magnitude, respectively.

If p and p_a are local static pressure and ambient back pressure, respectively, then the axial thrust F acting on the nozzle control surface DE is obtained as follows:

$$F={\displaystyle {\int }_D^E[(p-p_a)\tan \phi +\rho V^2\frac{\sin (\phi -\theta )\cos \theta }{\cos \phi }]dx}$$

(4)

The Lagrangian multiplier method is applied to solve the maximum value of thrust F under flow constraint ṁ and height constraint H. This problem is equivalent to solving the unconditional maximum value of the following integral, as follows:

$$I={\displaystyle {\int }_D^{E}f+{\lambda }_1{f}_1+{\lambda }_2{f}_{2}dx}$$

(5)

where

$$f=(p-p_a)\tan \phi +\rho V^2\frac{\sin (\phi -\theta )\cos \theta }{\cos \phi }$$

(6)

$$f_1=\rho V\frac{\sin (\phi -\theta )}{\cos \phi }$$

(7)

$$f_2=\tan \phi$$

(8)

The three conditions below are obtained by solving λ₁ and λ₂, as follows:

$$\phi =\theta +\mu$$

(9)

$${\lambda }_1=-V\frac{\cos (\theta -\mu )}{\cos \mu }$$

(10)

$${\lambda }_2=-\frac{1}{2}\rho V^2\sin 2\theta \tan \mu -p+p_a$$

(11)

The corner condition is satisfied by substituting (9)–(11) into (5), as follows:

$$\sin \theta =0$$

(12)

Therefore, the axial thrust is maximum when the air flow direction is along the axial direction everywhere on the last characteristic line. Moreover, the aerodynamics parameters are uniformly distributed on the last characteristic line when the axial thrust is maximum.

3.2. MOC

The governing equations of the steady supersonic flow field are hyperbolic; thus, they can be solved by the MOC. The MOC is used to solve steady, 2D, planar, irrotational, supersonic flow. The characteristic line equations along the left/right running characteristic line (C₊/C_-) and their corresponding compatibility equations are as follows:

$$\frac{dy}{dx}=\tan (\theta \pm \mu )$$

(13)

$$\frac{\sqrt{M{a}^2-1}}{\rho V^2}d{p}_{\pm }\pm d{\theta }_{\pm }+\delta \frac{\sin \theta }{yMa\cos (\theta \pm \mu )}d{x}_{\pm }=0$$

(14)

In the compatibility equation, δ = 0 is a 2D flow, and δ = 1 is an axisymmetric flow. The subscript “+” indicates the left-running characteristic line, and the subscript “−“ indicates the right-running characteristic line.

The characteristic line equation along the streamline (C₀) and its corresponding compatibility equations are as follows:

$$\frac{dy}{dx}=\tan \theta$$

(15)

$$\rho VdV+d\rho =0$$

(16)

$$dp-a^{2}d\rho =0$$

(17)

where a is the speed of sound.

The detailed unit process of MOC can be found in [15].

3.3. Design Process of SERN

The design process of the maximum thrust nozzle satisfying the height constraint based on the MOC is as follows:

(1) The ideal Mach number of outlet Ma_ideal can be obtained by (18) according to the geometry parameters of the inlet and outlet.

$$\frac{A_e}{A_{in}}=\frac{\frac{1}{M{a}_{ideal}}{[\frac{2}{\kappa +1}(1+\frac{\kappa -1}{2}M{a}_{ideal}{}^2)]}^{\frac{\kappa +1}{2(\kappa -1)}}}{\frac{1}{M{a}_{in}}{[\frac{2}{\kappa +1}(1+\frac{\kappa -1}{2}M{a}_{in}{}^2)]}^{\frac{\kappa +1}{2(\kappa -1)}}}$$

(18)

where A_e, A_in, Ma_in, and κ are outlet section area, inlet section area, inlet Mach number, and heat capacity ratio.

(2) The Mach region I_uOI_d in Figure 7 Schematic of MOC design procedure is calculated.

(3) The kernel region I_uCDAI_dO is calculated by considering θ_u and θ_d. When θ_u and θ_d are provided, the asymmetrical factor β can be obtained as follows:

$$\beta =\frac{{\theta }_d}{{\theta }_u}$$

(19)

(4) Point E on CDA is determined, satisfying the following two conditions:

$$M{a}_E=M{a}_{ideal}{\beta }_{Ma}$$

(20)

$$\theta =0^{\circ }$$

(21)

Here, β_Ma is a Mach number scale factor, which is a design parameter of the proposed design method. The mechanism of β_Ma’s influence on the nozzle is studied subsequently.

If point E satisfying the conditions cannot be found, then θ_u and θ_d should be provided again. The kernel region I_uCDAI_dO is calculated until the point E that satisfies the conditions is found.

(5) According to mass flow conservation, the last characteristic line EF of the upper turning region is derived, and the ramp profile CF is solved. Section 3 shows that if the nozzle has the maximum thrust under the height constraint, then the last characteristic line EF satisfies (20) and (21) everywhere. Therefore, the aerodynamic parameters on EF are evenly distributed.

(6) Point G on CDA is identified, and the last characteristic line GH of the lower turning region is solved according to mass flow conservation. Thus, the nozzle outlet HF can satisfy the nozzle height constraint. The aerodynamic parameters are evenly distributed on GH.

(7) The flap profile AH is calculated.

4. Parametric study of design parameters

4.1. Typical nozzle design

The design performance of the proposed nozzle can be evaluated based on an actual design case. The nozzle’s length is set to 800 mm. The flight Mach number Ma_∞ is 5.

Table 1 lists the parameters of the design condition and nozzle geometric constraints required by the aircraft. The aerodynamic parameters of the nozzle inlet and the geometric constraint of the nozzle outlet height are provided. The height of the inlet and outlet of the nozzle is 157 and 400 mm, respectively. The nozzle’s length is set to 800 mm. The flight Mach number Ma_∞ is 5.

In the design program, the initial θ₁, β_Ma, and β should be provided and are set to 30°, 1, and 0.75, respectively. The comparison of the normalized pressure contour between inviscid CFD and MOC is shown in Figure 8. The result shows the accuracy of the proposed design method.

4.2. Definition of nozzle performance parameters

The nozzle is the key component of the scramjet to generate thrust, which is an important performance parameter of the scramjet nozzle. The axial thrust and lift of the nozzle are the horizontal and vertical components of the integral of the nozzle wall pressure on the walls, respectively. The mathematical expressions are as follows:

$$F_x={\displaystyle \underset{wall}{\int }(p_w-p_a)dy}$$

(22)

$$F_y={\displaystyle \underset{wall}{\int }(p_w-p_a)dx}$$

(23)

$$F_{x,ideal}=\dot{m}\sqrt{\frac{2\kappa R}{\kappa -1}{T}_t[1-{(\frac{p_a}{p_t})}^{\frac{\kappa -1}{\kappa }}]}-I_{in}$$

(24)

$$I_{in}=\dot{m}{\overline{V}}_{in}+(p_{in}-p_a)A_{in}$$

(25)

$$C_{fx}=\frac{F_x}{F_{x,ideal}}$$

(26)

$$C_L=\frac{F_y}{F_{x,ideal}}$$

(27)

where p_w is the wall pressure, p_a is the ambient pressure, ṁ is mass flow rate, T_t is the nozzle inlet total temperature, p_t is the nozzle inlet total pressure, ${\overline{V}}_{in}$ is the average axial velocity of the nozzle inlet, p_in is the nozzle inlet static pressure, and R is the gas constant.

4.3. Influence of β

The asymmetrical factor β mainly changes the shape of the kernel region I_uCDAI_dO in Figure 7 and the nozzle geometry afterward. In the design process, θ_u remains constant when the inlet parameters, outlet height constraints, expansion arc radius, and β_Ma remain unchanged. The value of θ_d is determined when β is fixed.

Figure 9 shows the effect of different β values on the axial thrust coefficient and lift of the nozzle. In Figure 10, β changes only the length of the nozzle’s flap. However, it does not change the flap’s height and the ramp’s profile. The flap length decreases with the decrease in β, and the area change rate of the nozzle along the axial direction $\frac{dA}{dx}$ increases. Thus, when β is extremely small, the separation zone appears near the flap because of the large $\frac{dA}{dx}$. The effectiveness of the design method is lost due to the separation zone. Figure 9 shows that β is positively correlated with the axial thrust coefficient C_fx and negatively correlated with the lift force. Figure 10 shows that the streamline near the lower wall deviates from the axial direction downward with the decrease in β. This phenomenon results in decreased axial momentum on the outlet, decreased axial thrust, and increased lift.

4.4. Influence of β_Ma

The Mach number scale factor β_Ma can control the total nozzle length under the specified height constraint. The influence of β_Ma on the nozzle profile is shown in Figure 11. The figure shows that the increase in Mach number scale factor β_Ma increases the ramp length while decreasing the flap length to satisfy the height constraint and β.

Figure 12 shows the influence of different β_Ma values on the thrust coefficient and lift of the nozzle when β = 0.75. The figure shows that the increase in β_Ma increases the total length of the nozzle and decreases the length of the nozzle’s flap. Thus, the positive lift provided by the ramp increases, and the negative lift provided by the flap decreases simultaneously. The increase in nozzle geometric asymmetry increases the inhomogeneity of the aerodynamic parameters at the nozzle outlet, thereby decreasing the nozzle axial thrust coefficient. Moreover, the calculation results show that extremely large β_Ma causes serious nozzle deformation, further leading to serious deterioration of nozzle performance.

The influence of β_Ma in Figure 11 on the nozzle geometry and the nozzle design process in Figure 7 indicate that β_Ma can be considered a nozzle truncation design method. When β_Ma is unequal to 1, the upper wall I_uCF is designed according to Ma = Ma_E instead of Ma = Ma_ideal. Therefore, the nozzle design flow when β_Ma is larger than 1 is as follows: (1) The nozzle is designed according to the outlet Mach number Ma_c in the characteristic line design program, and the area ratio of the inlet and outlet of the nozzle conforms to (18). (2) The ramp remains motionless, and the flap is truncated to satisfy the specified height constraint.

4.5. Adjusting nozzle geometry by design parameters

Sections 0 and 4.4 state that in the proposed design method, the lengths of the flap and the ramp can be controlled using the asymmetrical factor β and the Mach number scale factor β_Ma under the given height constraint. Therefore, the two design parameters can be used to adjust nozzle geometry. The influences of β and β_Ma on the nozzle geometry are further illustrated by three design examples.

The three sets of geometric constraints are provided by maintaining the nozzle inlet conditions and height constraints in Section 4.1, as shown in Table 2. The nozzle in Section 4.1 is termed Nozzle B.

The design process is shown in Figure 13. The design should be initially carried out according to β_Ma = 1. If the total length constraint cannot be satisfied, then β_Ma should be adjusted and recalculated to satisfy the total length constraint. Then, β_Ma is maintained, and β is adjusted to satisfy the constraint of the flap length.

The nozzles satisfying the geometric constraints in Table 2 are designed using the de-sign program in Section 3.3, and they have the same height constraint as Nozzle B. The β and β_Ma of each nozzle are shown in Table 3, and their normalized pressure contours are shown in Figure 14.

5. Aerodynamic performance comparisons with traditional method

5.1. Performance comparisons

The proposed design method is compared with the traditional truncation design method to prove its superiority. The typical design method is used to design the fully expanded nozzle according to the inlet and outlet conditions. Then, this nozzle is cut proportionately to satisfy the geometric constraints. The specific design steps are shown in [16].

The nozzle designed by the truncation design method is referred to as Nozzle A. Figure 15 shows that the nonuniformity of pressure distribution on the outlet of Nozzle A is greater than that of Nozzle B when the nozzle outlet height and nozzle inlet parameters are the same. Table 4 lists the axial thrust coefficient and lift of Nozzle A and Nozzle B. The data in Table 4 are compared. The results indicate that the thrust and lift coefficients are significantly improved. The lift of Nozzle B is larger than that of Nozzle A, probably because the velocity angle equal to zero is satisfied on the last characteristic line when the ramp profile is designed. Thus, the ramp of Nozzle B produces less negative lift than that of Nozzle A.

5.2. Off-design Performance

The nozzle should provide sufficient thrust when the vehicle gains altitude and takes off horizontally. Thus, the aerodynamic performance under off-design operating points should be studied. The aerodynamic performances under the off-design operating points of Nozzle A and Nozzle B are compared to reflect the superiority of the proposed design method. Some typical off-design operating points are selected, as shown in Table 5. When NPR is 25.35, the nozzle is under the equal expansion condition at the inlet and outlet sizes according to (18), and NPR is 69.22 at the design point.

In Table 5, operating points Ma_∞ = 2.5 and Ma_∞ = 3 are in the overexpansion state. The operating points Ma_∞ = 3.5 and Ma_∞ = 4 are in the state of under expansion. However, the NPR is less than that at the design point. In addition, the transition point of the combined cycle engine’s working mode is usually near Ma_∞ = 2.5 in actual flight; thus, the typical off-design working conditions in Table 5 reflect the working performance of the nozzle in the full flight envelope under the ramjet mode.

Figure 16 shows the normalized pressure contours under various typical off-design conditions in Table 5. Table 6 shows the thrust and lift coefficient under various off-design conditions. The aerodynamics performance of Nozzle A under various off-design conditions in Table 5 is also calculated. The normalized pressure contours are shown in Figure 17, and the thrust and lift coefficients are shown in Table 7.

Figure 16 shows that the increase in NPR weakens the intensity of the oblique shock wave at the nozzle outlet, increases the downward angle of the streamline near the flap outlet, and makes the streamline near the outlet of the ramp relatively close to the horizontal direction. The change trend of the streamline near the flap in Figure 17 is the same as that in Figure 16. However, the streamline near the ramp is always upward. The positive lift generated by Nozzle B is larger than that by Nozzle A under the same working conditions because of the difference in the streamline on the ramp between Nozzle A and Nozzle B. After the thrust data comparison, the thrust coefficient of Nozzle B is slightly larger than that of Nozzle A. This finding proves that the aerodynamic performance of the nozzle by the proposed design method is better than that by the traditional proportional truncation method under off-design conditions.

6. Conclusion

In this study, the height constraint is directly introduced into the maximum thrust nozzle theory, and the maximum thrust nozzle satisfying the height constraint is designed using MOC. In this design method, the nozzle’s ramp and flap lengths are controlled by adjusting the asymmetrical factor β and the Mach number scale factor β_Ma under the specified outlet height constraint. The main conclusions are as follows:

1. β_Ma mainly controls the total length of the nozzle. The total length of the nozzle increases with the increase in β_Ma. The flap length decreases, and the positive lift of the nozzle increases to satisfy the height constraint. The increase in β_Ma increases the geometric asymmetry of the nozzle while decreasing C_fx and increasing C_L. In the actual design, extremely large β_Ma leads to a rapid deterioration of nozzle performance when β and the height constraint remain.

2. The influence of β on nozzle geometry is that only the flap length changes when the ramp profile and the flap height are unchanged. The flap length increases with the increase in β. Moreover, the streamline near the flap becomes close to the axial direction, thereby increasing the thrust coefficient and decreasing the lift. In the actual design, the nozzle is designed under the specified geometric constraints by adjusting β_Ma and β.

3. Compared with those of the nozzle (Nozzle A) designed using the traditional truncated design method, the thrust and lift coefficients of the nozzle (Nozzle B) designed using the proposed method are increased by 11.93% and 138.45%, respectively, at the design point. At the same off-design point, the thrust and lift coefficients of Nozzle B are greater than those of Nozzle A. Compared with the thrust coefficient of Nozzle A, that of Nozzle B is increased by 8.79% at most under various typical off-design working conditions. Moreover, the positive lift force is maintained in the entire working condition range.