International Journal of Aerospace Engineering
Volume 2021 (2021), Article ID 5574436, 20 pages
https://doi.org/10.1155/2021/5574436

Research Article

Design Optimization and Parameter Analysis of a Hybrid Rocket Motor-Powered Small LEO Launch Vehicle

Hao Zhu,1,2 Pengcheng Wang,1,2 Weile Xu,1,2 Yuanjun Zhang,1,2 Hui Tian,1,2 and Guobiao Cai1,2

¹School of Astronautics, Beihang University, Beijing 100191, China
²Key Laboratory of Spacecraft Design Optimization & Dynamic Simulation Technologies of Ministry of Education, Beihang University, Beijing 100191, China

Correspondence should be addressed to Yuanjun Zhang; zhangyuanjun@buaa.edu.cn

Received 19 January 2021; Revised 19 March 2021; Accepted 7 April 2021; Published 29 April 2021

Academic Editor: Joseph Majdalani

Copyright © 2021 Hao Zhu et al. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

Abstract.

In this paper, the effects of different grain shapes of a hybrid rocket motor (HRM) and different payload mass/orbit heights on the design of small launch vehicles (SLVs) are systematically discussed. An integrated overall design model for the hybrid rocket motor-powered small launch vehicle (HPSLV) is established, and two groups of three-stage SLVs capable of sending small payloads to the low earth orbit (LEO) are designed and optimized. In the first group, the SLVs with different grain shapes and different numbers of chambers in HRMs at the 1st and the 2nd stages are optimized and analyzed. In the second group, the SLVs capable of sending different payload mass to different orbit heights are optimized and analyzed. Pareto graphs of the design results show that the design of HRM at the 1st stage has the greatest impact on the take-off mass, total velocity increment, and maximum axial overload of the SLV. Self-organizing maps show that the take-off mass, maximum diameter, overall length, and velocity increment of the SLVs have the same variation tendency. For the 1-chamber HRM at the 1st stage, the wheel-shaped grain is better than circle-shaped and star-shaped grains in terms of reducing the total mass and length of the SLV, and the 4-chamber parallel HRM has more advantages over all 1-chamber designs for the same reason. The theoretical velocity increments are calculated by the Tsiolkovsky formula, and the actual velocity increments are obtained based on the trajectory simulation data. The results indicate that the HPSLV has a regular distribution in terms of the ratio of theoretical (actual) velocity increments at three different stages, and the estimated distribution ratio is around 1 : 1.55 : 1.69 (1 : 1.9 : 2.39), which can provide some reference for future development of HPSLV.

1. Introduction

With the increasing demands for low earth orbit (LEO) payloads, microsatellites are paid much attention in recent years. As a result, these satellites (especially formation-flying satellites [1]) spark a huge demand for transportation systems from Earth’s surface to LEO, which includes heavy-lift rockets capable of carrying multiple payloads (e.g., Falcon 9 in the U.S., CZ-5 in China, and PSLV in India) and small launch vehicles (SLV) capable of carrying one or more payloads (e.g., Fast Boat 1 in China, Pegasus in the U.S., and SS520 in Japan). Among these transportation systems, SLVs have the advantages of quick response, low cost, high reliability, and propellant safety and have therefore become a promising development direction today.

A good propulsion system for SLV can further reduce cost, improve reliability, and increase payload capacity at the same time. The hybrid rocket motor (HRM) is such an excellent propulsion system. Figure 1 illustrates the working principle of the HRM. Different from the traditional solid and liquid rocket motors which use an oxidizer and fuel in the same state of matter, the HRM uses a liquid oxidizer and solid fuel which are stored separately. This enables the advantages of high safety and high reliability of the HRM but also leads to the problems such as lower combustion efficiency and changes of the oxidizer-fuel ratio during work. Meanwhile, the more complex fuel regression rate laws also increase the challenge of design. Research worldwide is dedicated to addressing these problems and has made great progresses. Bernard et al. proposed a one-dimensional nonsteady flow predictive performance model for N₂O/paraffin HRMs and proved the efficiency by comparing the proposed model with hot-fire test data [2]. A regression rate model for a mixed oxidizer hybrid motor, whose grains are lightly oxidized, was proposed using a linear combination of traditional regression rate equations, and tests were conducted to determine the best coefficients for the presumed models [3]. A parametric model for the cost and performance of HRM was established and adopted in multiobjective HRM design optimization [4]. Combustion performance and flow of HRM with multisegmented grain are presented through three-dimensional numerical simulation [5]. These studies reveal the working mechanism of HRM and at the same time improve the combustion efficiency and specific impulse efficiency to a certain extent. Existing research provides a more accurate fuel regression rate model and some other parameters at the same time. Moreover, deeper research on pressure oscillations of a solid rocket motor [6] and on combustion instability of a liquid rocket motor [7] can be further introduced into HRM studies. All of these establish the foundation for the application of HRM.

Figure 1: Framework of the HRM working principle.

Various applications of HRM in transportation systems are also explored around the world. The application of HRMs in manned lunar landing has been studied, and the uncertain factors during the design phase are also considered [8]. The HRM was considered the upper stage of a three-stage launcher, and a means with a nested direct/indirect procedure was developed to optimize the design of the propulsion system and the trajectory at the same time [9]. Purdue University attempted to send a 4.5 kg payload into the 150 km LEO [10]. A 150-200 kg payload is planned to be launched into the 250-300 km orbit in Russia using a two-stage HRM-powered SLV [11]. 98% H₂O₂/paraffin is used for the air/ground-launched SLV to send a 20 kg payload into the 300 km circle orbit [12]. Rhee et al. have conducted research on application of the HRM with HTPB/LOX combination to the first stage of an air launch system [13]. A three-stage launch vehicle is designed using the HRM with H₂O₂/HTPB propellant [14].

However, most of the applications above focus on the overall design of HRM-powered LVs. In fact, the working process of HRM follows the law of diffusion combustion, and the internal ballistic performance of the motor depends on the mass flux of a liquid oxidizer, which leads to a strong coupling relationship between the fuel regression rate, oxidizer mass flux, and combustion chamber pressure during the operation time. This mechanism makes the design of HRM more complicated. Moreover, there are different shapes of fuel grains and chamber configurations, which lead to different flight characteristics of the SLV. To obtain a comprehensive understanding of the design criteria of the HRM used in SLV, the design parameters for different types of HRMs need to be compared. Thus, the HRM performances are simulated in this paper, and the overall scheme of the HRM-powered SLV (HPSLV) is designed and optimized. The impacts of design parameters on SLV overall performances, the influence of different grain shapes on the 1st- and 2nd-stage propulsion performances, and the relationship between the desired orbit height/payload mass and the velocity increment of HPSLV are studied. This paper is devoted to revealing the characteristics of the HRM used in SLV.

This paper is organized as follows. In Section 2, a detailed description of the HRM is given, and a parameterized model for the HPSLV is established. The corresponding mathematical model for overall design optimization of the HPSLV is given in Section 3. The specific design optimization for SLV with different fuel grain shapes/chamber configurations and different orbit heights/payload mass is carried out in Section 4, and the results and discussion of parameters are also given.

2. Detailed Description of the Problem

2.1. Overall Configuration of the LV

In this study, three-stage SLVs with tandem configuration are designed to send small payloads (from 100 kg to 200 kg) to the LEO (from 300 to 700 km circle orbit). The baseline configuration of SLV is shown in Figure 2, including the fairing, payload, payload adapter, electronic avionics device and attitude control equipment, two interstage cabins (housing stage separation mechanism), and three HRMs.

Figure 2: SLV configuration.

2.2. Calculation of Propulsion Performances

The main parts of the SLV are the three stages of HRM. Therefore, a preliminary design code of typical HRMs, which consist of an oxidizer feed system and several (one or four) thrusters with solid fuel grains, is developed. The main propulsion performances of each HRM, such as the propellant mass, specific impulse, and thrust, are computed by this code, as shown in Figure 3 [8]. The design variables include two parts: the initial values of the grain shapes and the working conditions of motor burning.

Figure 3: Framework of the propulsion system design module.

2.2.1. Grain Shape Design and Equal Thickness Grain Slice Cutting

In this study, three types of fuel grain are considered: the circle-shaped grain, star-shaped grain, and wheel-shaped grain, as shown in Figure 4. After the initial shape of fuel grain is known, the grain port area and the burning perimeter length , which will be used in the primary design of the motor and internal ballistic simulation, are both determined by grain size parameters and burning distance. The initial shape parameters and geometric formulas of different grains are shown in Table 1.

(a)

(b)

(c)

Figure 4: Parameters of different grain shapes: (a) circle-shaped grain; (b) star-shaped grain; (c) wheel-shaped grain.

Table 1: Definition of initial grain shape parameters.

Grain shape
Burning perimeter length
Grain port area

Circle-shaped grain

Star-shaped grain
Related to the outer diameter , number of star angle , grain thickness , star angle coefficient , angle of star root , radius of star tip arc , and radius of star root arc [15]

Wheel-shaped grain
Related to the outer diameter , number of wheel hole , grain thickness , center hole diameter , wheel channel inner diameter , wheel channel outer diameter , and radius of chamfer [16]

2.2.2. Propellant Selection and Thermodynamic Calculation

A 98% mass concentration of H₂O₂ and a solid fuel, which consists of 60% of HTPB, 28% of Al, 10% of Mg, and 2% of carbon, are selected as the propellant combination in this study. A thermal calculation code is used to obtain the performance parameters of the propellant combinations. These parameters, such as the specific impulse , adiabatic combustion temperature , specific heat ratio , and characteristic velocity , are exported by inputting the chamber pressure , oxidizer-fuel ratio , and nozzle expansion ratio . Two different processes, combustion in the chamber and flow in the nozzle, are both considered in HRM thermal calculation. Three assumptions are given in chamber thermodynamic calculation: (1) the adiabatic combustion process is considered, (2) the chemical equilibrium state is established for combustion products, and (3) the perfect gas hypothesis is adopted. Based on the above assumptions, the adiabatic-chemical model of the combustion process is established. Moreover, in the nozzle thermodynamic calculation, an assumption is considered that perfect gas with a uniform component has no irreversibility (such as friction, heat transfer, and other imbalances) during flowing. Thus, an equilibrium-frozen flow model is adopted to obtain thermal results. The - curves of the selected propellant combination under different chamber pressures are shown in Figure 5. It can be seen that the best specific impulse appears when the ratio is around 4.0.

Figure 5: - curves under different chamber pressures.

2.2.3. Primary Design

At the initial time of motor working, the flux of propellant, grain size, and area of nozzle throat can be all determined by the initial force , initial oxidizer-fuel ratio , and specific impulse (provided by thermodynamic calculation), using

where is the oxidizer mass flow rate, is the fuel mass flow rate, is the total propellant mass flow rate, is the length of fuel grain, is the density of fuel grain, is the initial burning perimeter length, is the regression rate of fuel grain, is the nozzle throat area, and is the initial pressure of the chamber. Other parameters, such as chamber length and nozzle exit area, can be further acquired.

As a critical parameter affecting the internal ballistic performance, the regression rate of the fuel grain depends on the oxidizer mass flux and the hydraulic diameter , as shown in

where is the regression rate coefficient and is the flux rate exponent. In the model of the regression rate, the scale effect of the grain is taken into account, and the hydraulic diameter in Equation (5) corresponds to the value in a single combustion chamber.

2.2.4. Internal Ballistic Simulation

Along with the process of grain shape design, thermodynamic calculation provides necessary values needed in the phase of internal ballistic simulation, which is conducted to calculate the internal ballistic parameters and main dimensions of the thruster [8].

The purpose of this phase is to obtain the HRM performance parameters, such as the - curve and - curve, and provide the basis for the design of HRM components. The combustion gas works in two ways: one is exhausted through the nozzle, while the other is used to improve the storage capacity of combustion gas, as shown in

where is the effective volume of chamber and is the combustion gas density in the chamber. Considering the HRM combustion equations and introducing the equilibrium pressure , the chamber pressure can be obtained as shown in Equation (7). Since the density of the combustion gas is much smaller than that of the propellant, so the third item on the right can be ignored.

where refers to the equilibrium pressure (corresponding to the zero-dimensional process in a solid rocket motor) [8], indicates the specific gas constant, and means the thickness of a unit fuel slice. After calculating the chamber pressure at the corresponding time, the thermal calculation results are used to calculate the current specific impulse and thrust coefficient , with combustion efficiency being 0.96 and nozzle efficiency being 0.93. Such values are slightly smaller than those in References [17, 18] for a conservative design. Then, the current thrust value is computed by , in which means the number of combustion chambers.

2.2.5. Validation of Simulation

A series of long-time working HRM firing tests is conducted [19], and the data of these tests are used to modify the regression rate model and evaluate the internal ballistic performance. As shown in Figure 6, based on the test data of an HRM with a circle-shaped grain, the test and simulation curves of the internal ballistic performance including chamber pressure and thrust are compared to verify the accuracy of the propulsion design model. They can be seen to have nearly the same behaviour during the working process, and the average error between the simulation and the test is less than 5%, which indicates that the propulsion design model can be used as the basis for the overall design of HRM.

(a) Chamber pressure

(b) Thrust

Figure 6: Experimental and predicted curves: (a) chamber pressure; (b) thrust.

2.3. Structural Design of Components

By the calculation of propulsion, HRM performance parameters are exported as outputs, which provides the parameters needed in calculation of the mass and dimension of HRM.

2.3.1. Overall Structural Description

The rocket take-off mass is a sum of the mass of all components. As shown in Figure 2, the fairing covers both the payload and the 3rd-stage motor as part of the skin of the 3rd-stage motor and improves the mass fraction of the 3rd-stage HRM. Components except for the motor, such as payload adapter, avionics device and attitude control equipment, interstages, and tail cabin, are set as constants based on expert experience. The Carbon Fiber- (CF-) reinforced composite material is used as the main shell material of the main components, including fairing, cover structure, nozzle, chamber, oxidizer tank, and gas bottle. Moreover, a semiellipsoid and a cylinder compose the whole fairing, and the thickness of fairing is set by experience.

2.3.2. HRM Description

The main components of an HRM are divided into two parts: solid part (solid fuel grain, chamber, nozzle, etc.) and liquid part (liquid oxidizer, oxidizer tank, gas bottle, valves and tubes, etc.), as shown in Table 2. H₂O₂ is stored in the oxidizer tank and transported into the chamber by a turbopump feed system or a gas pressure feed system. When it comes to the gas pressure feed system, the gas movement in the gas bottle is considered an adiabatic expansion process, and that in the oxidizer tank is considered a constant pressure process. The mass of pressurizing gas and the volumes of the oxidizer tank and gas bottle are obtained according to the first law of thermodynamics, the mass conservation law, and ideal gas state equation. The chamber is used to pack fuel grain, and the outer diameter of a single chamber consists of a diameter of fuel grain, thickness of heat insulating layer, and thickness of a chamber shell. A catalyst bed is also adopted in this study, which is used to decompose H₂O₂ into oxygen and water gas with high temperature to enhance combustion efficiency.

Table 2: Main structure of the SLV propulsion system.

Main structure
Material
Function

Liquid oxidizer
98% H₂O₂
Provide chemical energy

Solid fuel
HTPB-based grain
Provide chemical energy

Gas bottle
CF wound aluminum liner
Contain pressurization gas

Oxidizer tank
CF wound aluminum liner and oxidizer sac
Contain liquid oxidizer

Chamber and nozzle
CF wound shell and high-silica insulation
Contain fuel grain and high-pressure burn gas

Valves and tubes
Mainly aluminum alloy
Feed liquid propellant

Catalyst bed
Stainless-steel shell and silver mesh
Decompose H₂O₂

Turbopump
Details refer to Reference [20]
Increase the pressure of oxidizer

Thrust vector control actuator
Given by experience
Control the thrust vector direction

The mass and dimension of these components are calculated based on the results of the propulsion design module. All the thickness of the chamber, nozzle, oxidizer tank, and gas bottle is calculated based on the maximum stress intensity theory, while the mass and dimension of some components are given by experience, such as the catalyst bed, turbopump, actuator of thrust vector control, and valves and tubes.

In addition, the layout of these components in the HRM at different stages (1st/2nd/3rd stage) needs to be arranged properly.

(1) 1st and 2nd Stages of SLV. Two different HRM layouts are adopted in the 1st and the 2nd stage of SLV, one chamber and four identical paralleled chambers, as shown in Figure 7. Three different grain shapes (circle-shaped, star-shaped, and wheel-shaped) are adopted in the 1-chamber HRM, while the circle-shaped grain is adopted in each of the 4-chamber HRM. In the 1-chamber HRM, the chamber outer shell is also part of the rocket body. However, in the 4-chamber HRM, an extra frame-covering structure, whose mass and dimension can be deduced from the chamber dimension when the thickness and material are selected, is designed as part of the rocket body, as shown in Figure 7(b). The turbopump feed system is used in the 1st/2nd stage, and the main components are distributed in a tandem framework.

(a) 1-chamber HRM

(b) 4-chamber HRM

Figure 7: Different types of HRM in the 1st/2nd stage of SLV: (a) 1-chamber HRM; (b) 4-chamber HRM.

(2) 3rd Stage of SLV. Figure 8 describes the main structure of the 3rd stage of SLV. To reduce the total length of the stage and obtain a shape with a low length-diameter ratio for the whole SLV, the HRM at the 3rd stage adopts four thrusters in parallel. All the thrusters are distributed symmetrically around one oxidizer tank. As at the 3rd stage the SLV flights in a nearly vacuum space, a helium gas pressure feed system is used at this stage to simplify the structure. Four gas bottles are used and located in parallel around the oxidizer tank. The frame-covering structure is also used to integrate the main components of the motor. The circle-shaped grains are used at the 3rd stage, in which four chambers contain the same shape of grain.

Figure 8: Configuration of 3rd-stage HRM.

2.4. Trajectory Simulation

The mission of the SLV in this paper is to send a small payload to a LEO circle orbit, and the schematic diagram of the flight process is shown in Figure 9. In order to describe the motion of the rocket, several coordinates and their transformation need to be applied. The coordinates and angle relationship used in this study can be found in Reference [21].

Figure 9: Schematic diagram of the flight process.

2.4.1. Flight Plan Design

According to the flight process of SLV, different flight phases correspond to different control strategies. The pitch angle , attack angle , and trajectory tilt angle during the vertical take-off phase, program-turn phase, zero-attack-angle flight phase, and the 2nd- (or 3rd-) stage flight phase are decided by the following equations, respectively:

where is the maximum absolute value of attack angle during the program-turn phase, is a constant coefficient, is the pitch angular velocity, is the pitch angle at current time, and is time.

Besides, the above phases are powered by three HRMs, a free glide phase without power is considered to turn the velocity direction from radial direction to tangential direction of the required orbit, and the control equation is

Moreover, separation between stages occurs after shutdown of the 1st- (or 2nd-) stage HRM, as the 2nd- (or 3rd-) stage HRM fires several seconds later. The pitch angle during separation can be considered a constant, and the control equation is

2.4.2. Equations of Motion

The three-degree-of-freedom (3DOF) mass point trajectory considering the rotation of the earth is used in this paper. The dynamics equations and kinematics equations of the SLV flight are shown as Equations (14) and (15), respectively. Take the -direction as an example, is the position coordinate of the SLV, where is the initial position of the launch point. is the velocity component along the -direction in the launch coordinate, and the initial velocity of SLV . is the acceleration along the -direction in the launch coordinate. is the rocket mass at current time. is the thrust along the -direction in the launch coordinate (including controlling force) and can be calculated from Section 2.2 in the body coordinate without regarding thrust eccentricity and thrust deflection. is the aerodynamic force along the -direction in the launch coordinate and can be obtained according to the aerodynamic force in the body coordinate. , , and are the accelerations in the launch coordinate caused by gravity, centrifugal inertial force, and Coriolis inertial force, respectively. The trajectory parameters in the -direction and -direction are calculated using the same method as that in the -direction. Then, Equations (14) and (15) are solved with the fourth-order Runge-Kutta method to obtain the trajectory parameters.

2.4.3. Vehicle Aerodynamics

The equations for aerodynamic force are only a rough approximation based on small attack/sideslip angle hypothesis and are acquired as shown in Equation (16). The drag coefficient and the derivative of lift coefficient to attack angle are estimated using the data of the “Titan II” rocket [21], whose shape is similar to that of the SLV designed in this study. The specific values of and are shown in Table 3, and the 1976 U.S. Standard Atmosphere data are coded as the SLV working environment [22].

Table 3: Approximation aerodynamic coefficients of the “Titan II” rocket.

Aerodynamic coefficients
Velocity range
Formula

0.29

2.8

3.55

2.4.4. Calculation of Orbit Characteristic

Ignoring rarefaction gas and other negligible factors, the status at the moment of 3rd-stage shutdown directly determines the parameters of payload orbit. One method for acquiring these orbital characteristics is to use turnoff velocity, height, and local trajectory tilt angle, as shown in

where is the turnoff velocity of the 3rd stage, is the local trajectory tilt angle at the turnoff point, is the universal gravitational constant, is the earth mass, is the distance from the turnoff point to the earth’s core, is the orbital semimajor axis, is the orbital ellipticity, is the perigee altitude, and is the apogee altitude.

3. Design Optimization

3.1. Optimal Design Problem Statement

Based on difference of the grain shape, payload mass, and orbit height, two groups of different SLVs are designed and optimized. The SLVs in the first group are designed to deliver a 100 kg payload to a 300 km attitude orbit. Six cases with different 1st/2nd stages are considered in the first group to study how different grain shapes at different stages influence rocket performances. The SLVs in the second group are designed to deliver payloads with different mass (from 100 kg to 200 kg) to the orbits with different attitudes (from 300 km to 700 km). Five cases with the same SLV design scheme are considered. The purpose of the second group is to study how payload mass and orbit height influence HPSLV. All the cases mentioned above use the same solving method given in Section 2, and their design variables, optimization objectives, and boundary constraints are identical (but not of the same value) and are analyzed in the next section.

3.2. Variables, Target, and Constraints of SLV Optimal Design

3.2.1. Variables

Table 4 shows the range of design variables, where DV is the abbreviation of design variable. The lower bound and upper bound of the design variables are given according to engineering experience and previous works [14], which cover the general operating states of HPSLV.

Table 4: Description and ranges of variables.

Variable
Description
Unit
Lower bound
Upper bound

DV1
Diameter of the 1st stage
m
1
2.5

DV2
Diameter ratio of the 2nd stage to1st stage
—
0.7
1

DV3
Ratio of the envelope diameter of the 2nd HRM to the diameter of the 2nd stage
—
0.5
1

DV4
Diameter ratio of the 3rd tank to the 3rd stage
—
0.3
0.7

DV5
Ratio of the envelope diameter of the 3rd HRM to the diameter of the 3rd stage
—
0.3
1

DV6
Initial thickness of the 1st grain
m
0.05
0.4

DV7
Initial thickness of the 2nd grain
m
0.05
0.4

DV8
Initial thickness of the 3rd grain
m
0.05
0.4

DV9
Initial thrust at the 1st stage
kN
40
1250

DV10
Initial thrust at the 2nd stage
kN
2
650

DV11
Initial thrust at the 3rd stage
kN
0.25
85

DV12
Initial pressure at the 1st chamber
MPa
2
7

DV13
Initial pressure at the 2nd chamber
MPa
2
7

DV14
Initial pressure at the 3rd chamber
MPa
2
7

DV15
Initial oxidizer-fuel ratio at the 1st stage
—
2
7

DV16
Initial oxidizer-fuel ratio at the 2nd stage
—
2
7

DV17
Initial oxidizer-fuel ratio at the 3rd stage
—
2
7

DV18
Expansion ratio at the 1st stage
—
10
50

DV19
Expansion ratio at the 2nd stage
—
10
100

DV20
Expansion ratio at the 3rd stage
—
10
100

DV21
Maximum attack angle during the program-turn phase
°
0.1
3

DV22
Control parameter of attack angle during the program-turn phase
—
0.15
1

DV23
Glide time between 2nd and 3rd stages
s
0
300

As shown in Figure 10, refers to the envelope diameter of the 2nd HRM, while is that of the 3rd HRM. and correspond to the diameter of the 2nd and 3rd stages, respectively. The design variable in Table 4 is the ratio of to . According to the geometric relationship, the conversion between the diameter of the components and the diameter of the rocket can be calculated.

(a)

(b)

Figure 10: Envelope diameter of the 2nd HRM and 3rd HRM: (a) cross section of the 2nd stage (four-chamber); (b) cross section of the 3rd stage.

3.2.2. Constraints and Target of SLV Optimal Design

(1) Propulsion Constraints. Too large oxidizer mass flux will cause erosive burning. Therefore, the maximum oxidizer mass flux of each stage is set to 500 kg/(sm²).

(2) Structure Constraints. Considering that a too high is not good for the rocket structure strength, is set as smaller than 16. For protecting payload, the maximum axial overload is set below 16 g, while the maximum normal overload is set below 1 g.

(3) Trajectory Constraints. The maximum dynamic pressure during the flight needs to be smaller than 0.1 MPa. The final orbit height at the perigee, , must meet the target height of the orbit, , and thus, .

(4) Target. To acquire the optimal solution and a comprehensive consideration of cost and performance, the total take-off mass is set as the optimization target, while the payload mass is fixed in the optimization process.

3.3. Design Mathematical Model and Optimization Algorithm

According to Section 3.2, Equation (18) shows the mathematical model for the overall design optimization of HPSLV, including the variables, constraints, and target. For multidisciplinary design optimization involved in this paper, the differential evolution (DE) algorithm is selected as the unique optimization algorithm [23]. DE is a global optimization method and has been proven to be effective in solving complex optimization problems. Three parameters including the crossover ratio, mutation factor, and population size are used to control the DE process. New population members are generated by the operations of crossover and mutation from existing members according to certain rules. The operation of selection is finally carried out to save the members whose fitness values are better to the next generation.

4. Results and Discussion

4.1. Influence for Different Grain Shapes

This section analyzes the influence of different grain shapes on the main performance of SLV. All cases are aimed at sending a 100 kg payload into 300 km circle orbit. Table 5 shows six cases, A to F, with different grain shapes (circle-shaped, star-shaped, or wheel-shaped) and the number of chambers (1-chamber or 4-chamber) at the 1st and 2nd stages. The capital letters A~F in the “Case” column indicate different SLV schemes (the superscript/subscript indicates the HRM scheme at the 1st/2nd stage).

Table 5: SLV cases with different grain shapes and numbers of chambers.

Case
1st stage
2nd stage

Grain shape
Number of chambers
Grain shape
Number of chambers

Circle
1
Circle
4

Circle
1
Circle
1

Circle
1
Wheel
1

Circle
4
Circle
4

Wheel
1
Circle
4

Star
1
Circle
4

4.1.1. Optimal Results

As the mathematical model (as shown in Section 3.3) is used for design optimization in each case, the target gradually converges after about 500 iterations. Table 6 gives the optimal results of each SLV design case, including design variables, target, and constraints.

Table 6: Optimal results of SLV cases with different grain shapes.

Parameter
Description

DV1
(m)
1.466
1.477
1.458
1.824
1.687
1.499

DV2

0.994
0.851
0.999
0.718
0.773
0.983

DV3

0.936
0.891
0.580
0.968
0.978
0.943

DV4

0.498
0.402
0.465
0.415
0.428
0.418

DV5

0.989
0.843
0.928
0.943
0.960
0.903

DV6
(m)
0.058
0.062
0.059
0.242
0.055
0.087

DV7
(m)
0.254
0.058
0.310
0.276
0.235
0.253

DV8
(m)
0.283
0.329
0.308
0.291
0.333
0.284

DV9
(kN)
359.836
312.910
346.020
331.422
321.285
378.552

DV10
(kN)
92.776
69.271
92.970
83.086
85.685
102.660

DV11
(kN)
16.448
15.384
15.397
7.613
10.192
11.159

DV12
(MPa)
5.288
6.997
5.552
6.440
4.657
6.764

DV13
(MPa)
5.928
5.326
3.760
6.325
7.000
5.874

DV14
(MPa)
3.772
4.952
2.808
2.086
2.368
2.123

DV15

7.000
5.789
5.376
3.021
2.284
4.558

DV16

3.151
4.277
2.551
2.619
3.116
2.700

DV17

2.829
2.236
2.488
2.442
2.533
2.307

DV18

20.916
17.647
20.861
19.930
14.215
26.059

DV19

85.460
92.558
77.360
68.275
99.873
76.668

DV20

81.554
100.000
99.252
87.288
97.569
100.000

DV21
(°)
2.189
2.166
2.130
1.866
2.142
0.981

DV22

0.516
0.690
0.470
0.733
0.693
0.256

DV23
(s)
83.3
76.0
74.3
0
26.7
11.1

Target
(kg)
16489
16844
16257
15438
15705
15820

Constraint
(kg/(sm²))
98
83
89
265
77
86

Constraint
(kg/(sm²))
105
21
23
128
93
108

Constraint
(kg/(sm²))
59
109
68
28
60
33

Constraint

13.1
15.4
14.9
8.9
10.5
13.2

Constraint
(g)
15.8
15.1
16.0
15.9
15.9
15.9

Constraint
(g)
0.96
0.62
0.95
0.99
0.98
0.98

Constraint
(MPa)
0.066
0.052
0.065
0.054
0.057
0.069

Constraint
(km)
306
300
302
307
306
313

4.1.2. Discussion of the Main Performance Parameters

Figure 11 shows the overall structure of SLV in the six cases. Table 7 lists the take-off mass and overall size of each case. Compared with adopting the circle-shaped grain (in case A), adopting the wheel-shaped (in case E) and star-shaped grain (in case F) at the 1st stage can reduce the take-off mass by 4.75% and 4.06%, respectively, and the wheel-shaped grain has a better effect on reducing the take-off mass. Adopting the wheel-shaped grain can reduce the total length by 7.68% and increase the diameter by 15.08% at the same time. However, adopting the star-shaped grain has little influence on the overall size of the SLV.

Figure 11: Comparison of structures with different grain shapes.

Table 7: Take-off mass and overall dimension of different cases.

Case
Take-off mass (t)
Diameter (m)
Length (m)
Length-diameter ratio

16.489
1.466
19.163
13.071

16.844
1.477
22.800
15.436

16.257
1.458
21.662
14.859

15.438
1.824
16.241
8.906

15.705
1.687
17.692
10.486

15.820
1.499
19.731
13.161

The comparison between cases A, B, and D shows that the take-off mass and total length of the SLV are reduced when the 4-chamber configuration is adopted, but the diameter becomes larger. For instance, compared with those in case B, the take-off mass and total length in case D are reduced by 8.35% and 28.77%, respectively, and the diameter in case D increases by 23.49% correspondingly. This implies that the 4-chamber configuration makes the HRM more compact, which can effectively reduce the structural mass and shorten the overall length of the motor.

To evaluate the impact of different stages on the overall design of SLV, cases A, C, and E are selected for comparison. Compared with that in case A, the wheel-shaped grain at the 2nd stage is adopted in case C, and the take-off mass reduces by 1.41%. Meanwhile, the take-off mass of case E, in which the wheel-shaped grain is adopted at the 1st stage, reduces by 4.75% compared with case A. This means that compared with the situation when the wheel-shaped grain is adopted at the 1st stage, the impact of the 2nd grain is relatively small. Therefore, the mass of the 1st HRM is the main part affecting the overall mass of SLV. In general, the adoption of more complex grain and multicombustion chamber structure will help to achieve the key goal of reducing take-off mass but also increases the complexity of the propulsion system.

Figure 12 gives the curves of thrust and oxidizer-fuel ratio with time of the 1st stage in cases A, E, and F. It shows that in all the cases, the thrust of the HRM declines, while the oxidizer-fuel ratio increases during work. The value of the oxidizer-fuel ratio is related to the shape of grain. As shown in Figure 12(b), the circle-shaped grain has the highest oxidizer-fuel ratio, followed by the star-shaped grain and then the wheel-shaped grain. The main reason is that the wheel-shaped grain can provide larger burning area and can be designed according to the best oxidizer-fuel ratio. Therefore, compared with the circle-shaped (7.338) and star-shaped grain (5.768), the average ratio of the wheel-shaped grain (2.683) is closer to the best ratio of the propellant combination, which indicates that under the same requirement for total impulse and dimension, the mass of HRM with the wheel-shaped grain can be smaller than that in the other two cases, as shown in Table 7. The oxidizer-fuel ratio of circle-shaped grain shows large deviation compared with the best ratio, which will lead to loss of specific impulse.

(a)

(b)

Figure 12: Thrust and oxidizer-fuel ratio with time of the 1st stage in cases A, E, and F: (a) thrust; (b) oxidizer-fuel ratio.

Table 8 lists the flight time at each phase in different cases of SLV, and Figure 13 shows the corresponding height-time curves. It presents that the total flight time is less than 500 s (the shortest is 476.9 s) in cases A, B, and C, and the glide time between the 2nd and 3rd stages exceeds more than 74 s. Correspondingly, the working time of the 3rd HRM is shorter. On the contrary, in cases D, E, and F, more than 500 s (the longest is 546.8 s) is used to get into orbit, and the glide time is shorter or even nonexistent, while the 3rd HRM has smaller thrust and longer working time. Figure 13(b) regards the moment that the glide phase begins as the original time and shows the relation between the increased height of the rocket from the beginning of the glide phase to the injection of payload. It is shown that the increased heights in cases A, B, and C are smaller than those in the other three cases, which causes the adjustments of velocity vector in cases A (83.3 s), B (76.0 s), and C (74.3 s) needing a longer unpowered glide phase than those in other three cases (0 s, 26.7 s, and 11.1 s). On the contrary, larger increased heights are found in cases D, E, and F than in other three cases and allow for longer working time of the 3rd HRM with smaller thrust (314.9 s with 7.613 kN in case D, 261.7 s with 10.192 kN in case E, and 290.4 s with 11.159 kN in case F) than in other three cases (168.5 s with 16.448 kN in case A, 134.6 s with 15.384 kN in case B, and 188.1 s with 15.397 kN in case C), so that velocity vector can be adjusted and result in longer time from the beginning of the glide phase to the injection of payload. In addition, the maximum axial overloads and the maximum normal overloads are the smallest in case B, as shown in Table 6, which means that case B has the strongest ability to protect the payload and structure.

Table 8: Flight sequence of different cases of SLV.

Phase

1st motor work (s)
65.2
81.2
71.0
78.7
78.5
79.1

1st-stage separation (s)
8.0
8.0
8.0
8.0
8.0
8.0

2nd motor work (s)
145.9
182.1
127.5
137.2
127.0
135.2

2nd-stage separation (s)
8.0
8.0
8.0
8.0
8.0
8.0

Glide time (s)
83.3
76.0
74.3
0
26.7
11.1

3rd motor work (s)
168.5
134.6
188.1
314.9
261.7
290.4

Total time (s)
478.9
489.9
476.9
546.8
509.9
531.8

(a)

(b)

Figure 13: Height-time curve of different cases: (a) whole flight; (b) glide phase and 3rd-stage working phase.

4.1.3. Parameter Analysis: Sensitivity Analysis and SOM Analysis

(1) Sensitivity Analysis. To find out the design variables that have the greatest impact on the SLV performance, a sensitivity analysis of the main performance parameters is carried out based on case A. Near the optimal solution of case A, 1000 experiment points are selected by optimal Latin hypercube sampling. Then, a multiple quadratic regression model of the performance parameters is established based on these experiment points. Figure 14 shows the 10 factors that have the greatest impact on the main performance parameters (including take-off mass , maximum mass flux of each stage , maximum of axial overload ,maximum of normal overload , maximum dynamic pressure during the flight , height of the perigee , average specific impulse of each stage , effective propellant mass fraction at each stage , velocity increment of each stage , total velocity increment , and total velocity loss ). In the figure, the horizontal axis represents the contribution rate of the input variables to the response, where the blue colour represents a positive correlation and the red colour represents a negative correlation. The label “dv1” is the linear main effect, “dv6^2” means the quadratic main effect, and “dv1-dv6” corresponds to the interaction effect.

(a)

(b)

(c)

(d)

(e)

(f)

(g)

(h)

(i)

(j)

(k)

(l)

(m)

(n)

(o)

(p)

(q)

(r)

(s)

Figure 14: Pareto graphs of the main performance parameters in case A: (a) ; (b) ; (c) ; (d) ; (e) ; (f) ; (g) ; (h) ; (i) ; (j) ; (k) ; (l) ; (m) ; (n) ; (o) ; (p) ; (q) ; (r) ; (s) .

The result shows that the most affecting factors of are the diameter of the 1st stage (DV1) and the initial thickness of grain at the 1st stage (DV6). has the strongest relationship with the grain thickness at the corresponding stage. depends on 2 parameters most: the diameter of the grain at the 1st stage (DV1) and the initial thickness of grain at the 1st stage (DV6). is most affected by the initial thickness of grain at the 1st stage (DV6) and the control parameter of attack angle during the program-turn phase (DV22). The three factors that have the most significant impact on are the control parameter of attack angle (DV22), the initial thickness of grain at the 1st stage (DV6), and the initial thrust at the 1st stage (DV9). mainly depends on the diameter (DV1) and the initial thickness (DV6) of grain at the 1st stage. is determined by the grain thickness and initial oxidizer-fuel ratio at the corresponding stage. is closely related to the diameter of (DV1) and the initial thickness of the grain at the 1st stage (DV6), and the parameters which have great influence on and cover three size scale parameters (DV3, DV4, and DV5). shows similar tendency with , while and show definite association with the diameter (DV1) and the initial thickness of grain at the 1st stage (DV6).

(2) SOM Analysis. The self-organizing map (SOM) is a kind of visible image based on the artificial neural network and can perform unsupervised learning and clustering of data. The kernel of SOM is a neural network with only one input layer and one output layer. A node in the output layer represents the class that needs to be clustered, which can maintain the original topological structure of the input data and have a better visualization effect at the same time [24]. Figure 15 is the self-organizing map drawn according to the feasible solutions of case A. The points in the figure represent the clusters of feasible solutions which have a similar design state obtained through the self-organizing neural network. The level scale represents the value of the corresponding performance parameter. According to the figure, the take-off mass , 1st stage diameter , total length , and velocity increment show similar distribution laws, indicating that there is a positive correlation between these parameters, which is consistent with the results of sensitivity analysis.

(a)

(b)

(c)

(d)

Figure 15: Self-organizing map of case A: (a) ; (b) ; (c) ; (d) .

4.2. Analysis of Different Payloads and Different Orbit Heights

To analyze the influence of different payloads and different orbit heights on SLV performance, 5 SLV cases are set up for optimization design based on case A. The parameters of each case are shown in Table 9, in which the superscript corresponds to the payload mass and the subscript corresponds to the height of the orbit.

Table 9: Payload mass and orbit height of each case.

Case
Payload mass (kg)
Object height (km)

100
300

100
500

100
700

150
300

200
300

4.2.1. Optimal Results

Using the parameterized design model established in Section 2 and optimization methods described in Section 3, the SLVs in above five cases were designed and optimized. After about 500 iterations, the curve began to gradually converge. Table 10 gives the results of target, constraints, and design variables.

Table 10: Optimal results of SLV cases with different payload mass/orbit heights.

Parameter
Description

DV1
(m)
1.466
1.426
1.510
1.394
1.518

DV2

0.994
0.994
0.981
0.980
1.000

DV3

0.936
0.958
0.979
0.977
0.979

DV4

0.498
0.473
0.392
0.405
0.477

DV5

0.989
0.979
0.927
0.949
0.989

DV6
(m)
0.058
0.082
0.071
0.069
0.061

DV7
(m)
0.254
0.342
0.229
0.257
0.222

DV8
(m)
0.283
0.228
0.319
0.299
0.243

DV9
(kN)
359.836
364.718
465.276
350.354
464.752

DV10
(kN)
92.776
91.446
116.690
100.868
116.972

DV11
(kN)
16.448
10.369
11.016
10.958
19.369

DV12
(MPa)
5.288
5.529
6.994
6.902
6.081

DV13
(MPa)
5.928
7.000
6.288
6.952
6.501

DV14
(MPa)
3.772
5.656
2.037
5.579
2.095

DV15

7.000
4.343
4.890
3.701
4.444

DV16

3.151
3.014
2.923
3.289
3.046

DV17

2.829
2.072
2.351
2.755
2.250

DV18

20.916
22.300
30.020
17.212
32.685

DV19

85.460
62.732
88.940
84.065
93.787

DV20

81.554
100.000
97.231
91.644
99.179

DV21
(°)
2.189
1.682
0.958
2.734
2.168

DV22

0.516
0.563
0.269
0.524
0.437

DV23
(s)
83.285
131.649
131.352
27.845
84.972

Target
(kg)
16489
19675
22149
17256
20890

Constraint
(kg/(sm²))
98
105
115
100
105

Constraint
(kg/(sm²))
105
260
95
125
90

Constraint
(kg/(sm²))
59
25
40
40
45

Constraint

13.1
14.9
16.0
15.8
15.2

Constraint
(g)
15.8
15.7
16.0
15.8
15.2

Constraint
(g)
0.96
0.38
0.49
0.96
0.86

Constraint
(MPa)
0.066
0.048
0.057
0.066
0.066

Constraint
(km)
306
510
728
309
321

4.2.2. Discussion of the Main Performance Parameters

Obviously, Table 10 shows that when the payload mass is fixed, the take-off mass of SLV increases with the increase in the orbit height. Similarly, when the orbit height is fixed, the take-off mass increases with the increase in the payload mass. The performance of SLV is closely related to the velocity increments and distribution ratios of all stages. Table 11 lists the velocity increments and velocity distribution ratios of several typical three-stage SLVs. The influence of different payload mass and orbit heights on the velocity increment and velocity distribution ratio is discussed in the following.

Table 11: Velocity increments and velocity distribution ratios of several typical three-stage SLVs.

Name
Stage
Theoretical velocity increment (m/s)
Actual velocity increment (m/s)
Theoretical velocity distribution ratio
Actual velocity distribution ratio

Ро́кот (liquid) [25]
1st
3987
3139
44.12%
41.43%

2nd
2607
2327
28.85%
30.72%

3rd
2443
2110
27.03%
27.85%

Total
9037
7576
—
—

Pegasus (solid) [26]
Carrier
700
236
7.42%
2.90%

1st
2984
2307
31.64%
28.35%

2nd
2927
2787
31.04%
34.23%

3rd
2818
2810
29.89%
34.52%

Total
9429
8140
—
—

HPSLV of Purdue University (hybrid) [27]
1st
2773
—
26.76%
—

2nd
3734
—
36.03%
—

3rd
3851
—
37.20%
—

Total
10358
—
—
—

Figure 16 shows the velocity increment of each case and the corresponding velocity distribution ratio at each stage, in which the theoretical velocity increments (white) are calculated according to the Tsiolkovsky formula and the actual velocity increments (shaded) are calculated based on the trajectory simulation data. Compared with the SLV cases in Table 11, the velocity distribution ratio is closest to that of Purdue University’s HPSLV (the maximum deviation of the 1st, 2nd, and 3rd stages is smaller than 4.75%, 2.86%, and 5.85%, respectively).

(a) Theoretical and actual velocity increment

(b) Theoretical and actual velocity distribution ratio

Figure 16: Velocity increment and velocity distribution ratio: (a) theoretical and actual velocity increment; (b) theoretical and actual velocity distribution ratio.

A comparison of case , , and shows that with the increase in the orbit height, the SLV’s theoretical/actual velocity increments increase and the velocity loss also increases accordingly. Similarly, a comparison of case , , and shows that the theoretical/actual velocity increment and the proportion of the velocity loss continue to increase with the increase in payload mass.

From the perspective of the velocity distribution ratio at each stage, the theoretical/actual distribution ratio has slight change with the increase in payload mass and orbit height, and the estimated ratio is around 1 : 1.55 : 1.69/1 : 1.9 : 2.39. The actual velocity distribution ratio of the three-stage HPSLV is roughly 16% to 21% for the 1st stage, 32% to 38% for the 2nd stage, and 42% to 49% for the 3rd stage. This result can be used to estimate velocity increments of HPSLV in the future.

In all the cases, the velocity loss at each stage shows the same distribution law. The velocity at the 1st stage accounts for the largest proportion (23% to 35%), the 2nd stage follows (8% to 25%), and the 3rd stage is the smallest (1% to 10%). The loss of the 1st stage mainly comes from the atmospheric drag and gravity when the rocket flights through the dense atmosphere. At the 2nd stage, the SLV passes through the atmosphere during flight and is still greatly affected by atmospheric drag. At the 3rd stage, the SLV mainly works in the vacuum environment, and the loss mainly comes from the gravity when the altitude changes.

5. Conclusions

This paper carries out optimization design of HPSLVs with different grain shapes, chamber configurations, payload mass, and orbit heights. The results show that the star-shaped or wheel-shaped grains are more effective than circle-shaped grains in reducing the take-off mass and dimension of the SLV. The adoption of four combustion chambers makes the SLV design more compact, which has a significant effect on reducing the take-off mass and the total length. The analysis of SA and SOM shows that the optimal target of the take-off mass mainly depends on the design of the 1st stage. Furthermore, the influence of different payload mass and different orbit heights on the performance parameters of the HPSLV is investigated, and the velocity increments, velocity losses, and velocity distribution ratios at each stage are obtained. Compared with that in Purdue’s SLV, the velocity distribution ratio of the HPSLV in this paper shows the same regularity, which can provide some reference for future HPSLV design and development.

Nomenclature

: | Grain port area
: | Nozzle throat area
: | Regression rate coefficient, acceleration
: | Acceleration of Coriolis inertial force
: | Acceleration of centrifugal inertial force
: | Orbital semimajor axis
: | Drag coefficient
: | Thrust coefficient
: | Derivative of lift coefficient to attack angle
: | Characteristic velocity
: | Drag, diameter of the rocket
: | Outer diameter of grain
: | Centre hole diameter of wheel-shape grain
: | Wheel channel inner diameter of wheel-shaped grain
: | Wheel channel outer diameter of wheel-shaped grain
: | Hydraulic diameter
: | Grain thickness
: | Orbital ellipticity
: | Thrust
: | Universal gravitational constant, mass flux
: | Acceleration of gravity
: | Final orbit height at perigee
: | Specific impulse
: | Specific heat ratio
: | lift, length of the rocket
: | Length of the fuel grain
: | Take-off mass
: | The mass of the earth
: | The mass of the rocket
: | Mass flow rate
: | Number of combustion chambers, aerodynamic force
: | Axial overload
: | Normal overload
: | Number of star angle/wheel hole, flux rate exponent
: | Object height
: | Chamber pressure
: | Equilibrium pressure
: | Dynamic pressure
: | Specific gas constant
: | Radius of star tip arc, radius of chamfer of wheel-shaped grain, distance from the turnoff point to the earth’s core
: | Regression rate
: | Apogee altitude
: | Perigee altitude
: | Radius of star root arc
: | Burning perimeter length
: | Reference area
: | Glide time between 2nd and 3rd stages
: | Adiabatic combustion temperature
: | Time
: | Velocity increment
: | Effective volume of chamber
: | Velocity
: | Turnoff velocity of the third stage
: | Position in the -direction
: | Position in the -direction
: | Lateral force
: | Position in the -direction
: | Oxidizer to fuel ratio, attack angle
: | Control parameter of attack angle during the program-turn phase, sideslip angle
: | Star angle coefficient, nozzle expansion ratio
: | Angle of star root, trajectory tilt angle
: | Local trajectory tilt angle at the turnoff point
: | Effective propellant mass fraction
: | Diameter ratio coefficient in design variables
: | Pitch angle
: | Pitch angular velocity
: | Density
: | Velocity loss

Subscript

c: | Combustion gas
e: | Envelope
i: | Initial
o: | Oxidizer
f: | Fuel
max: | Maximum value
: | Current time
: | -direction
: | -direction
: | -direction
0: | Launch point
1: | First stage
2: | Second stage
3: | Third stage

Data Availability

The data used to support the findings of this study are included within the article.

Conflicts of Interest

The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this article.

References

B. Yi, D. Gu, x. Chang, and k. Shao, “Integrating BDS and GPS for precise relative orbit determination of LEO formation flying,” Chinese Journal of Aeronautics, vol. 31, no. 10, pp. 2013–2022, 2018.

B. Genevieve, M. Brooks, J. P. de la Beaujardiere, and L. Roberts, “Performance modeling of a paraffin wax/nitrous oxide hybrid rocket motor,” Tech. Rep., AIAA; Report No.: AIAA 2011-420, Reston, 2011.

M. Mascaro, D. A. Jones, D. M. Lineberry, R. A. Frederick, M. D. Moser, and K. Mahaffy, “Internal ballistics model for a mixed hybrid rocket motor,” Tech. Rep., AIAA; Report No.: AIAA 2015-3811, Reston, 2015.

P. C. Wang, H. Zhu, and G. B. Cai, “Optimization design for a cost-effective hybrid rocket motor,” Journal of Solid Rocket Technology, vol. 40, no. 5, pp. 537–544, 2017.

H. Tian, Y. Duan, and H. Zhu, “Three-dimensional numerical analysis on combustion performance and flow of hybrid rocket motor with multi-segmented grain,” Chinese Journal of Aeronautics, vol. 33, no. 4, pp. 1181–1191, 2020.

R. Taherinezhad and G. Zarepour, “Evaluation of pressure oscillations by a laboratory motor,” Chinese Journal of Aeronautics, vol. 33, no. 3, pp. 805–825, 2020.

A. Yang, B. Li, Y. Yan, S. Xue, and L. Zhou, “Experimental of combustion instability in NTO/MMH impinging combustion chambers,” Chinese Journal of Aeronautics, vol. 33, no. 5, pp. 1476–1485, 2020.

H. Zhu, H. Tian, and G. B. Cai, “Hybrid uncertainty-based design optimization and its application to hybrid rocket motors for manned lunar landing,” Chinese Journal of Aeronautics, vol. 30, no. 2, pp. 719–725, 2017.

L. Casalino and D. Pastrone, “Optimal design of hybrid rocket motors for launchers upper stages,” Journal of Propulsion and Power, vol. 26, no. 3, pp. 421–427, 2010.

J. Tsohas, L. J. Doppers, and S. D. Heister, “Sounding rocket technology demonstration for small satellite launch vehicle project,” Tech. Rep., AIAA; Report No.: AIAA RS4-2006-4004, Reston, 2006.

N. A. Davydenko, R. G. Gollender, A. M. Gubertov, V. V. Mironov, and N. N. Volkov, “Hybrid rocket engines: the benefits and prospects,” Aerospace Science and Technology, vol. 11, no. 1, pp. 55–60, 2007.

F. Costa, R. Contaifer, J. Albuquerque, S. Gabriel, and R. Marques, “Study of paraffin/H₂O₂ hybrid rockets for launching nanosats,” Tech. Rep., AIAA; Report No.: AIAA 2008-4542, Reston, 2008.

I. Rhee, C. Lee, and J. W. Lee, “Optimal design for hybrid rocket engine for air launch vehicle,” Journal of Mechanical Science and Technology, vol. 22, no. 8, pp. 1576–1585, 2008.

P. Wang, H. Tian, H. Zhu, and G. Cai, “Multi-disciplinary design optimization with fuzzy uncertainties and its application in hybrid rocket motor powered launch vehicle,” Chinese Journal of Aeronautics, vol. 33, no. 5, pp. 1454–1467, 2020.

C. X. Zhou, Y. T. Ju, and X. Chen, Rocket Design Theory, Beijing Institute of Technology Press, Beijing, 2014.

L. Xintian, C. Qiang, L. Yancheng, W. Xuekun, and D. Linpeng, “Research on wagon-wheel fuel grain parametric design and internal ballistics performance of hybrid rocket motor,” Journal of Beijing University of Aeronautics and Astronautics, vol. 46, no. 4, pp. 724–730, 2020.

Y. S. Chen and B. Wu, “Development of a small launch vehicle with hybrid rocket propulsion,” Tech. Rep., AIAA; Report No.: AIAA 2018-4835, Reston, 2018.

C. Carmicino and D. Pastrone, “Novel comprehensive technique for hybrid rocket experimental ballistic data reconstruction,” Journal of Propulsion and Power, vol. 34, no. 1, pp. 133–145, 2018.

Z. Hao, S. Xingliang, L. Chengen, T. Hui, and C. Guobiao, “Parameter analysis and transient-feature study of a long-time working hybrid rocket motor,” Journal of Aerospace Engineering, vol. 32, no. 6, article 04019086, 2019.

W. Werthman and C. Schroedea, “A preliminary design code for hybrid propellant rockets,” in 32nd Aerospace Sciences Meeting and Exhibit, Reno, NV, USA, 1994.

L. S. He and D. J. Xu, Solid Ballistic Missile and Launch Vehicle Conceptual Design, Beihang University Press, Beijing, 2017.

A. J. Krueger and R. A. Minzner, “A mid-latitude ozone model for the 1976 U.S. Standard Atmosphere,” Journal of Geophysical Research, vol. 81, no. 24, pp. 4477–4481, 1976.

R. Storn and K. Price, “Differential evolution: a simple and efficient adaptive scheme for global optimization over continuous spaces,” Journal of Global Optimization, vol. 23, no. 1, 1995.

T. Kohonen, “Self-organized formation of topologically correct feature maps,” Biological Cybernetics, vol. 43, no. 1, pp. 59–69, 1982.

Editorial Committee of Encyclopedia of Space Vehicles in the World, Encyclopedia of Space Vehicles in the World, China Aerospace Publishing House, Beijing, 2007.

Y. L. Kang, “Performance analysis of Pegasus,” Missile and Space Launch Technology, vol. 3, no. 3, pp. 27–32, 2002.

J. Tsohas, B. Appel, A. Rettenmaier, M. Walker, and S. Heister, “Development and launch of the Purdue hybrid rocket technology demonstrator,” Tech. Rep., AIAA; Report No.: AIAA 2009-4842, Reston, 2009.