Chinese Journal of Aeronautics 24 (2011) 150-163

Contents lists available at ScienceDirect

Chinese Journal of Aeronautics

journal homepage: www.elsevier.com/locate/cja

JOURNAL OF

AERONAUTICS

Hyper Heuristic Approach for Design and Optimization of Satellite

Launch Vehicle

Amer Farhan RAFIQUEa *, HE Linshua, Ali KAMRANb, Qasim ZEESHANa

aDepartment of Spacecraft Technology, School of Astronautics, Beihang University, Beijing 100191, China bDepartment of Space Propulsion, School of Astronautics, Beihang University, Beijing 100191, China Received 6 January 2010; revised 10 March 2010; accepted 25 September 2010

Abstract

Satellite launch vehicle lies at the cross-road of multiple challenging technologies and its design and optimization present a typical example of multidisciplinary design and optimization (MDO) process. The complexity of problem demands highly efficient and effective algorithm that can optimize the design. Hyper heuristic approach (HHA) based on meta-heuristics is applied to the optimization of air launched satellite launch vehicle (ASLV). A non-learning random function (NLRF) is proposed to control low-level meta-heuristics (LLMHs) that increases certainty of global solution, an essential ingredient required in product conceptual design phase of aerospace systems. Comprehensive empirical study is performed to evaluate the performance advantages of proposed approach over popular non-gradient based optimization methods. Design of ASLV encompasses aerodynamics, propulsion, structure, stages layout, mass distribution, and trajectory modules connected by multidisciplinary feasible design approach. This approach formulates explicit system-level goals and then forwards the design optimization process entirely over to optimizer. This distinctive approach for launch vehicle system design relieves engineers from tedious, iterative task and enables them to improve their component level models. Mass is an impetus on vehicle performance and cost, and so it is considered as the core of vehicle design process. Therefore, gross launch mass is to be minimized in HHA.

Keywords: multidisciplinary design optimization; satellite launch vehicle; heuristic optimization methods; hyper heuristic; air launched vehicles

1. Introduction

Multidisciplinary design optimization (MDO) is an emerging field in aerospace engineering that attempts to introduce a structured methodology to locate the best possible design in a multidisciplinary environment. In fact, MDO can be considered to be a discipline in and of itself with an intention of acting as an agent to bind the other disciplines together^. MDO methods can be from the gradient based class of methods, non-gradient based methods, heuristics, pa-

Corresponding author. Tel.: +86-13811123540. E-mail address: afrafique@yahoo.com

1000-9361/$ - see front matter © 2011Elsevier Ltd. All rights reserved. doi: 10.1016/S1000-9361(11)60019-8

rametric methods and so on. Each breed of methods has its strengths and weaknesses, and each is suitable for different types of problems[2].

In MDO the best choice of optimization method or combination of methods depends on the features of the given multidisciplinary problem. These features include heterogeneous mixes of analysis codes, discrete design parameter values; non-differentiable functions, large number of design variables and a large number of design constraints. In addition, it is dependent on the organization of data communication and execution paths of various interconnected disciplines. The subject of launch vehicle design and optimization is quite vast. Its application varies from small Scout launch vehicle to the current heavy launch system like Ariane. Many of the current state of the art approaches, in

launch vehicle design, employed basic heuristics, specially tailored meta-heuristics and hybrid heuris-tics[3-11]. MDO reduces the time required to execute the design process. In modern era, it is not only required to design a system which fits a customer's needs, but also required to deliver an optimized system. By using MDO methods, designers may quickly and efficiently handle alternative design points over a wide range of parameters[12].

MDO methods have evolved with computing power starting from conventional techniques to modern heuristics. Increase in computing power has limited the use of conventional methods as they are dominated by reformatting, transforming, and translating of data between design disciplines and analysis modules. This will lead to inferior solution. Meta-heuristics prove their superiority over conventional methods but cannot adopt universal settings for different design scenarios. The no free lunch theorem[13] ended the supremacy of individual heuristics. Present day demand of variable scale design necessitates the need of algorithms independent of the solution domain. Hyper heuristic approach (HHA) seems to be the way forward. Another motivating consideration for development of HHA is the involvement of the large number of interdisciplinary design variables in conceptual design of launch vehicles. As the number of variables increases, the volume of the design space increases exponentially. This increases the cost to optimize the design over even a small range for each variable. Furthermore, the feasible design space or the objective functions are non-convex and may contain multiple local optima that can trap optimizers and prevent them from locating the best design.

The approach presented in this study uses a non-learning random function (NLRF) applied to low-level meta-heuristics (LLMHs) in series. Genetic algorithm (GA), particle swarm optimization (PSO) and simulated annealing (SA) are the selected LLMHs. The proposed methodology attains the defined objective by providing diversity to population and forcing search direction by injecting feasible solutions. This, however, may lead to additional computational time but at the same time will serve as a tool to avoid local minima. The goal is not necessarily to "beat" previously employed methods but to obtain satisfactory results by employing a generalized method that can "solve" different problem scenarios with no or limited tuning of the algorithm.

This research effort is structured in three main sections. The first section covers HHA developed for launch vehicle system design. The second section presents assessment and validation of HHA through most widely used test functions. The third section explains multidisciplinary design problem of air launched satellite launch vehicle (ASLV). This section also covers subsystems of mass distribution, propulsion, aerodynamics and trajectory. ASLV design problem formulation and performance results are also presented in this

section. The last section covers brief discussion and conclusion.

2. Hyper Heuristic Optimization Approach

2.1. Need of HHA

Practical engineering problems cannot be tackled by exhaustive search, either to find an optimal solution or even to find a remarkably decent quality solution. The no free lunch theorem[13] proved that there is no one algorithm that could beat all other algorithms in all classes of problems. If an algorithm outperforms other algorithms on a specific class of problems, there must exist another class of problems on which this algorithm is worse than the others. Hence, a good way to raise the generality of heuristics is to apply different heuristics at different times of the search[14]. In this context, a generalized approach (termed hyper heuristics) is proposed[15],which broadly describes the process of using heuristics to choose heuristics to solve problems[16-17].

Although there has been a momentous research done in the field of launch vehicle design and optimization, it is still unclear what kind of search method should be employed that is good enough to be applied to the vast range of aerospace problems. There have been instances of applying individual heuristics, hybrid heuristics, probabilistic methods and even gradient methods in rocket based vehicle design and optimization. Therefore, in current research effort we aim to move a step further and apply HHA. HHA, presented in this research effort, selects optimizers of GA, swarm intelligence and SA to perform simultaneous search. Some key reasons for adopting HHA are as follows:

(1) The conceptual design of complex aerospace systems involves optimization with a large number of design variables involving multiple disciplines. This results in increase of volume of the design space exponentially and makes solution space critically sensitive of various mission scenarios.

(2) The feasible design space of large problems is often non-convex and may contain multiple local optima. This can trap optimizers and prevent them from locating the best design.

(3) The design of rocket based systems entails a large number of nonlinear constraints that also result in a non-convex feasible design space and numerous local optima.

(4) Heuristic optimization methods require significant parameter tuning of algorithms for the new problem or the new problem instance. Variation in design scenarios can lead to exhaustive search for the new problem instance.

(5) The aim is to devise an algorithm for solving a complex and multidisciplinary problem that is independent of problem scenario, reasonably comprehensible, trustable in terms of quality and repeatability and

■ 152 ■

with good behavior and perform exceptional search for global optima and avoid local minima.

2.2. Working of HHA

Working of HHA starts from setting the control function. The choice of control function can be critical in setting HHA. The bias of the control function can affect the performance in solution quality and computational cost, leading to selection of probabilistic approach. Statistically, random sampling has advantages; it produces unbiased estimates of the mean and the variance of the output variables. The choice of random function is attractive, from the perspective of increasing generality and avoiding bias of control function.

Mathematically, (Q, &, 6) is a probabilistic space and (9, fl) is a measurable observation, then random variable X is defined as

X: n^ 9 (1)

In the present case, the random variable is real value, thus

X: n^ R (2)

Such that

X(a< r)]e© Vr e R (3)

The approach presented in this study uses an NLRF applied to LLMHs in series. The working of the algorithm can be defined in the steps below and followed by the pseudo code (see Fig.1).

(1) A set of input defining objective, constraint function, design space specification (bounds of decision variables) and convergence criteria are established.

(2) Set total number of outer loop and inner loop iteration and total number of LLMHs available.

(3) Set bounds of operators in LLMHs; crossover and mutation for GA, inertia and acceleration factors for PSO, initial temperature and cooling schedule for SA.

(4) Outer iterations refer to the total number of calls to all LLMHs and inner iterations are the number of function evaluations done by individual LLMHs in one outer loop iteration.

(5) Random function divides the number of inner iterations and total population into three parts for three LLMHs used in this study. This division is done randomly while remaining within the limits of total number of inner iterations and total size of population. Based on this division each LLMH is called in series for set number of inner iterations to evaluate objective function and constraint violation.

(6) Call to each LLMH is stochastic (random) and not deterministic.

(7) Certain solutions obtained by each LLMH are stored as local and global best and worst.

(8) The above mentioned process is repeated till one outer loop iteration is done.

(9) After each outer loop iteration, new population is

created by each LLMH with values of mutation and crossover for GA, inertia and acceleration factors for PSO, initial temperature and cooling schedule for SA, are randomly chosen from bounds already provided. These values of operators are kept fixed for particular inner loop.

(10) Random function re-distributes the new population and inner iterations on random basis to each LLMH.

(11) Local and global best solutions are injected to new population after certain inner and outer iteration (on random basis) to diversify solution.

(12) The process is iterated till convergence is achieved.

• Set population sire for each LLMH

• Set lolal number or outer loop iterations

• Set total number of inner loop iterations

While (convergence NOT achieved)

• Create public-board to store information

• Generate sub populations for LLMH (random)

For i - I to outer loop iteration I OUTER LOOP] Call INDIVIDUAL HEURISTICS

CALL OFF-SPRINGS FUNCTIONS Send information to public-board

CALL POPULA TION HANDLING

Shufflem% of new population given to n optimizers Generate new sub population (NLRF) Generate new sub iterations (NLRF)

CALL POPULA TION DIVERSIFiCA TION (NLRF) Insert local best to sub population (random )

Insert global best to sub population (random )

Fig.1 Pseudo code of program.

The following is the brief description of the LLMHs selected for HHA. Ref.[18] developed GAs, which are capable of finding the global-optimal solution (or acutely near solutions) in complex, multidimensional search spaces. Theoretical details on GA can be found in Refs.[18]-[19] and is widely applied to aerospace system[3-4, 8 20] and subsystem[21] design. Swarm intelligence is a new realm of research in that the PSO technique takes inspiration from social behavior of insects and animals[22-23]. PSO has remarkable adaptability to continuous problems so it is applied to component level design and optimization of satellites[24], multidisciplinary optimization of ASLV[8, 11], and certain structural, multidisciplinary optimization prob-lems[25]. SA was originally proposed by Metropolis[ ] in the early 1950s as a model of the crystallization process. It was only in the 1980s that independent re-search[27-28] noted similarities between physical process of annealing and some combinatorial optimization problems. SA has proved to be a potent, stochastic search method and is applicable to a wide range of

engineering problems[29-32].

3. Experimental Setup

3.1. Benchmark test functions

Comprehensive empirical analysis is performed to investigate the performance of HHA. This analysis also evaluates performance of HHA with other popular non-gradient based optimization techniques. The following is the brief summary and mathematical description of selected test functions.

(1) The de Jong function[33] is continuous, convex, uni-modal and one of the simplest benchmark tests. Mathematical description is given by

f ( x) = £ xf

(-100 < x < 100, n = 30)

(2) The 2D six-hump camel back function[34] is a global optimization test function. There are six local minima and two global minima of this function. Mathematical description is given by

f (x) = x2(4 - 2.1x2 + x4/3) + xx2 + 4x4 ~ 4x2

(-3 < x < 3) (5)

(3) The Ackley problem[35] has several local minima but only one global minimum. It is a widely used multi-modal test function. Mathematical description is given by

-0.2 1¿xf 1cos(2xx; ) Vn.-0 _ en,.0

f ( x) = 20 + e - 20e

(-30 < x < 30, n = 10)

(4) Schwefel's function[36] is deceptive in that the global minimum is geometrically distant, over the parameter space, from the next best local minima. The surface of this function comprises a great number of peaks and valleys. Therefore, the search algorithms are potentially vulnerable to convergence in the wrong direction. This problem consists of numerous local minima and requires a large number of function evaluations to reach the global minimum. Mathematical description is given by

n , _,

f ( x)=S(-x )

(-500 < x < 500, n = 2)

(5) Rosenbrock function[37] is a classic optimization problem, also known as banana function or the second function of de Jong. The global minimum is present inside a long, narrow, parabolic-shaped flat valley. To locate the valley is trivial, however, convergence to the global optimum is difficult and hence it has been frequently used to assess the performance of optimization algorithms. Mathematical description is given by

f (x) = £ [(1 - x2) +100( xI+1 - x2)2] i=1

(-30 < x < 30, n = 20) (8)

(6) Rastrigin's function[38] is often used to evaluate global optimizers. This function is comparatively a difficult problem because of its large search space and the large number of local minima. The function is highly multi-modal, and the locations of the minima are regularly distributed. Mathematical description is given by

f (x) = £ [x2 -10 cos(2:rc xt) +10]

(-5.12 < x < 5.12,n = 20) (9)

(7) The Griewank function[39] has been widely employed as a test function for global optimization algorithms. The number of minima grows exponentially with increase in dimensions. Therefore, the global minimum becomes extremely difficult to detect. Mathematical description is given by

f ( x) =

—ixl-Г7 cos i^r

4 000 j! 1 H U (-600 < x < 600, n = 10)

3.2. Results of benchmark test functions

In experiments, population size is set at 100, and all test functions are executed for 100 times with the same algorithm parameter settings. In order to investigate the quality of the algorithm, a fair measure must be selected. The simulation is considered successful if the solution is less than |4 x10~41 of the true global minimum. Table 1 reveals the obtained results.

Table 1 Results of benchmark test functions obtained by using HHA

Function Success Mean Standard deviation Averag function

rate/% evaluation

de Jong (Sphere) 100 1.99x10~5 5.490 1> 10 -5 2 953

Six-hump camel back 99 -1.031 6 2.166 4> 10 -7 7 454

Ackley 98 7.890x 10~6 4.660 0> 10 -3 42 208

Schwefel 99 -4.189x102 3.390 0> 10 -10 7 551

Rosenbrock 100 1.180x10~6 4.880 0> 10 -6 11 912

Rastrigin 98 7.890x 10~6 1.400 0> 10 5 28 728

Griewank 98 2.210x 10~6 5.340 0> 10 5 32 028

Results reveal that HHA performed magnificently for solving these classical test functions. On the robustness of algorithms, HHA is significantly capable of finding the global optimum with outstanding success rate on all the test functions. HHA proves to be a highly stable algorithm, capable of obtaining consistent results in terms of accuracy and robustness.

3.3. Comparison with other methods

Performance of proposed HHA is compared to validate the superior performance over other methods available in literature. Table 2 presents the comparative summary of HHA with other methods with respect to success rate. There are other criteria according to which the performance of HHA is to be measured against widely used optimizers like GA, PSO and SA. Table 3 presents such a comparison.

Table 2 Comparison of HHA with other methods (success rate)

Function

de Jong (Sphere) Six-hump camel back Ackley Schwefel Rosenbrock

Rastrigin

Griewank

Success rate/%

HHA GA PSO SA

100 95 90 92

99 90 86 93

98 89 85 80

99 96 95 97

100 97 91 94

98 88 81 84

98 88 80 82

Table 3 Comparison of HHA with other methods (other criterion)

Criterion HHA GA PSO SA

Flexibility Best Good Worst Fair

Robustness Best Good Worst Fair

Ease of control Best Good Best Fair

Problem dependency Best Good Fair Fair

Simplicity Fair Fair Best Good

Computational cost Fair Fair Best Good

Adaptability Best Good Worst Fair

optimum. Table 3 shows that overall capability of HHA for solving complex problems is far better than widely applied heuristic optimization methods. Prime advantages of HHA include not being problem dependent, easy to control optimization parameters, robustness of solutions, flexible enough to apply to various kinds of problems and being adaptive in nature. The disadvantage of using HHA is the computational cost which can be neglected in favor of high quality solution and can also be reduced by improving computational resources.

4. Test Case (ASLV)

Recent trends within the global aerospace community require improvements to the satellite launch vehicle design process. The design of ASLV, capable of inserting specific payload to low Earth orbit (LEO), is a complex problem that must balance competing objectives and constraints. This segmented design process is iterative in nature and invariably leads to design compromises as system-level engineer's effort to make each component of the entire launch vehicle system compatible with each other. Therefore, the need arises for MDO and the use of artificial intelligence learning tool that can control the design of different components simultaneously. With MDO process, the design engineer can set explicit system-level goals and then turn the design optimization process entirely over to optimizer. This distinctive approach for launch vehicle system design relieves engineers from tedious, iterative task and enables them to improve their component level models.

4.1. Multidisciplinary design analysis of ASLV

ASLV design is a highly integrated process requiring synergistic compromise and tradeoffs of many parameters. The synthesis of an effective compromise requires balanced emphasis in subsystems, unbiased tradeoffs, and evaluation of many alternatives. Fig. 2 outlines the integrated disciplines considered in mul-

Table 3 reveals the exceptional superiority of HHA over GA, PSO and SA while searching for the global

Fig.2 Multidisciplinary design analysis.

tidisciplinary design analysis of ASLV and flow of information among these disciplines.

(1) Vehicle definition

The baseline design is launched from 40 000 ft (1 ft = 0.304 8 m) at Ma=0.8. This is intended as a representative number taken from launch conditions of similar vehicle[40]. The mission of ASLV is to deliver 400 kg payload (satellite) to LEO (500 km). Propulsion system is solid-fuelled solid rocket motor (SRM). In this research proposed ASLV is fixed as three-stage vehicle. The payload (satellite) is enclosed in fairing whose shape is known beforehand. The payload weight and volume requirements are specified in the problem definition before the optimization is computed. Payload fairing is 2 m long and with the same diameter as that of the first and second stages. Limits and bounds of certain flight performance parameters, for instance coast time, maximum angle of attack, axial overload coefficient (nx) and lateral overload coefficient (ny), are learnt from experience of similar and

other operational launch vehicles[40-41].

(2) Propulsion module

The design of next generation launch vehicles and space systems is dependent more than ever on high performance, increasingly efficient, reliable and affordable propulsion systems. The goal is to minimize the inert mass to impulse ratio. Ballistic and structural variables are strongly inter-linked and certain compromises are to be made to reach the optimal design configuration. Pressure, burn time, and nozzle area ratio dictate the inert mass of propulsion system and hence need to be optimized.

Propulsion analysis describes fundamental parameters like thrust, burn time, mass flow rate and nozzle parameters[42]. Chamber pressure is a vital design variable and has momentous influence on SRM specific impulse, burning rate of the propellant, size of expansion nozzle and thickness of casing materials to withstand the pressure stresses. In this analysis, we does not restrict to a particular shape of grain at the conceptual design level, rather a variable grain shape factor (ksi) is used to represent burning surface area (Sri) of grain as a function of length of grain (Li) and diameter of grain (Di) given by

ks i = Sr 1/D1 (11)

Mass of grain is calculated from design variables and propulsion analysis. Burn time (tbi), mass of ith stage grain (mgni) and mass flow rate ( mgni ) of ith

stage grain are calculated as Eqs.(12)-(14) respec-

tively[42-43].

nrlviDi

tu =■

grain, Agni fineness ratio of grain (grain length/diameter), and volumetric loading fraction. Nozzle throat area (At), expansion ratio (s) and nozzle exit area (Ae) are calculated as

PgnuAr

At =■

or c,max

1 r A l7 2r 1 - 7-1 ' f P V

l Pc J ]j r-1 l Pc J

Ae = AS

2 12ÔZÏÏ

where Ab is burning area of grain, Rc the gas constant, Tc the temperature in the combustion chamber, pe the exit pressure, pc the chamber pressure, and y the specific heat ratio of gas. Refs. [42]-[43] provided the calculation of vacuum specific impulse (/svpac) and thrust (F) as

F = C mgn, - Pa Ae , (i = 1,2, •••, N)

rvac _ 7-a I

r R;T;

gl Isap

where N is the number of stages, pa atmospheric pressure, /sap the average specific impulse, and g0 the acceleration due to gravity. (3) Mass module

Using a combination of physics-based methods and empirical data, the mass of the foremost components for the solid stages is determined by Ref.[43]. The total mass of multistage ASLV includes the masses of pro-pellants and their tanks, related structures and payload mass. The gross launch mass (GLM) (mM) for a multistage ASLV is expressed as

m01 = mPAY + X (mgni + msti + mai + mii) (21) i=1

where mPAY is mass of the payload, msti mass of the ith stage SRM structure, mai total mass of the ith stage aft skirt, mfi total mass of the ith stage forward skirt. Design mission dictates mass of the payload, 400 kg in this case. Total mass of the ith stage aft and forward skirt are simplified as

m = m„„. + m„

mfi = mf„ + m,

mai + mfi = N,m0,

mgn, = - Pgmks,\n,Df mgrn = PgnU,Sr, = PgmU1ks1\mD!

"gni /^gni"^ ri /^gni"T "sr "gm^i (14)

where u is burning rate of propellant, pgni density of

where msvi is mass of control system, safety self-destruction system, servo, and cables inside the ith stage aft skirt; masi mass of the ith stage aft skirt including shell structure, equipment rack, heating protector structure, and direct subordinate parts for integration;

mfei mass of equipment and cables inside the ith stage forward skirt; mfsi mass of the ith stage forward skirt including shell structure, equipment rack, and direct subordinate parts for integration; m0i mass of the ith stage SRM. Skirt mass ratio or structural coefficient N has small dispersions which can be selected from the statistical data[43].

The simplified form of N-stage launch vehicle mass equation comes to be

"0; - N

n[1 - N -Kgn,^,(1 + «.„■)]

where Kgni is propellant reserve coefficient and its values can be selected from statistical data. Relative mass coefficients of effective grain and structure mass fraction asti are given by

M» = (26)

®.ti = — (27)

where mgeni is the mass of effective grain.

The main parameter for designing a multistage ASLV is structure mass fraction. It is dependent upon structural material, grain shape, as well as the parameters of internal ballistics of SRM.

The mass of the ith stage SRM structure (see Eq.(28)), comprises mass of the motor cylinder (mcyi), motor dome ends (mcii and mc2i), forward and aft skirts (mqi), forward and aft attachments (mj1i and mj2i), forward and aft insulation liners (min,c1 i and min,c2i), cylindrical section insulation liners (min,cyi), nozzle expansion cone (mnoz,eci), nozzle spherical head (mnoz,shi), nozzle insulation (mnoz,ini), igniter (migi), thrust vector control (mTVCi), cables (mcabi) and attachment parts (mapi). These component masses are presented in Eqs.(29)-(47). Mass of structure is used to calculate osti and so GLM of ASLV

msti = mcy, + mc1i + mc2i + mqi + mj1i + mj2i + min,c1, + min,c2i + min,cy, +

m ■=■

«TVC. + mcab, + map, 3Kcy * ffDl pX

where Kcy is ratio of the cylindrical length to grain length which can be obtained from statistics, f factor of safety, fp relative pressure in chamber (f,=pc,max /pc),

Dch diameter of chamber, and crb ultimate strength of chamber material.

M2 Dc3h ffpPc

8(2e2 - 1)acos e2

1 + (2e2 - 1)sin26>2

where Âe is ellipsoid ratio, a strength ratio (a=ob/pcl,

where pcl is the density of closure material), and 0 is normally 60°-65°.

M2 D3h ffp Pc

mc2; =-2-—

c2; 8(le2 - 1)acos2ft

42Dc2h - (42 - 1)d2 42 Dc2h

1 + (le2 - 1)sin2<92

where d is diameter of rear end opening for nozzle.

r,2 c lqi + lq2 mq i = ^DchAA

where pq is density of skirt material, ¿>q thickness of skirts, and lqi and lq2 are lengths of forward and aft skirts.

mj1, = ■

-Kp, {rfSXl - yf(ôlt +^2, ) + 1 [(2Ä02 + r2)V b2 - r2 - (2bl + Rf)4 bl - R2 ]} (33)

mj2; = mj1;

where i = 1 for forward metal joint, and i = 2 for aft metal joint, pj is density of joint material, Ri radius of case cylindrical section; S1i, S2i, yi and r are explained in Fig.3, and

bo =-f -^cy (35)

where ¿>cy thickness of case cylindrical section, and Dcy diameter of case cylindrical section.

Fig.3 Forward and aft metal joint.

mrn,c1, = TAnDRatb

An ^DcyRa tbl

ln,c2; 4(2e2 -1)

42Dc2y - (^2 - 1)^2

where pm is density of insulation material, and Ra rate of ablation.

(^in CmPm +agi

' cy^cy f-'m.

min,cy; = nDcy LcyPm

«gCin Pie

L 2ein (ccypcy3cy ln + agltb )a

1 (^inCinAn + agl

OT_, =

where Lcy is length of case cylindrical section, sin heat transfer coefficient of insulation material, cm specific heat capacity of insulation, agi heat transfer coefficient from gas to insulation, ccy specific heat capacity of cylindrical section, pcy density of case cylindrical section material, Tg temperature of gas, Tcy allowable temperature of cylindrical section and T initial temperature of cylindrical section.

mnoz,ec = Pec , ■ „ (~T ~ 1) '

4sin^n At 0.67 f jd- -1(1 - S^Ëf

where pec is material density of the expansion cone, pu nozzle expansion half angle, dt diameter of throat, de nozzle exit diameter, S submerged coefficient of nozzle, crec ultimate strength of expansion cone material and Eec elastic modulus of expansion cone material.

!,shi = 3.656Ah dt3

where psh is material density of spherical head of nozzle.

noz,ini Pin,nz Snz^in,nz

where pin,nz is material density of the nozzle insulation, Snz surface area of nozzle and 8,„„z thickness of the

nozzle insulation.

mig i =1.454

f 7< ^

mTVCi = (1.3-1.7) mn

= (1 + 3 Lcy, ) Pei (46)

where Lcyi is the ith stage SRM length and pci the ith stage cable density.

mapi = (6.13 x10-7)D,2 L^48 (47)

(4) Aerodynamics module

Speedy and cost-effective estimations of aerodynamic stability and control characteristics are frequently required in preliminary design phase of ASLV. Thus, a need arises for use of time-efficient computer software that can predict aerodynamic properties over a range of flight conditions. For this purpose, U.S. Air Force Missile DATCOM 1997 (digital)[44] has been widely used in the aerospace industry. DATCOM is capable of quickly and economically estimating the aerodynamics of a wide variety of design configuration and in different flow field regions that the ASLV encounters during atmospheric flight. DATCOM is used in the present research for estimation of coefficient of lift (CL) and coefficient of drag (CD).

(5) Trajectory module

Three degree-of-freedom (DOF) approach is time efficient and diminishes the requirements of a large number of input variables calculated with lesser accuracy due to unavailability of detailed shape of the

structure and other aerodynamic parameters in conceptual design phase. Therefore, in this research effort 3DOF trajectory analysis is preferred over 6DOF mod-

el[45-46]

3DOF model is developed and simulated in SIMU-LINK to analyze the flight path. Trajectory simulation obtained from 3DOF model is computationally efficient and serves the purpose at the conceptual design level. The trajectory analysis depends on inputs from the aerodynamic, mass and propulsion modules. The flight program and results obtained from other disciplines compute the flight trajectory. This investigation treats ASLV as a point-mass and flight in 2D over a spherical and non-rotating earth is assumed. This implies that the Coriolis and centrifugal pseudo forces are negligible. Fig.4 shows forces acting on ASLV and Eq.(47) presents the set of governing equations of motion.

Fig.4 Forces acting on ASLV. dv F cosa- D

--go sin$

d5 _ Fsina+ L g0cos5+f v ^

dt dh dt dl

mv = v sin5

V Re + h y

V Re + h y

v cos 3

a = r/ + p-3

1= — Re

a = «pro(t ) L = Cl 2 pv2 Sref

D = CD 2 Pv Sref

where v is velocity, D the drag, L the lift, m mass of

vehicle, a angle of attack (AOA), cp flight path angle, 6 trajectory angle, r) range angle, 3 trajectory angle in local frame, apro(t) programmed angle of attack, Re radius of Earth, h height about ground, l range, and Sref surface area.

The powered flight trajectory of ASLV is composed of several distinct phases. After a short, horizontal launch phase, the ASLV undergoes a launch maneuver by introducing a carefully selected AOA profile. The flight program is optimized under the constraints of the maximum allowed AOA (amax), lateral overloads and axial overloads. The lateral and axial overloads are limited to ensure the structural integrity.

Each stage shuts down one after another and separates to shed inert mass. Because the equations of motion change discontinuously at the shutdown points, the trajectory must be divided into intervals, the number of which corresponds to number of stages. The phases are explained as follows.

1) Horizontal launch phase

The launch sequence of proposed ASLV begins with the release from the carrier aircraft at an altitude of approximately 40 000 ft (12 195 m) at Ma=0.80. Stage 1 ignition occurs approximately 5 s after drop when the proposed ASLV has cleared carrier aircraft. An initial condition is that the flight path angle is 0°. Aerodynamic lifting surfaces (wings, horizontal and vertical tail) provide stability and lift after launch and further within atmospheric flight.

2) Powered ascent phase

This phase starts from the end of transonic phase (0.8<Ma<1.3) up to the shutdown of second-stage motor. AOA must be 0° when vehicle is in transonic regime. This is because the aerodynamic forces vary acutely in transonic speeds and large AOA would yield strong disturbance to control system and servo overload to the body of ASLV. Furthermore, drag loss will be more in the transonic region. AOA must approach to 0° at the time of the first-stage separation. However, there is no such restriction on separation of upper stages.

Flight program in rocketry refers to the prescribed variation of vehicle pitching angle during the flight. This has substantial influence on reachable altitude, orbit injection accuracy, aerodynamic loading and heating etc. Flight AOA is programmed using the following relations:

apro(t ) = a

f (t ) =

sin2 f (t )

rc (t ~ ti) k (t2 -1) + (t -11) L -1

t2 - tm

where am is launch maneuver variable, t time of flight, t1 time of start of power ascent phase (turning), t2 time at end of turning and tm time corresponding to maximum AOA. Fig.5 explains the powered ascent phase of ASLV.

Fig.5 Powered ascent phase.

Axial overload coefficient should not exceed the allowable range. Otherwise, the control system, the strength of grain and the bond strength between grain and motor case will be destroyed. Lateral overload coefficient should also not exceed the allowable range. Constraints on these overloads are enforced in the optimization problem. The axial and lateral overload coefficients are calculated as

F + L sin a- D cos a n =-< «

ny =■

mgo F sina + L

■< n,,

3) Coasting phase

The final stage carrying the payload separates and flies up with no thrust in an elliptical orbit at the end of the 2nd stage.

4) Kick phase

Finally, an apogee kick puts the payload (satellite) into a circular LEO. The constraints in this phase are applied to flight path angle, dynamic pressure, normal acceleration and body axial acceleration.

4.2. ASLV optimization problem formulation

ASLV design and optimization problem is formulated by defining design objective, variables and design constraints.

(1) Design objective

ASLV design can pose different objective functions for the optimization problem depending upon mission requirements. For example, one could minimize cost, maximize payload for a fixed launch weight, maximize injection accuracy in orbit, and minimize launch weight for placing a fixed payload in a particular orbit. Traditionally, minimum take-off weight concepts have been sought in launch vehicle design. This is because weight (or mass) is a strong driver on vehicle performance and cost, and so plays a central role in vehicle design process. Design objective for ASLV optimization problem is to minimize the GLM of the entire vehicle to inject a payload (400 kg) into LEO. Mathematical description of design objective is provided by min GLM(X) s.t.

gj (X) < 0

hk (X) = 0

Xlb ^ X, < Xub (54)

where X is design variables, Xib lower bound of design variables and Xub upper bound of design variables. From Eq.(25), GLM of ASLV can be written as

m0i = m0, (N, mPAY , Kgn, , Mki ,«sa) (55)

By exploring asti and introducing propulsion variables (pci, pei and u,), trajectory variables (omax and am) and removing independent variables (N, mPAY and Kgni), we can write the set of design variables (X) as

X = f , D,, Pc,, Pe,, U,, ks,, «max , On ) (56)

(2) Design variables

Table 4 lists the system design variables for each stage with respective discipline. There are a total of 19 design variables that govern the integrated design and optimization problem of ASLV.

Table 4 Discipline-wise design variables

Discipline Design variable Symbol

Structure Propulsion Relative mass coefficient of grain (Stages 1-3) №

Structure Propulsion Aerodynamics Body diameter (Stages 1-3) Di/m

Structure Propulsion Chamber pressure (Stages 1-3) Pci/bar

Structure Propulsion Trajectory Exit pressure (Stages 1-3) Pei/bar

Structure Propulsion Coefficient of grain shape (Stages 1-3) ksi

Propulsion Trajectory Grain burning rate (Stages 1-3) ui/(mm-s_1)

Aerodynamics Trajectory Max AOA (Stage 1) «max/(°)

Aerodynamics Trajectory Launch maneuver variable (Stage 1) am

(3) System constraints

Final mission velocity (vf) and corresponding altitude (Hf) are trajectory constraints. The overall structure of the system should be extremely strong to survive the high g-loads. Therefore, an axial overload constraint is implemented to restrict it below 12g for the 1st and 2nd stages, and lateral overload to restrict below 2g for the 1st stage. During launch maneuver, the maximum AOA is constrained to be below 22° and ensured that it is 0° during transonic phase and at first

stage separation. Flight path angle is constrained to be 0° at orbit insertion (the 3rd stage). Nozzle exit diameters are constrained to be less than stage diameters. The diameters of the 1st and 2nd stages are con-

strained to be equal. If any of these conflicts occur, the program is set to send back extremely poor performance values in each goal area so that it will learn not to evaluate these designs in the future. Constraints are formulated as

Ct > 0 (i = 1,2, — ,9)1 Ct = 0 (i = 10,11) J

where C{ is given as

Cl = vf > 7 600 m/s

C2 = Hf > 450 km

C3 = nxi * 12 g

C4 = nx 2 ^ 12 g

C5 = ny1 < 2 g

C6 = «max ^ 22°

C7 = De2 < D

C8 = De3 < D1

C9 = D1 = D2

C = a= 0° (0.8 <Ma < 1.3)

C11 = ф = 0° (At orbit insertion)

Constraints are included in the objective function and are handled by dynamic penalty function. Eq.(58) represents the symbolic problem statement.

min GLM(X) = f (X) + h(k)£max{0, gt (X)} (58)

where f (X) is the fitness function, h(k) modified penalty value and k is the current iteration number of the algorithm. The function gi(X) is a relative violated function of the constraints[47].

4.3. Optimization parameters

Table 5 explains the optimization parameters for optimization problem of ASLV using HHA.

Table 5 Optimization parameters for HHA

Entity Control parameter Value

Population size 100

Outer loop iteration 50

Inner loop iteration 120

NLRF Sub-population size for LLMH Random

Sub-iterations for LLMH Random

Maximum function evaluation 100 000

Function tolerance 10~6

GA Crossover rate Two point (0.6-0.8) random

Mutation rate Uniform (0.1-0.3) random

Acceleration constant, c1 (1.8-2.2) random

PSO Acceleration constant, c2 (1.8-2.2) random ЦС1,С2) = 4

Initial inertia weight (0.7-0.9) random

Final inertia weight (0.3-0.5) random

SA Initial temperature/°C 90-100

Cooling schedule r. /-v(Iteration) т 0.9

Note: T0 is the start time for SA.

4.4. Results of ASLV test case

Multidisciplinary design and optimization using HHA is successfully implemented for complex design problem of ASLV under stringent mission objectives and performance constraints. MDO and HHA prove to be more effective in terms of exploring the design space and fulfilling design objectives of ASLV design problem.

Fig.6 illustrates the performance graphs of optimized configurations of ASLV. The required orbit injection velocity (> 7 600 m/s) and altitude (= 500 km) are perfectly achieved. Axial overloads are within the allowable ranges and required thrust is achieved to complete the mission of ASLV.

Fig.6 Performance graphs of ASLV

Table 6 contains optimal values for all the design variables. Optimized design variables of design space of ASLV lie between respective upper and lower bounds. Results reveal that ASLV is capable of completing the specified mission through hyper heuristic-based MDO approach. Table 7 demonstrates the stage layout and propulsion parameters of optimal configurations of ASLV. Table 8 presents the characteristics of launch vehicle optimization problem.

The computational efficiency can never be ignored in conceptual design and optimization of any aerospace system. Conceptual design phase requires that various configurations can be tested in a remarkably short span of time. The proposed 3DOF simulations give ASLV optimization routine the computational

Table 6 Optimal results of design variables

Design variable Symbol Lower bound Upper bound Optimal solution

m 0.60 0.80 0.70

Relative mass coefficient of grain №2 /«il 0.90 1.00 0.95

m /rn 0.90 1.00 0.95

Body diameter D1/m D^m 1.25 0.80 1.40 0.90 1.30 0.82

pc1/bar 65.0 80.0 70.00

Chamber pressure pc2/bar 65.0 80.0 72.29

pC3/bar 40.0 60.0 44.84

pe1/bar 0.12 0.16 0.15

Exit pressure pe2/bar 0.08 0.16 0.13

Pe3/bar 0.05 0.10 0.08

u1/(mm-s_ 1) 5.00 8.00 5.10

Grain burning rate u2/(mms 1) 5.00 8.00 6.50

u3/(mm-s 1) 4.00 6.00 4.86

is1 1.50 2.30 2.05

Coefficient of grain shape is2 1.50 2.30 2.12

is3 1.50 2.30 1.90

Max AOA «max/(°) 1.00 22.0 21.98

Launch maneuver variable aj(°) 0.01 0.10 0.012

Table 7 Optimal configuration

Optimized configuration

GLM: 22 310 kg Payload: 400 kg

Stage 3

Gross mass: 1 108.27 kg Propellant mass: 952.86 kg Thrust: 50.46 kN; Burn time: 53.21 s Stage 2

Gross mass: 4 222.92 kg Propellant mass: 3 811.24 kg Thrust: 179.00 KN; Burn time: 60.00 s Stage 1

Gross mass: 16 577.74 kg Propellant mass: 15 616.25 kg Thrust: 521.79 kN Burn time: 83.83 s

efficiency required for the conceptual design. 5. Discussion and Conclusions

HHA, which integrates GA, PSO and SA through a random function, is proposed in this paper. To illustrate the efficiency and effectiveness of the proposed method, standard test functions in the field of optimization are solved. The results indicate better performance in terms of both diversity and convergence. The proposed methodology attains the defined objective by providing diversity in the population and alters search direction. Hence, it is proved to be significantly efficient in finding global optima and avoiding being trapped in local minima.

The proposed approach is successfully applied to practical and complex design case study of ASLV This methodology can handle optimization problem in the area of launch vehicle design, where the optimization of several multidisciplinary design variables and considering conflicting constraints is required. The proposed approach tackles a complex and multidiscipli-nary problem that is independent of problem scenario, reasonably comprehensible, trustable in terms of quality and repeatability and with good behavior. It performs exceptional search for global optima and avoids local minima. It can be concluded that the proposed approach has shown enormous potential and can be applied to complex optimization problems in design of launch vehicles and other aerospace vehicles.

The results can be used as a basis for detailed design. The optimization results and performance are to be considered as preliminary (proof-of-concept) only, comparable to the existing systems.

Acknowledgement

Amer Farhan RAFIQUE wishes to thank Higher Education Commission (HEC) of Pakistan for award of scholarship for Ph.D. studies. He would also like to thank Mr. Sadique AHMED for the valuable help in programming.

References

Table 8 Optimization characteristics of launch vehicle [2]

Vehicle description Velocity Minimum launch mass Altitude Objective [3]

Critical constraint Insertion velocity Flight path angle Axial, lateral overloads AOA Chamber pressure Trajectory Trajectory Structure Aerodynamics Propulsion [4]

Optimum design All constraints satisfied GLM = 22 310 kg

Sobieszczanski-Sobieski J, Haftka R T. Multidiscipli-nary aerospace design optimization: survey of recent developments. Journal of Structural Optimization 1997; 14(1): 1-23.

Olds J R. The suitability of selected multidisciplinary design techniques to conceptual aerospace vehicle design. AIAA-1992-4791, 1992.

Bayley D J, Hartfield R J, Burkhalter J E, et al. Design optimization of a space launch vehicle using a genetic algorithm. Journal of Spacecraft and Rockets 2008; 45(4): 733-740.

Rafique A F, He L S, Zeeshan Q, et al. Multidiscipli-nary design and optimization of an air launched satellite launch vehicle using a hybrid heuristic search algorithm. Journal of Engineering Optimization 2011;

43(3): 305-328.

[5] Rafique A F, He L S, Zeeshan Q, et al. Integrated system design of air launched small space launch vehicle using genetic algorithm. AIAA-2009-5506, 2009.

[6] Rafique A F, Zeeshan Q, He L S. Conceptual design of a small satellite launch vehicle using hybrid optimization. In: Sivasundaram S. Seventh International Conference on Mathematical Problems in Engineering, Aerospace and Sciences. Cambridge: Cambridge Scientific Publishers Ltd., 2008; 681-697.

[7] LeMoyne R. Multidisciplinary cost and performance optimization of a two-stage liquid propulsion based launch vehicle. AIAA-2008-2642, 2008.

[8] Rafique A F, He L S, Kamran A, et al. Multidiscipli-nary design of air launched satellite launch vehicle: performance comparison of heuristic optimization methods. Acta Astronautica 2010; 67(7-8): 826-844.

[9] Takeshi T, Takashige M. Optimal conceptual design of two-stage reusable rocket vehicles including trajectory optimization. Journal of Spacecraft and Rockets 2004; 41(5): 770-778.

[10] Jahangir J, Masoud E, Jafar R. Multidisciplinary design optimization of a small solid propellant launch vehicle using system sensitivity analysis. Journal of Structural and Multidisciplinary Optimization 2009; 38(1): 93-100.

[11] Rafique A F, He L S, Zeeshan Q, et al. Multidiscipli-nary design of air-launched satellite launch vehicle using particle swarm optimization. In: Anderssen R S, Braddock R D, Newham L T H. 18th World IMACS Congress and MODSIM09 International Congress on Modelling and Simulation. Canberra: Modelling and Simulation Society of Australia and New Zealand and International Association for Mathematics and Computers in Simulation. 2009; 418-424.

[12] Rowell L F, Korte J J. Launch vehicle design and optimization methods and priority for the advanced engineering environment. NASA TM-2003-212654, 2003.

[13] Wolpert D H, Macready W G. No free lunch theorems for optimization. IEEE Transactions on Evolutionary Computation 1997; 1(1): 67-82.

[14] Burke E, Hart E, Kendall G, et al. Hyper-heuristics: an emerging direction in modern search technology. In: Glover F, Kochenberger G. Handbook of Meta-heuris-tics. Dordrecht: Kluwer, 2003; 457-474.

[15] Denzinger J, Fuchs M, Fuchs M. High performance ATP systems by combining several AI methods. 15th International Joint Conference on Artificial Intelligence. 1997; 102-107.

[16] Cowling P, Soubeiga E. Neighborhood structures for personnel scheduling: a summit meeting scheduling problem. 3rd International Conference on the Practice and Theory of Automated Timetabling. 2000.

[17] Ross P. Hyper-heuristics. In: Burke E, Kendall G. Search Methodologies: Introductory Tutorials in Optimization and Decision Support Techniques. Amsterdam: Springer, 2005; 529-556.

[18] Goldberg D. Genetic algorithms in search, optimization and machine learning. Reading: Addison-Wesley Longman, 1989.

[19] Holland J H. Adaptation in natural and artificial systems: an introductory analysis with applications to biology, control, and artificial intelligence. Cambridge, MA: MIT Press, 1992.

[20] Zeeshan Q, Dong Y F, Nisar K, et al. Multidisciplinary

design and optimization of multistage ground-launched boost phase interceptor using hybrid search algorithm. Chinese Journal of Aeronautics 2010; 23(2): 170-178.

[21] Kamran A, Liang G Z. Design and optimization of 3D radial slot grain configuration. Chinese Journal of Aeronautics 2010; 23(4): 409-414.

[22] Kennedy J, Eberhart R. Particle swarm optimization. IEEE Proceedings International Conference on Neural Networks. 1995; 1942-1945.

[23] Kennedy J, Eberhart R C. Swarm intelligence. San Mateo, CA: Morgan Kaufmann, 2001.

[24] Hassan R, Cohanim B, Weck O D, et al. A comparison of particle swarm optimization and the genetic algorithm. AIAA-2005-1897, 2005.

[25] Venter G, Sobieszczanski-Sobieski J. Multidisciplinary optimization of a transport aircraft wing using particle swarm. Journal of Structural and Multidisciplinary Optimization 2004; 26(1): 121-131.

[26] Metropolis N, Rosenbluth A W, Rosenbluth M N, et al. Equations of state calculations by fast computing machines. Journal of Chemical Physics 1953; 21(6): 1087-1092.

[27] Kirkpatrick S, Gelatt C D, Vecchi M. Optimization by simulated annealing. Journal of Science 1983; 220 (4598): 498-516.

[28] Cerny V. Thermodynamical approach to the travelling salesman problem: an efficient simulation algorithm. Journal of Optimal Theory and Applications 1985; 45(1): 41-51.

[29] Chattopadhyay A, Seeley C E. Simulated annealing technique for multiobjective optimization of intelligent structures. Journal of Smart Materials and Structures 1994; 3(2): 98-106.

[30] Rafique A F, He L S, Zeeshan Q, et al. Multidiscipli-nary design of air-launched space launch vehicle using simulated annealing. In: Mertsching B, Hund M, Aziz Z. Lecture Notes in Computer Science 5803 LNAI. Berlin Heidelberg: Springer-Verlag, 2009; 719-726.

[31] Tekinalp O, Utalay S. Simulated annealing for missile trajectory planning and multidisciplinary missile design optimization. AIAA-2000-684, 2000.

[32] Tekinalp O, Bingol M. Simulated annealing for missile optimization: developing method and formulation techniques. Journal of Guidance, Control, and Dynamics 2004; 27(4): 616-626.

[33] de Jong K A. An analysis of the behavior of a class of genetic adaptive systems. PhD thesis, Michigan: University of Michigan, 1975.

[34] Dixon L C W, Szego G P. The optimization problem: an introduction. In: Dixon L C W, Szego G P. Towards Global Optimization II. New York: North Holland, 1978.

[35] Ackley D H. A connectionist machine for genetic hill-climbing. Boston: Kluwer, 1987.

[36] Schwefel H P. Numerical optimization of computer models. Chichester: Wiley & Sons, 1981.

[37] Dixon L C W, Mills D J. Effect of rounding errors on the variable metric method. Journal of Optimization Theory and Applications 1994; 80(1): 175-179.

[38] Torn A, Zilinskas A. Global optimization. Lecture Notes in Computer Science. Berlin: Springer-Verlag, 1989.

[39] Griewank A O. Generalized descent for global optimization. Journal of Optimization Theory and Applications 1981; 34(1): 11-39.

[40] Isakowitz S J. International reference guide to space launch systems. 3rd ed. Reston: American Institute of Aeronautics and Astronautics, 1999.

[41] Turner M J L. Rocket and spacecraft propulsion. 3rd ed. Chichester: Springer Praxis, 2009.

[42] Sutton G P, Biblarz O. Rocket propulsion elements. 7th ed. New York: Wiley-Interscience, 2001.

[43] He L S. Launch vehicle design. Beijing: Beihang University Press, 2004.

[44] Blake W B. Missile DATCOM: user's manual-1997 FORTRAN 90 revision. Oklahoma: Wright-Patterson AFB, 1998.

[45] Zipfel P H. Modeling and simulation of aerospace vehicle dynamics. Reston: AIAA, 2007.

[46] Xiao Y L. Rocket ballistics and dynamics. Beijing: Beihang University Press, 2005.

[47] Crossley W A, Williams E A. A study of adaptive penalty functions for constrained genetic algorithm-based optimization. AIAA-1997-0083, 1997.

Biography:

Amer Farhan RAFIQUE Born in 1980, he received B.E. (Aerospace) degree from College of Aeronautical Engineering, National University of Science and Technology, (NUST), Pakistan in 2002, M.S. degree from University of Engineering and Technology, Taxila, Pakistan in 2006, and Ph.D. degree from Beihang University, China in 2010. His main research interest includes multidisciplinary design and optimization of aerospace systems. E-mail: afrafique@yahoo.com; afrafique@sa.buaa.edu.cn