1. Introduction
A Launch Vehicle (LV) design requires complex decisions to be made from the start of the design phase. Mission designers must decide on inter-related factors, which include the most feasible engine design and propellants to reach the mission’s requirements, environmental impact, manufacturing processes, overall cost and the availability and skillset of human resource required to manufacture, test and deploy a launch vehicle [
1]. These factors are all inter-connected, highlighting the intricacies involved in making major decisions at an early stage of the design process.
A typical rocket engine uses solid, liquid or hybrid propellants. Compared to solid or hybrid propellants, liquid propellants have the following advantages: (1) relatively high energy density, which leads to a higher Specific Impulse (
), resulting in increased engine efficiency; (2) propellant storage tanks require a reasonably low volume; (3) engine throttling (if required) is easily achieved by varying the propellant flow rate in order to control the thrust; (4) allowance for multiple burns as the engine can be repeatedly shut down and restarted [
2]. Liquid Propellant Rocket Engines (LPREs) are classified by their propellant feed system: either pressure fed or pump fed. Mission designers prefer pump-fed LPREs due to the higher combustion chamber pressure (
) produced. Higher
results in a relatively smaller engine; therefore, the payload carrying capacity of the overall LV increases. However, pump-fed engines, also known as turbopump engines, come with varying complexities. There are three types of turbopump LPRE cycles: Gas Generator (GG), Staged Combustion (SC) and Expander. GG and SC cycles contain a secondary combustion chamber where tapped-off fuel and oxidizers are ignited. The resulting exhaust gases are used to power the turbine. This turbine powers the pumps, which raise the propellant pressure to the required
. In SC, the exhaust gases are redirected back to the main chamber for complete combustion. Expander cycles do not contain a secondary combustion chamber. Such engines are only used for cryogenic propellants as their heat of vaporization is greater; therefore, more energy can be extracted and used in mechanical work. Thus, Expander cycle LPREs are only used as upper stage engines in LVs. Hence, for single-stage LVs or first-stage engines, there is often a tradeoff between GG and SC engines: SC engines typically achieve higher
; therefore, the LV can carry greater payloads. However, this increase in payload comes at the expense of greater design complexity.
The design phase for an LV consists of three stages: conceptual, preliminary and detailed design. The propellant pair and engine cycle are finalized in the conceptual design phase. The mission designer’s main aim is to maximize the LPRE performance while simultaneously minimizing the overall LV weight in order to increase the payload [
3]. According to Hammond [
4], the Life-Cycle Cost (LCC) of an LV follows the Pareto Principle: 80% of a vehicle’s LCC is finalized at the conceptual design stage.Thus, any changes in the preliminary or detailed design stages do not result in exponential improvements. Therefore, to minimize the risk of major revisions at a later stage in the design cycle, research recommends a thorough analysis of all major sub-systems and components at the conceptual design phase [
5]. The use of multiple trade-off studies during the conceptual design stage enables mission designers to quantitatively determine the effect of mission requirements on the engine’s performance and overall LV size. Such studies aid mission designers in understanding the effect of interconnected variables such as propellant mass, LV weight,
, thrust and manufacturing complexity [
6].
Trade-off studies rely on mathematical tools to model the propulsion sub-systems and provide LV mass estimates. These tools are based on Mass Estimate Relationships (MERs), which are empirical models derived from parametric equations between thrust and engine mass [
7]. Therefore, the accuracy of the model is proportional to the number of data points used. A review of MERs determined that there is no single set of relationships that is more accurate [
8]. Wildvank [
9], therefore, recently developed a structural dry mass model for an LV, with the aim of moving away from empirical-based mass estimates.
Optimization techniques initially used on LPRE cycles included non-gradient- and gradient-based methods such as Genetic Algorithm (GA) and Sequential Quadratic Programming (SQP) for both entire cycles as well as specific sub-systems [
10,
11,
12]. Akhtar et al. and Bayley et al. [
13,
14] used GA-based techniques to optimize a single-stage and a four-stage LV, respectively. Da mota [
15] utilized SQP to optimize two LPRE cycles with the objective of maximizing payloads. He concluded that algorithms such as Particle Swarm Optimization (PSO) should also be studied. Zhou et al. [
16] modelled the optimal starting and regulating characteristics of an electric pump for LPRE using the GA technique. A fuzzy Multi-Objective GA (MOGA) technique has been applied to a GG cycle [
17] as well as on a hybrid rocket [
18] with the objectives of maximizing altitude and minimizing Gross Lift-Off Weight (GLOW). Orgeira-Crespo et al. [
19] compiled a database of existing LPREs and examined conceptual multistage LVs based on those existing engines. A GA technique was used on these conceptual LVs to optimize the GLOW based on different trajectories. Castellini et al. [
20] compared different techniques to optimize the LV ascent trajectory and concluded that GA and PSO techniques were the most feasible techniques to use.
A surrogate model-based two-step method was used by Mack et al. [
21] to optimize a radial turbine in an Expander cycle LPRE. First, Response Surface Approximation was used to search the solution space to develop an approximate model. Secondly, GA was used to refine the search area and obtain a global optimum point. An alternative to this method is CFD, which is computationally intensive and time-consuming. Cui et al. [
22] also used a surrogate-based model to optimize the trajectory of an air-launched solid rocket. Neural Network (NN)-based models were applied by Dresia et al. [
23] to LPRE problems such as estimating the fatigue life of reusable thrust chambers and estimating the maximum temperature in nozzle cooling channels. In both cases, NN models were created as large sets of training data were available: 120,000 data points to estimate fatigue life and 20,000 CFD simulations to calculate the maximum temperature. NN-based controllers are proposed by Yu et al. [
24] to ensure steady-state spacecraft formation when in orbit. Therefore, similarly to [
25], surrogate and NN models are beneficial when problems are computationally intensive and already have numerous data sets available. For the purposes of this research, the numerical calculations are not computationally intensive; therefore, the benefits of such techniques are minimized and non-gradient- and gradient-based algorithms can be directly applied.
When determining GLOW for an LV, according to the best of our investigations, researchers, apart from [
19], used mass models based on MERs to provide estimates for the masses. While minimizing GLOW and maximizing payload have been used as Objective Functions (OFs), to the best of our knowledge, other OFs such as maximizing
or Thrust-to-Weight ratio (
) have not been investigated. In addition, the optimization models mentioned previously require Thrust (
F) and
as inputs. This results in running the optimizer numerous times at different user defined
F and
values in order to obtain a Pareto front of optimal points. Furthermore, although newer techniques such as the Cuckoo Search Algorithm (CSA) and Pigeon-Inspired Optimization (PIO) have successfully been used in various aerospace applications [
26], to our knowledge, these techniques have not been used in LPREs.
This research has three aims. First, to estimate the mass and dimensions of a pump-fed cycle LPRE using mathematical model and not via statistical or empirical estimates. Second, to compare seven techniques, including a novel hybrid GA-PSO technique (GA-Swarm), and four OFs on GG and SC cycles and determine the most suitable technique and function for each cycle. Thirdly, a comparison is made between the two cycles for a given mission.
The structure of this paper is as follows: a description of the modelling and optimization process used in this research is provided in
Section 2; a validation of the mass model and results of the optimization algorithms and OFs are discussed in
Section 3; conclusions and future work are provided in
Section 4.
2. Methodology
A flow chart describing the process used to estimate the mass and dimensions of a pump fed LPRE is shown in
Figure 1. The framework developed in this research has a modular format. Therefore, the propulsion sub-systems are ’assembled’ in various combinations to model different cycles.
The required inputs are propellant pair, the LPRE cycle to be modelled (SC or GG) and the mission. The mission as defined in this process includes the required change in velocity (
) and the required altitude (
H). In this research,
,
F and burn time (
) are calculated during the optimization process, unlike previously mentioned models in
Section 1 where
and
F were required inputs.
2.1. LV Model
Once the propellant pair is entered, Chemical Propellant Equilibrium Program (CPropEP), an external software based on NASA calculations [
27], is used to calculate propellant thermodynamic properties. Detailed information on CPropEP is given in [
28,
29]. This software calculates a range of output variables, from which the characteristic velocity (
), thrust coefficient (
), combustion chamber temperature (
), isentropic expansion factor (
), molar mass (M) and the expansion ratio (
) are required for further calculations in this research.
As is not fixed and as optimum is required, the above variables would have to be recalculated in every iteration during the optimization process. In order to reduce computing times, an additional code was written to calculate the optimum Mixture Ratio (MR) for the propellant combination at different values of . The data were fitted using non-linear curve-fitting. The equations generated are then used to recalculate the thermodynamic variables in each iteration during optimization instead of repeatedly running CPropEP again. Therefore, the advantages of this additional step are as follows:
Time saved during the optimization process as CPropEP is not repeatedly run for each iteration. Instead, it only runs once to generate equations for each of the required variables.
The optimum MR for each is calculated to achieve the maximum at that specific . Maximizing is important as it indicates a higher combustion efficiency and an increase in engine performance.
The thermodynamic properties are then used in the engine model developed for the purpose of this research.
This method also minimizes the optimization processing time as otherwise CPropEP would load and run for each value of
.
MR is calculated via the equation below:
The output variables
,
,
,
M and
are calculated in succession using the following equation:
where
y alternately represents
,
,
, M and
.
Equation (
3) is used to calculate
.
These equations and the values of the variables a, b and c are unique to the combination of a given propellant pair at a given exit pressure. If the performance of an LPRE cycle needs to be analyzed with a different propellant pair or the performance of the propellant pair needs to be evaluated at a different value of , data need to be generated from the start to obtain the required values for each equation.
2.1.1. Propellant Mass and Tank Calculations
A forward calculation is used in the optimization techniques, which requires Thrust
F,
and burn time (
) to be specified. Therefore, initial values are given; however, the final values are obtained as an optimized output. Using the output values calculated in
Section 2.1 and standard equations as defined in [
30], specific impulse,
, propellant mass flow rate,
, and overall propellant mass,
, are calculated using the following equations:
where
is the standard value of Earth’s gravity defined as 9.81 m/s
. The multiplication by
is for the purpose of unit conversion only as
is defined as total impulse per unit weight of propellant. Using the optimum MR calculated in the previous section, the mass of fuel,
, is calculated via the following.
As fuel has been predefined, its density is already known; therefore, the volume of fuel and, hence, volume of fuel tank can be calculated as per Equation (
9). An additional 10% volume is added to account for ullage volume.
Assuming the tank is cylindrical in shape, with a Length:Diameter (L/D) of 7, rearranging the standard equation for volume of a cylinder, the radius and length are calculated in the following equations.
The selected propellant tank material for both propellants is Stainless Steel 316 [
30]. Safety factors of 1.25 and 1.5 are, respectively, applied to the yield strength and ultimate strength to calculate the maximum allowable operating stress,
.
The total operating pressure is calculated as given below where the pressure drops due to the injector head assumed to be 10 bar, and an additional 10% is added due to losses in the fittings.
As the tank material and allowable operating stress is known, the equation for Hoop Stress as defined in [
30] is rearranged to obtain the tank’s thickness via the following equation.
Finally, the propellant fuel tank mass is determined using Equation (
14):
where
is already known as the tank material is pre-defined. Equations are then repeated to calculate the oxidizer mass,
, and oxidizer tank mass,
. Instead of Equation (
7), the following equation is used to calculate the oxidizer mass.
In this framework, it is assumed both tank radii are the same. Therefore, once the radii of both tanks are calculated, the maximum value is taken as the radii for both, and the tank lengths and thicknesses are recalculated accordingly.
2.1.2. Pressurant Mass and Tank Calculations
The selected pressurant is Helium, He, and the pressurant tank material is Ti 6Al-4V [
30]. The initial tank pressure,
, is set at 450 bar. The number of moles of He required to pressurize the fuel and oxidizer tanks is calculated using the Ideal Gas Equation:
where
R is the molar gas constant, 8.314 J/Kmol, and temperature,
T, is assumed at room temperature, 298.16K. The volumes for fuel (
) and oxidizer, (
) were previously calculated using Equation (
8).
The ratio of He left in the tank,
, is assumed to be 2. The pressure of He left in the tank,
is given by the following.
The minimum number of moles of
He, (
nHe), therefore, required to pressurize both propellant tanks is calculated as follows.
The mass of He,
, is obtained via converting the number of moles of He to kg.
The Ideal Gas Equation is again used to calculate the volume of He,
, in the tank. This is assumed to be the inner tank’s volume. The pressurant tank is assumed to be spherical; therefore, using the standard equation for volume of a sphere, the inner radius is calculated via the following.
Given the allowable stress,
S (assumed 1100 MPa), and joint efficiency,
E (0.9), of the material and the equation from [
30], the tank’s thickness,
, is determined as follows.
The tank’s outer radius,
, is simply the following.
The He tank outer,
, and inner
volumes are calculated using the volume of a sphere and the respective inner and outer tank radii given in equations. The mass of the pressurant tank,
, is therefore calculated below with an additional 20% mass as the safety margin.
2.1.3. Combustion Chamber Mass Calculations
The combustion chamber’s shape is assumed to be cylindrical with a converging-diverging nozzle attached at the lower end. The assumptions, based on those given in [
30], used are chamber area (
) equivalent to three times the throat area (
), nozzle characteristic length (
) of 1.143 m, nozzle converging section is a truncated cone and the chamber and nozzle thickness is 3 mm. A 45°angle is assumed for the converging nozzle section and 15°angle for the diverging section. The nozzle and chamber’s material comprises stainless steel. Therefore, once
is calculated as in Equation (
25), and as material
is known, standard equations for volume of a cylinder and truncated cone are used to calculate the remaining dimensions and, consequently, the nozzle’s mass,
, and combustion chamber’s mass,
.
2.1.4. Pre-Burner and Gas Generator Dimensions and Masses
The pre-burner and gas generator dimensions and mass (
and
are calculated in a similar manner to the main thrust chamber given in the previous section). The assumptions about the material, shape, converging angle and nozzle thickness are as given in
Section 2.1.3.
The mixture ratio for both GG and PB is fixed at 0.5. CPropEP is used to calculate the subsequent thermodynamic parameters required for mass and dimension calculations.
2.1.5. GG Parameters
For the GG cycle, the combustion chamber pressure of the GG, , is assumed to be the same as the main chamber, .
The propellant mass flow rate for GG,
, is assumed to be 2% of the main chamber flow rate, as given in [
30].
2.1.6. SC Parameters
For the SC cycle, the combustion chamber pressure of the PB,
, is calculated via the critical pressure ratio, as given in [
30]:
The propellant mass flow rate for
PB is assumed to be fuel rich [
30]. It is assumed that the entire fuel volume flows through the
PB. The total fuel mass flow rate through the
PB is assumed at
. As the mixture ratio is already known, as is
, the
PB oxidizer mass flow rate is then calculated.
2.1.7. Engine Balancing
The design of a turbopump system depends on the power balance between the pumps and turbine. The mission designer must ensure that the power generated by the turbine matches the power required by the pumps.
The power required by each pump is a function of the volumetric flow rate,
Q, change in pressure,
, and the pump efficiency,
, as shown in the following equation.
The power required by each pump is calculated separately and summed to provide the total power required by the pumps, .
The actual power produced by the turbine is calculated via the following:
where
is the turbine efficiency,
is the turbine mass flow rate,
is the specific heat capacity at constant pressure of the gas entering the turbine,
is the Turbine Inlet Temperature,
is the turbine pressure ratio (
), and
is the ratio of the specific heats.
The equations given above undergo several iterations until the power produced by the turbine matches the power required by the pumps.
2.1.8. GG Power Balance
GG is an open cycle LPRE. The turbine pressure ratio across the turbine is kept fixed at 20 [
30]. In addition, the turbine mass flow rate is assumed equal to the GG mass flow rate. As mentioned in
Section 2.1.5, the mixture ratio is fixed as is the GG flow rate. Therefore, the oxidizer and fuel mass flow rates through the GG can be calculated as per
Section 2.1.5.
2.1.9. SC Power Balance
For closed cycles, such as SC and Expander, the LPRE tool keeps the turbine pressure ratio as a variable to be optimized. According to [
30], this ratio is approximately 1.5 for SC and Expander cycles. The turbine must generate enough power to drive two oxidizer pumps and one fuel pump. The PB mass flow rate and individual propellant mass flow rates through the PB are calculated as mentioned in
Section 2.1.6. The volumetric flow rate through the secondary oxidizer pump,
is taken as half of the total fuel volumetric flow rate,
. The turbine mass flow rate is assumed to be equal to the PB mass flow rate; however, this flow rate is greater in magnitude compared to the GG flow rate as the entire fuel is directed to the PB.
2.1.10. System Mass and Dimensions
The overall system mass,
, is then a sum of the combustion chamber, nozzle, tanks, propellant and pressurant masses, turbopump assembly, Pre-Burner (or Gas Generator) and the payload, as shown in Equation (
29). This is also known as Gross Lift-Off Weight or GLOW. Calculations for system mass and dimensions are as per the model previously described in [
31].
For the SC cycle, would replace and the secondary oxidizer pump, , would also be included.
2.2. Flight
This section describes the models used to calculate the forces acting on the LPRE (and hence the amount of Thrust, F, required) as it ascends through the atmosphere to reach the required mission altitude. Changes in gravity, air density and atmospheric pressure and temperature are taken into account. The framework assumes vertical lift-off and flight.
Once engine parameters are calculated, F and Drag (D) are calculated at pre-determined intervals of time to determine the performance of the LPRE from takeoff until reaching the required H.
2.2.1. Atmospheric Model
Atmospheric temperature, (
), and atmospheric pressure, (
), vary with altitude,
H. This variation influences the air density (
), which in turn affects the drag on any vehicle passing through the Troposphere, lower Stratosphere and upper Stratosphere. In this research, the US Standard Atmosphere Model, 1976 [
32], is used to calculate the pressure, temperature and density with increasing altitudes as per the following equations.
For the Troposphere with an altitude below 11,000 m,
and
are calculated via the following:
When the LPRE is travelling through the lower Stratosphere where the altitude is between 11,000 m and 25,000 m, the following set of equations are used to obtain
and
:
The upper Stratosphere lies above 25,000 m. For these higher altitudes, (
) and (
) are determined as per the following equations:
The three equations are given above for calculating return values in °C. Therefore, an additional step is required to add 273.1 to each value to convert to Kelvin. Similarly, the values of are in kPa and must be converted to Pa in order to ensure all units are uniform throughout the framework.
The air density,
, is a function of
and
and is calculated via the followin.
2.2.2. Gravitational Model
Newton’s Law of Universal Gravitation is used to determine the change in acceleration due to gravity,
, at each change in altitude. The values of the constants
G, Gravitational constant,
, Earth’s mass and
, Earth’s radius, are given in
Table 1. These values are used in the following equation to determine
at a particular altitude.
2.2.3. Flight Model
The values calculated in
Section 2.2.1 and
Section 2.2.2 are used to determine the Thrust,
F, and Drag,
D, during the LPRE’s flight from takeoff to the required altitude,
H. For this framework, a time interval,
, of 0.1 s is used when calculating the velocity and height during the flight. The Drag Equation (assuming a fixed drag coefficient (
) of 0.2) and General Thrust Equation are, respectively, used to calculate the Drag and Thrust, as given in the following equations:
where
and
are defined as per [
30] using the universal gas constant
R assumed to be 8.314 J/Kmol.
,
M and
are calculated via
Section 2.1 and used as follows.
The exit pressure (
) is fixed and assumed at 0.7 bar for this research study. Therefore, for each interval of time,
is recalculated as per the equations given in
Section 2.2.1 and then used in the Thrust equation. The second term in the Thrust Equation is the pressure component, which varies as
changes as the LPRE gains altitude.
F, therefore, increases when the nozzle exit velocity,
, is at a maximum point (as propellant mass flow rate,
is fixed).
In the Drag Equation, is assumed as the frontal surface area. This is a fixed value, and as is kept at a fixed value, changes in D are only due to changes in and velocity, v.
The system mass at each interval of time and
D,
F,
,
,
and gravity,
, calculated earlier are used to determine the system mass,
, velocity,
v, and subsequently the altitude,
H, at each time interval as per the following equations.
2.3. Optimization
There is a need to quantitatively compare different optimization techniques to evaluate which methods are most suitable for a particular problem. There may not be a single method that consistently performs better than others for different objective functions. This is why the following techniques are evaluated according to their performance: Simulated Annealing (SA), Nelder Mead (NM), Cuckoo Search Algorithm (CSA), Pigeon-Inspired Optimization (PIO), Particle Swarm Optimization (PSO), Genetic Algorithm (GA) and a novel hybrid GA-PSO technique named GA-Swarm, which is described in further detail in
Section 2.3.1.
The importance of selecting a relevant OF for use in multidisciplinary design studies cannot be overemphasized. Da Mota et al. [
15] tried to minimize the Gross Lift-Off Weight (GLOW), the structural mass (SM) and the expected total cost to investigate the effect of changing the OF using a single optimization technique. Dresia et al. [
33] investigated the effect of minimizing GLOW or Structural Mass (SM) on a GG cycle using GA techniques and Akhtar et al. [
13] studied the effect of GLOW using GA on a pressure-fed cycle.
These seven optimization techniques and four Objective Functions (OFs) are used to optimize each LPRE cycle and to determine the most suitable OF and technique. A list of the four OFs, design variables and constraints is provided in
Table 2.
2.3.1. Genetic Algorithm-Swarm (GA-Swarm)
GA-Swarm is a novel technique developed in this research by incorporating the ’swarm’ behaviour of PSO within GA. GA-PSO hybrid techniques have been used previously [
34,
35,
36]. However, these techniques either incorporate mutation within PSO to ensure the technique can escape local optima or apply the GA ’fitness’ function to swarm particles in order to divide them and then apply crossover and mutation functions on the ’better’ half of the population.
In GA-Swarm, as shown in
Figure 2, PSO is incorporated into the GA algorithm by ensuring that after crossover, the offspring ’swarm’ towards the best solution in order to ensure a more efficient search within the solution space.
The offspring is given a bias toward the best solution if a random number is less than a predefined probability, pSwarm, which is usually equal to 0.2.
For the strength of bias, sSwarm is taken to be equal to 0.2.
2.3.2. Constraints and Boundaries
There are two constraints for all techniques and OF:
and
H. In order to ensure that constraints are met every time the algorithm is run, a penalty-based method is applied. This weightage factor is provided in Equation (
46):
where there is a penalty for each constraint, as shown in Equation (
47):
where the Objective Function is as given in
Table 2: i.e., GLOW,
,
or maximizing both
and
. In every iteration, there is a check whether the velocity and height calculated match the required values of
and
H, respectively. Both constraints are multiplied by a penalty coefficient of 100,000 for
and 10,000 for
H. If the constraints are violated, then the measure of violation will be non-zero. If the constraints are not violated then the measure of violation will be zero as the penalty coefficient multiplication has no effect.
The upper and lower bounds for the three variables:
F,
and
are given in
Table 3. These variables are the same for both GG and SC cycles.
The number of iterations per run for all techniques was kept at 1000, except for PSO and GA-Swarm where each run had 150 iterations. For PSO, the number of particles in the ‘swarm’ was 50; the number of ‘pigeons’ in PIO was 128; CSA algorithm contained 125 ‘nests’, and both GA and GA-Swarm algorithms had 150 ’chromosomes’.
3. Results
3.1. Model Validation
Shahab-3 is an Iranian-designed single-stage GG rocket. The first documented flight was in 1998; however, this rocket is currently still in use.
Table 4 lists the approximate engine performance parameters as sourced from the open literature [
37].
The values of the variables used to calculate
,
,
,
,
M and
as per Equations (
1)–(
3) are provided in
Table 5. A graphical visualization of this data is included in
Figure A1.
Two optimized results are included in
Table 6 using PSO as the optimization technique. PSO is an established optimization technique and has provided consistent results for all four OFs used in this research. Hence, PSO is used to optimize the given model. The first result gives the lowest GLOW as per the OF used in the optimization technique. The second result also used GLOW as OF; however, the
returned is close to the claimed
value of Shahab-3 as per the literature. The reason behind including this second set of results is to highlight that optimization techniques can return a range of values within the solution space. Therefore, according to this particular optimization technique and the OF used, this engine could have a lower GLOW but with a higher
for the same propellant pair and cycle. The result with a similar
(48 bar) to the actual engine also calculated a similar GLOW compared to the reported engine values.
Differences in results compared to the actual engine are due to a range of factors. The choice of material for tanks and engine, L/D used, turbine speed and actual turbine power produced in the actual engine are possible causes that account for the difference in values.
A visual summary of the optimized parameters is given in
Figure 3.
3.2. Gas Generator Cycle: Comparison of Techniques and Objective Functions
Four different objective functions and seven optimization techniques are tested on a RFNA-RP1 GG cycle with mission constraints
3800 m/s and an
H of 120 km, as shown in
Figure 4, with a payload of 1000 kg. A ’box’ is included in each graph to depict the operating range of
generally ascribed to GG cycle engines, as compiled by [
11]. The results of three OF (minimum GLOW, maximum
and maximum
) are provided in
Figure 4, plotting the optimum values of each OF against
. The results of the fourth OF are displayed separately in
Figure 5.
The purpose of testing multiple objective functions and techniques is to determine the most suitable function and technique for the GG LPRE cycle.
3.2.1. Nelder Mead
With reference to
Figure 4a, the NM technique returned results that met the mission’s constraints for two out of three OFs. Therefore, the results for maximizing the Thrust-to-Weight ratio have been excluded. A probable reason for this inability to meet the constraints is due to NM’s susceptibility to the initial guess. This technique searches for the local maximum within the vicinity of the initial guess. Therefore, if the initial guess is close to the global maximum, the probability of NM converging on the local maxima is higher than when the initial guess is further away from the global maximum.
Out of these results for NM, the best-performing OF is , as the optimum result both clearly lies within the accepted operating range of for a GG cycle, and occurs at a of 71 bar, which is lower than the minima obtained using GLOW OF. For GLOW OF, there are distinct bands of results at 20, 30 and 40 bar , with the minima lying at 97 bar. This minima lies within the operating range for a GG cycle; however, it is relatively higher than the obtained using as the OF. Using as an OF provided a relatively larger cluster of results between 11 and 18 bar , which is an infeasible operating range for GG cycle. Another smaller band of results between 40 and 45 bar coincided with a similar band of results when using GLOW as the OF. Results of both OFs demonstrate the presence of multiple local minima in the solution space.
3.2.2. Simulated Annealing
Figure 4b shows the results for three OFs using the Simulated Annealing technique. Using GLOW as the OF provides a band of results between 70 and 74 bar
. The results also follow a fairly linear path within the GG cycle
operating range. However, this technique and OF return results between 12 and 13 bar and one result at 479 bar, all outside the accepted
’s operating range for this LPRE cycle.
Similarly to the NM technique, using as the OF provides a relatively greater number of results below the acceptable range for this cycle as compared to results within the accpetable range. In contrast, the T/W function follows the same trend as GLOW by providing a greater range of results within the acceptable range. The T/W function provides a greater number of results between 43 and 60 bar. These bands demonstrate SA’s strength in avoiding local maxima.
These technique required multiple runs in order to generate a sufficient data set for analysis as the technique did not meet the constraints on every run. Therefore, although results are promising, this technique is not recommended due to the additional computational time involved in generating numerous runs and filtering out erroneous data.
3.2.3. Cuckoo Search Algorithm
Looking at
Figure 4c, the majority of CSA results lie between
of 29 and 86 bar for all three OF. There is relatively less variation between the minimum and maximum
for all three OFs, as compared to NM and SA. The maximum
achieved with any of the three OF using CSA is 84 bar, unlike NM and SA where a
of greater than 90 bar was noted. Similarly to SA, the three OFs provide an optimum result within a relatively close
range: minimum GLOW at 86 bar, maximum T/W at 84 bar and maximum
at 77 bar. This narrow range of results highlights CSA’s immunity to the objective function.
The values for GLOW lie on a fairly linear curve with distinct clusters of results: 39–41 bar, 56–63 bar, 70–72 bar and 78–83 bar. also provides a cluster of results in a similar range of 69–72 bar. The curve for T/W highlights two clusters: 44–55 bar and 64–75 bar. Therefore, unlike NM and SA, all three OFs return clusters of results near the global optimum result, highlighting the presence of multiple solutions near the optima within the solution space. Apart from GLOW OF, CSA failed to meet the constraints for the other two OFs on numerous runs; hence, multiple runs were required to generate a sufficient data set of results.
3.2.4. Particle Swarm Optimization
In
Figure 4d, the trend for all three OFs using the PSO technique follows that of SA and NM with a wide range of values returned. Although the majority of results for all OFs lie within the GG cycle
’s operating range, both
and T/W OFs return results below and above this range. In particular, the global maxima for both OFs lies above the operating
range, with a maximum
at 192 bar and maximum T/W at 175 bar, whereas GLOW OF returns a global minima at 112 bar.
The results for GLOW OF are contained in two distinct clusters: 30–33 bar and 57–65 bar. also returns a cluster of results at a similar range of 54–56 bar, underlining the existence of multiple local minima in the solution space. T/W also returns results in three clusters: 23–26 bar, 52–59 bar and 81–89 bar.
It is pertinent to mention that this technique returned results for T/W (111 bar) and (114 bar) functions at a value of near 112 bar (the that resulted in minimum GLOW). Therefore, this indicates the presence of an optimal value for each OF in the search space. CSA also meets the mission constraints for each optimization run.
3.2.5. Pigeon-Inspired Optimization
In
Figure 4e, similarly to CSA, the results for all three OFs using the PIO technique vary between a relatively short range: 23 bar to 95 bar. Despite multiple runs, the technique failed to meet both constraints on numerous occasions; therefore, only a limited number of values that met both constraints are included in this research study. However, the majority of these limited values lie within the operating
range for a GG cycle. GLOW OF returns a minima at 84.5 bar
, whereas a maximum
is attained at 95.7 bar and maximum T/W is attained at 97 bar.
The data points for each OF are in small clusters, demonstrating multiple local minima. With GLOW OF, there is a cluster near the global minima between 79 and 84 bar. The T/W function returns three groups of values at the following: 30–25 bar, 55–61 bar and 93–97 bar. There is one cluster of values for ranging from 48 to 61 bar, which overlaps with a similar (but smaller) cluster using T/W as the OF.
3.2.6. Genetic Algorithm
The results for all OF using the GA technique, as shown in
Figure 4f, are similar to CSA as all results lie within a relatively narrow range of
(38–89 bar) compared to other techniques. However, unlike CSA, all results for GA lie within the operating
range for this cycle. The minimum GLOW that meets the constraints is at 96 bar
, the maximum
occurs at 89.9 bar and the maximum T/W occurs at a lower
value of 75 bar.
For as the OF, the results are in three distinct clusters, 38–48 bar, followed by an equal number in the 58–67 bar and then 80–89 bar range. T/W also returns a dense group of results between 62 and 68 bar, as does GLOW (59–68 bar), reinforcing the presence of multiple local minima. This clustering shows that GA does not become stuck in multiple local minima.
3.2.7. GA-Swarm
The results of the final technique, GA-Swarm, are displayed in
Figure 4g. Due to the hybrid nature of this technique, it not only meets the constraints for each run but also provides results with minimal variance and within a relatively narrow range of
, similarly to GA. The maximum
that meets the constraints is at 86 bar, the minimum GLOW is at 87 bar
, and the maximum T/W occurrs at a lower
value of 79 bar. These values are similar to those obtained using the GA technique.
Unlike GA, the values obtained using as OF are not grouped as distinct clusters, however, a dense group of results is observed between 63–66 bar and 71 bar. Using T/W as an OF returns three clusters at: 38–45 bar, 63–68 bar and 71–77 bar. The results for GLOW as OF follow a similar trend as the GA technique, with a relatively large cluster between 50 and 59 bar.
3.2.8. Objective Function: T/W vs.
The results of a fourth OF used on GG cycle are shown in
Figure 5 where maximizing
and T/W are both given equal weightage. Values of the global optima for each technique are given in
Table 7. Five techniques return results with a maximum
between 262 and 272 s; however, only CSA and SA provide theoretically higher performance values for this Objective Function.
F is directly proportional to
; therefore, a linear trend is expected if
is plotted against
, which is verified by the results as shown in
Figure 5. This Objective Function is included in this research to validate the seven techniques. All techniques returned a similar trend line, as shown in
Figure 5. If a technique had not returned values other than those shown in
Figure 5, the Objective Function would, therefore, have invalidated that particular technique.
3.3. Evaluating Objective Functions and Techniques
Results for all techniques are shown in
Table 8,
Table 9 and
Table 10. When comparing the OF results from each technique, some techniques such as CSA, PIO and GA-Swarm are immune to the influence of one technique over another as the results returned are within a narrow range regardless of the OF. However, PIO (along with NM and SA) and CSA will still be regarded as techniques that fail to meet the mission constraints on numerous occasions. GA-Swarm meets the mission constraints on each run.
Using GLOW as an OF provides relatively few results which lie near the global minimum (as clusters are observed in local minima). However, when comparing the results for each technique, results from T/W provided a maxima at a relatively lower when using GA-Swarm compared to the other OF. Based on experience, the end user prefers an OF with emphasis on low . Therefore, using T/W as an OF may be more feasible as compared to using GLOW.
3.4. Staged Combustion Cycle: Comparison of Techniques and Objective Functions
Similarly to
Section 3.2, the same objective functions and optimization techniques are tested on an RFNA–RP1 SC cycle. Mission constraints and payload remain the same and the results for three OF are shown in
Figure 6. A ‘box’ is included in each graph to depict the operating range of
generally ascribed to SC cycle engines [
11].
This given mission has returned values (for all techniques and objective functions) that lie below the operating range for SC cycles. However, for the purposes of this research, as the SC cycle is compared to GG cycle for a given mission, the overall trends and results for each technique will be discussed.
3.4.1. Nelder Mead
With reference to
Figure 6a, unlike for GG cycle, the NM technique returned results that met the mission constraints for all functions. The minimum GLOW is obtained at a
of 122 bar, while with
as the OF, the maximum
is determined at 99 bar and the maximum T/W occurs at 92 bar.
For OF, there are distinct bands of results from 26 to 33 bar, 47 to 52 bar and 65 to 82 bar. These clusters coincide with similar groups of results when T/W is used as OF, highlighting the presence of local minima in the solution space.The results for GLOW as an OF follow a fairly linear trend without clusters of results. There is also relatively less variance in values as compared to the other two OFs.
3.4.2. Simulated Annealing
Figure 6b shows the results for the three OFs using the Simulated Annealing technique. All three OFs return the optimum result within a relatively close
range: minimum GLOW at 94 bar, maximum T/W at 96 bar and maximum
at 96 bar.
The results for all three OFs using the SA technique display visible clusters of results, once again highlighting the effectiveness of this technique in evading local minima. However, this technique is unable to meet the mission constraints for all OFs on each run, requiring additional runs in order to obtain a minimum data set of points.
3.4.3. Cuckoo Search Algorithm
The results from CSA technique are shown in
Figure 6c. This technique follows the same trend as when run for GG cycle, with results for the three OFs lying within a relatively narrow range of
, between 46 bar and 92 bar.
Another similarity between the results of both cycles using this technique is the presence of multiple clusters using all three OFs. Once again, the clusters are close to the global optimum’s result, which demonstrates CSA’s superior search ability compared to NM and SA.
3.4.4. Particle Swarm Optimization
The trend for all three OFs using PSO technique is shown in
Figure 6d. The graph shows similar results for GG cycle using this technique. Although each run meets exactly the mission’s constraints, however, a wide range of values are returned, which eventually converge on the global optimum solution.
This technique also returns optimum results for each function at a relatively higher compared to other techniques: minimum GLOW at 133 bar, maximum T/W at 111 bar and maximum at 110 bar.
3.4.5. Pigeon-Inspired Optimization
The results for PIO, as shown in
Figure 6e, are similar to both CSA and PSO. The similarity with CSA lies in the fact that the results for all three OFs using the PIO technique vary between a relatively short range: 30 bar to 104 bar. Similarly to PSO, there are no distinct clusters of results (as shown in CSA, SA and NM); instead, results show a gradual convergence to the optimum solution for each technique. This is also a contrast to the GG cycle where data were shown in clusters. Moreover, PIO managed to meet the criteria in each run, unlike for the GG cycle.
3.4.6. Genetic Algorithm
For SC cycle, the results of all OF with GA techniques are shown in
Figure 6f. For GLOW and
functions, the technique failed to meet the constraints on numerous runs. There are only a limited number of data points shown. For T/W, the results show clusters, with a relatively dense cluster between 47 and 70 bar. This cluster coincides with two and one data points using GLOW and
functions, respectively, demonstrating that for this particular cycle, T/W is a relatively better function for the GA technique.
3.4.7. GA-Swarm
For the hybrid technique, GA-Swarm results are shown in
Figure 6g. There is a marked improvement compared to GA results as this technique met the mission constraints for each run, which is a characteristic of the PSO technique. In addition, unlike PSO, GA-Swarm provided results in a relatively narrow range, demonstrating qualities similar to GA. The maximum
that meets the constraints is at 93 bar, and the minimum GLOW and maximum T/W both occur at a lower
value of 81 bar and 86 bar, respectively.
Unlike GA, the values obtained using as OF are not grouped as distinct clusters; however, a density of results is observed between 63 and 66 bar and 71 bar. Using T/W as an OF returns three clusters at 38–45 bar, 63–68 bar and 71–77 bar. The results for GLOW as OF follow a similar trend as the GA technique, with a relatively large cluster between 50 and 59 bar.
3.4.8. Objective Function: T/W vs
In
Figure 7, the results of all seven techniques using the fourth OF on SC cycle are displayed.
Table 11 displays values for each variable at max
and T/W. For this cycle, NM returns the maximum result for
and T/W; however, due to NM’s susceptability to the initial guess, these results are disregarded. SA, CSA and GA-Swarm returned similar values of
while SA and GA-Swarm also returned same values of T/W. As mentioned in
Section 3.2.8, this Objective Function was specifically included in order to further validate each technique. Similarly to
Figure 5, all techniques follow a similar trend with
linearly increasing with
. As the same mission was optimized using both cycles, the magnitude of results for GG and SC using this Objective Function is also similar.
3.5. Evaluating Objective Functions and Techniques
The results for all techniques are shown in
Table 12,
Table 13 and
Table 14. The most promising techniques for SC cycle are CSA, PSO, PIO and GA-Swarm. All four techniques consistently met the mission’s constraints and provided a wide range of results. Apart from PSO, the use of GLOW as an OF returns the lowest value of
compared to the other two OFs, with
F and
values remaining nearly equal. Therefore, if the mission designer requires a low
that still provides high engine performance, then GLOW should be used as an OF for the SC cycle.
3.6. GG vs. SC Cycle for Given Mission
The primary aim of comparing GG cycle and SC cycle with the same mission requirements is to aid the mission designer in decision making when faced with the task of selecting a particular engine cycle for a mission.
For GG cycles, all techniques using GLOW OF returned a minima between 52 and 112 bar and an average F of 259 kN. The same techniques on SC, however, resulted in the minima at a relatively higher of 98 to 102 bar with a resulting average increase in F to 272 kN. This higher is still within the operating range for an SC cycle engine.
Therefore, for this mission, a GG cycle would be more beneficial as the same GLOW can be achieved at a relatively lower without a significant loss of F. In addition, the additional design, manufacture and operational complexity of a SC cycle engine do not need to meet the requirements of this mission.
4. Conclusions
This research includes a preliminary design of a conceptual pump-fed LPRE cycle, including mass and dimension estimates of all sub-systems using mathematical models and not empirical models. Unlike previous research, the only inputs required for this model were the mission requirements (H and ) and the propellant pair. The model then calculates the optimum set of engine performance parameters. In addition, seven optimization techniques and four Objective Functions were analyzed in order to provide the maximum performance of the cycle for a given mission, propellant pair and payload. An Objective Function simultaneously maximizing both and was used to validate all techniques. One of the optimization techniques, GA-Swarm, used was a novel technique. This research concludes that the best-performing OF for the GG cycle was T/W as the optimized engine performance values for all variables were better compared to using GLOW as an Objective Function. The best-performing Objective Function for SC was GLOW. However, when comparing two cycles, the benchmark OF to consider for both is GLOW.
Looking at the optimization techniques, Nelder Mead is extremely susceptible to initial guess across both cycles and, therefore, cannot be used for an initial optimization. Simulated Annealing shows variance, provides a wide range of values and also returns results in clusters, highlighting that the results are immune to local minima. However, Simulated Annealing does not match the constraints on each run and, therefore, requires additional computational time to generate a data set of values for analysis. Cuckoo Search Algorithm technique works well for SC cycle only. PSO again shows variance by providing a wide range of results for all cycles and missions; however, this technique meets the mission constraints for each run. Pigeon-Inspired Optimization provided optimum performances for the SC cycle only. Similarly to CSA, for the GG cycle, this technique struggled to meet constraints and provided a greater scatter effect in the results. GA provides results within a wider range of acceptability. However, the constraints were not met on multiple instances. Therefore, although this is a robust technique, the additional computational time negates the use of this technique. GA-Swarm, due to its hybrid nature, combines the search ability of GA by efficiently exploring the solution space with the fast convergence of PSO for all cycles and OFs. In conclusion, GA-Swarm performs better than GA for all cycles and both PSO and GA-Swarm provide a set of optimum results for each objective function for all cycles.
Suggestions for future work would include running the optimizer using different propellant pairs to check if GG is still the more feasible cycle for this mission. This model can be modified to model a pressure-fed LPRE cycle and to incorporate multiple stage LVs. A database can be developed by comparing various LPRE cycles against each other for different missions. In addition, the application of GA-Swarm for a variety of engineering problems is also recommended in order to further estabish its robustness.