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

Authors: Hao Zhu; Pengcheng Wang; Weile Xu; Yuanjun Zhang; Hui Tian; Guobiao Cai
Journal: International Journal of Aerospace Engineering
 (2021)

---

## 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.


---

## Body

## 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. Figure1 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 N2O/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% H2O2/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 H2O2/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 Section2, 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 Figure2, 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 Figure3 [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 Figure4. After the initial shape of fuel grain is known, the grain port area Ap and the burning perimeter length Sc, 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.Figure 4
Parameters of different grain shapes: (a) circle-shaped grain; (b) star-shaped grain; (c) wheel-shaped grain.
(a)(b)(c)Table 1
Definition of initial grain shape parameters.
Grain shapeBurning perimeter lengthGrain port areaCircle-shaped grainπDp‐2eπDp‐2e2/4Star-shaped grainRelated to the outer diameterDp, number of star angle n, grain thickness e, star angle coefficient ε, angle of star root ϑ, radius of star tip arc r, and radius of star root arc r1 [15]Wheel-shaped grainRelated to the outer diameterDp, number of wheel hole n, grain thickness e, center hole diameter Di, wheel channel inner diameter Dpi, wheel channel outer diameter Dpo, and radius of chamfer r [16]


#### 2.2.2. Propellant Selection and Thermodynamic Calculation

A 98% mass concentration of H2O2 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 Is, adiabatic combustion temperature Tf, specific heat ratio k, and characteristic velocity c∗, are exported by inputting the chamber pressure Pc, 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 Is-α 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
Is-α 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 forceFi, initial oxidizer-fuel ratio αi, and specific impulse Is (provided by thermodynamic calculation), using
(1)ṁoi=ṁiαiαi+1=FiαiIsαi+1,(2)ṁfi=ṁiαi+1=FiIsαi+1,(3)Lp=ṁfiρfSciṙ,(4)At=ṁic∗Pci,where ṁo is the oxidizer mass flow rate, ṁf is the fuel mass flow rate, ṁi is the total propellant mass flow rate, Lp is the length of fuel grain, ρf is the density of fuel grain, Sci is the initial burning perimeter length, ṙ is the regression rate of fuel grain, At is the nozzle throat area, and Pci 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 fluxGo and the hydraulic diameter D0, as shown in
(5)ṙ=aGonD0−0.2=aṁoApn4ApSc−0.2,where a is the regression rate coefficient and n 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 thePc-t curve and F-t 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
(6)ṁo+ṁf=PcAtc∗+Vcdρcdt+ρcdVcdt,where Vc is the effective volume of chamber and ρc is the combustion gas density in the chamber. Considering the HRM combustion equations and introducing the equilibrium pressure Peq, the chamber pressure Pc 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.
(7)Pc=Peq−ṙc∗VcAtRTfdPcde=ṁf+ṁoc∗At−ṙc∗VcAtRTfdPcde,where Peq refers to the equilibrium pressure (corresponding to the zero-dimensional process in a solid rocket motor) [8], R indicates the specific gas constant, and de 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 Is and thrust coefficient CF, 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 F=NCFPcAt, in which N 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.Figure 6
Experimental and predicted curves: (a) chamber pressure; (b) thrust.
(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 Figure2, 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 Table2. H2O2 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 H2O2 into oxygen and water gas with high temperature to enhance combustion efficiency.Table 2
Main structure of the SLV propulsion system.
Main structureMaterialFunctionLiquid oxidizer98% H2O2Provide chemical energySolid fuelHTPB-based grainProvide chemical energyGas bottleCF wound aluminum linerContain pressurization gasOxidizer tankCF wound aluminum liner and oxidizer sacContain liquid oxidizerChamber and nozzleCF wound shell and high-silica insulationContain fuel grain and high-pressure burn gasValves and tubesMainly aluminum alloyFeed liquid propellantCatalyst bedStainless-steel shell and silver meshDecompose H2O2TurbopumpDetails refer to Reference [20]Increase the pressure of oxidizerThrust vector control actuatorGiven by experienceControl the thrust vector directionThe 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.Figure 7
Different types of HRM in the 1st/2nd stage of SLV: (a) 1-chamber HRM; (b) 4-chamber HRM.
(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 Figure9. 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:
(8)α=0,φ=π2,(9)αt=4αmaxeβt1−teβt1−t−1,φ=ϑ+α,(10)α=0,φ=ϑ,(11)ωφ=dφtdt=const,where αmax is the maximum absolute value of attack angle during the program-turn phase, β is a constant coefficient, ωφ is the pitch angular velocity, φt is the pitch angle at current time, and t 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(12)φ=ϑ.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(13)ωφ=dφtdt=0.


#### 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 x-direction as an example, xt is the position coordinate of the SLV, where x0 is the initial position of the launch point. vx is the velocity component along the x-direction in the launch coordinate, and the initial velocity of SLV vx0=0. ax is the acceleration along the x-direction in the launch coordinate. m is the rocket mass at current time. Fx is the thrust along the x-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. Nx is the aerodynamic force along the x-direction in the launch coordinate and can be obtained according to the aerodynamic force in the body coordinate. gx, aex, and acx are the accelerations in the launch coordinate caused by gravity, centrifugal inertial force, and Coriolis inertial force, respectively. The trajectory parameters in the y-direction and z-direction are calculated using the same method as that in the x-direction. Then, Equations (14) and (15) are solved with the fourth-order Runge-Kutta method to obtain the trajectory parameters.
(14)axayaz=1mFxFyFz+NxNyNz+gx−aex−acxgy−aey−acygz−aez−acz,(15)xtytzt=∫0tvxdt∫0tvydt∫0tvzdt+x0y0z0.


#### 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 CD and the derivative of lift coefficient to attack angle CLα 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 CD and CLα are shown in Table 3, and the 1976 U.S. Standard Atmosphere data are coded as the SLV working environment [22].
(16)DLZ=−CDqSM−αCLαqSM−βCLαqSM.Table 3
Approximation aerodynamic coefficients of the “Titan II” rocket.
Aerodynamic coefficientsVelocity rangeFormulaCD0≤Ma≤0.80.290.8≤Ma≤1.068Ma−0.51Ma>1.0680.091+0.5Ma−1CLα0≤Ma≤0.252.80.25≤Ma≤1.12.8+0.447Ma−0.251.1≤Ma≤1.62.8+0.660Ma−1.11.6≤Ma≤3.62.8+0.350Ma−1.6Ma>3.63.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(17)aorbit=r2−vk2r/GMearth,eorbit=1−vk2rGMearth2−vk2rGMearthcos2ϑk,rporbit=aorbit1−eorbit,raorbit=aorbit1+eorbit,where vk is the turnoff velocity of the 3rd stage, ϑk is the local trajectory tilt angle at the turnoff point, G is the universal gravitational constant, Mearth is the earth mass, r is the distance from the turnoff point to the earth’s core, aorbit is the orbital semimajor axis, eorbit is the orbital ellipticity, rporbit is the perigee altitude, and raorbit is the apogee altitude.


## 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 Figure2, 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 Figure3 [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 Figure4. After the initial shape of fuel grain is known, the grain port area Ap and the burning perimeter length Sc, 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.Figure 4
Parameters of different grain shapes: (a) circle-shaped grain; (b) star-shaped grain; (c) wheel-shaped grain.
(a)(b)(c)Table 1
Definition of initial grain shape parameters.
Grain shapeBurning perimeter lengthGrain port areaCircle-shaped grainπDp‐2eπDp‐2e2/4Star-shaped grainRelated to the outer diameterDp, number of star angle n, grain thickness e, star angle coefficient ε, angle of star root ϑ, radius of star tip arc r, and radius of star root arc r1 [15]Wheel-shaped grainRelated to the outer diameterDp, number of wheel hole n, grain thickness e, center hole diameter Di, wheel channel inner diameter Dpi, wheel channel outer diameter Dpo, and radius of chamfer r [16]


### 2.2.2. Propellant Selection and Thermodynamic Calculation

A 98% mass concentration of H2O2 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 Is, adiabatic combustion temperature Tf, specific heat ratio k, and characteristic velocity c∗, are exported by inputting the chamber pressure Pc, 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 Is-α 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
Is-α 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 forceFi, initial oxidizer-fuel ratio αi, and specific impulse Is (provided by thermodynamic calculation), using
(1)ṁoi=ṁiαiαi+1=FiαiIsαi+1,(2)ṁfi=ṁiαi+1=FiIsαi+1,(3)Lp=ṁfiρfSciṙ,(4)At=ṁic∗Pci,where ṁo is the oxidizer mass flow rate, ṁf is the fuel mass flow rate, ṁi is the total propellant mass flow rate, Lp is the length of fuel grain, ρf is the density of fuel grain, Sci is the initial burning perimeter length, ṙ is the regression rate of fuel grain, At is the nozzle throat area, and Pci 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 fluxGo and the hydraulic diameter D0, as shown in
(5)ṙ=aGonD0−0.2=aṁoApn4ApSc−0.2,where a is the regression rate coefficient and n 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 thePc-t curve and F-t 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
(6)ṁo+ṁf=PcAtc∗+Vcdρcdt+ρcdVcdt,where Vc is the effective volume of chamber and ρc is the combustion gas density in the chamber. Considering the HRM combustion equations and introducing the equilibrium pressure Peq, the chamber pressure Pc 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.
(7)Pc=Peq−ṙc∗VcAtRTfdPcde=ṁf+ṁoc∗At−ṙc∗VcAtRTfdPcde,where Peq refers to the equilibrium pressure (corresponding to the zero-dimensional process in a solid rocket motor) [8], R indicates the specific gas constant, and de 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 Is and thrust coefficient CF, 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 F=NCFPcAt, in which N 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.Figure 6
Experimental and predicted curves: (a) chamber pressure; (b) thrust.
(a)
Chamber pressure(b)
Thrust


## 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 Figure4. After the initial shape of fuel grain is known, the grain port area Ap and the burning perimeter length Sc, 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.Figure 4
Parameters of different grain shapes: (a) circle-shaped grain; (b) star-shaped grain; (c) wheel-shaped grain.
(a)(b)(c)Table 1
Definition of initial grain shape parameters.
Grain shapeBurning perimeter lengthGrain port areaCircle-shaped grainπDp‐2eπDp‐2e2/4Star-shaped grainRelated to the outer diameterDp, number of star angle n, grain thickness e, star angle coefficient ε, angle of star root ϑ, radius of star tip arc r, and radius of star root arc r1 [15]Wheel-shaped grainRelated to the outer diameterDp, number of wheel hole n, grain thickness e, center hole diameter Di, wheel channel inner diameter Dpi, wheel channel outer diameter Dpo, and radius of chamfer r [16]


## 2.2.2. Propellant Selection and Thermodynamic Calculation

A 98% mass concentration of H2O2 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 Is, adiabatic combustion temperature Tf, specific heat ratio k, and characteristic velocity c∗, are exported by inputting the chamber pressure Pc, 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 Is-α 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
Is-α 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 forceFi, initial oxidizer-fuel ratio αi, and specific impulse Is (provided by thermodynamic calculation), using
(1)ṁoi=ṁiαiαi+1=FiαiIsαi+1,(2)ṁfi=ṁiαi+1=FiIsαi+1,(3)Lp=ṁfiρfSciṙ,(4)At=ṁic∗Pci,where ṁo is the oxidizer mass flow rate, ṁf is the fuel mass flow rate, ṁi is the total propellant mass flow rate, Lp is the length of fuel grain, ρf is the density of fuel grain, Sci is the initial burning perimeter length, ṙ is the regression rate of fuel grain, At is the nozzle throat area, and Pci 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 fluxGo and the hydraulic diameter D0, as shown in
(5)ṙ=aGonD0−0.2=aṁoApn4ApSc−0.2,where a is the regression rate coefficient and n 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 thePc-t curve and F-t 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
(6)ṁo+ṁf=PcAtc∗+Vcdρcdt+ρcdVcdt,where Vc is the effective volume of chamber and ρc is the combustion gas density in the chamber. Considering the HRM combustion equations and introducing the equilibrium pressure Peq, the chamber pressure Pc 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.
(7)Pc=Peq−ṙc∗VcAtRTfdPcde=ṁf+ṁoc∗At−ṙc∗VcAtRTfdPcde,where Peq refers to the equilibrium pressure (corresponding to the zero-dimensional process in a solid rocket motor) [8], R indicates the specific gas constant, and de 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 Is and thrust coefficient CF, 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 F=NCFPcAt, in which N 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.Figure 6
Experimental and predicted curves: (a) chamber pressure; (b) thrust.
(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 Figure2, 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 Table2. H2O2 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 H2O2 into oxygen and water gas with high temperature to enhance combustion efficiency.Table 2
Main structure of the SLV propulsion system.
Main structureMaterialFunctionLiquid oxidizer98% H2O2Provide chemical energySolid fuelHTPB-based grainProvide chemical energyGas bottleCF wound aluminum linerContain pressurization gasOxidizer tankCF wound aluminum liner and oxidizer sacContain liquid oxidizerChamber and nozzleCF wound shell and high-silica insulationContain fuel grain and high-pressure burn gasValves and tubesMainly aluminum alloyFeed liquid propellantCatalyst bedStainless-steel shell and silver meshDecompose H2O2TurbopumpDetails refer to Reference [20]Increase the pressure of oxidizerThrust vector control actuatorGiven by experienceControl the thrust vector directionThe 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.Figure 7
Different types of HRM in the 1st/2nd stage of SLV: (a) 1-chamber HRM; (b) 4-chamber HRM.
(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.3.1. Overall Structural Description

The rocket take-off mass is a sum of the mass of all components. As shown in Figure2, 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 Table2. H2O2 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 H2O2 into oxygen and water gas with high temperature to enhance combustion efficiency.Table 2
Main structure of the SLV propulsion system.
Main structureMaterialFunctionLiquid oxidizer98% H2O2Provide chemical energySolid fuelHTPB-based grainProvide chemical energyGas bottleCF wound aluminum linerContain pressurization gasOxidizer tankCF wound aluminum liner and oxidizer sacContain liquid oxidizerChamber and nozzleCF wound shell and high-silica insulationContain fuel grain and high-pressure burn gasValves and tubesMainly aluminum alloyFeed liquid propellantCatalyst bedStainless-steel shell and silver meshDecompose H2O2TurbopumpDetails refer to Reference [20]Increase the pressure of oxidizerThrust vector control actuatorGiven by experienceControl the thrust vector directionThe 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.Figure 7
Different types of HRM in the 1st/2nd stage of SLV: (a) 1-chamber HRM; (b) 4-chamber HRM.
(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 Figure9. 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:
(8)α=0,φ=π2,(9)αt=4αmaxeβt1−teβt1−t−1,φ=ϑ+α,(10)α=0,φ=ϑ,(11)ωφ=dφtdt=const,where αmax is the maximum absolute value of attack angle during the program-turn phase, β is a constant coefficient, ωφ is the pitch angular velocity, φt is the pitch angle at current time, and t 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(12)φ=ϑ.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(13)ωφ=dφtdt=0.


### 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 x-direction as an example, xt is the position coordinate of the SLV, where x0 is the initial position of the launch point. vx is the velocity component along the x-direction in the launch coordinate, and the initial velocity of SLV vx0=0. ax is the acceleration along the x-direction in the launch coordinate. m is the rocket mass at current time. Fx is the thrust along the x-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. Nx is the aerodynamic force along the x-direction in the launch coordinate and can be obtained according to the aerodynamic force in the body coordinate. gx, aex, and acx are the accelerations in the launch coordinate caused by gravity, centrifugal inertial force, and Coriolis inertial force, respectively. The trajectory parameters in the y-direction and z-direction are calculated using the same method as that in the x-direction. Then, Equations (14) and (15) are solved with the fourth-order Runge-Kutta method to obtain the trajectory parameters.
(14)axayaz=1mFxFyFz+NxNyNz+gx−aex−acxgy−aey−acygz−aez−acz,(15)xtytzt=∫0tvxdt∫0tvydt∫0tvzdt+x0y0z0.


### 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 CD and the derivative of lift coefficient to attack angle CLα 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 CD and CLα are shown in Table 3, and the 1976 U.S. Standard Atmosphere data are coded as the SLV working environment [22].
(16)DLZ=−CDqSM−αCLαqSM−βCLαqSM.Table 3
Approximation aerodynamic coefficients of the “Titan II” rocket.
Aerodynamic coefficientsVelocity rangeFormulaCD0≤Ma≤0.80.290.8≤Ma≤1.068Ma−0.51Ma>1.0680.091+0.5Ma−1CLα0≤Ma≤0.252.80.25≤Ma≤1.12.8+0.447Ma−0.251.1≤Ma≤1.62.8+0.660Ma−1.11.6≤Ma≤3.62.8+0.350Ma−1.6Ma>3.63.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(17)aorbit=r2−vk2r/GMearth,eorbit=1−vk2rGMearth2−vk2rGMearthcos2ϑk,rporbit=aorbit1−eorbit,raorbit=aorbit1+eorbit,where vk is the turnoff velocity of the 3rd stage, ϑk is the local trajectory tilt angle at the turnoff point, G is the universal gravitational constant, Mearth is the earth mass, r is the distance from the turnoff point to the earth’s core, aorbit is the orbital semimajor axis, eorbit is the orbital ellipticity, rporbit is the perigee altitude, and raorbit is the apogee altitude.


## 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:
(8)α=0,φ=π2,(9)αt=4αmaxeβt1−teβt1−t−1,φ=ϑ+α,(10)α=0,φ=ϑ,(11)ωφ=dφtdt=const,where αmax is the maximum absolute value of attack angle during the program-turn phase, β is a constant coefficient, ωφ is the pitch angular velocity, φt is the pitch angle at current time, and t 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(12)φ=ϑ.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(13)ωφ=dφtdt=0.


## 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 x-direction as an example, xt is the position coordinate of the SLV, where x0 is the initial position of the launch point. vx is the velocity component along the x-direction in the launch coordinate, and the initial velocity of SLV vx0=0. ax is the acceleration along the x-direction in the launch coordinate. m is the rocket mass at current time. Fx is the thrust along the x-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. Nx is the aerodynamic force along the x-direction in the launch coordinate and can be obtained according to the aerodynamic force in the body coordinate. gx, aex, and acx are the accelerations in the launch coordinate caused by gravity, centrifugal inertial force, and Coriolis inertial force, respectively. The trajectory parameters in the y-direction and z-direction are calculated using the same method as that in the x-direction. Then, Equations (14) and (15) are solved with the fourth-order Runge-Kutta method to obtain the trajectory parameters.
(14)axayaz=1mFxFyFz+NxNyNz+gx−aex−acxgy−aey−acygz−aez−acz,(15)xtytzt=∫0tvxdt∫0tvydt∫0tvzdt+x0y0z0.


## 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 CD and the derivative of lift coefficient to attack angle CLα 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 CD and CLα are shown in Table 3, and the 1976 U.S. Standard Atmosphere data are coded as the SLV working environment [22].
(16)DLZ=−CDqSM−αCLαqSM−βCLαqSM.Table 3
Approximation aerodynamic coefficients of the “Titan II” rocket.
Aerodynamic coefficientsVelocity rangeFormulaCD0≤Ma≤0.80.290.8≤Ma≤1.068Ma−0.51Ma>1.0680.091+0.5Ma−1CLα0≤Ma≤0.252.80.25≤Ma≤1.12.8+0.447Ma−0.251.1≤Ma≤1.62.8+0.660Ma−1.11.6≤Ma≤3.62.8+0.350Ma−1.6Ma>3.63.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(17)aorbit=r2−vk2r/GMearth,eorbit=1−vk2rGMearth2−vk2rGMearthcos2ϑk,rporbit=aorbit1−eorbit,raorbit=aorbit1+eorbit,where vk is the turnoff velocity of the 3rd stage, ϑk is the local trajectory tilt angle at the turnoff point, G is the universal gravitational constant, Mearth is the earth mass, r is the distance from the turnoff point to the earth’s core, aorbit is the orbital semimajor axis, eorbit is the orbital ellipticity, rporbit is the perigee altitude, and raorbit 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 Section2, 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

Table4 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.
VariableDescriptionUnitLower boundUpper boundDV1Diameter of the 1st stageD1m12.5DV2Diameter ratio of the 2nd stage to1st stageξ1—0.71DV3Ratio of the envelope diameter of the 2nd HRM to the diameter of the 2nd stageξ2—0.51DV4Diameter ratio of the 3rd tank to the 3rd stageξ3—0.30.7DV5Ratio of the envelope diameter of the 3rd HRM to the diameter of the 3rd stageξ4—0.31DV6Initial thickness of the 1st graine1m0.050.4DV7Initial thickness of the 2nd graine2m0.050.4DV8Initial thickness of the 3rd graine3m0.050.4DV9Initial thrust at the 1st stageF1kN401250DV10Initial thrust at the 2nd stageF2kN2650DV11Initial thrust at the 3rd stageF3kN0.2585DV12Initial pressure at the 1st chamberPc1MPa27DV13Initial pressure at the 2nd chamberPc2MPa27DV14Initial pressure at the 3rd chamberPc3MPa27DV15Initial oxidizer-fuel ratio at the 1st stageα1—27DV16Initial oxidizer-fuel ratio at the 2nd stageα2—27DV17Initial oxidizer-fuel ratio at the 3rd stageα3—27DV18Expansion ratio at the 1st stageε1—1050DV19Expansion ratio at the 2nd stageε2—10100DV20Expansion ratio at the 3rd stageε3—10100DV21Maximum attack angle during the program-turn phaseαmax°0.13DV22Control parameter of attack angle during the program-turn phaseβ—0.151DV23Glide time between 2nd and 3rd stagesTGlides0300As shown in Figure10, D2e refers to the envelope diameter of the 2nd HRM, while D3e is that of the 3rd HRM. D2 and D3 correspond to the diameter of the 2nd and 3rd stages, respectively. The design variable ξ2ξ4 in Table 4 is the ratio of D2eD3e to D2D3. According to the geometric relationship, the conversion between the diameter of the components and the diameter of the rocket can be calculated.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.
(a)(b)


#### 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 Gomaxi is set to 500 kg/(s∗m2).(2) Structure Constraints. Considering that a too high L/D is not good for the rocket structure strength, L/D is set as smaller than 16. For protecting payload, the maximum axial overload Nxmax is set below 16 g, while the maximum normal overload Nymax is set below 1 g.(3) Trajectory Constraints. The maximum dynamic pressure during the flight qmax needs to be smaller than 0.1 MPa. The final orbit height at the perigee, Hp, must meet the target height of the orbit, ObjectH, and thus, Hp/ObjectH>1.(4) Target. To acquire the optimal solution and a comprehensive consideration of cost and performance, the total take-off mass M0 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 Section3.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.
(18)findxminM0=fxs.t.g1x=Gomax1500−1<0g2x=Gomax2500−1<0g3x=Gomax3500−1<0g4x=L16D−1<0g5x=Nxmax16−1<0g6x=Nymax−1<0g7x=qmax0.1−1<0g8x=1−HpObjectH<0xL≤x≤xU.


## 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 Section2, 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

Table4 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.
VariableDescriptionUnitLower boundUpper boundDV1Diameter of the 1st stageD1m12.5DV2Diameter ratio of the 2nd stage to1st stageξ1—0.71DV3Ratio of the envelope diameter of the 2nd HRM to the diameter of the 2nd stageξ2—0.51DV4Diameter ratio of the 3rd tank to the 3rd stageξ3—0.30.7DV5Ratio of the envelope diameter of the 3rd HRM to the diameter of the 3rd stageξ4—0.31DV6Initial thickness of the 1st graine1m0.050.4DV7Initial thickness of the 2nd graine2m0.050.4DV8Initial thickness of the 3rd graine3m0.050.4DV9Initial thrust at the 1st stageF1kN401250DV10Initial thrust at the 2nd stageF2kN2650DV11Initial thrust at the 3rd stageF3kN0.2585DV12Initial pressure at the 1st chamberPc1MPa27DV13Initial pressure at the 2nd chamberPc2MPa27DV14Initial pressure at the 3rd chamberPc3MPa27DV15Initial oxidizer-fuel ratio at the 1st stageα1—27DV16Initial oxidizer-fuel ratio at the 2nd stageα2—27DV17Initial oxidizer-fuel ratio at the 3rd stageα3—27DV18Expansion ratio at the 1st stageε1—1050DV19Expansion ratio at the 2nd stageε2—10100DV20Expansion ratio at the 3rd stageε3—10100DV21Maximum attack angle during the program-turn phaseαmax°0.13DV22Control parameter of attack angle during the program-turn phaseβ—0.151DV23Glide time between 2nd and 3rd stagesTGlides0300As shown in Figure10, D2e refers to the envelope diameter of the 2nd HRM, while D3e is that of the 3rd HRM. D2 and D3 correspond to the diameter of the 2nd and 3rd stages, respectively. The design variable ξ2ξ4 in Table 4 is the ratio of D2eD3e to D2D3. According to the geometric relationship, the conversion between the diameter of the components and the diameter of the rocket can be calculated.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.
(a)(b)


### 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 Gomaxi is set to 500 kg/(s∗m2).(2) Structure Constraints. Considering that a too high L/D is not good for the rocket structure strength, L/D is set as smaller than 16. For protecting payload, the maximum axial overload Nxmax is set below 16 g, while the maximum normal overload Nymax is set below 1 g.(3) Trajectory Constraints. The maximum dynamic pressure during the flight qmax needs to be smaller than 0.1 MPa. The final orbit height at the perigee, Hp, must meet the target height of the orbit, ObjectH, and thus, Hp/ObjectH>1.(4) Target. To acquire the optimal solution and a comprehensive consideration of cost and performance, the total take-off mass M0 is set as the optimization target, while the payload mass is fixed in the optimization process.


## 3.2.1. Variables

Table4 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.
VariableDescriptionUnitLower boundUpper boundDV1Diameter of the 1st stageD1m12.5DV2Diameter ratio of the 2nd stage to1st stageξ1—0.71DV3Ratio of the envelope diameter of the 2nd HRM to the diameter of the 2nd stageξ2—0.51DV4Diameter ratio of the 3rd tank to the 3rd stageξ3—0.30.7DV5Ratio of the envelope diameter of the 3rd HRM to the diameter of the 3rd stageξ4—0.31DV6Initial thickness of the 1st graine1m0.050.4DV7Initial thickness of the 2nd graine2m0.050.4DV8Initial thickness of the 3rd graine3m0.050.4DV9Initial thrust at the 1st stageF1kN401250DV10Initial thrust at the 2nd stageF2kN2650DV11Initial thrust at the 3rd stageF3kN0.2585DV12Initial pressure at the 1st chamberPc1MPa27DV13Initial pressure at the 2nd chamberPc2MPa27DV14Initial pressure at the 3rd chamberPc3MPa27DV15Initial oxidizer-fuel ratio at the 1st stageα1—27DV16Initial oxidizer-fuel ratio at the 2nd stageα2—27DV17Initial oxidizer-fuel ratio at the 3rd stageα3—27DV18Expansion ratio at the 1st stageε1—1050DV19Expansion ratio at the 2nd stageε2—10100DV20Expansion ratio at the 3rd stageε3—10100DV21Maximum attack angle during the program-turn phaseαmax°0.13DV22Control parameter of attack angle during the program-turn phaseβ—0.151DV23Glide time between 2nd and 3rd stagesTGlides0300As shown in Figure10, D2e refers to the envelope diameter of the 2nd HRM, while D3e is that of the 3rd HRM. D2 and D3 correspond to the diameter of the 2nd and 3rd stages, respectively. The design variable ξ2ξ4 in Table 4 is the ratio of D2eD3e to D2D3. According to the geometric relationship, the conversion between the diameter of the components and the diameter of the rocket can be calculated.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.
(a)(b)


## 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 Gomaxi is set to 500 kg/(s∗m2).(2) Structure Constraints. Considering that a too high L/D is not good for the rocket structure strength, L/D is set as smaller than 16. For protecting payload, the maximum axial overload Nxmax is set below 16 g, while the maximum normal overload Nymax is set below 1 g.(3) Trajectory Constraints. The maximum dynamic pressure during the flight qmax needs to be smaller than 0.1 MPa. The final orbit height at the perigee, Hp, must meet the target height of the orbit, ObjectH, and thus, Hp/ObjectH>1.(4) Target. To acquire the optimal solution and a comprehensive consideration of cost and performance, the total take-off mass M0 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 Section3.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.
(18)findxminM0=fxs.t.g1x=Gomax1500−1<0g2x=Gomax2500−1<0g3x=Gomax3500−1<0g4x=L16D−1<0g5x=Nxmax16−1<0g6x=Nymax−1<0g7x=qmax0.1−1<0g8x=1−HpObjectH<0xL≤x≤xU.


## 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. Table5 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.
Case1st stage2nd stageGrain shapeNumber of chambersGrain shapeNumber of chambersA4CircleCircleCircle1Circle4BCircleCircleCircle1Circle1CWheelCircleCircle1Wheel1D4Circle4CircleCircle4Circle4E4CircleWheelWheel1Circle4F4CircleStarStar1Circle4


#### 4.1.1. Optimal Results

As the mathematical model (as shown in Section3.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.
ParameterDescriptionA4CircleCircleBCircleCircleCWheelCircleD4Circle4CircleE4CircleWheelF4CircleStarDV1D1 (m)1.4661.4771.4581.8241.6871.499DV2ξ10.9940.8510.9990.7180.7730.983DV3ξ20.9360.8910.5800.9680.9780.943DV4ξ30.4980.4020.4650.4150.4280.418DV5ξ40.9890.8430.9280.9430.9600.903DV6e1 (m)0.0580.0620.0590.2420.0550.087DV7e2 (m)0.2540.0580.3100.2760.2350.253DV8e3 (m)0.2830.3290.3080.2910.3330.284DV9F1 (kN)359.836312.910346.020331.422321.285378.552DV10F2 (kN)92.77669.27192.97083.08685.685102.660DV11F3 (kN)16.44815.38415.3977.61310.19211.159DV12Pc1 (MPa)5.2886.9975.5526.4404.6576.764DV13Pc2 (MPa)5.9285.3263.7606.3257.0005.874DV14Pc3 (MPa)3.7724.9522.8082.0862.3682.123DV15α17.0005.7895.3763.0212.2844.558DV16α23.1514.2772.5512.6193.1162.700DV17α32.8292.2362.4882.4422.5332.307DV18ε120.91617.64720.86119.93014.21526.059DV19ε285.46092.55877.36068.27599.87376.668DV20ε381.554100.00099.25287.28897.569100.000DV21αmax (°)2.1892.1662.1301.8662.1420.981DV22β0.5160.6900.4700.7330.6930.256DV23TGlide (s)83.376.074.3026.711.1TargetM0 (kg)164891684416257154381570515820ConstraintGomax1 (kg/(s∗m2))9883892657786ConstraintGomax2 (kg/(s∗m2))105212312893108ConstraintGomax3 (kg/(s∗m2))5910968286033ConstraintL/D13.115.414.98.910.513.2ConstraintNxmax (g)15.815.116.015.915.915.9ConstraintNymax (g)0.960.620.950.990.980.98Constraintqmax (MPa)0.0660.0520.0650.0540.0570.069ConstraintHp (km)306300302307306313


#### 4.1.2. Discussion of the Main Performance Parameters

Figure11 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.
CaseTake-off mass (t)Diameter (m)Length (m)Length-diameter ratioA4CircleCircle16.4891.46619.16313.071BCircleCircle16.8441.47722.80015.436CWheelCircle16.2571.45821.66214.859D4Circle4Circle15.4381.82416.2418.906E4CircleWheel15.7051.68717.69210.486F4CircleStar15.8201.49919.73113.161The 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.Figure12 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.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.
(a)(b)Table8 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 Nxmax and the maximum normal overloads Nymax 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.
PhaseA4CircleCircleBCircleCircleCWheelCircleD4Circle4CircleE4CircleWheelF4CircleStar1st motor work (s)65.281.271.078.778.579.11st-stage separation (s)8.08.08.08.08.08.02nd motor work (s)145.9182.1127.5137.2127.0135.22nd-stage separation (s)8.08.08.08.08.08.0Glide time (s)83.376.074.3026.711.13rd motor work (s)168.5134.6188.1314.9261.7290.4Total time (s)478.9489.9476.9546.8509.9531.8Figure 13
Height-time curve of different cases: (a) whole flight; (b) glide phase and 3rd-stage working phase.
(a)(b)


#### 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 M0, maximum mass flux of each stage Gomaxi, maximum of axial overload Nxmax,maximum of normal overload Nymax, maximum dynamic pressure during the flight qmax, height of the perigee Hp, average specific impulse of each stage Isi, effective propellant mass fraction at each stage μi, velocity increment of each stage Vi, total velocity increment V, and total velocity loss ΔV). 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.Figure 14
Pareto graphs of the main performance parameters in case A: (a)M0; (b) Gomax1; (c) Gomax2; (d) Gomax3; (e) Nxmax; (f) Nymax; (g) qmax; (h) Hp; (i) Is1; (j) Is2; (k) Is3; (l) μ1; (m) μ2; (n) μ3; (o) V1; (p) V2; (q) V3; (r) V; (s) ΔV.
(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 ofM0 are the diameter of the 1st stage (DV1) and the initial thickness of grain at the 1st stage (DV6). Gomaxi has the strongest relationship with the grain thickness at the corresponding stage. Nxmax 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). Nymax 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 qmax 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). Hp mainly depends on the diameter (DV1) and the initial thickness (DV6) of grain at the 1st stage. Isi is determined by the grain thickness and initial oxidizer-fuel ratio at the corresponding stage. μ1 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 μ2 and μ3 cover three size scale parameters (DV3, DV4, and DV5). Vi shows similar tendency with μi, while V and ΔV 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 M0, 1st stage diameter D1, total length L, and velocity increment V show similar distribution laws, indicating that there is a positive correlation between these parameters, which is consistent with the results of sensitivity analysis.Figure 15
Self-organizing map of case A: (a)M0; (b) L; (c) D1; (d) V.
(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 Table9, 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.
CasePayload mass (kg)Object height (km)A300100100300A500100100500A700100100700A300150150300A300200200300


#### 4.2.1. Optimal Results

Using the parameterized design model established in Section2 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.
ParameterDescriptionA300100A500100A700100A300150A300200DV1D1 (m)1.4661.4261.5101.3941.518DV2ξ10.9940.9940.9810.9801.000DV3ξ20.9360.9580.9790.9770.979DV4ξ30.4980.4730.3920.4050.477DV5ξ40.9890.9790.9270.9490.989DV6e1 (m)0.0580.0820.0710.0690.061DV7e2 (m)0.2540.3420.2290.2570.222DV8e3 (m)0.2830.2280.3190.2990.243DV9F1 (kN)359.836364.718465.276350.354464.752DV10F2 (kN)92.77691.446116.690100.868116.972DV11F3 (kN)16.44810.36911.01610.95819.369DV12Pc1 (MPa)5.2885.5296.9946.9026.081DV13Pc2 (MPa)5.9287.0006.2886.9526.501DV14Pc3 (MPa)3.7725.6562.0375.5792.095DV15α17.0004.3434.8903.7014.444DV16α23.1513.0142.9233.2893.046DV17α32.8292.0722.3512.7552.250DV18ε120.91622.30030.02017.21232.685DV19ε285.46062.73288.94084.06593.787DV20ε381.554100.00097.23191.64499.179DV21αmax (°)2.1891.6820.9582.7342.168DV22β0.5160.5630.2690.5240.437DV23TGlide (s)83.285131.649131.35227.84584.972TargetM0 (kg)1648919675221491725620890ConstraintGomax1 (kg/(s∗m2))98105115100105ConstraintGomax2 (kg/(s∗m2))1052609512590ConstraintGomax3 (kg/(s∗m2))5925404045ConstraintL/D13.114.916.015.815.2ConstraintNxmax (g)15.815.716.015.815.2ConstraintNymax (g)0.960.380.490.960.86Constraintqmax (MPa)0.0660.0480.0570.0660.066ConstraintHp (km)306510728309321


#### 4.2.2. Discussion of the Main Performance Parameters

Obviously, Table10 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.
NameStageTheoretical velocity increment (m/s)Actual velocity increment (m/s)Theoretical velocity distribution ratioActual velocity distribution ratioРо́кот (liquid) [25]1st3987313944.12%41.43%2nd2607232728.85%30.72%3rd2443211027.03%27.85%Total90377576——Pegasus (solid) [26]Carrier7002367.42%2.90%1st2984230731.64%28.35%2nd2927278731.04%34.23%3rd2818281029.89%34.52%Total94298140——HPSLV of Purdue University (hybrid) [27]1st2773—26.76%—2nd3734—36.03%—3rd3851—37.20%—Total10358———Figure16 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).Figure 16
Velocity increment and velocity distribution ratio: (a) theoretical and actual velocity increment; (b) theoretical and actual velocity distribution ratio.
(a)
Theoretical and actual velocity increment(b)
Theoretical and actual velocity distribution ratioA comparison of caseA300100, A500100, and A700100 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 A300100, A300150, and A300200 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.


## 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. Table5 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.
Case1st stage2nd stageGrain shapeNumber of chambersGrain shapeNumber of chambersA4CircleCircleCircle1Circle4BCircleCircleCircle1Circle1CWheelCircleCircle1Wheel1D4Circle4CircleCircle4Circle4E4CircleWheelWheel1Circle4F4CircleStarStar1Circle4


### 4.1.1. Optimal Results

As the mathematical model (as shown in Section3.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.
ParameterDescriptionA4CircleCircleBCircleCircleCWheelCircleD4Circle4CircleE4CircleWheelF4CircleStarDV1D1 (m)1.4661.4771.4581.8241.6871.499DV2ξ10.9940.8510.9990.7180.7730.983DV3ξ20.9360.8910.5800.9680.9780.943DV4ξ30.4980.4020.4650.4150.4280.418DV5ξ40.9890.8430.9280.9430.9600.903DV6e1 (m)0.0580.0620.0590.2420.0550.087DV7e2 (m)0.2540.0580.3100.2760.2350.253DV8e3 (m)0.2830.3290.3080.2910.3330.284DV9F1 (kN)359.836312.910346.020331.422321.285378.552DV10F2 (kN)92.77669.27192.97083.08685.685102.660DV11F3 (kN)16.44815.38415.3977.61310.19211.159DV12Pc1 (MPa)5.2886.9975.5526.4404.6576.764DV13Pc2 (MPa)5.9285.3263.7606.3257.0005.874DV14Pc3 (MPa)3.7724.9522.8082.0862.3682.123DV15α17.0005.7895.3763.0212.2844.558DV16α23.1514.2772.5512.6193.1162.700DV17α32.8292.2362.4882.4422.5332.307DV18ε120.91617.64720.86119.93014.21526.059DV19ε285.46092.55877.36068.27599.87376.668DV20ε381.554100.00099.25287.28897.569100.000DV21αmax (°)2.1892.1662.1301.8662.1420.981DV22β0.5160.6900.4700.7330.6930.256DV23TGlide (s)83.376.074.3026.711.1TargetM0 (kg)164891684416257154381570515820ConstraintGomax1 (kg/(s∗m2))9883892657786ConstraintGomax2 (kg/(s∗m2))105212312893108ConstraintGomax3 (kg/(s∗m2))5910968286033ConstraintL/D13.115.414.98.910.513.2ConstraintNxmax (g)15.815.116.015.915.915.9ConstraintNymax (g)0.960.620.950.990.980.98Constraintqmax (MPa)0.0660.0520.0650.0540.0570.069ConstraintHp (km)306300302307306313


### 4.1.2. Discussion of the Main Performance Parameters

Figure11 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.
CaseTake-off mass (t)Diameter (m)Length (m)Length-diameter ratioA4CircleCircle16.4891.46619.16313.071BCircleCircle16.8441.47722.80015.436CWheelCircle16.2571.45821.66214.859D4Circle4Circle15.4381.82416.2418.906E4CircleWheel15.7051.68717.69210.486F4CircleStar15.8201.49919.73113.161The 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.Figure12 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.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.
(a)(b)Table8 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 Nxmax and the maximum normal overloads Nymax 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.
PhaseA4CircleCircleBCircleCircleCWheelCircleD4Circle4CircleE4CircleWheelF4CircleStar1st motor work (s)65.281.271.078.778.579.11st-stage separation (s)8.08.08.08.08.08.02nd motor work (s)145.9182.1127.5137.2127.0135.22nd-stage separation (s)8.08.08.08.08.08.0Glide time (s)83.376.074.3026.711.13rd motor work (s)168.5134.6188.1314.9261.7290.4Total time (s)478.9489.9476.9546.8509.9531.8Figure 13
Height-time curve of different cases: (a) whole flight; (b) glide phase and 3rd-stage working phase.
(a)(b)


### 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 M0, maximum mass flux of each stage Gomaxi, maximum of axial overload Nxmax,maximum of normal overload Nymax, maximum dynamic pressure during the flight qmax, height of the perigee Hp, average specific impulse of each stage Isi, effective propellant mass fraction at each stage μi, velocity increment of each stage Vi, total velocity increment V, and total velocity loss ΔV). 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.Figure 14
Pareto graphs of the main performance parameters in case A: (a)M0; (b) Gomax1; (c) Gomax2; (d) Gomax3; (e) Nxmax; (f) Nymax; (g) qmax; (h) Hp; (i) Is1; (j) Is2; (k) Is3; (l) μ1; (m) μ2; (n) μ3; (o) V1; (p) V2; (q) V3; (r) V; (s) ΔV.
(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 ofM0 are the diameter of the 1st stage (DV1) and the initial thickness of grain at the 1st stage (DV6). Gomaxi has the strongest relationship with the grain thickness at the corresponding stage. Nxmax 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). Nymax 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 qmax 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). Hp mainly depends on the diameter (DV1) and the initial thickness (DV6) of grain at the 1st stage. Isi is determined by the grain thickness and initial oxidizer-fuel ratio at the corresponding stage. μ1 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 μ2 and μ3 cover three size scale parameters (DV3, DV4, and DV5). Vi shows similar tendency with μi, while V and ΔV 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 M0, 1st stage diameter D1, total length L, and velocity increment V show similar distribution laws, indicating that there is a positive correlation between these parameters, which is consistent with the results of sensitivity analysis.Figure 15
Self-organizing map of case A: (a)M0; (b) L; (c) D1; (d) V.
(a)(b)(c)(d)


## 4.1.1. Optimal Results

As the mathematical model (as shown in Section3.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.
ParameterDescriptionA4CircleCircleBCircleCircleCWheelCircleD4Circle4CircleE4CircleWheelF4CircleStarDV1D1 (m)1.4661.4771.4581.8241.6871.499DV2ξ10.9940.8510.9990.7180.7730.983DV3ξ20.9360.8910.5800.9680.9780.943DV4ξ30.4980.4020.4650.4150.4280.418DV5ξ40.9890.8430.9280.9430.9600.903DV6e1 (m)0.0580.0620.0590.2420.0550.087DV7e2 (m)0.2540.0580.3100.2760.2350.253DV8e3 (m)0.2830.3290.3080.2910.3330.284DV9F1 (kN)359.836312.910346.020331.422321.285378.552DV10F2 (kN)92.77669.27192.97083.08685.685102.660DV11F3 (kN)16.44815.38415.3977.61310.19211.159DV12Pc1 (MPa)5.2886.9975.5526.4404.6576.764DV13Pc2 (MPa)5.9285.3263.7606.3257.0005.874DV14Pc3 (MPa)3.7724.9522.8082.0862.3682.123DV15α17.0005.7895.3763.0212.2844.558DV16α23.1514.2772.5512.6193.1162.700DV17α32.8292.2362.4882.4422.5332.307DV18ε120.91617.64720.86119.93014.21526.059DV19ε285.46092.55877.36068.27599.87376.668DV20ε381.554100.00099.25287.28897.569100.000DV21αmax (°)2.1892.1662.1301.8662.1420.981DV22β0.5160.6900.4700.7330.6930.256DV23TGlide (s)83.376.074.3026.711.1TargetM0 (kg)164891684416257154381570515820ConstraintGomax1 (kg/(s∗m2))9883892657786ConstraintGomax2 (kg/(s∗m2))105212312893108ConstraintGomax3 (kg/(s∗m2))5910968286033ConstraintL/D13.115.414.98.910.513.2ConstraintNxmax (g)15.815.116.015.915.915.9ConstraintNymax (g)0.960.620.950.990.980.98Constraintqmax (MPa)0.0660.0520.0650.0540.0570.069ConstraintHp (km)306300302307306313


## 4.1.2. Discussion of the Main Performance Parameters

Figure11 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.
CaseTake-off mass (t)Diameter (m)Length (m)Length-diameter ratioA4CircleCircle16.4891.46619.16313.071BCircleCircle16.8441.47722.80015.436CWheelCircle16.2571.45821.66214.859D4Circle4Circle15.4381.82416.2418.906E4CircleWheel15.7051.68717.69210.486F4CircleStar15.8201.49919.73113.161The 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.Figure12 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.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.
(a)(b)Table8 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 Nxmax and the maximum normal overloads Nymax 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.
PhaseA4CircleCircleBCircleCircleCWheelCircleD4Circle4CircleE4CircleWheelF4CircleStar1st motor work (s)65.281.271.078.778.579.11st-stage separation (s)8.08.08.08.08.08.02nd motor work (s)145.9182.1127.5137.2127.0135.22nd-stage separation (s)8.08.08.08.08.08.0Glide time (s)83.376.074.3026.711.13rd motor work (s)168.5134.6188.1314.9261.7290.4Total time (s)478.9489.9476.9546.8509.9531.8Figure 13
Height-time curve of different cases: (a) whole flight; (b) glide phase and 3rd-stage working phase.
(a)(b)


## 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 M0, maximum mass flux of each stage Gomaxi, maximum of axial overload Nxmax,maximum of normal overload Nymax, maximum dynamic pressure during the flight qmax, height of the perigee Hp, average specific impulse of each stage Isi, effective propellant mass fraction at each stage μi, velocity increment of each stage Vi, total velocity increment V, and total velocity loss ΔV). 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.Figure 14
Pareto graphs of the main performance parameters in case A: (a)M0; (b) Gomax1; (c) Gomax2; (d) Gomax3; (e) Nxmax; (f) Nymax; (g) qmax; (h) Hp; (i) Is1; (j) Is2; (k) Is3; (l) μ1; (m) μ2; (n) μ3; (o) V1; (p) V2; (q) V3; (r) V; (s) ΔV.
(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 ofM0 are the diameter of the 1st stage (DV1) and the initial thickness of grain at the 1st stage (DV6). Gomaxi has the strongest relationship with the grain thickness at the corresponding stage. Nxmax 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). Nymax 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 qmax 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). Hp mainly depends on the diameter (DV1) and the initial thickness (DV6) of grain at the 1st stage. Isi is determined by the grain thickness and initial oxidizer-fuel ratio at the corresponding stage. μ1 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 μ2 and μ3 cover three size scale parameters (DV3, DV4, and DV5). Vi shows similar tendency with μi, while V and ΔV 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 M0, 1st stage diameter D1, total length L, and velocity increment V show similar distribution laws, indicating that there is a positive correlation between these parameters, which is consistent with the results of sensitivity analysis.Figure 15
Self-organizing map of case A: (a)M0; (b) L; (c) D1; (d) V.
(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 Table9, 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.
CasePayload mass (kg)Object height (km)A300100100300A500100100500A700100100700A300150150300A300200200300


### 4.2.1. Optimal Results

Using the parameterized design model established in Section2 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.
ParameterDescriptionA300100A500100A700100A300150A300200DV1D1 (m)1.4661.4261.5101.3941.518DV2ξ10.9940.9940.9810.9801.000DV3ξ20.9360.9580.9790.9770.979DV4ξ30.4980.4730.3920.4050.477DV5ξ40.9890.9790.9270.9490.989DV6e1 (m)0.0580.0820.0710.0690.061DV7e2 (m)0.2540.3420.2290.2570.222DV8e3 (m)0.2830.2280.3190.2990.243DV9F1 (kN)359.836364.718465.276350.354464.752DV10F2 (kN)92.77691.446116.690100.868116.972DV11F3 (kN)16.44810.36911.01610.95819.369DV12Pc1 (MPa)5.2885.5296.9946.9026.081DV13Pc2 (MPa)5.9287.0006.2886.9526.501DV14Pc3 (MPa)3.7725.6562.0375.5792.095DV15α17.0004.3434.8903.7014.444DV16α23.1513.0142.9233.2893.046DV17α32.8292.0722.3512.7552.250DV18ε120.91622.30030.02017.21232.685DV19ε285.46062.73288.94084.06593.787DV20ε381.554100.00097.23191.64499.179DV21αmax (°)2.1891.6820.9582.7342.168DV22β0.5160.5630.2690.5240.437DV23TGlide (s)83.285131.649131.35227.84584.972TargetM0 (kg)1648919675221491725620890ConstraintGomax1 (kg/(s∗m2))98105115100105ConstraintGomax2 (kg/(s∗m2))1052609512590ConstraintGomax3 (kg/(s∗m2))5925404045ConstraintL/D13.114.916.015.815.2ConstraintNxmax (g)15.815.716.015.815.2ConstraintNymax (g)0.960.380.490.960.86Constraintqmax (MPa)0.0660.0480.0570.0660.066ConstraintHp (km)306510728309321


### 4.2.2. Discussion of the Main Performance Parameters

Obviously, Table10 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.
NameStageTheoretical velocity increment (m/s)Actual velocity increment (m/s)Theoretical velocity distribution ratioActual velocity distribution ratioРо́кот (liquid) [25]1st3987313944.12%41.43%2nd2607232728.85%30.72%3rd2443211027.03%27.85%Total90377576——Pegasus (solid) [26]Carrier7002367.42%2.90%1st2984230731.64%28.35%2nd2927278731.04%34.23%3rd2818281029.89%34.52%Total94298140——HPSLV of Purdue University (hybrid) [27]1st2773—26.76%—2nd3734—36.03%—3rd3851—37.20%—Total10358———Figure16 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).Figure 16
Velocity increment and velocity distribution ratio: (a) theoretical and actual velocity increment; (b) theoretical and actual velocity distribution ratio.
(a)
Theoretical and actual velocity increment(b)
Theoretical and actual velocity distribution ratioA comparison of caseA300100, A500100, and A700100 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 A300100, A300150, and A300200 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.


## 4.2.1. Optimal Results

Using the parameterized design model established in Section2 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.
ParameterDescriptionA300100A500100A700100A300150A300200DV1D1 (m)1.4661.4261.5101.3941.518DV2ξ10.9940.9940.9810.9801.000DV3ξ20.9360.9580.9790.9770.979DV4ξ30.4980.4730.3920.4050.477DV5ξ40.9890.9790.9270.9490.989DV6e1 (m)0.0580.0820.0710.0690.061DV7e2 (m)0.2540.3420.2290.2570.222DV8e3 (m)0.2830.2280.3190.2990.243DV9F1 (kN)359.836364.718465.276350.354464.752DV10F2 (kN)92.77691.446116.690100.868116.972DV11F3 (kN)16.44810.36911.01610.95819.369DV12Pc1 (MPa)5.2885.5296.9946.9026.081DV13Pc2 (MPa)5.9287.0006.2886.9526.501DV14Pc3 (MPa)3.7725.6562.0375.5792.095DV15α17.0004.3434.8903.7014.444DV16α23.1513.0142.9233.2893.046DV17α32.8292.0722.3512.7552.250DV18ε120.91622.30030.02017.21232.685DV19ε285.46062.73288.94084.06593.787DV20ε381.554100.00097.23191.64499.179DV21αmax (°)2.1891.6820.9582.7342.168DV22β0.5160.5630.2690.5240.437DV23TGlide (s)83.285131.649131.35227.84584.972TargetM0 (kg)1648919675221491725620890ConstraintGomax1 (kg/(s∗m2))98105115100105ConstraintGomax2 (kg/(s∗m2))1052609512590ConstraintGomax3 (kg/(s∗m2))5925404045ConstraintL/D13.114.916.015.815.2ConstraintNxmax (g)15.815.716.015.815.2ConstraintNymax (g)0.960.380.490.960.86Constraintqmax (MPa)0.0660.0480.0570.0660.066ConstraintHp (km)306510728309321


## 4.2.2. Discussion of the Main Performance Parameters

Obviously, Table10 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.
NameStageTheoretical velocity increment (m/s)Actual velocity increment (m/s)Theoretical velocity distribution ratioActual velocity distribution ratioРо́кот (liquid) [25]1st3987313944.12%41.43%2nd2607232728.85%30.72%3rd2443211027.03%27.85%Total90377576——Pegasus (solid) [26]Carrier7002367.42%2.90%1st2984230731.64%28.35%2nd2927278731.04%34.23%3rd2818281029.89%34.52%Total94298140——HPSLV of Purdue University (hybrid) [27]1st2773—26.76%—2nd3734—36.03%—3rd3851—37.20%—Total10358———Figure16 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).Figure 16
Velocity increment and velocity distribution ratio: (a) theoretical and actual velocity increment; (b) theoretical and actual velocity distribution ratio.
(a)
Theoretical and actual velocity increment(b)
Theoretical and actual velocity distribution ratioA comparison of caseA300100, A500100, and A700100 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 A300100, A300150, and A300200 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.



---

*Source: 5574436-2021-04-29.xml*