Loïc Brevault
email: loic.brevault@onera.fr

Mathieu Balesdent
email: mathieu.balesdent@onera.fr

Nicolas Bérend
email: nicolas.berend@onera.fr

Rodolphe Le Riche
email: leriche@emse.fr

Decoupled UMDO formulation for interdisciplinary coupling satisfaction under uncertainty

Keywords: Multidisciplinary Design Optimization, Uncertainty, Interdisciplinary couplings, Polynomial Chaos Expansion

Abstract

At early design phases, taking into account uncertainty in the optimization of a multidisciplinary system is essential to establish its optimal characteristics and performances. Uncertainty Multidisciplinary Design Optimization (UMDO) methods aim at efficiently organizing not only the different disciplinary analyses, the uncertainty propagation, the optimization, but also the handling of interdisciplinary couplings under uncertainty. A new decoupled UMDO formulation (named Individual Discipline Feasible -Polynomial Chaos Expansion) ensuring the coupling satisfaction for all the instantiations of the uncertain variables is presented in this paper. Coupling satisfaction in instantiations is essential to guarantee the equivalence between the coupled and decoupled UMDO formulations and therefore to ensure the physical relevance of the obtained designs. The proposed approach relies on the iterative construction of Polynomial Chaos Expansions in order to represent, at the convergence of the optimization problem, the coupling functional relations between the disciplines. The performances of the proposed formulation is assessed on an analytic test case and on the design of a new Vega launch vehicle upper stage.

I. Introduction

Multidisciplinary Design Optimization (MDO) is a set of engineering methodologies to optimize systems that are modeled as a set of coupled disciplinary analyses. For example, a launch vehicle is customarily decomposed into interacting submodels for propulsion, aerodynamics, trajectory, mass and structure. Taking into account the different disciplines requires to model and manage the interactions between them all along the optimization process. Martins et. al. [START_REF] Martins | Multidisciplinary design optimization: Survey of architectures[END_REF] explained that the designers can improve system design and decrease design cycle cost and time by using MDO in the early design phases of a system. In these phases, the system architecture is not defined and the exploration of the design space requires repeated discipline simulations. In order to find the optimal conceptual design at an affordable computational cost, low fidelity disciplinary analyses are mostly employed, introducing model (epistemic) uncertainties. Handling the uncertainties at early design phases is essential to efficiently characterize the optimal system design and its performances because it can reduce time and cost of the next design phases. [START_REF] Zang | Needs and opportunities for uncertainty-based multidisciplinary design methods for aerospace vehicles[END_REF] Uncertainty Multidisciplinary Design Optimization (UMDO) aims at solving MDO problems in the presence of uncertainty. This induces several new challenges compared to deterministic MDO: uncertainty modeling, uncertainty propagation, optimization under uncertainty and the interdisciplinary coupling satisfaction under uncertainty. In this paper we focus on the last issue as it is essential to ensure the system physical consistency. In deterministic MDO, different strategies have been developed to manage disciplinary couplings. The interactions are handled through coupling variables and the system consistency is described as a set of interdisciplinary equations that have to be satisfied. Two types of problem formulations can be distinguished [START_REF] Balling | Optimization of coupled systems-a critical overview of approaches[END_REF] : the coupled versus the decoupled formulations. The coupled formulations [START_REF] Balling | Optimization of coupled systems-a critical overview of approaches[END_REF] ensure system consistency by repeating interdisciplinary equation solving (MultiDisciplinary Analyses or MDA) all along the optimization process. In the coupled formulations, the decision variables handled by the system optimizer are the design variables. The decoupled formulations [START_REF] Allison | On the use of analytical target cascading and collaborative optimization for complex system design[END_REF][START_REF] Balling | Optimization of coupled systems-a critical overview of approaches[END_REF][START_REF] Braun | Use of the collaborative optimization architecture for launch vehicle design[END_REF][START_REF] Sobieszczanski-Sobieski | Optimization by decomposition: a step from hierarchic to non-hierarchic systems[END_REF] relax interdisciplinary couplings during optimization and satisfy the system consistency only at the optimal design. The system optimizer handles the design variables and the coupling variables. The decoupled approaches have been proposed to avoid repeated time consuming MDA. UMDO formulations have to efficiently organize the different disciplinary analyses, the uncertainty propagation, the optimization, but also the handling of interdisciplinary couplings under uncertainty. In the design of a launch vehicle for example, uncertainties can arise from low fidelity models to estimate the stage dry mass or from a wind gust during the vehicle launch. The introduction of uncertainty in MDO problems increases the number of required discipline evaluations due to the uncertainty propagation. Among the different formalisms (probability theory [START_REF] Kolmogorov | F oundations of the Theory of Probability[END_REF] , interval analysis [START_REF] Moore | Introduction to interval analysis[END_REF] , evidence theory [START_REF] Dempster | Upper and lower probabilities induced by a multivalued mapping[END_REF] , possibility theory [START_REF] Zadeh | Fuzzy sets as a basis for a theory of possibility[END_REF] , etc.) that can model the uncertainty, we adopt here the probability formalism. Within this framework, coupled UMDO formulations have been proposed [START_REF] Koch | Facilitating probabilistic multidisciplinary design optimization using kriging approximation models[END_REF][START_REF] Oakley | Performance optimization of multidisciplinary mechanical systems subject to uncertainties[END_REF] in which Monte Carlo (MC) method combined with MDA is used to propagate uncertainty. For each instantiation of the uncertain variables, the interdisciplinary equations are solved (MDA). Whereas in MDO there is only one value for the coupling variables that satisfies the interdisciplinary coupling equations for a given design, in UMDO, the coupling variables have to satisfy the interdisciplinary equations for each instance of the uncertain parameters. The computational cost of the coupled approaches becomes prohibitive due to the required number of discipline evaluations (in a straightforward implementation, it is the number of samples in MC multiplied by the number of discipline evaluations needed for an MDA). In the literature, decoupled UMDO formulations have therefore been investigated. The existing decoupled UMDO formulations ensure the system consistency for the statistical moments of the coupling variables [START_REF] Liu | Probabilistic analytical target cascading: a moment matching formulation for multilevel optimization under uncertainty[END_REF][START_REF] Mcallister | Multidisciplinary robust design optimization of an internal combustion engine[END_REF] , at the Most Probable failure Point (MPP) [START_REF] Du | Sequential optimization and reliability assessment for multidisciplinary systems design[END_REF] or by constructing the coupling variable probability densities [START_REF] Sankararaman | Likelihood-based approach to multidisciplinary analysis under uncertainty[END_REF] . The objective of this paper is to propose a new decoupled UMDO formulation to satisfy the interdisciplinary couplings under uncertainty at the optimal design in instantiations. By instantiations we mean that the system couplings must be consistent for all uncertain variable instances and not just for some particular values (e.g., statistical moments, MPP, etc). The approach relies on an iterative construction of surrogate models of the coupling functional relations. At the optimum, the surrogate models represent these mappings as would do an MDA under uncertainty, ensuring the equivalence between the proposed decoupled and coupled formulations. The rest of the paper is organized as follows. Section 2 introduces MDO and UMDO problem formulations, notations and key concepts. Section 3 provides an overview of the existing methods to satisfy interdisciplinary couplings under uncertainty, highlighting advantages and drawbacks. Section 4 presents the proposed decoupled UMDO formulation, named IDF-PCE (Individual Discipline Feasible -Polynomial Chaos Expansion). In section 5, the proposed method is illustrated with an analytic test case and with the design of a new Vega launch vehicle upper stage.

II. Mathematical formulation of MDO and UMDO problems

A. MDO problem formulation

A general MDO problem can be formulated as follows:

Minimize F(z,y) (1)
with respect to z, y

subject to g(z, y) ≤ 0 (2) h(z, y) = 0 (3) ∀i = j, y ij = c ij (z i , y .i ) (4) z min ≤ z ≤ z max (5)
z is the design variable vector. We note z i the input design variable vector of the discipline i ∈ {1, ..., N } and z = N i=1 z i without duplication. F is the multiobjective function (also called performance) to be optimized, g is an inequality constraint vector and h is an equality constraint vector. In a multidisciplinary environment, the disciplines exchange coupling variables, y. The latters link the different disciplines to model the interactions enabling to evaluate the multidisciplinary system consistency. c ij (z i , y .i ) are coupling functions to calculate the output coupling variable vector which is calculated by discipline i and input to discipline j. y .i refers to all the input coupling variable vector of discipline i and y ij is the input coupling variable vector which is input to discipline i and output from discipline j. We note y = N i=1 y .i = N i=1 y i. without duplication. From the design variables and the input coupling variables of the discipline i, the output coupling variables are computed with the coupling function: c i. (z i , y .i ). The couplings between the disciplines i and j are said to be satisfied (also called feasible or consistent) when the following interdisciplinary system of equations is verified:

y ij = c ij (z i , y .i ) y ji = c ji (z j , y .j ) (6)
When all the couplings are satisfied, i.e. when Eqs.( 6) are satisfied ∀i ∈ {1, ..., N }, ∀j ∈ {1, ..., N }, j = i, the system is said to be multidisciplinary feasible. The satisfaction of the interdisciplinary couplings is essential as it is a necessary condition for the system model to be physically realistic.

B. MDO coupling satisfaction

In MDO, we can distinguish two categories of methods to satisfy the interdisciplinary couplings. [START_REF] Coelho | Bi-level model reduction for coupled problems[END_REF] • Coupled formulations perform a MultiDisciplinary Analysis (MDA) (Figure 1) to ensure the interdisciplinary couplings at each iteration of the system level optimization. In the rest of the paper, MDO coupled formulations are defined by the optimization problem (Eqs.(7-10)): (10) with y(z) the coupling variable vector satisfying the interdisciplinary equations. • Decoupled formulations impose equality constraints on the coupling variables in the MDO formulation at the system level (Eq. 4) to ensure the interdisciplinary couplings at the optimal design (Figure 2). In the rest of the paper, MDO decoupled formulations are defined by the optimization problem described by Eq. (1-5).

Minimize F (z,y(z)) (7) with respect to z subject to g (z,y(z)) ≤ 0 (8) h (z,y(z)) = 0 (9) z min ≤ z ≤ z max
An overview of the two classes of methods is provided in the next paragraphs. MultiDisciplinary Analysis (coupled formulations). MDA is an auxiliary analysis aiming to find an equilibrium between the disciplines by solving the system of interdisciplinary equations (Eqs.( 6)). Two classical MDA methods are distinguished: either the Fixed Point Iteration (FPI) or an auxiliary optimization process minimizing the residual of the interdisciplinary equations. [START_REF] Coelho | Bi-level model reduction for coupled problems[END_REF] FPI is an iterative procedure involving a loop between the disciplines with no control on the coupling variables (excepted for the initialization) which directly result from the discipline simulations. If the interdisciplinary set of equations defines a contraction mapping, FPI converges to consistent couplings [START_REF] Ortega | Stability of difference equations and convergence of iterative processes[END_REF] . Alternatively, MDA can be solved by minimizing the discrepancy between the input coupling variables and the different coupling outputs [START_REF] Coelho | Bi-level model reduction for coupled problems[END_REF] :

Minimize y 1. -c 1. (z, y .1 ) 2 +...+ y N. -c N. (z, y .N ) 2 (11)
with respect to y with y i. the input coupling variables of all the disciplines linked to discipline i. An efficient auxiliary optimization algorithm requires often fewer calls to the discipline i than FPI, as the optimization process chooses the steps more freely than FPI [START_REF] Sankararaman | Likelihood-based approach to multidisciplinary analysis under uncertainty[END_REF] . Newton-Raphson or staggered solution approach [START_REF] Felippa | Partitioned analysis of coupled mechanical systems[END_REF] are examples of optimization algorithms applied to MDA. More details on MDA can be found in [START_REF] Keane | Computational approaches for aerospace design: the pursuit of excellence[END_REF] .

Decoupled formulations. Instead of performing the MDA optimization (Eq. ( 11)) at each MDO process iteration in z, equality constraints can be imposed between the input and the output coupling variables in the MDO formulation at the same level as the constraints g and h: ∀i, ∀j = i, y ij = c ij (z i , y .i ) (Eq.( 4)).

In the equality constraint approach, the system level optimizer handles both the design variables and the input coupling variables. The equality constraints on coupling variables may not be satisfied at each iteration but they participate to the convergence.

These two methods of interdisciplinary coupling satisfaction have been incorporated in various MDO formulations that can be classified in four categories:

• Single level approaches by application of MDA: Multi Discipline Feasible (MDF), [START_REF] Balling | Optimization of coupled systems-a critical overview of approaches[END_REF] • Multi level approaches by application of MDA: Concurrent SubSpace Optimization (CSSO), [START_REF] Sobieszczanski-Sobieski | Optimization by decomposition: a step from hierarchic to non-hierarchic systems[END_REF] Bi-Level Integrated System Synthesis (BLISS), 51 • Single level approaches with equality constraints on the coupling variables: Individual Discipline Feasible (IDF), [START_REF] Balling | Optimization of coupled systems-a critical overview of approaches[END_REF] All At Once (AAO), 9 • Multi level approaches with equality constraints on the coupling variables: Collaborative Optimization (CO), [START_REF] Braun | Use of the collaborative optimization architecture for launch vehicle design[END_REF] Analytical Target Cascading (ATC). [START_REF] Allison | On the use of analytical target cascading and collaborative optimization for complex system design[END_REF] The multi level approaches have optimization processes at the system level and at the disciplinary level. MDO formulations satisfying the interdisciplinary equations with MDA ensure the system feasibility at each system level optimization iteration. Among the formulations relying on MDA, MDF is the most usual. [START_REF] Balesdent | A survey of multidisciplinary design optimization methods in launch vehicle design[END_REF] MDF is a single level optimization formulation in which the system performance is evaluated with a disciplinary iterative process. CSSO and BLISS formulations use MDA to ensure interdisciplinary couplings but enable decoupled discipline evaluations. IDF, CO, ATC and AAO are fully decoupled formulations with satisfaction of the couplings by additional variables and equality constraints in the formulations. The decoupled MDO formulations offer several advantages compared to MDF [START_REF] Balesdent | A survey of multidisciplinary design optimization methods in launch vehicle design[END_REF][START_REF] Martins | Multidisciplinary design optimization: Survey of architectures[END_REF] :

• The optimization process allows parallel analyses and/or disciplinary level optimizations,

• The multi-level methods facilitate the system optimization by distributing the problem complexity over the different dedicated discipline optimizations, • In the multi-level approaches, the discipline optimizers handle local design variables, decreasing the system level design space size which handles the shared design variables between several disciplines, • The number of calls to the computationally expensive discipline codes can be notably decreased by avoiding expensive MDA calculations.

However, compared to MDF, the decoupled MDO formulations require an appropriate interdisciplinary coupling handling and involve an optimization problem with more variables in total (the design variables plus the coupling variables that can be distributed among the system and the local disciplinary optimizers in the case of multi-level approaches) and more constraints. If uncertain variables are considered in a MDO problem, new formulations are required.

C. UMDO problem formulation

The introduction of uncertainty in a MDO problem leads to the new general UMDO problem [START_REF] Yao | Review of uncertainty-based multidisciplinary design optimization methods for aerospace vehicles[END_REF] :

Minimize Ξ [F(z, Θ [Y] , U)] (12)
with respect to

z, Θ [Y] subject to K [g(z, Θ [Y] , U)] ≤ 0 (13) ∀i = j, Y ij = c ij (z i , Y .i (Θ [Y]), U i ) (14) z min ≤ z ≤ z max (15)
Besides the coupling equations ( 14), we consider only inequality constraints. Important differences exist between the UMDO formulation and the classical MDO formulation:

• U is the uncertain variable vector. We note U i , the input uncertain variable vector of the discipline i and U = N i=1 U i without duplication. In this paper, it is assumed that the uncertain variables are modeled with the probability theory, and that the input variable distributions are known. Note that other uncertainty modelings exist such as evidence theory, 17 possibility theory 57 or interval analysis. [START_REF] Moore | Introduction to interval analysis[END_REF] Moreover, it is assumed that the design variables are deterministic variables, and all the uncertainties are represented by U. We note (Ω, σ Ω , P Ω ) be a probability space with Ω the sample space for U, a sigma-algebra σ Ω , and a probability measure P Ω . We note f U the probability density function of the uncertain variable vector U. A real-valued random variable U is a function:

U : Ω → R (16) w → U (w) (17)
To simplify the notation in the rest of the paper, for all the uncertain variables, the instantiation U (w) is noted u.

• Because of the presence of the uncertain variable vector U, the disciplinary coupling variable vector Y (Eq.( 14)) is also an uncertain variable vector. In the decoupled formulations, coupling variables have to be handled by the optimizer, however in presence of uncertainty the optimizer cannot directly handle the coupling variables. Θ [Y] are parameters characterizing the uncertain coupling vector Y. These parameters can be instantiations of the uncertain variables, the statistical moments, the parameters of the probability density function, etc.

• Ξ denotes the objective function measure vector. A performance measure reflects the uncertainty in the objective function to be optimized due to the presence of the uncertain variables. Within the probability formalism, the expected value E [F(z, Y, U)] or an aggregation of the expected value and the standard deviation are commonly used to quantify the uncertainty in the objective function [START_REF] Baudoui | Optimisation robuste multiobjectifs par modèles de substitution[END_REF] .

• Two classical measures of uncertainty exist for the constraints in a UMDO problem and are formulated as follows:

- )] and σ[g(.)] are the expected value and the standard deviation of the constraint functions. The robust formulation is based on the statistical moments of the inequality functions to ensure that despite the uncertainty, the system will stay feasible. k indicates the restriction of the feasible region to k standard deviations away from the mean values of the constraint functions.

K [g(z, Y, U)] = E [g(z, Y, U)]+kσ [g(z, Y, U)]. E[g(.
- )] stands for the measures of uncertainty for the inequality constraint functions. The uncertainty measures of the constraints have to be at most equal to Λ reqd . [START_REF] Agarwal | Uncertainty quantification using evidence theory in multidisciplinary design optimization[END_REF] It reflects the requirement for the optimized system to lie in the feasible region with a given reliability despite the uncertainty. As we assume that the uncertain variables are modeled within the probability theory,

K [g(z, Y, U)] = Λ [g(z, Y, U) ≥ 0] -Λ reqd . Λ[g(.
K i [g i (z, Y, U)] = P [gi(z,Y,U)≥0] -Λ reqdi = Fi f U (u)du -Λ reqdi ,
Ξ[F(z,Θ[Y],U)] κ[g(z,Θ[Y],U)] Y..=c..(z,Θ[Y],U)
Calculation(of(F,g Uncertainty(simulation Incorporating uncertainty in MDO raises a number of additional challenges that need to be addressed (uncertainty modeling, uncertainty propagation, optimization under uncertainty) [START_REF] Yao | Review of uncertainty-based multidisciplinary design optimization methods for aerospace vehicles[END_REF] . In this paper, we focus on the interdisciplinary coupling satisfaction for the decoupled UMDO formulations. Similarly to the MDO formulations, two approaches can be distinguished for the satisfaction of the interdisciplinary couplings in UMDO: coupled approaches which perform MDA under uncertainty (Figure 3), or decoupled approaches which impose equality constraints on the coupling variables (Figure 4). MDA under uncertainty has been investigated to solve the UMDO problem [START_REF] Mahadevan | Probability, reliability, and statistical methods in engineering design[END_REF] and is briefly discussed in the following.

D. UMDO coupled formulations (by performing MDA)

The most straightforward approach to ensure the coupling satisfaction in UMDO is to use Monte Carlo simulations (MC) to propagate uncertainty in MDA [START_REF] Koch | Facilitating probabilistic multidisciplinary design optimization using kriging approximation models[END_REF][START_REF] Oakley | Performance optimization of multidisciplinary mechanical systems subject to uncertainties[END_REF] (Figure 3). The following formulation is a coupled UMDO formulation:

Minimize Ξ [F(z, Y(z, U), U)] (18)
with respect to z

subject to K [g(z, Y(z, U), U)] ≤ 0 (19) z min ≤ z ≤ z max (20)
In MC, one sample corresponds to one instantiation of the uncertain variables. For a given design variable vector z 0 , to evaluate the performance measures Ξ [F(z 0 , Y(z 0 , U), U)], it is necessary to propagate the uncertainty from the input discipline uncertainties characterized by U to the objective function. In the coupled formulations, Y is the coupling variable vector satisfying the following system of interdisciplinary equations:

y ij = c ij (z i , y .i , u) y ji = c ji (z j , y .j , u) (21)
∀ u ∈ Ω and ∀i ∈ {1, ..., N }, ∀j = i ∈ {1, ..., N }, with N the number of disciplines. We assume that for a given instantiation of the uncertain variables u 0 , there exists a unique set of coupling variables such that the couplings satisfy: ∀i ∈ {1, ..., N }, ∀j ∈ {1, ..., N }, i = j, ∃!y ij | y ij = c ij (z 0 , y .i , u 0 ). To compute the uncertainty measure of the performance Ξ[F(z 0 , Y(z 0 , U), U)], repeated MDA are performed for a set of uncertain variable instantiations sampled by MC. This may be computationally expensive due to the repeated evaluations of the disciplines. The computational cost of MDA under uncertainty corresponds to that of one MDA multiplied by the number of uncertain variable instantiations [START_REF] Mahadevan | Probability, reliability, and statistical methods in engineering design[END_REF] . To overcome the computational cost barrier introduced by the repetitive MDAs, Du et al. [START_REF] Du | Efficient uncertainty analysis methods for multidisciplinary robust design[END_REF] propose to replace the MDA by a surrogate model. The surrogate model is obtained by a first order Taylor series expansion and serves to estimate the first two statistical moments (the mean and the variance). This approach allows to model the coupling variables as a function of the uncertain variables. The method enables to find the optimal design while ensuring interdisciplinary couplings for all the uncertain variable instantiations. However, the method has several limits: first order Taylor approximation is only valid for functions that can be locally approximated as linear functions and the method requires to perform a MDA to locally build the surrogate model.

An alternative way to avoid to perform MDA under uncertainty is to use the same approaches as in deterministic MDO and to perform UMDO on the decoupled multidisciplinary system (Figure 4). It would allow to benefit from the same advantages as the deterministic decoupled MDO formulations highlighted in the previous section. To guarantee the interdisciplinary coupling satisfaction at the UMDO problem optimum, three approaches are distinguished in the literature:

• Coupling satisfaction by moment-matching for the uncertain coupling variables [START_REF] Liu | Probabilistic analytical target cascading: a moment matching formulation for multilevel optimization under uncertainty[END_REF][START_REF] Mcallister | Multidisciplinary robust design optimization of an internal combustion engine[END_REF] ,

• Coupling satisfaction at the Most Probable Point of failure (MPP) for the uncertain coupling variables 20, 21 , • Coupling satisfaction in instantiations by constructing the coupling variable probability distribution densities [START_REF] Sankararaman | Likelihood-based approach to multidisciplinary analysis under uncertainty[END_REF] .

These three approaches are further discussed in the next section, highlighting the main characteristics, advantages and drawbacks.

III. Existing approaches for interdisciplinary coupling satisfaction in decoupled UMDO formulation

In the presence of uncertainty, the coupling variables become uncertain variables and several methods have been introduced in order to enable the optimizer to handle these uncertain variables. One of the earliest proposition was done by Gu et al. [START_REF] Gu | An investigation of multidisciplinary design subject to uncertainty[END_REF] propagating worst case linearized uncertainty estimates of coupling variables. The worst case approach provides conservative results [START_REF] Gu | An investigation of multidisciplinary design subject to uncertainty[END_REF] and consequently degrades the objective function due to an overestimation of the uncertainty effects. The moment-matching approach has been introduced to evaluate the uncertainty encountered in the design process of multidisciplinary systems [START_REF] Liu | Probabilistic analytical target cascading: a moment matching formulation for multilevel optimization under uncertainty[END_REF][START_REF] Mcallister | Multidisciplinary robust design optimization of an internal combustion engine[END_REF] . This method is detailed in the next paragraph.

A. Moment-matching approach for interdisciplinary coupling satisfaction

In reference [START_REF] Mcallister | Multidisciplinary robust design optimization of an internal combustion engine[END_REF] , the authors investigated hierarchical frameworks for the uncertainty propagation in MDO. The authors developed a Robust Collaborative Optimization (RCO) formulation in which the mean and the variance of the random coupling variables are used to establish the subsystem-level objectives and the interdisciplinary compatibility constraints. By handling the statistical moments of the coupling variables Y, the optimizer controls a finite number of parameters that describe an infinite number of variable values. This approach has been used in reference [START_REF] Xiong | A moment-matching robust collaborative optimization method[END_REF] and adapted in other UMDO formulations such as Probabilistic ATC [START_REF] Liu | Probabilistic analytical target cascading: a moment matching formulation for multilevel optimization under uncertainty[END_REF] . The general UMDO formulation with a matching-moment approach is:

Minimize Ξ [F(z, M 1 [Y] , ..., M q [Y] , U)] (22)
with respect to

z, M 1 [Y] , ..., M q [Y] subject to K [g(z, M 1 [Y] , ..., M q [Y] , U)] ≤ 0 (23) ∀i, ∀j = i, ∀l ∈ {1, ..., q}, M l [Y ij ] = M l [c ij (z i , M 1 [Y .i ] , ..., M q [Y .i ] , U i )] (24) z min ≤ z ≤ z max (25)
with M l the l th statistical moment vector. Most of the formulations 6, 40 consider the first two moments of the coupling variables. Eq.( 24) imposes equality constraints between the statistical moments of the input coupling variables Y ij and the output coupling variables computed with the discipline simulations c ij (z i , Y .i , U i ). These constraints ensure the interdisciplinary compatibility for the statistical moments of the uncertain coupling variables. The statistical moments can be computed by different methods: Taylor series expansion, Monte Carlo simulation, quadrature rules. These moment-matching UMDO formulations are closed to the deterministic MDO formulations and present some advantages. First, they allow to decouple the discipline simulations in the MDO way. Second, at the optimum the statistical moments of the couplings are satisfied. Third, if Taylor approximations are used to compute the first two statistical moments, the computational cost is reduced compared to MDA with MC uncertainty propagation. Fourth, this formulation introduces a limited number of new coupling variables handled by the optimizer. The principal drawback of these formulations is that they do not ensure the interdisciplinary coupling satisfaction for all possible instantiations of the uncertain variables, i.e., they do not match the probability density function of the coupling variables. To be equivalent to MDF under uncertainty in terms of coupling satisfaction, a decoupled approach has to ensure the coupling satisfaction for each uncertainty variable instantiation.

B. Coupling satisfaction based on the Most Probable failure Point in UMDO

In the context of UMDO formulation involving probability constraints, several approaches have been proposed to allow the evaluation of the UMDO objectives and constraints. They require Multidisciplinary Reliability Analysis (MRA) to estimate the inequality constraint Λ g [g(z, Y, U) ≥ 0] -Λ reqdg ≤ 0. Again, an appropriate incorporation of MRA in the UMDO formulation is essential as the MRA is the most computationally expensive step [START_REF] Ahn | Sequential approach to reliability analysis of multidisciplinary analysis systems[END_REF][START_REF] Du | Sequential optimization and reliability assessment for multidisciplinary systems design[END_REF] . Several approaches UMDO formulations exist 7 :

• Nested loop approach 32 : at each iteration of the system level optimizer, MRA is performed,

• Single loop approach [START_REF] Agarwal | An inverse-measure-based unilevel architecture for reliability-based design optimization[END_REF][START_REF] Dormohammadi | Reliability-based design optimization within analytical target cascading framework[END_REF] : the nested loop is replaced by a single loop by transforming the reliability constraints with the Karush-Kuhn-Tucker optimality conditions, • Sequential approach [START_REF] Du | Sequential optimization and reliability assessment method for efficient probabilistic design[END_REF][START_REF] Du | Sequential optimization and reliability assessment for multidisciplinary systems design[END_REF] : the system level optimization and the evaluation of the constraints with MRA are performed sequentially.

We focus in the next paragraphs on the sequential approach as it is the most used approach to solve decoupled UMDO problems in the literature [START_REF] Ahn | An efficient strategy for reliability-based multidisciplinary design optimization using bliss[END_REF][START_REF] Chiralaksanakul | Decoupled approach to multidisciplinary design optimization under uncertainty[END_REF][START_REF] Du | Sequential optimization and reliability assessment for multidisciplinary systems design[END_REF][START_REF] Li | A hierarchical hybrid strategy for reliability analysis of multidisciplinary design optimization[END_REF][START_REF] Li | An efficient strategy for multidisciplinary reliability design and optimization based on csso and pma in sora framework[END_REF][START_REF] Zhang | A collaborative approach for multidisciplinary systems reliability design and optimization[END_REF] .

Coupling satisfaction at the Most Probable failure Point (MPP) with a sequential approach Du et al. [START_REF] Du | Sequential optimization and reliability assessment for multidisciplinary systems design[END_REF] proposed the Sequential Optimization and Reliability Assessment (SORA) for UMDO. In this approach, the optimization of the system level objective Ξ[F(.)] and the MRA are performed sequentially.

The UMDO problem is decomposed in a sequence of deterministic MDO problems. SORA 20 replaces the probabilistic reliability constraints by a deterministic approximation of the reliability constraints evaluated at the MPP. MRA is performed by First Order Reliability Method (FORM) [START_REF] Rackwitz | Reliability analysisa review and some perspectives[END_REF] to find the MPP (noted u * ).

The MPP is the closest failure point to the origin in the standard normalized space in which the variable U are transformed into Ȗ. The MPP is noted ȗ * in the standard space. To distinguish the original and the transformed inequality constraint vector, we denote by g the transformed constraint vector such that with In SORA [START_REF] Du | Sequential optimization and reliability assessment for multidisciplinary systems design[END_REF] , four steps are distinguished (Figure 5):

FORM 1 : g(z, y, ȗ * ) = Λ g [g(z, y, u) ≥ 0] -Λ reqd (26
• Step 1: the first deterministic MDO problem is solved with the uncertain variables fixed at their mean value. At the k th SORA iteration, the deterministic MDO problem is solved with the uncertain variables fixed at the MPP found at the (k -1) th iteration.

• Step 2: MRA is performed to identify the MPP u * (k) of all the inequality constraints by FORM with the design variables fixed at the optimal design z (k)

opt found at the k th iteration of Step 1. The objective function is computed: F z (k) opt , y c z (k) opt , u * (k) , u * (k) ,
U = u * (k) .
The deterministic MDO problem of Step 1 can be solved with the classical decoupled MDO methods (IDF, AAO, BLISS, ATC, etc). With the IDF method, the deterministic MDO problem at the SORA k th -cycle, (k ≥ 2), is formulated as follows [START_REF] Du | Sequential optimization and reliability assessment for multidisciplinary systems design[END_REF] :

Minimize F(z, y, u * (k-1) ) ( 27
)
with respect to z, y subject to g(z, y, u * (k-1) ) ≤ 0 (28)

∀i = j, y ij = c ij (z i , y .i , u * (k-1) i ) (29) z min ≤ z ≤ z max (30)
Once converged, the deterministic IDF by SORA ensures the system feasibility for the coupling variables at the MPP. After the design variables z

(k)
opt have been found, MRA is performed based on FORM in the standard uncertain space [START_REF] Chiralaksanakul | Decoupled approach to multidisciplinary design optimization under uncertainty[END_REF][START_REF] Du | Sequential optimization and reliability assessment for multidisciplinary systems design[END_REF] :

Maximize g(z (k) opt , y, ȗ) (31)
with respect to ȗ, y subject to ȗT ȗ

1 2 = -Φ -1 (Λ reqd ) (32) ∀i = j, y ij = c ij (z (k) opt , y .i , ȗi ) (33)
with Φ(.) the standard Gaussian cumulative distribution function. This optimization provides the MPP value ȗ * (k) for the uncertain variables at the SORA k th -cycle. MRA is performed on a decoupled multidisciplinary system and the interdisciplinary couplings are satisfied at the MPP in the standard normalized space (Eq. ( 33)). By decoupling MRA from the deterministic MDO, SORA tends to decrease the number of calls to the disciplinary functions compared to the nested approach [START_REF] Du | Sequential optimization and reliability assessment for multidisciplinary systems design[END_REF] . SORA with CO [START_REF] Li | A hierarchical hybrid strategy for reliability analysis of multidisciplinary design optimization[END_REF][START_REF] Zhang | A collaborative approach for multidisciplinary systems reliability design and optimization[END_REF] , CSSO [START_REF] Li | An efficient strategy for multidisciplinary reliability design and optimization based on csso and pma in sora framework[END_REF][START_REF] Zhang | A collaborative strategy for reliability-based multidisciplinary design optimization[END_REF] or BLISS 4 has been implemented, but the coupling satisfaction relies on the same approach: satisfaction at the MPP of the coupling variables.

The interdisciplinary coupling satisfaction within SORA presents several advantages:

• Possibility to perform disciplinary analysis in parallel,

• Satisfaction of the interdisciplinary couplings at the MPP value for the uncertain and coupling variables at the optimum, • Reduction of the computational cost with a sequence of deterministic MDO and MRA compared to the nested approach involving MDA with MC.

However, SORA has also several limitations. MRA is performed by FORM which locally linearizes in the standard space the inequality constraints that can lead to inaccurate estimation of the probability of failure. FORM also assumes the uniqueness of the MPP which is a limiting hypothesis in practical applications [START_REF] Dubourg | Metamodel-based importance sampling for structural reliability analysis[END_REF] . Furthermore, in terms of interdisciplinary coupling satisfaction, the couplings are ensured only at the MPP. The method does not ensure that the couplings will be satisfied for any other instantiation of the uncertain variables. The coupling is verified only at one point, the MPP.

C. Coupling satisfaction in instantiations with probability density estimation

Instead of satisfying the couplings at the MPP, Sankararaman et al. proposed to satisfy the couplings in instantiations by constructing the coupling variable probability densities [START_REF] Sankararaman | Likelihood-based approach to multidisciplinary analysis under uncertainty[END_REF] . The proposed approach has not been applied to a UMDO problem although it could. Let consider the figure 6 representing a partially decoupled multidisciplinary system. The couplings from the discipline i to the discipline j are removed. For each instantiation of the uncertain variables U and the coupling variables Y ij , it is possible to compute:

c ij (z i , u i , c ji (z j , u j , y ij )) and e ij (z, u, y ij ) = y ij -c ij (z i , u i , c ji (z j , u j , y ij )
) . E ij is the error between the input coupling variables and the corresponding output coupling variables. In this partially decoupled approach, y ji = c ji (z j , u j , y ij ) results directly from the simulation of the discipline j. To simplify, we denote c ij (z i , U i , c ji (z j , U j , Y ij )) by c ij (.). For a particular input instantiation: (u 0 , y 0 ij ), the couplings between the disciplines i and j are satisfied if

e 0 ij = y 0 ij -c ij (z i , u 0 i , y 0 ji ) = 0 with y 0 ji = c ji (z j , u 0 j , y 0 ij ). Discipline i Discipline j Yji cij Eij=Yij-cij Ξ[F(z,Y,U)] κ[g(z,Y,U)] z U Yij Yij
Y ij |E ij (.) = 0, 47 f ij (.) = f [cij (.)|Eij (.)=0] (.) (34)
To satisfy the couplings in instantiations, the error E ij must be equal to 0 for any instantiation of the uncertain variables. To construct f ij (.), the authors consider a particular instantiation y 0 ij . The uncertainties due to U are propagated in the partially coupled system and the authors assessed the probability of satisfying the couplings conditioned to the particular instantiation y 0 ij : P (c ij(.) = y 0 ij |y 0 ij ). The authors highlighted that this conditional probability is proportional to the likelihood L(y 0 ij ) ∝ P (c ij (.) = y 0 ij |y 0 ij ). However, the probability of having c ij (.) = y 0 ij is theoretically equal to zero. To overcome this issue, an infinitesimal interval [START_REF] Pawitan | In all likelihood: statistical modelling and inference using likelihood[END_REF] Finally, the conditional input density probability is given by Bayes' theorem [START_REF] Sankararaman | Likelihood-based approach to multidisciplinary analysis under uncertainty[END_REF] :

[y 0 ij -2 , y 0 ij + 2 ] is considered, leading to: L(y 0 ij ) ∝ f cij (.) (c ij (.) = y 0 ij |y 0 ij ).
f ij (y ij |y ij = c ij (.)) = f cij (.) (c ij (.) = y ij |y ij )f Yij (y ij ) f cij (.) (c ij (.) = y ij |y ij )f Yij (y ij )dy ij (35) = L(y ij )f Yij (y ij ) L(y ij )f Yij (y ij )dy ij (36)
Sankararaman et al. assumed a uniform prior distribution for y ij leading to:

f ij (y ij |y ij = c ij (.)) = L(y ij ) L(y ij )dy ij ( 37
)
This approach is interesting for the following reasons:

• The couplings are satisfied in instantiations and not just at the MPP or at the mean of Y, • The multidisciplinary system is partially decoupled allowing to decrease the number of calls to the disciplines.

However, this approach has limits in terms of interdisciplinary coupling satisfaction:

• The functional dependency between the uncertain variables U and the input coupling variables Y ij is not considered in the construction of the input coupling probability density functions. Normally, Y depends on z and U. This functional dependency is essential because it ensures that for a particular instantiation of the uncertain variables there exists a unique converged value for the coupling variables as with MDA under uncertainty, • In Bayes' theorem, a uniform prior distribution for Y ij is assumed and can be a non valid assumption.

To overcome these limitations, it would be necessary to model and parametrize the joint probability distribution between U and Y based on their respective marginals and a copula that would model the dependency between these variables. However, parametrizing a copula modeling a functional dependency arises numerous challenges as the copula must not model a statistical dependency but a functional one.

In the next section, we propose an interdisciplinary coupling satisfaction method for decoupled UMDO formulations in which the couplings are satisfied for all the instantiations of the uncertain variables at the UMDO optimal design. In this approach, the coupling functional relations are modeled with surrogates in order to incorporate the dependency between Y and U at a given z.

IV. Decoupled and consistent UMDO formulations

In the proposed approach, interdisciplinary coupling satisfaction is ensured for all the instantiations of the uncertain variables to ensure that the system is multidisciplinary feasible. Moreover, it does not require any MDA, allowing to reduce the number of calls to the computationally expensive disciplines. The proposed method relies on the iterative construction of surrogate models of the coupling functional relations along with the system level UMDO optimization. At the UMDO optimum, the metamodels of the coupling functional relations simulate these mappings as would MDA under uncertainty does. The surrogate models, written ŷ(U, α), provide a functional representation of the dependency between the uncertain variables U and the input coupling variables with α the vector of the metamodel parameters. This approach ensures that the functional dependency between the uncertain variables and the coupling variables is taken into account. Note that, in order to keep ŷ simple, the dependency between ŷ and z is not described here: ŷ is not a function of z, it is learned for the specific z where the optimization converges. We propose to model the coupling functional relations with Polynomial Chaos Expansion (PCE) because this surrogate model presents advantages in terms of uncertainty analysis and propagation [START_REF] Eldred | Recent advances in non-intrusive polynomial chaos and stochastic collocation methods for uncertainty analysis and design[END_REF] . An overview of PCE is provided in the next paragraph.

A. Polynomial Chaos Expansion (PCE) as a functional representation of the uncertain coupling variables

Let consider a function c : Ω → R depending on U ∈ Ω, a vector of uncertain variables in the sample space Ω. PCE is an approximation of the function c on its entire domain of definition [START_REF] Eldred | Recent advances in non-intrusive polynomial chaos and stochastic collocation methods for uncertainty analysis and design[END_REF] . PCE consists of an expansion of the function c over a polynomial orthogonal basis [START_REF] Hosder | Stochastic response surfaces based on non-intrusive polynomial chaos for uncertainty quantification[END_REF] :

c(U) = a 0 + ∞ i1=1 a i1 P 1 (U i1 ) + ∞ i1=1 i1 i2=1 a i1,i2 P 2 (U i1 , U i2 ) + ... (38)
with {P 1 , P 2 , ..., P r , ...} a basis of orthogonal polynomials, with P r of degree r. The choice of the polynomial basis is made consistently with the distribution of the input random variables. The polynomial basis is orthogonal to the weighting function [START_REF] Eldred | Recent advances in non-intrusive polynomial chaos and stochastic collocation methods for uncertainty analysis and design[END_REF] of the input uncertain variable distributions. In practice, the expansion (Eq.38) is truncated to a degree d and is reorganized (Eq.39) to have a one-to-one correspondence between the coefficients and the polynomials:

c(U) d j=0 α j Ψ j (U) = ŷ(U, α) (39)
where α j and Ψ j correspond to a i1,i2,...,ir and P r (U i1 , U i2 , ..., U ir ). Two types of truncation can be distinguished. The total order expansion includes a complete basis of polynomials up to a total order specification p. The number of PCE coefficients is given by: d + 1 = (n+p)! n!p! with n the number of uncertain variables and p the total order specification. The other approach, the tensor product expansion, does not bound the total order expansion but only truncates on a per-dimension basis. [START_REF] Eldred | Recent advances in non-intrusive polynomial chaos and stochastic collocation methods for uncertainty analysis and design[END_REF] It allows to have different truncation orders p i and therefore enables anisotropy in the polynomial order for each dimension. The number of PCE coefficients is given by: d+1 = n i=1 (p i +1). The difficulty in PCE is the estimation of the polynomial coefficients [START_REF] Poles | A polynomial chaos approach to robust multiobjective optimization[END_REF] . Intrusive and non intrusive approaches exist to compute the coefficients. The intrusive approaches require to modify the function c(.) in order to compute the PCE coefficients. In this paper we focus on non intrusive methods as the disciplines are considered as black box functions. Two main techniques can then be employed: the orthogonal spectral projection or the regression [START_REF] Xiong | Weighted stochastic response surface method considering sample weights[END_REF] . The orthogonal spectral projection consists in projecting the output c(.) on each polynomial basis function using the orthogonality [START_REF] Eldred | Recent advances in non-intrusive polynomial chaos and stochastic collocation methods for uncertainty analysis and design[END_REF] :

α j = < c, Ψ j > < Ψ 2 j > = 1 < Ψ 2 j > Ω c(u)Ψ j (u)f U (u)du (40)
with < ., . > the inner product on functions, < . [START_REF] Agarwal | Uncertainty quantification using evidence theory in multidisciplinary design optimization[END_REF] > the norm squared. The multivariate integral can be estimated by sampling or other the numerical integration methods. [START_REF] Eldred | Recent advances in non-intrusive polynomial chaos and stochastic collocation methods for uncertainty analysis and design[END_REF] The regression method relies on a least square method. Given M sample points {u [START_REF] Agarwal | An inverse-measure-based unilevel architecture for reliability-based design optimization[END_REF] , ..., u (M ) }, the polynomial chaos coefficients α = [α 0 , ..., α d ] are determined by [START_REF] Eldred | Comparison of non-intrusive polynomial chaos and stochastic collocation methods for uncertainty quantification[END_REF] :

min α M i=1 c(u (i) ) -ŷ(u (i) , α) 2 (41)
Xiong et al. [START_REF] Xiong | Weighted stochastic response surface method considering sample weights[END_REF] proposed to use the roots of the orthogonal polynomial basis as the sample points and to use a weighted least square regression to represent the higher contribution of sample points in the higher frequency region of the input random variables. PCE metamodeling is interesting in the uncertainty context as it provides analytic statistical moment formula linking the PCE coefficients to the output statistical moments. For instance:

µ ŷ d j=0 α j < Ψ j > (42)
σ 2 ŷ d j=0 α 2 j < Ψ 2 j > (43)
These analytical statistical moments converge to the real moments as the truncation degree increases. Higher statistical moments also have analytical expressions [START_REF] Eldred | Comparison of non-intrusive polynomial chaos and stochastic collocation methods for uncertainty quantification[END_REF] . The assessment of the fitting property of a surrogate model is measured with a loss function [START_REF] Vapnik | The nature of statistical learning theory[END_REF] L, which quantifies the error between the function modeled and the surrogate model for a particular u 0 . Different loss functions exist among which the commonly used square loss function [START_REF] Vapnik | The nature of statistical learning theory[END_REF] :

L(u 0 , α) = [c(u 0 ) -ŷ(u 0 , α)] 2 (44)
Statistical learning theory defines the generalization error J to quantify the error of the surrogate model as the expectation of the loss function [START_REF] Vapnik | The nature of statistical learning theory[END_REF] :

J = E (c(U) -ŷ(U, α)) 2 = Ω [c(u) -ŷ(u, α)] 2 f U (u)du (45)
In the proposed approach (IDF-PCE), IDF is adopted to decouple the multidisciplinary problem. The input coupling variables are modeled through PCE decomposed on an orthogonal polynomial chaos basis according to the input distribution of the uncertain variables. The generalization error is used to ensure the interdisciplinary coupling satisfaction in instantiations.

B. Proposed formulation: IDF-PCE

The IDF-PCE is a decoupled UMDO formulation for a single objective, relying on a functional representation of the coupling variables. IDF-PCE is formulated as follows:

Minimize Ξ [F (z, α, U)] (46)
with respect to z, α (ij) , ∀i = j

subject to K [g(z, α, U)] = E [g(z, α, U)] + kσ [g(z, α, U)] ≤ 0 (47)
J (mij ) = Ω c ij z i , u i , ŷ.i u, α (.i) -ŷij u, α (ij) 2 (mij )
f U (u)du = 0, ∀i = j (48)

z min ≤ z ≤ z max (49)
with J (mij ) the interdisciplinary constraint for the m th coupling from the discipline i to the discipline j and ŷ.i u, α (.i) the PCEs of all the input coupling variables of the discipline i. In this formulation, the system level optimizer handles the design variables and the PCE coefficients of the coupling variables α (ij) . The handling of the PCE coefficients at the system level allows to decouple the disciplines and to evaluate them in parallel (Figure 7). In comparison to the coupled formulations, the dimension of the design space is therefore increased by the number of parameters α (ij) . To ensure the multidisciplinary feasibility at the optimum, equality constraints involving the generalization error are imposed (Eq.( 48)). These constraints involve the input coupling variables modeled by PCE and the output coupling variables resulting from the discipline simulations. The constraints have an integral form to ensure the coupling satisfaction for all the possible instantiations of the uncertain variables. If we have: ∀i = j, J (mij ) (.) = 0, then the couplings are satisfied for all the instantiations u ∈ Ω, excepted perhaps for a finite number of instantiations. Indeed, due to the integral form of the generalization error, the loss function could be different from 0 only for a countable number of instantiations u ∈ Ω if the disciplines present some discontinuities. The vector J stands for the distances with respect to the MDA coupling satisfaction. Indeed, in the MDA approach, J(z, U) = 0, ∀(z, U) ∈ [z min , z max ] × Ω. In the IDF-PCE, J(z, U) = 0, ∀U ∈ Ω, has to be satisfied only at the UMDO optimum z = z opt . The interdisciplinary feasibility is not ensured all along the optimization. In practice, the multidimensional integrals associated to expectations, standard deviations and J (Eqs. ( 46), ( 47) and ( 48)) are difficult to compute. We propose three techniques to estimate these integrals, by Monte Carlo, by quadrature rules and by decomposition of the output coupling variables over a PCE. These three approaches are detailed in the next paragraphs.

Monte Carlo. MC approximation of Eq.( 48) relies on repeated sampling in the input uncertain space Ω. The integral in Eq.( 48) is approximated by:

J (mij ) J MC (mij ) = 1 M M k=1 c ij z i , u ik , ŷ.i u k , α (.i) -ŷij u k , α (ij) 2 (mij ) (50)
with M the number of samples. Similar estimations are carried out for the expectations and standard deviations in Eqs.( 46) and (47). J MC (mij ) is called the empirical error [START_REF] Vapnik | The nature of statistical learning theory[END_REF] and is related to the Root Mean Squared Error (RMSE) metric by:

RM SE = J MC (mij ) (51)
The smaller the RMSE, the better the interdisciplinary couplings are satisfied. The uncertain variables are sampled at the system level and are propagated in the decoupled subsystems (ŷ and c). Figure 8 illustrates the IDF-PCE method. MC is easy to implement and it can reach any level of accuracy if enough samples are calculated. The convergence of MC to the integral value is in order of 1 √ M . The convergence is slow and if the disciplines are computationally expansive, MC becomes intractable. Optimizer

α U12] α PCEyof Y21 PCEyof Y12 y21 y12 Subsystemy1 Subsystemy2 J12=yyy||y12-c12|| 1 M J21=yyy||y21-c21|| 1 M 2 2 Ξ[FUz,α,U]] κ[gUz,α,U]] c21 c12 z u ^Û ncertaintyysimulation y21 ^y12 ÛMC] UMC]
Figure 8: IDF-PCE for a coupled system of two disciplines with the generalization errors evaluated by Monte Carlo Quadrature rules. Quadrature rules approximate the integrals as a weighted sum of function values at specified points within the domain of integration [START_REF] Davis | Methods of numerical integration[END_REF] . Instead of sampling randomly as in MC simulation, for a function: c : Ω → R, with dim(Ω) = D, a set of specific points is used to approximate the multivariate integral:

Ω c(u)f U (u)du M1 i1=1 M2 i2=1 ... M D i D =1 (w i1 ⊗ w i2 ... ⊗ w i D )c(u i1 , u i2 , ..., u i D ) (52)
where w are weights and ⊗ is the tensor product operator. The simplest approximation to the integral is done through a tensor product quadrature with M 1 , ..., M D the number of specific points in each dimension. The points and the weights are computed according to the input variable uncertainty distributions : the points are the roots of the polynomials which are orthogonal (in the inner product sense) to the weighting function of the input uncertain variable distributions [START_REF] Davis | Methods of numerical integration[END_REF] . Although the proposed formulation stands for any Ξ that is an integral, we continue the article with the expected value of the performances Ξ [F (z, α,

U)] = E [F (z, α, U)].
The performances and the constraints are all multivariate integrals and can be computed by quadrature rules as follows:

E [F (z, α, U)] = Ω F (z, α, u) f U (u)du (53) M1 i1=1
...

M D i D =1 (w i1 ⊗ ... ⊗ w i D )F (z, ŷ ([u i1 , ..., u i D ], α) , [u i1 , ..., u i D ]) (54) E [g(z, α, U)] = µ g = Ω g (z, α, u) f U (u)du (55) M1 i1=1
...

M D i D =1 (w i1 ⊗ ... ⊗ w i D )g (z, ŷ ([u i1 , ..., u i D ], α) , [u i1 , ..., u i D ]) (56) σ [g(z, α, U)] 2 = Ω (g (z, α, u) -µ g ) 2 f U (u)du (57) M1 i1=1
...

M D i D =1 (w i1 ⊗ ... ⊗ w i D ) [g (z, ŷ ([u i1 , ..., u i D ], α) , [u i1 , ..., u i D ]) -µ g ] 2 (58)
J (mij ) = Ω c ij z i , u i , ŷ.i u, α (.i) -ŷij u, α (ij) 2 (mij ) f (u)du (59) M1 k1=1
...

M D k D =1 (w k1 ⊗ ... ⊗ w k D ) c ij z i , [u k1 , ..., u k D ] i , ŷ.i [u k1 , ..., u k D ], α (.i) -ŷij [u k1 , ..., u k D ], α (ij) 2 (mij ) (60)
The quadrature rule requires D i=1 M i discipline evaluations to propagate the uncertainty. Compared to MC, the approximation based on tensor product is efficient for a small number of input uncertain variables, but the method suffers from the curse of dimensionality. [START_REF] Eldred | Recent advances in non-intrusive polynomial chaos and stochastic collocation methods for uncertainty analysis and design[END_REF] Sparse grid approaches can be implemented in order to decrease the number of function evaluations while preserving the accuracy for high dimensional integrals. More details on the quadrature rules can be found in [START_REF] Davis | Methods of numerical integration[END_REF] .

Polynomial Chaos Expansion. Another way to impose coupling constraints and to benefit from the PCE uncertainty propagation is to decompose the output coupling variables with another PCE with respect to U. Instead of propagating the uncertainty into the subsystems by MC, an output PCE can be constructed to model the output coupling variables. We denote by system level PCEs, ŷ.i U, α (.i) , the PCEs of the coupling variables whose coefficients are handled by the system level optimizer. We also denote by discipline level PCEs, the PCEs of the coupling variables c ij z i , ŷ.i U, α (.i) , U i . Two further assumptions allow to simplify the expression of the discipline level PCE: like the system level PCE, the effect of the design variables z is not described (so the discipline level PCE is valid only at the convergenced design variables); Y .i is approximated by the system level PCE, ŷ.i U, α (.i) . Therefore, once α are chosen, the discipline level PCE is only a function of the uncertainties and its own coefficients, cij U, α(ij) . For the m th component of the disciplinary couplings between the disciplines i and j:

c ij z i , ŷ.i U, α (.i) , U i (mij ) d k=0 αk (m ij ) Ψ k (U) = cij U, α(ij) (mij ) (61)
The output PCE coefficients α(ij) correspond to the coefficients of the output coupling functions c ij (.) of the subsystem i (input PCEs + the discipline i ). They are calculated using the techniques detailed in the previous section. With the orthogonal spectral projection, we have:

αk (m ij ) = < c ij (m ij ) , Ψ k > < Ψ 2 k > = 1 < Ψ 2 k > Ω c ij z i , u i , ŷ.i u, α (.i) (mij ) Ψ k (u)f U (u)du (62)
We can approximate Eq.( 48) as:

J (mij ) = Ω c ij z i , u i , ŷ.i u, α (.i) -ŷij u, α (ij) 2 (mij ) f U (u)du (63) Ω cij u, α(ij) -ŷij u, α (ij) 2 (mij )
f U (u)du ( 64)

Ω d k=0 α(ij) k -α (ij) k Ψ k (u) 2 (mij ) f U (u)du (65)
To ensure the coupling satisfaction at the optimum between the disciplines i and j for all the uncertain variable instantiations based on PCE coefficients, and avoid the calculation of the integral, Eq.( 65) is replaced by: ∀i = j, α (ij) -α(ij) 2 = 0 (66)

If Eq.( 66) is satisfied then the integral ( 65) is null and Eq.( 48) is approximately satisfied. Indeed, if the system level and discipline level PCE coefficients are equals, therefore the sum in the integral of Eq.( 65) is composed of terms equal to zero. The constraints on the input and output PCE coefficients ensure the coupling satisfaction between the disciplines i and j for each uncertainty variable instantiation at the UMDO optimum under two conditions:

• ∀ u ∈ Ω, c ij z opt , u i , ŷ.i u, α (.i) -cij u, α(ij) < ,
with a tolerance and u an instantiation of the uncertain variables (good approximation of the coupling variables by the output PCE), • Eq. ( 66) is verified.

In this case, as the two uncertain functions cij (.) and ŷij (.) have the same polynomial chaos expansion, for any uncertain instantiation: cij (.) ŷij (.) c ij (.). The statistical moments of the objective function and of the constraint functions are calculated by quadrature rules as in the previous paragraph. The important difference with IDF-PCE (quadrature) lies in the quadratic coupling constraints Eq.(66) in IDF-PCE (PCE) which only involve the PCE coefficients and could facilitate the optimizer convergence. This approach allows to avoid the propagation of the uncertainty for a large number of uncertain variable instantiations as in MC, and to decrease the number of calls to the disciplines.

IDF-PCE presents several advantages:

• It ensures the coupling satisfaction at the optimum for each uncertain variable instantiation,

• It models the functional dependency between the uncertain variables, the design variables and the coupling variables, which may be useful beyond optimization, e.g. for a late sensitivity analysis, • It allows to decouple the UMDO problem and to simulate in parallel the disciplines, • It does not require any complete MDA, • The quadrature and the PCE approaches allow to propagate uncertainty more efficiently (using less calls to the disciplines) compared to MC.

However, this formulation presents some drawbacks:

• It increases the number of variables handled by the performance optimizer (design variables + PCE coupling variable coefficients), • And it increases the number of equality constraints at the system level in the optimization formulation.

V. Applications to analytic and launch vehicle design examples

A. Analytic test case

Numerical comparisons between the reference MDF (using MDA and Monte Carlo) and the proposed formulations are provided for an analytic test problem. First the mathematical problem is presented, then the four UMDO formulations are described and finally the results are analyzed. + 100 calculation of F and g:

F = 1 5 (z sh -5) 2 + (z 1 -3) 2 + (z 2 -7) 2 + (c 21 + z 1 * z 2 ) 0.6 + (U sh + 9) 2 g = 150 + exp(-U 2 1 * 0.01) * z sh * z 1 -0.02 * z 3 2 * U 5 2 + 0.01 * c 2.5 12 * z 2 * exp(-U sh * 0.1) • 3 design variables: z 1 ∈ [0, 1], z 2 ∈ [0, 1]
and the shared variable z sh ∈ [0, 1], we note:

z = [z sh , z 1 , z 2 ],
• 3 uncertain variables: U 1 = Uniform(-1,1), U 2 = Normal(0,1) and the shared uncertain variable U sh = Normal(0,1), we note:

U = [U sh , U 1 , U 2 ],
• 2 coupling variables:

Y = [Y 12 , Y 21 ],
•

1 objective function: Ξ[F (z, Y, U)] = E[F (z, Y, U)],
the expected value of the function F,

• 1 constraint function: E[g(z, Y, U))] + 3σ[g(z, Y, U))] ≤ 0
This analytic UMDO problem has been implemented for the following reasons:

• It has been empirically observed for all z and u values tried that this problem is such that MDA converges (it is a contraction mapping by FPI), and it converges to a unique coupling value.

• The disciplines are non linear and provide non gaussian coupling variable distributions.

MDF under uncertainty approach

The MDF under uncertainty approach is formulated as follows:

Minimize E[F (z, Y(z, U), U)] (67)
with respect to

z sh , z 1 , z 2 subject to E[g(z, Y(z, U), U))] + 3σ[g(z, Y(z, U), U)] ≤ g (68) z sh , z 1 , z 2 ∈ [0, 1] 3 ( 69
)
The dimension of the design space is 3. Uncertainties are propagated with Monte Carlo method on MDA for each instantiation of the uncertain variables with a sample size M = 150000 in order to have an error lower than 10 -3 on the estimation of the statistical moments. Because the means and the standard deviation are estimated by Monte Carlo, the objective function is noisy, therefore gradient based optimizers are not appropriate for this test case. Diverse derivative-free algorithms have been proposed in the literature to handle such problems 34 . To compare the three proposed methods to MDF with MDA we choose an Ant Colony optimizer (ACOmi) from the Matlab DOTcvp toolbox [START_REF] Hirmajer | Dotcvp: Dynamic optimization toolbox with cvp approach for handling continuous and mixed-integer do problems[END_REF] , as it is a derivative free optimization method that handles constraints with a penalization method. [START_REF] Schlüter | The oracle penalty method[END_REF] Optimizations are stopped if there are no progress in 500 consecutive objective function evaluations with a tolerance of 10 -3 on the objective function and the constraint. The MDA convergence criteria for the Fixed Point Iterations has been set to 10 -4 as it corresponds to a variation in the objective and constraint functions smaller than 10 -3 . g = -0.004 is a conservative tolerance due to the estimation of the mean and the standard deviation of the constraint by MC to ensure that the measure of the constraint is inferior or equal to 0.

IDF-PCE formulation

The general proposed IDF-PCE formulation is given by Minimize

E [F (z, α, U)] (70)
with respect to z sh , z 1 , z 2 , α [START_REF] Braun | Use of the collaborative optimization architecture for launch vehicle design[END_REF] , α [START_REF] Du | Sequential optimization and reliability assessment for multidisciplinary systems design[END_REF] subject to E[g(z, α, U)] + 3σ[g(z, α, U)] ≤ g (71)

J 12 = Ω c 12 z sh , z 1 , u sh , u 1 , ŷ21 u, α (21) -ŷ12 u, α (12) 2 f U (u)du ≤ (72)
J 21 = Ω c 21 z sh , z 2 , u sh , u 2 , ŷ12 u, α (12) -ŷ21 u, α (21) 2 f U (u)du ≤ (73) z sh , z 1 , z 2 ∈ [0, 1] 3 (74)
The methods to compute the multivariate integrals are detailed in the next paragraph.

IDF-PCE-MC formulation

In IDF-PCE with Monte Carlo the interdisciplinary constraints are computed with:

J 12 1 M M k=1 c 12 z sh , z 1 , u sh k , u 1 k , ŷ21 u k , α (21) -ŷ12 u k , α (12) 2 ≤ (75) J 21 1 M M k=1 c 21 z sh , z 2 , u sh k , u 2 k , ŷ12 u k , α (12) -ŷ21 u k , α (12) 2 ≤ (76)
The mean of the objective function and the mean and the standard deviation of the constraint g are computed by Monte Carlo. Moreover, the system level coupling variables are decomposed according to:

ŷij (U, α (ij) ) = d k=0 α (ij)
k Ψ k (U), with Ψ k the product of Hermite and Legendre polynomials of total degree k as these polynomial bases are orthogonal to the input density distributions (Gaussian and uniform). The PCE are decomposed with a total order expansion of degree 3. As there are three uncertain variables for the decomposition, dim α [START_REF] Braun | Use of the collaborative optimization architecture for launch vehicle design[END_REF] = dim α [START_REF] Du | Sequential optimization and reliability assessment for multidisciplinary systems design[END_REF] = (3+3)! 3!3! = 20. The design space is of dimension 43. The uncertainties are propagated with Monte Carlo on the decoupled system with a sample size of M = 150000. The interdisciplinary constraints on the couplings are such that = 10 -3 in order to have in average a coupling error under .

IDF-PCE-quadrature formulation

In the proposed decoupled formulation IDF-PCE with quadrature rules the constraints are computed as follows: [START_REF] Du | Sequential optimization and reliability assessment for multidisciplinary systems design[END_REF] (77)

J 12 = M sh i=1 M1 j=1 M2 k=1 (w shi ⊗ w 1j ⊗ w 2 k ) c 12 z sh , z 1 , u shi , u 1j , ŷ21 u shi , u 1j , u 2 k , α
-ŷ 12 u shi , u 1j , u 2 k , α (12) 2
J 21 = M sh i=1 M1 j=1 M2 k=1 (w shi ⊗ w 1j ⊗ w 2 k ) c 21 z sh , z 2 , u shi , u 2 k , ŷ12 u shi , u 1j , u 2 k , α (12)
(78)

-ŷ 21 u shi , u 1j , u 2 k , α (21) 2
The expected value of F is also approximated by the quadrature rules:

E[F (z, α, U)] M sh i=1 M1 j=1 M2 k=1 (w shi ⊗ w 1j ⊗ w 2 k )F z, ŷ u shi , u 1j , u 2 k , α , u shi , u 1j , u 2 k (79)
The mean and the standard deviation of the constraint g are computed as explained in the quadrature rules paragraph. The quadrature rules used to compute the multidimensional integrals correspond to the tensor product of the one dimensional Gauss-Hermite and Gauss-Legendre quadratures. The number of sampling points in each direction is: M sh = M 2 = 8 and M 1 = 10, resulting in a tensor product of 640 discipline evaluations to compute the multivariate integrals. This number of quadrature points is selected in order to have an error less than 10 -3 compared to a Monte Carlo computation of the integrals with 1 000 000 points. The decomposition of the coupling variables is the same as in the IDF-PCE with Monte Carlo formulation.

IDF-PCE (PCE) formulation

In this approach, the discipline level PCE coefficients α(ij) are computed by orthogonal spectral projection based on Eq.( 62) in which the multivariate integrals are computed by quadrature rules:

α(12) k = < c 12 , Ψ k > < Ψ 2 k > = 1 < Ψ 2 k > Ω c 12 z sh , z 1 , u sh , u 1 , ŷ21 u, α (21) Ψ k (u)f U (u)du (80) 1 < Ψ 2 k > M sh l=1 M1 m=1 M2 k=1 (w sh l ⊗ w 1m ⊗ w 2 k )c 12 z sh , z 1 , u sh l , u 1m , ŷ21 u sh l , u 1m , u 2 k , α (21) Ψ k (u sh l , u 1m , u 2 k ) (81)
The interdisciplinary constraints J 12 and J 21 are replaced by:

α (12) -α(12) 2 ≤ α (82) α (21) -α(21) 2 ≤ α (83)
To compute the discipline level PCE coefficients, we use the same quadrature rules as in the IDF-PCE (quadrature) formulation: M sh = M 2 = 8 and M 1 = 10, resulting in a tensor product of 640 discipline evaluations. The constraints on the couplings are such that α = 0.5 as it generates an error on the average smaller than 10 -3 . The main difference with the IDF-PCE quadrature formulation is in the constraints to ensure the interdisciplinary couplings in instantiations. In the IDF-PCE (PCE) the quadratic constraints only involve the PCE coefficients and could facilitate the optimizer convergence.

Results

Due to the presence of uncertainty and the use of a global optimizer, each optimization is repeated 5 times and the results given in Table 1 are the average of the 5 optimizations. The ratio of the standard deviation over the expected value of the results is added between bracket in order to quantify the robustness of the results.

Table 1: Analytic test case results for the different proposed IDF-PCE formulations. Between brackets, standard deviation over average of each result.

Results

MDF-MDA (Ref)

IDF-PCE (MC)

IDF-PCE (quadrature)

IDF-PCE (PCE)

Objective µF = 0.928 (0.64%) µF = 0.926 (0.65%) µF = 0.926 (0.70%) µF = 0.914 (0.49%) Design z sh = 0.520 (0.63%) z sh = 0.511 (0.86%) z sh = 0.514 (1.34%) In order to highlight the importance of incorporating uncertainty in MDO problems, a deterministic MDF optimization with the uncertain variables set to their mean values is performed. The optimal objective value is 0.466 and the set of optimal design variables is: z sh = 0.504, z 1 = 0.452, z 2 = 0.682. For this analytical test case, the presence of uncertainty modifies both the optimal objective and the set of optimal design variables.

z sh = 0.523 (1.03%) variables z1 = 0.340 (1.13%) z1 = 0.339 (1.11%) z1 = 0.340 (1.27%) z1 = 0.349 (1.13%) z2 = 0.658 (1.55%) z2 = 0.661 (1.30%) z2 = 0.661 (1.68%) z2 = 0.649 (0.95%) Coupling |c12 -y12| 2 ≤ 0.0001 J12 = 0.00067 (1.23%) J12 = 0.00054 (1.08%) α (21) -α(21) 2 = 0.48 (1.56%) constraints |c21 -y21| 2 ≤ 0.0001 J21 = 0.00057 (1.15%) J21 = 0.00074 (1.12)% α (21) -α(21) 2 = 0.3 (2.13%) Constraint K value -0.
The MDF under uncertainty formulation is considered as the reference for interdisciplinary coupling satisfaction as the couplings are satisfied for each uncertain variable instantiation. The convergence criteria for the Fixed Point Iterations is chosen in order to generate an error smaller than 10 -3 in both the objective and the constraints as the ACOMI convergence criteria is set to 10 -3 for the objective and constraint functions. MDF and the proposed formulations converge to the same design variables z with errors inferior to 2.04% in terms of distance to MDF results. IDF-PCE (quadrature) provides the smaller error. The higher error stems from IDF-PCE (PCE) due to the approximations introduced by the discipline level PCEs. In terms of objective values, the relative error compared to MDF is of 0.21% for IDF-PCE (MC) and IDF-PCE (quadrature) and 1.5% for IDF-PCE (PCE). The design space dimension is 3 in the MDF approach whereas it is 43 in the proposed formulations. The number of constraints is 1 in MDF whereas it is 3 in the proposed formulations. The increase in number of design variables and constraints leads to more optimization iterations. The number of calls to each discipline is 1512 * 10 6 for MDF. It is divided by 1.80 for IDF-PCE (MC) formulation. IDF-PCE (quadrature) and IDF-PCE (PCE) allow to decrease respectively to 3.52 * 10 6 and to 3.37 * 10 6 the number of calls to each discipline. Compared to MDF, it corresponds to a decrease by a factor of 430 in the number of calls to each discipline in IDF-PCE (quadrature) and by 449 in IDF-PCE (PCE). The reduction of the number of calls to the disciplines is due to the absence of complete MDA and the uncertainty propagation technique (quadrature and PCE) instead of Monte Carlo. While the absence of MDA decreases the number of calls to the disciplines, it generates higher errors in the interdisciplinary coupling satisfaction : the couplings are satisfied with a precision of 10 -4 in MDF for all the instantiations of the uncertain variables and with a precision of 6.7 * 10 -4 on the average in IDF-PCE (MC). The replacement of Monte Carlo in IDF-PCE enables to decrease the number of calls to the disciplines while ensuring coupling satisfaction with a precision of 7.4 * 10 -4 on the average. The distribution of the performance value given by MDF and by the proposed approaches are similar (Figure 11). The distributions of the coupling variable Y 12 have the same tendencies for MDF and the decoupled approaches. Moreover, it can be noted that the proposed approaches succeed to handle multimodal probability density for the coupling variables. However, differences can be noted in the distribution tails (Figure 11). This is due to the error introduced by the PCE metamodeling of the coupling relations. The distributions of the coupling errors for the proposed formulation are given Figure 12. The distributions have a small dispersion around 0, and the higher dispersion arises for IDF-PCE (PCE).

B. Launch vehicle design test case

The launch vehicle test case consists of the design of a new upper stage for the Vega launch vehicle. Vega is currently composed of three solid rocket motor stages (P80, zefiro-23 and zefiro-9A) and a liquid propulsion (UDMH/NTO) upper stage AVUM (Attitude Vernier Upper Module). This architecture enables to launch 1500kg into a 700km circular polar orbit. To comply with new demands in payload mass, several studies have been performed [START_REF] Sippel | Study trade-offs on future european expendable launchers[END_REF][START_REF] Wisse | Venus-conceptual design for vega new upper stage[END_REF] to enhance the vehicle capabilities. In reference [START_REF] Wisse | Venus-conceptual design for vega new upper stage[END_REF] , a sensitivity analysis has been conducted to take into account the main uncertainties at this early design phase.

In this test case, we are interested in designing a new upper stage for Vega with UMDO methodologies and in comparing coupled and decoupled formulations. The targeted mission is the launch of 1800kg into a 700km polar circular orbit. Five disciplines are involved: the propulsion, the mass budget and geometry design, the aerodynamics, the trajectory and the stage fallout.

Design and uncertain variables

The expected upper stage mass of Vega is minimized. The problem has six design variables: The problem incorporates also three uncertain variables:

• 1 st stage specific impulse: Isp 1 ∼ Uniform(174.5, 275.5) and MMH has higher density and slightly higher performance than UDMH (Unsymmetrical DiMethylHydrazine) [START_REF] Wisse | Venus-conceptual design for vega new upper stage[END_REF] . The chamber pressure (Pc), the mixture ratio (OF), the throat cross-sectional area (At) and the nozzle expansion ratio ( ) are inputs to the propulsion discipline. CEA computes the characteristic velocity (c* ) and the coefficient of thrust (c τ ) with shifting equilibrium until the nozzle throat and frozen expansion after. The outputs of the propulsion discipline are the upper stage thrust (T ), the upper stage specific impulse (Isp), the upper stage mass flow rate ( ṁ) assumed constant and the exit nozzle area (Ae).

Mass budget and geometry design

The mass budget and geometry design discipline aims at estimating the mass of the upper stage and its geometry. The 1 st stage is a P100 booster which replaces the P80 (in order to provide enough thrust to take off), the 2 nd stage is the Zefiro-23 and the 3 rd stage is the Zefiro-9 with their geometries and masses fixed to the existing baseline. Only the new upper stage is designed. The dry mass of the upper stage (M d ) is the sum of the masses of the tanks, turbopumps, combustion chamber, nozzle, pressurization system and structural masses. The mass and geometry models used are derived from the engineering models for the conceptual design of launch vehicle developed by Castellini [START_REF] Castellini | Multidisciplinary design optimization for expendable launch vehicles[END_REF] . Interdisciplinary couplings between the mass budget and geometry design discipline and the trajectory is required to model the dependencies between the dry mass and the loads undergone during the flight. Two load parameters constitute the couplings: the maximal axial load factor N ax max and the maximal dynamic pressure P dyn max .

Aerodynamics

The aerodynamics discipline consists in computing the aerodynamics coefficients such as the drag and lift coefficients required to compute the aerodynamics loads during the rocket atmospheric flight. A zero-lift model is used. The calculations of the drag coefficients are based on the US AirForce computer program MissileDATCOM [START_REF] Blake | Missile datcom: User's manual-1997 fortran 90 revision[END_REF] which relies on an experimental data base to determine the aerodynamic forces and coefficients of complex rocket geometries. It takes as inputs the geometry of the launcher, the angle of attack, the Mach and the altitude. As the payload fairing is fixed and only the upper stage is modified, the drag coefficient (C d ) is interpolated from the angle of attack (a) and the Mach (M ) based on MissileDATCOM simulations. The drag coefficient table is directly given to the trajectory discipline enabling to remove the feedback loop between the trajectory and the aerodynamics. This model is generally sufficient in the early design studies.

Trajectory

A three dimensional model with rotating round Earth is used. The flight control corrects all along the flight eventual non nominal propulsion events or non nominal environmental flight conditions in order to reach the mission injection point. A simplified approach is adopted here to take into account the uncertainties, a trajectory optimization is performed for each uncertain variable instantiation. First, an optimization yields a nominal feasible trajectory (by minimizing the distance to the target injection point) where the uncertain variables are set to their expected values. Then, for the uncertainty propagation, local trajectory optimizations are performed for each instantiation of the uncertain variables. The trajectory optimization consists in defining crossing points for the pitch angle and modifying them to satisfy the specifications of the mission. The pitch angle is calculated by piecewise linear function. The trajectory enables ballistic flight for the upper stage. The trajectory discipline computes the loads (N ax max and P dyn max ) required to simulate the mass budget and geometry design discipline. Moreover, to simulate the critical 3 rd stage fallout impact point in the fallout discipline, the state of the 3 rd stage at the moment of the separation (Sep states ) is computed.

Stage fallout

To protect life and property on Earth, it is necessary to control the fall back zone of the stages. On the Vega launcher, the 3 rd stage fallout is critical as it might fall back on the ground whereas the two first stages fall off at the north of Kourou. A three dimensional ballistic flight model with rotating round Earth is used to simulate the fallout.

UMDO formulations

Due to computational time constraints, only one run of MDF and IDF-PCE (quadrature) have been tested for this test case (3 weeks of calculation per formulation on a computer, Intel Xeon CPU E5-16030, 2.80GHz). The IDF-PCE (quadrature) formulation has been selected because among the proposed formulations it is the one with the most promising results highlighted in the analytic test case.

MDF formulation

Minimize

E[M 4 th (z, Y(z, U), U)] (84)
with respect to

z = [D, A t , P c , OF, M ox , ] subject to E[g(z, Y(z, U), U)] + 3σ[g(z, Y(z, U), U)] ≤ 0 (85)
where M 4 th is the mass of the upper stage (dry mass plus propellant masses) and g is the success mission constraint. The success mission constraint takes into account the constraints to reach the target injection point (in altitude 700km, in velocity 7005m/s and in flight path angle 0 • ) but also the 3 rd stage impact point that must be in the ocean.

Only the feedback couplings between the trajectory and the mass budget and geometry disciplines are decoupled. It corresponds to two couplings: the maximal axial load and the maximal dynamic pressure during the flight (red arrow in Figure 13). These two couplings correspond to feedback couplings in the Design Structure Matrix (Figure 13). We denote by N axmax (.) the maximal axial load factor which is an output of the trajectory discipline. In the same way, we denote by P dyn max (.) the maximal dynamic pressure. Naxmax (.) and P dyn max (.) are the system level PCEs and their PCE coefficients (α (Nax max ) , α (P dynmax ) ) are handled by the optimizer.

IDF-PCE formulation

Minimize

E [M 4 th (z, α, U)] (86)
with respect to z = [D, A t , P c , OF, M ox , ], α (Nax max ) , α (P dynmax )

subject to E[g(z, α, U)] + 3σ[g(z, α, U)] ≤ 0 (87)

J Nax max = Ω N axmax (z, y T raj. , u) -Naxmax u, α (Nax max ) 2 f U (u)du ≤ (88)
J P dynmax = Ω P dynmax (z, y T raj. , u) -Pdynmax u, α (P dynmax ) 2 f U (u)du ≤ (89)

Results

As in the analytic case, ACOmi is used in MDF and IDF-PCE (quadrature). In MDF, uncertainties are propagated with Monte Carlo with a sample size of M = 700 in order to have an error lower than 1kg in the estimation of the objective function. Due to computational expensive disciplines (especially the trajectory analysis which involves an optimization problem), a limited budget for the UMDO optimization of 15 * 10 4 calls to the disciplines has been chosen as the stopping criteria. The MDF and IDF-PCE (quadrature) optimizations are initialized with the same set of design variables (identified by deterministic MDO optimization) in order to be able to compare the results. In MDF, the design space is of dimension 6. The MDA is solved by Fixed Point Iterations and the convergence criteria is set to 0.1g for the maximal axial load factor and to 2.24Pa for the maximal dynamic pressure in order to have an error on the upper stage dry mass lower than 1kg. The FPI is initialized in order to converge on the average in two iterations.

The PCE are decomposed with a total order expansion of degree 2. As there are three uncertain variables, we have dim(α (Nax max ) ) = dim(α (P dynmax ) ) = (2+3)! 2!3! = 10. The design space is of dimension 26. To compute the multivariate integral by quadrature the tensor product, it requires 216 discipline evaluations. This number of quadrature points is selected in order to also have an error lower than 1kg in the estimation of the objective function compared to a Monte Carlo computation of the integral with 10 000 points.

Starting from the same initialization and with the same number of calls to the disciplines, IDF-PCE (quadrature) formulation finds an upper stage mass (2896.2kg) lower than the MDF formulation (2983.0kg) (Table 2). To ensure that the approximations introduced in IDF-PCE (quadrature) are acceptable, a Monte Carlo simulation (1000 samples) on a MDA has been performed at the design found by IDF-PCE (quadrature). The objective function found with MDA is E[M 4 th ] = 2895.7kg and the constraint value K = -0.0364. Figures 18 and 19 present the distribution of the coupling error with IDF-PCE (quadrature). J Nax max and J P dynmax are in the same order as the MDA convergence criteria, however compared to the coupled approach, the error in the couplings with IDF-PCE (quadrature) is higher. It can be explained by the choice of a PCE decomposition order of 2, and the UMDO optimization that does not have reached complete convergence. The coupling distributions are represented in Figures 22, 23, 24, 25 with the coupled and decoupled approaches at the optimal solution found by IDF-PCE (quadrature). The stage 3 impact points are presented in Figure 26. All the instantiations of the impact point for the 3 rd stage fall back in the ocean (Figure 26).

VI. Conclusions

This paper provides a new decoupled UMDO formulation that ensures the system multidisciplinary feasibility for all the instantiations of the uncertain variables. The satisfaction of the interdisciplinary couplings in instantiations is important to ensure the equivalence between the coupled and decoupled UMDO formulations and therefore the physical relevance of the obtained designs. An overview of existing decoupled UMDO formulations in the literature has been presented, highlighting the advantages and the drawbacks of each approach. The proposed Individual Discipline Feasible -Polynomial Chaos Expansion formulation is based on the iterative construction of surrogate models of the functional coupling relations. The PCE coefficients are handled by the optimizer and integral coupling constraints are imposed to ensure the coupling satisfaction. Three techniques have been used to compute the multidimensional integrals involved in the statistical moment calculations and the coupling constraints: Monte Carlo, quadrature rules and discipline level PCE decomposition. Numerical comparisons between MDF formulation and the proposed formulations have been performed on an analytic and a launch vehicle design test cases. The proposed formulations enable to decouple the multidisciplinary system and to decrease the number of calls to the disciplines while ensuring coupling satisfaction in instantiations. Further tests, in particular for launch vehicle design, should be performed to better qualify IDF-PCE formulation. Additional research is necessary to propose a decoupled UMDO formulation incorporating reliability analysis (to compute the probability of a constraint function to be above a certain threshold) instead of a robust formulation of the constraints (based on the statistical

Figure 2 :

 2 Figure 2: Multidisciplinary Design Optimization, decoupled approach

Figure 4 :

 4 Figure 4: Uncertainty Multidisciplinary Design Optimization decoupled process

Figure 5 :

 5 Figure 5: SORA procedure for UMDO 21

Figure 6 :

 6 Figure 6: Partially decoupled multidisciplinary system 47 and related notations

Figure 7 :

 7 Figure 7: Decoupled UMDO with the surrogate models of the coupling functional relations

Figure 9 :

 9 Figure 9: IDF-PCE with the numerical evaluation of the generalization errors by PCE

Figure 10 : 1 * c 0.47 21 -

 10121 Figure 10: Analytic test case of a multidisciplinary coupled system

4 Figure 12 :

 412 Figure 12: Distribution of the coupling errors. Left column, estimations of J 12 , right column, estimations of J 21 , from 150000 U samples.

•

  D: upper stage diameter (m) ∈ [1.8, 3.2] • At: upper stage throat cross-sectional area (cm 2 ) ∈ [10, 50] • P c: upper stage chamber pressure: (bar) ∈ [30, 75] • OF : upper stage mixture ratio ∈ [1.9, 3.5] • M ox : upper stage oxydizer mass (kg) ∈ [1100, 2000] • : upper stage nozzle expansion ratio ∈ [230, 300]

  to each discipline N d = 15 * 10 4 N d = 15 * 10 4

  to the disciplines Expected value of the upper stage mass (kg) Convergence curve of the MDF and IDF-PCE (quadrature) objective function MDF IDF-PCE (quadrature) 1 st feasible points

Figure 17 :Figure 18 :Figure 19 :

 171819 Figure 17: Convergence curves for MDF and IDF-PCE (quadrature) obtained for the feasible points

Figure 20 :Figure 21 :

 2021 Figure 20: Distribution of the upper stage mass (1000 samples) with IDF-PCE formulation

Figure 22 :

 22 Figure 22: Distribution of the Nax max coupling (1000 samples) with IDF-PCE formulation

  load factor (g)

Figure 23 :

 23 Figure 23: Distribution of the Nax max coupling (1000 samples) with MDA

Figure 24 :Figure 25 :Figure 26 :

 242526 Figure 24: Distribution of the P dynmax coupling (1000 samples) with IDF-PCE formulation

  with g i the i th component of the inequality constraint vector, F i = {U ∈ Ω|g i (z, Y, U) ≥ 0} and f U (.) the probability density function of U.

	Design,variables:,z Optimizer Uncertainty,simulation			Optimizer
					Uncertain,variables:,U		Design(variables:(z	Coupling(variables:(Θ[Y]
									Uncertain(variables:(U
	Discipline 1 c1.	Discipline Coupling,variables:,Y		Discipline 1	Discipline 2	Discipline N
	Coupling,variables:,Y	2	c2.	Discipline	Ξ[F(z,Y(z,U),U)] κ[g(z,Y(z,U),U)]	c1.	c2.	cN.
					N
					cN.
	Multidisciplinary,Design,Analysis
	Figure 3: Uncertainty Multidisciplinary Design Op-
	timization coupled process

  with y c the coupling variable vector satisfying interdisciplinary couplings (see Eq.33). • Step 3: the convergence is checked. If the inequality constraints (Eq.26) are verified and the objective function becomes stable 21 , the solution is found. • Step 4: if the convergence is not reached, or the inequality constraints are violated, a new deterministic MDO problem is formulated for

Table 2 :

 2 Vega upper stage design results for MDF and IDF-PCE (quadrature) formulations Mox = 1473.0kg Mox = 1339.7kg Coupling |Nax max -yN axmax | 2 ≤ 0.01g 2 JN axmax = 0.0147g 2 constraints |P dynmax -yP dynmax | 2 ≤ 5P a 2 JP dynmax = 3.32P a 2

	Results		MDF-MDA	IDF-PCE (quadrature)
	E[M 4 th ]		2983.0kg	2896.2kg
	Design		OF = 2.71	OF = 3.49
	variables		Pc = 70bar	Pc = 50.3bar
			= 300	= 288.9
			D = 1.82m	D = 2.16m
			At = 28cm 2	At = 33cm 2
	Constraint K value	-0.0234	-0.0358
	Dimension	design
	space

of

American Institute of Aeronautics and Astronautics

of 32 American Institute of Aeronautics and Astronautics

of 32 American Institute of Aeronautics and Astronautics

of 32 American Institute of Aeronautics and Astronautics

Acknowledgments : The work presented in this paper is part of a CNES/ONERA PhD thesis.

Propulsion

The propulsion discipline is based on the NASA computer program CEA (Chemical Equilibrium with Applications) [START_REF] Gordon | Computer program for calculation of complex chemical equilibrium compositions and applications[END_REF] which computes chemical equilibrium compositions and properties of complex mixtures but also the flow dynamics of the expansion of gases through the engine nozzle. The chosen propellants are the MMH/NTO (MonoMethyl Hydrazine / Nitrogen TetrOxide) because they are storable liquid propellants moments of the constraint function) as adopted in this paper.