Probabilistic Multidisciplinary Design
Optimization on a high-pressure sandwich
wall in a rocket engine application
Author: Dennis WAHLSTROM¨ Supervisors: S ¨oren KNUTS Robert TANO Jan ¨OSTLUND
A thesis submitted in fulfillment of the requirements for the degree of Master of Science
Declaration of Authorship
I, Dennis WAHLSTROM¨ , declare that this thesis titled, “Probabilistic Multidisciplinary
Design Optimization on a high-pressure sandwich wall in a rocket engine applica-tion” and the work presented in it are my own. I confirm that:
• This work was done wholly or mainly while in candidature for a master degree at Ume˚a University for GKN Aerospace.
• Confidential information has been removed or modified for publication. • Where I have consulted the published work of others, this is always clearly
attributed.
• Where I have quoted from the work of others, the source is always given. With the exception of such quotations, this thesis is entirely my own work.
• I have acknowledged all main sources of help.
“The important thing is not to stop questioning. Curiosity has its own reason for existence. One cannot help but be in awe when he contemplates the mysteries of eternity, of life, of the marvelous structure of reality. It is enough if one tries merely to comprehend a little of this mystery each day.”
Ume˚a University
Abstract
Department of Physics Master of Science
Probabilistic Multidisciplinary Design Optimization on a high-pressure sandwich wall in a rocket engine application
by Dennis WAHLSTROM¨
A need to find better achievement has always been required in the space industry through time. Advanced technologies are provided to accomplish goals for human-ity for space explorer and space missions, to apprehend answers and widen knowl-edges. These are the goals of improvement, and in this thesis, is to strive and de-mand to understand and improve the mass of a space nozzle, utilized in an upper stage of space mission, with an expander cycle engine.
The study is carried out by creating design of experiment using Latin Hypercube Sampling (LHS) with a consideration to number of design and simulation expense. A surrogate model based optimization with Multidisciplinary Design Optimization (MDO) method for two different approaches, Analytical Target Cascading (ATC) and Multidisciplinary Feasible (MDF) are used for comparison and emend the con-clusion. In the optimization, three different limitations are being investigated, de-sign space limit, industrial limit and industrial limit with tolerance.
Optimized results have shown an incompatibility between two optimization ap-proaches, ATC and MDF which are expected to be similar, but for the two limita-tions, design space limit and industrial limit appear to be less agreeable. The ATC formalist in this case dictates by the main objective, where the children/subproblems only focus to find a solution that satisfies the main objective and its constraint. For the MDF, the main objective function is described as a single function and solved subject to all the constraints. Furthermore, the problem is not divided into subprob-lems as in the ATC.
Ume˚a University
Sammanfattning
Department of Physics Master of Science
Probabilistic Multidisciplinary Design Optimization on a high-pressure sandwich wall in a rocket engine application
by Dennis WAHLSTROM¨
Ett behov av att ¨oka prestandan har alltid varit eftertraktat inom rymdindustrin. Avancerade teknologier skapar nya m ¨ojligheter f ¨or m¨anniskan att genomf ¨ora rym-dresor f ¨or att kunna besvara fr˚agor och f ¨orb¨attra b˚ade kunskap och f ¨orst˚aelse om rymden. I den h¨ar uppsatsen handlar om att f ¨orst˚a och minimera massan p˚a ett rymdmunstycke som anv¨ands i ett ¨ovre raketstegsuppdrag med en expandercykel motor.
Studien ¨ar genomf ¨ord genom att skapa f ¨ors ¨oksplan med Latin Hypercube Sam-pling (LHS) som tar h¨ansyn till behovet av ett stort antal design och simuleringstid. Surrogat modeller ¨ar framtagna och optimering ¨ar genomf ¨ord med Multidisciplin¨ar Design Optimering (MDO) metod f ¨or tv˚a olika ansatser, Analytical Target Cascad-ing (ATC) och Multidisciplinary Feasible (MDF) som har anv¨ants f ¨or att j¨amf ¨ora och verifiera slutsatsen. I optimeringen, har tre olika begr¨aningar unders ¨okts, design-srum gr¨ans, industriell gr¨ans, och industriell gr¨ans med tolerans.
Resultat fr˚an optimeringen har visat p˚a en inkompatibilitet mellan de tv˚a olika ansatserna, d¨ar ATC och MDF f ¨orv¨antades vara lika f ¨or de tv˚a begr¨ansningarna, designsrum och industriell gr¨ans. I det h¨ar fallet, styrs ATC formuleringen av hu-vudproblemet, d¨ar underproblemen fokuserar enbart p˚a att hitta en l ¨osning som uppfyller huvudproblemet och dess tv˚angsvillkor. F ¨or MDF formuleringen, ¨ar hu-vudproblemet definierad som ett singel problem som l ¨oses med avseende p˚a alla tv˚angsvillkor, dessutom har problemet inte delats in i underproblemen som i ATC formuleringen.
Valet av surrogatmodell har stor p˚averkan p˚a optimeringens exakthet, och robus-theten i optimeringen har unders ¨okts med en annan typ av f ¨ors ¨oksplan. I den h¨ar studien skapas en f ¨ors ¨oksplan med full faktoriellt metod som ¨ar centrerad n¨ara den optimala l ¨osningen. Inom en s˚adan region har resultatet visat en ¨overenst¨ammelse f ¨or de flesta surrogatmodellerna f ¨orutom max temperatur, livsl¨angd och t ¨ojning p˚a den varmaste delen och har visat sig ha st ¨orst effekt p˚a innerv¨aggstjockleken p˚a rymdmunstycket.
Jan ¨Ostlund, who have been great support to this thesis and given me this unique opportunity to work with a very interesting project and what I like to do. Especially, S ¨oren Knuts for his high effort for supporting and giving me a clear figure about this project and how to solve it, this has been an important step to keep solving the problems and moving forward. During the thesis at GKN I have learned many things such as how to cooperate in a group and about space industry, things that I did not know before.
My advisors, Raja Visakha at GKN and Eddie Wadbro at Ume˚a University for suggesting and providing me an idea and understanding about Multidisciplinary Design Optimization method. Magnus Wir for supporting on how to describe the production cost and manufacturing limits. Lars and Ingegerd Ljungkrona for help-ing me with ModeFrontier Software. Arne Boman for his aerodynamics and thermo-dynamics explanations, Martin Carlsson for the help for understanding FEI space nozzle better, Oscar Linde for helping with linux OS, and many thanks to everyone at space nozzle’s team that was very supportive and able to provide me a happy work place. Alexander Hall, another thesis workers and trainees for all the enter-tainments.
Contents
Declaration of Authorship iii
Abstract vii Sammanfattning ix Acknowledgements xi 1 Introduction 1 1.1 GKN Aerospace . . . 1 1.2 Background . . . 2
1.3 Main component of a space rocket . . . 2
1.4 Space Nozzle. . . 3
1.5 Objective of the thesis. . . 6
2 Theoretical Background 7 2.1 Physics of Space Nozzle . . . 7
2.1.1 Thermodynamics and Chemical Reaction . . . 8
2.1.2 Aerodynamics. . . 9
2.1.3 Solid Mechanics. . . 10
Fatigue crack initiation. . . 10
Ductile Rupture . . . 11
2.1.4 Physics in a sandwich wall . . . 12
2.2 Statistic Background and Design of Experiment. . . 13
2.2.1 Latin Hypercube Sampling . . . 13
2.2.2 Full Factorial Analysis . . . 14
2.2.3 Analysis of Variance . . . 15
2.3 Multidisciplinary Design Optimization Architecture . . . 15
2.3.1 Single Objective, Multidisciplinary Design Optimization Ar-chitecture . . . 16
2.3.2 Multi Objective, Multidisciplinary Design Optimization Ar-chitecture . . . 17
2.3.3 Individual Disciplinary Feasible . . . 18
Analytical Target Cascading . . . 18
Augmented Lagrangian Relaxation . . . 19
2.3.4 Multidisciplinary Feasible . . . 20
3 Implementation and Multidisciplinary optimization 21 3.1 Main study methodology . . . 21
3.2 Decomposition of cooling system . . . 23
3.3 Implementation . . . 23
3.3.1 Data Implementation. . . 24
xiv
3.3.3 Solid mechanics simulation and finite element method . . . 27
3.4 Prediction of surrogate model . . . 27
3.4.1 Natural behavior of the surrogate model . . . 28
3.5 Multidisciplinary Optimization . . . 28
3.5.1 Optimization constraints . . . 28
3.5.2 Optimization problem formulations . . . 30
3.5.3 Local and global optimization . . . 33
3.5.4 Different limitations . . . 33
3.6 Certainty of Surrogate model . . . 34
4 Data analysis and Results with Discussion 35 4.1 Data analysis. . . 35
4.1.1 Fitting of the surrogate model. . . 35
4.1.2 Quality of response surface and surrogate model . . . 36
4.1.3 The limitations used to optimize . . . 36
4.2 Results of optimization . . . 37
4.3 Probabilistic analysis of uncertainty . . . 41
4.3.1 Accuracy of Surrogate Model . . . 41
4.3.2 Source of uncertainties . . . 43
4.3.3 Impact of uncertainties. . . 44
5 Conclusion and Recommendation 47 5.1 Conclusion . . . 47
5.2 Recommendation . . . 48
Bibliography 49 A Derivation 51 A.1 Main jet expansion . . . 51
B Figures and Tables from Surrogate model 53
List of Figures
1.1 Schematic picture of Ariane 6 main components. (ESA) . . . 3
1.2 Schemetic system picture over an expander cycle engine. . . 4
1.3 Picture of ETID (left) and FEI (right) [5]. . . 5
1.4 Two identical channels in a sandwich wall, observed from top view point. . . 5
2.1 Sections of space nozzle. . . 8
2.2 Cooling system in sandwich wall.. . . 12
2.3 Full factorial experiment for 22. . . . 14
2.4 Single-Objective, MDO architecture. . . 16
2.5 An illustration of Hierarchical in the ATC case. . . 18
2.6 An illustration Nof on-Hierarchical in the ATC case. . . 19
3.1 MDO Methodology.. . . 22
3.2 Half of cooling channel. Left side is the flame side and right side is toward the outside, with meshes. . . 23
3.3 LHS for 7 design factors for two different design spaces 0.5 − 1.5 and 0.2 − 1.7, with 50 and 200 as number of designs. . . 25
3.4 Condensation operating point with wall temperature toward and in the flame side, plotted with the lowest temperature different ∆Tcondense and the saturation line. . . 26
3.5 Channels in cooling system [5]. . . 29
3.6 System Level. . . 30
3.7 Hieararchical ATC problem formulation. . . 30
3.8 ATC Algorithm. . . 32
3.9 An example concept and illustration on the different limits, Design Space (D.S.), Industrial (Ind.) and Industrial with Tolerance, are de-fined with function’s minimum.. . . 34
4.1 Optimized geometry from different limits and methods. . . 40
4.2 Plots of ANOVA, Actual vs. Surrogate model with its residual. . . 41
4.3 Plots of ANOVA, Actual vs. Surrogate model with its residual. . . 42
4.4 Effect influence for each design variable. . . 44
B.1 Correlation matrix. . . 54
B.2 Plot of mass over design variables. . . 54
B.3 Plot of ∆Tcover design variables.. . . 55
B.4 Plot of ∆PC in logarithmic scale over design variables. . . 55
B.5 Plot of Tmaxover design variables. . . 56
B.6 Plot of damage in logarithmic scale over design variables. . . 56
B.7 Plot of heat pick-up over design variables.. . . 57
B.8 Plot of ∆1 over design variables. . . 57
xvi
B.10 Plot of mach number over design variables. . . 58
B.11 Plot of temperature different, ∆Tcondenseover design variables.. . . 59
B.12 Plot of ∆Pcover design variables.. . . 59
B.13 Plot of Damage over design variables. . . 60
C.1 Different strain peaks. . . 61
C.2 Optimized with mass as objective and no physical constraints. . . 61
List of Tables
1.1 Tolerance cost for weight reductions 1 kg or 1% of the structural weight. [4] . . . 2 3.1 Description of design variables. . . 23 3.2 List of constraints and each of its bound with description . . . 29 4.1 Number of data left to create surrogate model for 0.5 to 1.5 and 0.2 to
1.7interval. . . 35 4.2 Fitting quality to the data for each variable using 200 designs between
a LHS interval of 0.2 to 1.7. . . 36 4.3 Different limitations for each design variable. . . 36 4.4 Total number of local solver evaluated and converged in the MDF and
ATC for different limits. . . 37 4.5 Result from the optimization for MDF and ATC, for different limits.
The results are relative to reference values. . . 37 4.6 Result from the optimization for MDF and ATC, with a comparison
between surrogate models and actual results without industrial lim-its. The results are relative to reference values. . . 38 4.7 Result from the optimization for MDF and ATC, with a comparison
between surrogate models and actual results with industrial limits. The results are relative to reference values. . . 38 4.8 Result from the optimization for MDF and ATC, with a comparison
between surrogate models and actual results with industrial limits and tolerances. The results are relative to reference values. . . . 39 B.1 Fitting quality to the data for each variables using 50 designs between
an interval of 0.5 to 1.5. . . 53 B.2 Fitting quality to the data for each variables using 50 designs between
List of Abbreviations
ATC Analytical Target Cascading
ASL Airbus Safran Launchers
DS Design Space
DoE Design of Experiment
ESA European Space Agency
ETID Expander Technology Integrated Demonstrator
FEI Flight Engine Image
FEM Finite Element Method
HCF High Cycle Fagtigue
IDF Individual Discipline Feasible
LCF Low Cycle Fagtigue
LHS Latin Hypercube Sampling
MDO Multidisciplinary Design Optimization
MDF MultiDisciplinary Feasible
MOMDOA Multi Objective Multidisciplinary Design Optimization Architecture
NE Nozzle Extension
SOMDOA Single Objective Multidisciplinary Design Optimization Architecture
SQP Sequential Quadratic Programming
List of Symbols
Ae nozzle exhaust area m2 E modulus of elasticity Pa
hg heat transfer coefficient W m−2K−1
H enthalpy J
Isp specific impulse m/s2 l area moment of inertia m4
m mass kg
˙
m mass change kg/s
np particle density mol−1m−3 Nc number of cycle Unit-less
P pressure Pa
Pa ambient pressure Pa Pc coolant pressure Pa Pe exhaust pressure Pa
q heat flux W/m2
qg heat flux of gas W/m2
Q heat W
Qtot heat pick-up W r radius of gyration m
R ideal gas constant J mol−1 K−1 R strain ratio Unit-less Rσ stress ratio Unit-less
T temperature K Tc coolant temperature K u velocity m/s U internal energy J V volume m3 ve exhaust velocity m/s W work J
∆ physical different Unit-less
strain rate Unit-less
γ heat capacity ratio Unit-less
ρ density kg/m3
Chapter 1
Introduction
Space is a subject that has been curious to us for many centuries and even before Christian. Aristotle ideal about cosmology and heliocentric model was grounded by Aristarchus of Samos had stunned many individuals, but at that time it was lacked mathematics and unable to prove, and was barely accepted as a fundamental con-tent [1]. These doctrines had become more underwritten many centuries later when telescope had been invented and able to observe and prove this theory with both ob-servational and mathematical. Copernicus presented a heliocentric model indepen-dently of Aristarchus based on his astronomical observations. Later in seventeenth century, Galilei could prove the heliocentric model with a scientific method, during the time in his house arrest caused by the heliocentric ideology instead of geocen-tric ideology that was approved by the church. The heliocengeocen-tric model describes planetary motion, and in fact was one of the starting point for Kepler. He advanced this ideology and provided laws that planets unavoidably must obey and known as Kepler’s law. The result was essential to physicist at that time, and its derivation of Kepler’s law was utilized by Newton to generalize motion laws and universal gravi-tation, Newton’s law [2] [3].
Curiosity is one of the strongest ”sense” in human nature, and it is a center of evolution to strive to answer on these inquires, what, why and how. Without these perspectives, the evolution and development of technology would have been differ-ent from today and this comes to our cause of matter.
In this section, a short introduction about GKN Aerospace will be reviewed as well for background and objective of this thesis, together with an explanation of space rocket mission, and briefly introduce the importance of space nozzle that will be focusing.
1.1
GKN Aerospace
2 Chapter 1. Introduction
1.2
Background
The development of space rocket is essential for us to provide a possibility to travel outer space and explore the outer world. This conception requires a more advanced technology that we do not have today, the components for each space rocket part need to be more effective, which is challenging for engineers and scientists. Space rocket and space shuttle have been manufactured since middle of twentieth cen-tury, when Yuri Gagarin went outer space and Neil Armstrong went to the moon, and when the first artificial satellite Sputnik in-landed in the Earth’s orbit, and ever since the space industry keep expanding. New methods and technologies provide possibility in opportunity to make software more effective and complex with more functions able to interpret and handle issues and bugs cause from the methods. The issues and causes that today stand as a center is that the solutions are not robust and need to be investigated. Last but not least, cost of all actions need to be reduced to their minimum for maintenance. Table1.1show costs of weight reduction between vehicle, aircraft, large-aircraft and space rocket. An optimized design solution for space rocket has a greater impact for cost standpoint comparing with aircraft, there-fore this is an essential task for engineers and scientists to examine.
TABLE 1.1: Tolerance cost for weight reductions 1 kg or 1% of the structural weight. [4] Te/kg Te/wt.% Vehicle 0.01 0.1 Aircraft 0.5 100 Large-aircraft 1 1000 Space rocket 10 10000
One project that GKN participates in is Ariane 6 program, originally proposed by European Space Agency (ESA) and sub-divided to Airbus Safran Launcher (ASL). GKN is enrolling in and responsible for the development and production of turbine engine for upper stage and first stage, an upgrade engine from Vulcain II use in Ar-iane 5. For the upcoming ArAr-iane 6, the nozzle is called, Sandwich Advanced Nozzle (SWAN). In the future, with great technology GKN is processing and wanting to advance to design and manufacture even for the upper stage nozzles.
1.3
Main component of a space rocket
FIGURE1.1: Schematic picture of Ariane 6 main components. (ESA)
1.4
Space Nozzle
4 Chapter 1. Introduction
FIGURE1.2: Schemetic system picture over an expander cycle engine.
There are different types of space nozzles suited to specific engines. In Ariane 6, an upgrade Vulcain II is designed to use a gas generator cycle, differently from the upper stage engine Vinci utilizes another principle. That is a gas expander cycle, for reusing the fuel in a different way. The gas expander cycle is shown in figure1.2 which will be focused on in this thesis.
Gas expander cycle utilizes heat/temperature generated from the engine in a way that fuel travels down and absorbs the heat, and later is being reused to drive the turbine. Another concept is a gas generator cycle, fuel separates and divides into two ways, first, to a pre-burner and second to the nozzle acts as coolant. Further-more, the fuel that has cycled around the nozzle will be disposed in the combustion chamber without being reused. The pre-burner is connected to oxidizer pump and to fuel pump, which some of the oxidizer is divided to the pre-burner and the com-bustion chamber for decreasing the temperature. Fuel and oxidizer that accumulates in the pre-burner is used to drive the turbine, and discarded after being used.
FIGURE1.3: Picture of ETID (left) and FEI (right) [5].
There is another important thing concerning this thesis, is the sandwich wall of cooling channels, and it is shown in figure 1.4. The observing point is from top view, where the considered system is for two channels with a complete symmetric between them when it comes to size and this two channels system is repeated into the whole size space nozzle. In the figure, an idea is to illustrate the channels that will be studied and a fully and more rigorously explanation will be provided in2.1.4 and3.3.1sections.
6 Chapter 1. Introduction
1.5
Objective of the thesis
The objective of this thesis is to find the ”optimal” structural design of sandwich wall of FEI space nozzle with mass as the main optimization objective, using Multidis-ciplinary Design Optimization (MDO) method. By creating Design of Experiment and forming a surrogate model of aerodynamics, thermodynamics and solid me-chanics, and analyze their interactions and addiction on the design variables. This is for widening knowledge and understanding of the space nozzle. A compari-son between two different MDO problem formulations, Analytical Target Cascading (ATC) and Multidisciplinary Feasible (MDF) will be investigated and concluded for strengthening the reliability of the solutions. The optimal solutions of the space noz-zle need to be considered to the manufacture and production perspective, that GKN has high value and expectation.
Chapter 2
Theoretical Background
In chapter 1, a short introduction about space nozzle and objective of this thesis have been introduced, where the space nozzle that will be used in this study is FEI which is an upper stage engine with an expander cycle system. The need to understand such system is important, with a goal in chapter 2 to give a concise and wide intro-duction about physics theory and statistic theory that will be needed to understand this thesis.
The study can be done in many ways, with a conception to find surrogate models that describes physics in term of geometry. Besides, the analysis is consisting of fi-nite element method (FEM), with that cause the simulations are time expensive and needed to construct and formulate the problem in such way that reduces the num-ber of simulations. Therefore, the decomposition of the problem becomes essential. One of the ways to engage the problem is by creating a Design of Experiment (DoE) scheme and benefiting its feasibility of probabilistic spread of the data and simu-late the physics, aerodynamics, thermodynamics and solid mechanics, that will be previewed in this chapter.
2.1
Physics of Space Nozzle
Space nozzle experiences high thermal and pressure loads created from turbine and chemical reaction between fuel and oxidizer. This energy production is used to ac-celerate space rocket, higher velocity results in higher thrust force, and at the same time the material needs to sustain such extreme conditions. A need to understand the energy travels and acts in the space nozzle is a key to visualize the thermal and pressure loads affecting the material, and with that knowledge, the cooling system can be designed designate to these physical behaviors and properties.
8 Chapter 2. Theoretical Background
FIGURE2.1: Sections of space nozzle.
2.1.1 Thermodynamics and Chemical Reaction
Space nozzle regions can be seen in figure2.1. In the first region, combustion cham-ber, an interaction between fuel and oxidizer occurs, that is exchanging and releasing energy in form of thermal energy. This phenomenon is called exothermic reaction:
reactant → product + energy (2.1) For a reaction between liquid hydrogen and liquid oxygen:
2H2+ O2 → 2H2O (2.2)
Resulting in bonding energy described as an enthalpy ∆H and has its physical for-mula:
∆H = HCreation− HAnnihilation (2.3) where HCreationis the enthalpy used to create a system and HAnnihilationis the energy use to destroy a system, with the definition of enthalpy:
H ≡ U + P V (2.4)
An amount of thermal energy U created from the substances obeys the first law of thermodynamics and the ideal gas law to the first order. The first law of thermody-namics is defined as:
P = npkbT (2.6) with npas particle density, kb Boltzmann constant and T temperature. This concep-tion has not yet violated any physical law, but in an experimental concept, the ideal gas law is not the fundamental way of describing the pressure in a space nozzle, ex-cept it contains more terms to equation (2.6), by using of the perturbation theory and getting the correction terms. This also applies to the diverging section. In addition to the diverging section, its physical properties of thermodynamics and aerodynam-ics, their character is behaving as an isentropic, meaning it is both adiabatic and reversible [6].
When the wall of a space nozzle has low heat resistance the temperature tends to decrease and causes the gas to change its property due to the phase transition. The condensation is determined by temperature and pressure of the hot gas, with a different in density between the two phases, liquid and gas. This phase transition is decided by a saturation line, which is a line that separates different states and phases of the gas, solid and liquid. A well-known theory of condensation on a solid interface is resolved by Nusselt, which in this thesis will be referred to the original Nusselt’s work [7].
2.1.2 Aerodynamics
In the combustion chamber the temperature and chemical reaction are dominating, but in converging section, throat and diverging section the aerodynamics is the focal problem. In addition to the earlier statement, the ideology of the space nozzle’s sec-tions is to convert the thermal energy into kinetic energy, with the thrust force of the space rocket depends on the exit velocity. This can be illustrated as the momentum law. This also implies that each section has its specific velocity, and in the converging section the fluid has a subsonic pattern, while at the throat the velocity increases and becomes a sonic. Conceptional, the desired outcome can be regulated by decreasing area and control the mass flow of the fluid. The mass flow over a surface Ω1is given by: dm dt = Z Z Ω1 ρ1u1dA1 (2.7)
Imaginably, the entry of the converging section as 1 and at the throat as 2, by using the mass flow ratio between two sub-systems assumed to be equal, a relation can be determined:
10 Chapter 2. Theoretical Background The exit velocity from a main jet expansion is:
ve= v u u t2 γ γ − 1 RTc M 1 − Pe Pc γ−1 γ ! (2.9) where the derivation of this can be seen in AppendixA. The exit velocity veis essen-tial to the space rocket’s thrust force and is given by:
F = ( ˙mve+ PeAe) − PaAe (2.10) where ˙mis mass change, Peis pressure at the exit, Aeis area of the nozzle at the exit and Pa, ambient pressure. Same thrust force can be derived for impulse:
F = ˙mIsp (2.11)
with Ispas the specific impulse. The thrust force is depending on the change of mass ˙
m. When a fluid undergoes phase transition its density is then changed affecting the specific impulse and the momentum of the space rocket rocket. A well-known equation of motion for a space rocket is described by Tsiolkovsky’s rocket equation, with the conservation of momentum derived:
∆v = veln( m0 m1
) (2.12)
where m0is initial mass and m1is final mass.
2.1.3 Solid Mechanics
The material that is used to manufacture a space nozzle is important whereas it must be able to keep its properties due to loads from aerodynamics and thermodynamics, such as pressure and heat. Pursuant to earlier sections2.1.1and 2.1.2, the largest temperature production is held in the combustion chamber and can cause material to break or fail, due to fatigue or ductile rupture.
Fatigue crack initiation
Fatigue is a way to describe metal properties, and the limit is studied by testing many cycles. There are two different types of fatigue, Low Cycle Fatigue (LCF) and High Cycle Fatigue (HCF), these are depending on various factors.
LCF is used to describe life under high temperature and low number of cycles, which the cause of this is predominant of thermal loads and cause the material to fail [9]. But for HCF, the material can withstand many numbers of cycles and the loading is caused mainly from the elasticity [10]. The stress ratio is given by:
Rσ = σmin/σmax (2.13)
where σminand σmaxare the minimum stress and maximum stress respectively. The strain ratio is:
R = min/max (2.14)
This equation has excluded the mean stress, which in some case the solution can be better approximated by taking account of its mean value, with the mean stress correction gives by the Morrow hypothesis [11]:
∆ 2 = (σf0 − σmean) E (2Nc) b+ 0 f(2Nc)c (2.16) where E is the modulus of elasticity. The value of b and c can be fitted by creating a curve (Strain range versus number of cycles to failure). The mean stress may have a significant effect on life, especially for HCF and for large strain ranges the effect of the mean stress is negligible, independently on LCF or HCF. Another important physical property is Damage, and the damage caused by a cycle is defined as:
Damage ≡ 1 Nc
(2.17) For several cycles, the damage is added linearly:
Damagetot = D1+ D2+ ... (2.18) Observing equation (2.18), the Damage is defined as an inverse proportional to the number of cycles, with few numbers of cycles the Damage is high, vice versa for many numbers of cycles the Damage is low.
Ductile Rupture
12 Chapter 2. Theoretical Background
2.1.4 Physics in a sandwich wall
FIGURE2.2: Cooling system in sandwich wall.
In sections2.1.1,2.1.2and2.1.3the main source of heat and pressure from the rocket engine has been presented. In this section, is to utilize previous knowledge and fo-cus on the sandwich wall channels. For a bi-channel system in the sandwich wall a 2D illustration can be seen in figure2.2, where the inner wall to the gas flame is positioned toward us and opposite for the outer wall. Moreover, this bi-channel is repeated into the fully size space nozzle, figure1.3. Starting from top, fuel is divided into each channel equivalently and under its way down, it experiences pressure loss due to frictions from the wall surfaces. Channel split is used to separate fuel into two evenly channels, which has an importance for the space nozzle when the area and the number of channels have increased, shown in figure2.2. Such structural design is necessary to the inner wall’s prospect of damage and temperature. Continuously, this fuel travels to the second channel and out to the outlet manifold with a hydro-dynamics equilibrium as the physical requirement at the outlet manifold. Whereas, this requirement is the source for the velocity at the inlet of the second channel to be higher than the first channel.
The decisive distant of fuel is between the inlet manifold and outlet manifold, which has a significant physics to the turbine. Energy the fuel successfully picking up during the transportation can directly be contributed to the turbine. Pressure difference of the fluid/coolant is given by:
∆Pc= Pcin− Pcout (2.19)
Temperature difference:
Qtot= qdA (2.21) Where q is the heat flux produces from the thermodynamics reactions in the com-bustion chamber. Observing that the temperature of the fuel increases due to the heat, equation (2.20) both with the heat pick-up equation (2.21) are having a com-mon relation. The heat transfer from hot flame side to a solid wall is given by:
qg = hg(Twall− Tgas) (2.22) This description shows the heat from the gas transfers to the inner wall with a heat transfer factor hg.
2.2
Statistic Background and Design of Experiment
There are many available methods for designing an experiment, all having advan-tages and disadvanadvan-tages. Factorial analysis is one of the common methods used to design an experiment, for two levels the number of design grows exponential in term of factors, 2k. Unavoidable, the designs become large for many factors, but this can be useful when the problem has decomposed and reduced into fewer pa-rameters. A way to decompose the problem is by the Plackett-Burman screening, where this method objectively finds the relations and the main effects between fac-tors, and reduce the number of variables by neglecting those with low correlation and main effect. In this study, one wants to collect information even those with low main effect, therefore, a more suitable method for many factors problem is the Latin Hypercube Sampling (LHS).
2.2.1 Latin Hypercube Sampling
LHS assembles and generates the needed values for each variable Xj by statistical span the hyperspace/hypercube. Each variable is generated by multidimensional distribution Dj with equal probable/uniformly and disjoint. For N numbers, the probability for each sample is 1/N , this repeats for each variable and stratifies into a matrix Xij, where i = 1, ..., N and j = 1, ..., M . To illustrate this, assume that X is continuous with a cumulative distribution function FX(x).
X ∼ P (X ≤ x) ≡ FX(x) (2.23)
Let also φi be sampling probability with equation (2.23) describes random variable, X. The sample xican be described by the inverse transformation of the sampling φi:
xi = FX−1(φi) (2.24)
Formulation of equation (2.24) fulfills The inversion principle[12].
14 Chapter 2. Theoretical Background For Yi= g(Xi), its mean value is given by:
¯ Y = 1 N N X i=1 Yi (2.25)
and variance for the random sample: VRS( ¯YRS) =
1
NV (Y ) (2.26)
For stratify case:
VS( ¯YS) = VRS( ¯YRS) − 1 N2 N X i=1 (µi− µ)2 (2.27) The variance of LHS: VLS( ¯YLS) = VRS( ¯YRS) + 1 (NM +1(N − 1)M −1) X RS (µi− µ)(µj− µ) (2.28)
with µi, µj as the expectation values in cell i or cell j in a stratify or LHS case, and µ = E(Y ). This has evidently shown a reduction of the variance of the mean value of a factor 1/N2, equations (2.27). However, in case of LHS, equation (2.28), the variance is a combination between the random sampling and stratify sampling, in such a way the method gives higher stability to the probabilistic scheme [14].
2.2.2 Full Factorial Analysis
Full factorial analysis is one of the DoE methods used to study the variations be-tween factors, for two levels (low and high) the number of studies increase by 2kfor kfactors mentioned earlier. The conception of this is to create points for every edge (low and high) in a k number of variable, and for the 2 dimensional its full factorial can be illustrated in figure2.3[15].
MA= i=1
(¯yA=i− ¯y)2 (2.29)
where ¯yis average value of the responses and yA=iis response of the study A.
2.2.3 Analysis of Variance
Analysis of Variance (ANOVA) is a way to investigate variations between variables and determine quality of the data. The coefficient of determine is given by:
R2= 1 −SSE SST
(2.30) This coefficient is used to evaluate the deviation between data and predicted func-tion, where SSE and SST are error sum square and total sum of square.
SSE = n X i=1 (yi− ˆyi)2= n X i=1 e2i (2.31) SST = n X i=1 (yi− ¯y)2 (2.32)
Sum square of error is used to describe the residual between estimated points and data. Differently from the total sum of square, which is a sum of total error that has possibility to occur and the sum square of error.
There is other type of coefficient of determine with a consideration to the data points used to interpolate and predict a function, this is expressed as R2adj:
Radj2 = 1 − SSE/dfe SST/dft
(2.33) with degree of freedom dft= n − 1and dfe= n − p − 1for n number of designs and pnumber of variables used.
2.3
Multidisciplinary Design Optimization Architecture
16 Chapter 2. Theoretical Background and its mathematics expression:
Minimize : F (X) (2.34)
with respect to : X
subject to : g(X) ≤ 0 ψ(X) = 0 R(X) = 0
with g = g1, g2, ..., gm represents unequal constraints, ψ = ψ1, ψ2, ..., ψl equal con-straint/consistency constraints and R residual constraints, where X = X1, X2, ..., Xn are the design variables.
2.3.1 Single Objective, Multidisciplinary Design Optimization Architec-ture
Single Objective, Multidisciplinary Design Optimization Architecture (SOMDOA) can be seen in figure2.4.
FIGURE2.4: Single-Objective, MDO architecture.
its environment is based on the SOMDOA except with multiple objectives. These objectives have possibility to link and depend on each other and its formulation is shown below, equation (2.35).
Minimize F(X) (2.35)
with respect to X
subject to : g(X) ≤ 0 ψ(X) = 0 R(X) = 0 where F = F1, F2, ..., Fpis a vector of p objectives.
There are four different recipes in the MDO method for the monolithic formula-tion. All-at-Once (AaO) recipe is the general problem formulation of a MDO prob-lem, whereas no assumption has applied, equation (2.35). In turn, this can create problem in simulation expense due to its difficulty in converging depend on the original problem matrix, and therefore no guarantee on a feasible solution. Nor-mally, the AaO formulation is rare to use because of in a real-world problem, con-sistency of the constraints or residual of the model is often needed to be assumed and in some occasion negligible [17]. This is for keeping the problem in a range that is solvable and estimable for understanding the study outcomes. The second recipe is Simultaneous Analysis and Design (SAND), which is a simplification of AaO problem by eliminating the equal consistency constraints ψ of equation (2.35). This assumption can result in a less problematic in convergence comparing with the AaO [18].
18 Chapter 2. Theoretical Background
2.3.3 Individual Disciplinary Feasible
Elimination of residual constraints R(X) from the general MDO formulation equa-tion (2.35) results in an IDF formulaequa-tion.
Minimize : F(X) (2.36)
with respect to : X subject to : g(X) ≤ 0
ψ(X) = 0
IDF solves problems individually and for a large problem the method decomposes in such way, that increases the feasibility of the problem and easier to find the conver-gence, which is a great benefit for a large problem [17]. There are plenty of methods in the IDF formulation, but in this the focus will be on Analytical Target Cascading (ATC).
Analytical Target Cascading
ATC is one of the methods that divides problem into subproblems and increases the convergence ability in comparison with other MDO methods[19] [20]. This allows the user to optimize subproblems individually, and keep consistency between all the subproblems [21]. The feature of this is to utilize penalty function that each optimize loop condition is violated, will result in an addition penalty value in term of consistency constraint differences between subproblems and objective.
There are two different ways to solve and describe the problem in the ATC for-mulation, hierarchical and non-hierarchical shown in figures2.5and2.6.
FIGURE2.6: An illustration Nof on-Hierarchical in the ATC case.
Independently on how the main objective problem is described, in hierarchical ATC, the problems level i − 1 share variables with level i, inconsiderate to the sub-problems i − 1 that might have the same dependent variables. Instead, they strictly need to share their variables through the upper level problem i. Non-hierarchical despites the earlier formulation and allows all the problems and subproblems to share variables with each other freely, figure2.6. The number of sharing variables increase with the number of problem and subproblem, which can be difficult to keep tracking of all shared variables when the problems and the systems are large.
Augmented Lagrangian Relaxation
Individual Discipline Feasible operates by keeping track on sharing variables in each optimize function, holding consistency within the defined system. This consistency is formulated by Lagrange Multiplier and one conventional penalty function is the Augmented Lagrangian (AL) [23] [22].
The AL relaxation function is described:
π(c) = vTc + wT(c ◦ c) (2.37) where c is constraint, and in this case is the response and the target value from sys-tem i, ci = Ti− Ri. The AL can be applied to the ATC formulation function for an objective F :
Minimize: F (xj, gh(xj)) + πk(cj) (2.38) subject to: xj, gh(xj)
20 Chapter 2. Theoretical Background The update of v and w in equation (2.37) is formulated as:
vk+1= vk+ 2wk◦ wk◦ ck (2.39) wk+1i = ( wki if|cki| ≤ γ|ck−1i | βwki if|ck i| > γ|ck−1i | For constants β > 1 and 0 < γ < 1.
2.3.4 Multidisciplinary Feasible
Another way to decompose the AaO is Multidisciplinary Feasible (MDF) by elim-inating consistency constraints ψ(X) and R(X) residual constraint from equation (2.35). The MDF problem formulation becomes:
Minimize : F(X) (2.40)
with respect to : X
subject to : g(X) ≤ 0
Chapter 3
Implementation and
Multidisciplinary optimization
The theories to understand this thesis have been presented in the previous chapter and now we are able to apply these theories and deepen into our main objective.
In this chapter, the goal is to give readers a clear figure about implementation and optimization, including assumptions that are needed to be made with a consid-eration to the simulation expense. The implementation will carefully be introduced from construction of DoE scheme and to surrogate model. In the optimization part, the theory of focus to our problem will be formulated as well for the definition of different study optimize limitations.
3.1
Main study methodology
Work frame in this thesis is divided into six main steps and are presented in figure 3.1.
First step: This is to find and define what is it that needs to be studied in the cooling system. Importantly, understand the physics of how the engine contributes the main jet and the pressure inside the cooling channels and how these affect the material.
Second step: When the physical properties are known, the design variables that influence the physics can be created by making a DoE scheme for probabilistic vari-ation.
Third step: The created DoE is putting into black boxes to get physical relations influence by the design variables. The black boxes in this case are Tccool and Ansys, which will be discussed more in section3.3.1.
Fourth step: Interpret the outputs from the black boxes, including understand if the received results are logic and using to create surrogate models that describes the physical properties, such as, aerodynamics, thermodynamics and solid mechanics.
Fifth step: Optimization, using the MDO method includes formulation of objec-tive problem and creation of surrogate models application to the previous step. The MDO method will be used with two different approaches for two distinguish solv-ing ways. Before optimizsolv-ing the problem, three different limitations will be defined are implied to the MDO method, Design Space limit, industrial limit and industrial limit with tolerance.
22 Chapter 3. Implementation and Multidisciplinary optimization
as explained in the theory part. Further decomposition of the problem is to neglect manifold at the bottom, with an assumption of low impact to the results.
3.3
Implementation
This section is about the main set-up of the study, with further and specific assump-tions will be presented.
FIGURE3.2: Half of cooling channel. Left side is the flame side and right side is toward the outside, with meshes.
TABLE3.1: Description of design variables.
Design variable Description
x1 Number of channel
x2 Inner wall thickness
x3 Channel height
x4 Mid wall thickness x5 Outer wall thickness
x6 Mill radius
24 Chapter 3. Implementation and Multidisciplinary optimization
3.3.1 Data Implementation
The methodology uses to solve the problem is shown in figure3.1. Firstly, in the implementation step a DoE scheme is created using LHS method, which is necessity needed to create surrogate models that describes and suits to the physical problems. LHS scheme is built on seven design variables (mainly are geometries) for fifty and two hundred designs, for getting a wide range of study as possible by taking simu-lation time into account.
These design variables are shown in figure 3.2 and table 3.1. Observing, the design variables x4 and x1 are divided by a factor 2, because the presented piece is only half of the original channel and it is positioned at the top of the space nozzle in the combustion chamber. The notation a2 represents a constant, preventing the ”buckling” phenomenon in the space nozzle for different load of pressure between inside and outside [25], this constant only exist at the top of the space nozzle. The stress acts on the space nozzle’s wall is given by equation (3.1) which is originally derived by Euler. Differently, at bottom of the space nozzle the channel height is twice larger than the inlet, and the constant a2which is shown in figure3.2goes to zero. σ = F A = π2E l r 2 (3.1)
where F is force, A area, E modulus of elasticity, l area moment of inertia and r radius of gyration. This also shows for a constant stress σ, then F ∝ A.
There is one more important thing that must be considered when creating a new design variable for the channel height x3. When the design space has large varia-tions, this can result channel height to be smaller than the sum of the radii, x6 and x7, and cause problem when creating mesh elements and geometry for the finite element analysis. The problem is manipulated by defining a constant a1as:
a1 ≡ x3ref − (x6ref + x7ref) (3.2)
and make variations on this constant a1, instead of the channel height x3. Thereby a new channel height is taking form:
x3i = αi3a1+ αi6x6ref + αi7x7ref (3.3)
(A) Design space 0.5 − 1.5, with 50 designs. (B) Design space 0.2 − 1.7, with 50 designs.
(C) Design space 0.2 − 1.7, with 200 designs.
FIGURE3.3: LHS for 7 design factors for two different design spaces 0.5 − 1.5and 0.2 − 1.7, with 50 and 200 as number of designs.
The first study is to contribute LHS scheme value between 0.5 to 1.5, figure3.3 with a probability density 1/N described in section2.2.1. Wide range of the LHS can cause trouble in simulation, e.g. thin walls and large channel give a high hydro-dynamics pressure load on the wall and thin walls cannot withstand such pressure. In some exceptional case, if the walls are thin and the mill and the weld radii are large, this could compensate the fact that the walls are thin. These are depending on several factors and difficult to determine and face the problem directly. Therefore, by creating a DoE scheme with wide range is important to the study for finding the limitation of the material, and possibly designing space nozzles that today unable to manufacture. The factors generated from the LHS 0.5 to 1.5 will be multiplied into a reference value.
26 Chapter 3. Implementation and Multidisciplinary optimization
3.3.2 Aerodynamics and thermodynamics simulation
For given designs, the varies of the geometries are putting in black boxes, such as Tccool and Ansys for physics study, aerodynamics, thermodynamics and solid me-chanics, ∆Tc, ∆Pc, Qtot, Mass, Tmax etc. these are mentioned in section 2.1.4. The aerodynamics and thermodynamics are the results from the Tccool script, which is used and developed by GKN with the goal to simulate the space nozzle’s channels for given geometries of cooling system, figure2.2. It also neglects effects from the bottom manifold and the study can be chosen between one dimensional or two di-mensional along axial axis. In this study one didi-mensional will be selected. The physical result returns from the Tccool is a mean value of aerodynamics and ther-modynamics, where the bi-channel is assumed to be repeated into the whole space nozzle. Conceptional, is to minimize the study problem without losing any impor-tant physical properties. The Mach number uses in this study is the largest value in the first channel, figure2.1, where another physical has been selected depending on their maximum or significant for their specific position in the bi-channel.
There are two operating points that needs to be considered, hot operating point and cold operating point, which is also called condensation point. The differences of these points are physical initial values of the fuel and the engine. The condensation phenomenon has been mentioned in section2.1.1and this ∆Tcondensein figure3.4is representing the minimum of difference between the wall temperature toward the flame and the saturation temperature.
FIGURE3.4: Condensation operating point with wall temperature to-ward and in the flame side, plotted with the lowest temperature
and simulated in Ansys, to get the LCF for different geometrical designs. A three-dimensional simulation is time expensive, and the study problem formulation is decomposed even further into a two dimensional by giving a simulation position in axial axis, specific for the LCF study. This selection is analyzed from the Tccool script of the aerodynamics and thermodynamics with the theoretical support. Accordingly, LCF is expected to have the largest impact influences by the thermal load in the combustion chamber due to the thermodynamics reactions, and such axial position will be selected in this study. At the bottom of the space nozzle, pressure drop of fuel is low relative to the inlet. Besides, the pressure force acts on the wall is two times larger than at the inlet due to equation (3.1), and can results in an increasing of maxfor cycles Nc. The study will first be simulated in the hottest region and second at the bottom, where the data will be analyzed and selected to create a safety factor that each geometry is forbidden to exceed.
A script CUMFAT is developed by GKN and used to evaluate the solid mechan-ical properties, such as strain and damage with the theoretmechan-ical support described in the section2.1.3. Normally, the stress is described in a three-dimensional tensor with the CUMFAT script specializes in search and select the highest value in a specific axis and return as result.
3.4
Prediction of surrogate model
A CAE program ModeFrontier is used to create response surfaces and predict surro-gate models, and to analyze behavior of the data. In this study, the surrosurro-gate models have chosen to be approximated and described as a second order polynomial in R2. The model Y of these for n dimensional is given by:
Y = β0+ n X i=1 βixi+ n X i=1 K X k≥i βikxkxi+ n X i=1 βiix2i (3.4)
For a well-defined polynomial in an interval a and b, the function says to be contin-uous if there exist a value |f (a) − f (b)| ≥ L|a − b| for an arbitrary Lipschitz constant L. The proof of continuity of a polynomial will not be performed, except the concept of understanding the surrogate model is essential for usage in the optimization [26]. If discrete data points can be interpolated into a surrogate model and able to form a polynomial that does not contradict the Lipschitz condition, then the discrete data points can be transformed into an unobservable continuous function and is differ-entiable.
The number of design is essential when the computational time is limited, with the need to reduce the simulation time and still be able to create a well-behaved surrogate model suits to the physics. An approximate number of design that needs for creating a polynomial is given by [27]:
28 Chapter 3. Implementation and Multidisciplinary optimization where mcis the number of coefficients, n is the number of variables and k, degree of polynomial. For 7 design variables and third order polynomial, number of designs required is 120 and for a second order polynomial, 36 designs.
In some case, the data is unable to describe as a second order polynomial, another surrogate model approximation will be applied. This approximation is experimental described the data as an exponential function, which is objectively linearized the data by natural logarithm before fitting the new data into a second order polynomial.
3.4.1 Natural behavior of the surrogate model
The difficulty in this study is to decide if the surrogate models are well behaved for describing the physical properties. Taking mass as an example, by describing as a function of geometries x, it is expected to have the third order polynomial, but in some case in a specific interval the function can vaguely be assumed and described as a second order polynomial. This is inconsistent with the theory, but as far the function has revealed the ”original” behavior of its physical, then this violation can be ignored. For instant, the Mass,
Mass ∝ x3[m3] (3.6)
This can be discussed even further, at this design space range it is not necessary for such physical function to fulfil the theoretical approach, which often is a global definition for the mass, x ∈ [0, inf). But if one has more points, then it is possible to create such surrogate model that fits to the third order polynomial and fulfills the theoretical point of view.
3.5
Multidisciplinary Optimization
The optimization problem will be engaged with two different multidisciplinary ap-proaches, ATC and MDF, for comparison between two different methods. A surro-gate model based optimization will be optimized using a non-linear optimization method and in this case the Sequential Quadratic Programming (SQP), with an im-plementation in MATLAB. Following from earlier, when optimizing, a restriction of search interval is needed with lower and upper bound which is this case defined with three different limitations considerable to the manufacturing aspect.
3.5.1 Optimization constraints
straints on the cooling system, listed in table3.2.
TABLE3.2: List of constraints and each of its bound with description
Constraint Description
g1 = ∆Pc− C1 ≤ 0 Pressure drop g2 = −∆Tc+ C2 ≤ 0 Temperature drop g3 = Tmax− C3≤ 0 Max temperature g4 =Damage − C4 ≤ 0 Damage
g5 = M − C5≤ 0 Mach number g6 =Producible − C6≤ 0 Producible g7 = Qtot− C7 ≤ 0 Heat pick-up
g8 = ∆1− C8 ≤ 0 Strain range top edge g9 = ∆3− C9 ≤ 0 Strain range bottom edge
For the strain range ∆1 and ∆3, their values are based on the strain for each cycle and not the max strain, a picture of them can be seen in figureC.1where the max strain max is more complicated to create a surrogate model. The notations 1 and 3 are the positions shown in figure3.2and max temperature, Tmaxis usually the temperature at point 1 where the temperature is the largest.
In this, a consideration to the manufacture cost will be considered with an esti-mation of the cost, Producible describes as:
Producible = vp1 x1 x1ref + vp2 x2+ x3 x2ref + x3ref + vp3 x3 x3ref + vp4 x5 x5ref (3.7) whereP
ivpi = 1 are weighted coefficients describe the weight of each variable’s
property in the producible, equation (3.7). In this case the weighted coefficients are assumed to be vp1 = 0.60, vp2 = 0.20, vp3 = 0.15and vp4 = 0.05, where the number
of channels x1 has the highest impact comparing with the rest of the factors. The vp2 and vp4 represent the weight of the sheets for both inner and outer sheet, where
the last vp3 is for the channel height. A concept of how to manufacture a sandwich
wall is shown in figure3.5 together with equation (3.7) has shown an incomplete model of the producible, where the times that requires and the fact that it becomes more difficult to manufacture when the size of geometries are small is pretermitted. Except the concept provides a relation of how the producible can be estimated.
30 Chapter 3. Implementation and Multidisciplinary optimization
3.5.2 Optimization problem formulations
The main goal is to optimize mass of the space nozzle with respect to design vari-ables and subject to constraints, with the problem formulations for both ATC and MDF are shown below.
FIGURE3.6: System Level.
Hierarchical ATC formulation is used to formulate the problem with constraints g1, g2and g3 reserve in the Mass main objective function, as for the remaining con-straints are selected to include in the penalty functions πk(c). This selection lies with an argument that the cooling system is based on the ∆Pc, ∆Tc and Tmax, that are the main constraints of the fuel and the inner wall toward the flame side, thereby, reduce number of the subproblems.
with respect to : X
subject to: g1(X) ≤ 0 g2(X) ≤ 0 g3(X) ≤ 0
where k = 4, 5, ..., 9 are equivalent to the number of constraints and c is a vector of sharing variables. An additional need is to minimize the penalty function AL, with the requirement to satisfy the constraints gk:
Minimize: πk(c) (3.9)
with respect to : X
subject to: gk(X) ≤ 0
Optimization tolerance between the design variables for each system is numerically chosen: T ol ≤ 0.001, where T ol = |XT arget− XResponse| with update factors γ = 0.8 and β = 2.1 [22]. These factors are presented in the theory section.
32 Chapter 3. Implementation and Multidisciplinary optimization
FIGURE3.8: ATC Algorithm.
The number of solvers become large for many constraints defined for the ATC formulation. On contrary, the MDF formulation its way of describing the problem is easier to formulate independent to the number of constraints which are gathered into a single solver, equation (3.10).
Minimize M ass(X) (3.10) with respect to : X subject to : g1(X) ≤ 0 g2(X) ≤ 0 .. . g9(X) ≤ 0
ficult to receive, when the solver is programmed only for finding a solution regard-less to the global solution. One of the ways to handle this is to evaluate many initial points, but this can be problematic for many design variables and constraints. Al-ternatively, such feasibility of global evaluation exists in the MATLAB’s Global Op-timization toolbox and will be utilized in this study. The MATLAB’s ”GlobalSearch” constructs the problem and initialize search for the best objective value.
Sequential Quadratic Programming (SQP) method will be applied into the opti-mization, where this requires the objective function and the constraints to be twice differentiable and continuous, and in our case the functions are described in second order polynomial or exponential function, which both are differentiable and contin-uous.
3.5.4 Different limitations
34 Chapter 3. Implementation and Multidisciplinary optimization
FIGURE3.9: An example concept and illustration on the different lim-its, Design Space (D.S.), Industrial (Ind.) and Industrial with
Toler-ance, are defined with function’s minimum.
3.6
Certainty of Surrogate model
Chapter 4
Data analysis and Results with
Discussion
Based on chapter 2 and implemented according to chapter 3, the main result of this thesis is revealed. This chapter contains three parts, data analysis, main optimization result and probabilistic analysis of uncertainty.
The results that are presented in the thesis have been normalized with reference values and modified in such a way that they are able to publish.
4.1
Data analysis
Results from creating the surrogate models using the ModeFrontier program are presented in this section with different optimization limitations.
4.1.1 Fitting of the surrogate model
The surrogate model described in equation (3.4) was created in ModeFrontier, where the original data contained ”bad data/design”. To improve the accuracy of the surro-gate model, designs had been manually excluded where the selection was depend-ing on the main behavior of the data points. The number of useful results accorddepend-ing to the three different studies with different LHS design spaces and design numbers, after manual exclusion are shown in table4.1.
TABLE4.1: Number of data left to create surrogate model for 0.5 to 1.5and 0.2 to 1.7 interval.
Nr. of design Nr. of removed data Nr. of remaining data Design range
50 3 47 0.5 − 1.5
50 8 42 0.2 − 1.7
200 90 110 0.2 − 1.7
36 Chapter 4. Data analysis and Results with Discussion
4.1.2 Quality of response surface and surrogate model
TABLE4.2: Fitting quality to the data for each variable using 200 de-signs between a LHS interval of 0.2 to 1.7.
Property Polynomial order R2 R2 adjust Tmax 2 0.9975 0.9963 ∆Tc 2 0.9989 0.9905 ∆Pc 2 0.9869 0.9808 ln(∆Pc) 2 0.9977 0.9966 Mass 2 0.9996 0.9994 ln(Damage) 2 0.9869 0.9808 ∆1 2 0.9941 0.9914 ∆3 2 0.9941 0.9914 Heat pick-up 2 0.9996 0.9994 Mach number 2 0.9982 0.9974
The fitting quality for 50 designs for 0.5 − 1.5 and 0.2 − 1.7 intervals are shown in tablesB.1 andB.2, where all the surrogate models have high number fitting qual-ity. With a good fitting quality, the effect on the physical outputs from the design variables are better known than with poor fitting quality. Pressure drop ∆Pc and Damage are described in a natural logarithmic scale due to their complexity. For the ∆Pcthe improvement of R2 for defining as the natural logarithm ≈ 1%, table4.2. The choice to logarithm the data points before creating the polynomial is based on the data’s character. One reason is to avoid negative value which is unphysical and not to be expected when defining the function as a second order polynomial.
4.1.3 The limitations used to optimize
TABLE4.3: Different limitations for each design variable.
4.2
Results of optimization
The number of local solvers evaluated and converged are shown in table4.4 and are increasing as the study range decreasing for the ATC approach. The explanation is that when the design range is decreasing, there is less number of local minima to converge on, and by using the same number of points in the optimizer, more points will then converge at the same local minimum. For the MDF, the number of local solvers have shown to be compatible between the limitations, where most of the results successfully converge at the same solution. But for the industrial limit and industrial limit with tolerance the outcomes have shown, that the number of local solver evaluated is inconsistent with the number of local solver converged to the solution, by the cause of difficulty of finding physical allowance in the particu-lar region. Adding into the statement, all the functions are described in one single problem and by suffocating its search interval the solver has also trouble to find the solution, due to the size of information the optimizer needs to handle.
An essential remark to the ATC solutions are that the results have converged in the first iteration which means that there is no penalty update applied to the solu-tion. This kind of result is not unexpected, since the optimizer evaluates many points and succeeds to find the ”best” solution with respect to all the conditions. If one needs to understand how the convergence behaves for the two MDO approaches, the ”GlobalSearch” can be removed when trying to find the solution. In other way is to suffocate the constraints even more while solving, in order to make the optimizer works harder to find the common solution.
TABLE4.4: Total number of local solver evaluated and converged in the MDF and ATC for different limits.
D.S. limit Ind. limit Ind. limit with Tol.
ATC MDF ATC MDF ATC MDF
Local Solver Converged
53 39/39 258 51/53 307 31/43
TABLE4.5: Result from the optimization for MDF and ATC, for dif-ferent limits. The results are relative to reference values.
D.S. limit Ind. limit Ind. limit with Tol. Design
Variable
ATC MDF ATC MDF ATC MDF
38 Chapter 4. Data analysis and Results with Discussion
TABLE 4.6: Result from the optimization for MDF and ATC, with a comparison between surrogate models and actual results without
in-dustrial limits. The results are relative to reference values.
ATC MDF
Physical and Pro-ducibility Surro. Model Actual Surro. Model Actual Mass 0.48 0.49 0.49 0.51 ∆Pc 0.86 0.81 0.86 0.86 ∆Tc 1.01 1.01 1.01 1.01 Tmax 0.94 0.95 0.93 0.95 ∆1 0.95 0.95 0.95 1.00 ∆3 1.00 1.00 1.00 1.00 Damage 0.25 0.33 0.21 0.34 Producibility 0.65 0.65 0.92 0.92 Mach number 0.96 0.92 0.90 0.87 Heat pick-up 1.00 1.01 1.01 1.01 max 1.55 1.15 ∆Tcondense 0.28 0.29
TABLE 4.7: Result from the optimization for MDF and ATC, with a
comparison between surrogate models and actual results with indus-trial limits. The results are relative to reference values.
ATC MDF
ATC MDF Physical and
Pro-ducibility Surro. Model Actual Surro. Model Actual Mass 1.01 1.04 1.04 1.04 ∆Pc 1.01 0.99 1.01 0.99 ∆Tc 1.01 1.01 1.01 1.01 Tmax 1.05 1.05 1.05 1.05 ∆1 1.10 1.10 1.10 1.10 ∆3 1.00 1.00 1.00 1.00 Damage 0.59 0.55 0.59 0.54 Producibility 0.74 0.74 0.73 0.73 Mach number 1.01 1.10 1.02 1.01 Heat pick-up 1.01 1.01 1.01 1.01 max 1.10 1.10 ∆Tcondense 0.57 0.56
Optimized designs can be seen in table4.5. The designs are different depend-ing on the limits defined in the optimization. For the D.S. limit solution, the result has shown an incompatibility between ATC and MDF approaches, but both lead to the similar objective function Mass. Using figureB.2, the Mass is less affected by increasing x1 compared to increasing x5. Differences in the results can be origi-nated from the problem formulation of the ATC and utilization of the ”GlobalSearch”. Combining these two facts, the ATC problem formulation is dictated by the main objective function and all initial conditions defined at the beginning. Where the children/subproblems equation (3.9), only requires to satisfy the responses from the main objective and its own constraints.
The industrial limit has shown the same behavior of x1 for ATC and MDF as earlier which can be directed to the above argument. In the opposite, the indus-trial limit with tolerance gives similar results between the two different MDO ap-proaches, which is not observed in the D.S and the industrial limit. Meaning that, in the defined lower and upper bounds there are fewer remaining solutions, and therefore, the two disparate formulations of optimization tend to seek and converge in the same optimum.
40 Chapter 4. Data analysis and Results with Discussion higher damage, see figures in Appendix.B. For high temperature in the inner wall at a constant pressure, ∆Tcondensewill increase inversely with the risk for gas phase change/condensation that can affect the space rocket’s specific impulse.
Another important measure to each solution is the max strain, max. The indus-trial limit with tolerance has the lowest max value which is also in an agreement with ∆1, meaning that the max strain is independent to the number of cycles, or it has found its steady state, such solution can be considered as preservable. The man-ufacture tolerance limit itself has limited the solution in the way that it is realizable. On the other hand, the industrial limit and the D.S. limit indicate an increase maxi-mum strain cycle by cycle. For many cycles, this may lead to material failure much earlier than for the case with a stable strain. Meaning for many cycles, the damage and strain are accumulated for each cycle and when it exceeds the limit it will crack.
(A) ATC without industrial limit. (B) MDF without industrial limit.
(C) ATC with industrial limit. (D) MDF with industrial limit.
(E) ATC with industrial limit and tolerance. (F) MDF with industrial limit and tolerance. FIGURE4.1: Optimized geometry from different limits and methods.
The optimized geometries in the combustion chamber of the cooling system are shown in figure4.1. These geometries have shown an incompatibility between meth-ods, which can indicate to that, there exist many local minima and to the MDO prob-lem formulation.
