Daniela Domenico

Keywords: Robust design optimization, uncertainty quantification, organic Rankine cycles, turbomachinery

I want to thank Mr. Benoit Obert, for having followed

Résumé

L'optimisation robuste (RDO) est un outil important pour la conception de produits industriels sous incertitude. Elle combine des algorithmes d'optimisation et des techniques de quantification de l'incertitude (UQ). La quantification d'incertitudes est généralement trop coûteuse pour des modèles numériques complexes de systèmes en ingénierie. Dans le but de développer des stratégies de RDO efficaces conçues pour des applications industrielles, le couplage de techniques UQ parcimonieuses avec un algorithme génétique multi-objectif basé sur des modèles substituts (SMOGA) a été étudié. A cet égard, une technique RDO prometteuse a été utilisée, basée sur le couplage de deux modèles substituts imbriqués: le premier est utilisé pour l'UQ, tandis que la surface de réponse du second est utilisée pour accélérer l'optimisation; un critère d'enrichissement est utilisé pour actualiser le modèle substitut pendant la convergence de l'optimiseur. Plusieurs méthodes d'UQ utilisant des informations sur les gradients de la solution par rapport aux variables incertaines ont été mises en oeuvre et comparées en termes de précision et coût de calcul. Nous avons ensuite sélectionné une méthode UQ dite «basse fidélité», c'est-à-dire peu coûteuse mais pas très précise, et une méthode «haute fidélité» afin de construire un modèle substitut multi-fidélité pour l'optimisation robuste. Ce modèle permet d'avoir une précision proche du modèle haute fidélité pour un coût de calcul bien moindre. Les méthodes étudiées ont été appliquées à la RDO de cycles thermodynamiques de

Introduction

The Organic Rankine Cycle (ORC) is a viable technology for the exploitation of renewable energies like concentrated solar power, geothermal power, biomass or waste heat recovery. In these applications, it usually outperforms classic steam cycles for its simplicity, the lower operational costs and the higher thermodynamic efficiency [START_REF] Macchi | Organic Rankine Cycle (ORC) Power Systems[END_REF].

ORCs are Rankine Cycles employing as working fluid complex organic compounds (hydrocarbons, silicon oils or refrigerants), instead of steam: they are closed cycles involving at least a pump, that compresses the working fluid, a group of hot heat exchangers, usually composed by one or more pre-heaters, an evaporator and sometimes also a superheater, where the working fluid is heated by an external heat source to become vapour; afterwards, there is a turbine, that converts the thermodynamic power of the fluid in mechanical one, and a condenser, where the residual heat is released at the environment allowing the fluid to come back to the liquid status. The mechanical power at the turbine shaft is converted in electricity by a generator.

During the last ten years, Organic Rankine Cycles (ORCs) have become a competitive technical solution for the exploitation of low-medium temperature heat sources of limited capacity [START_REF] Macchi | Organic Rankine Cycle (ORC) Power Systems[END_REF], bringing about an extraordinary growth of their market, in particular for geothermal, biomass and waste heat recovery (WHR) applications. Data plotted in Fig. 1.1 and Fig. 1.2 show an explosion in the number of the ORCs installed since the early 2000s simultaneously combined with a rapid increase in the average size of these plants.

From Fig. 1.1 it arises that ORCs main field of application is geothermal power, even if this technology is characterized by a low number of installations: therefore, geothermal plants generally Figure 1.2: Installation of ORC per year [2] have a considerable size. In such applications, ORCs are mainly employed for the conversion of liquid-dominated geothermal reservoirs at temperatures of 120-150 degC [START_REF] Colonna | Organic Rankine Cycle Power Systems: From the Concept to Current Technology, Applications, and an Outlook to the Future[END_REF]. In case of two-phase geothermal fluid, brine and steam are separated: the first is employed to warm the liquid working fluid up in the pre-heater(s), while the second one is used to evaporate the organic working fluid.

In any case, geothermal applications use generally saturated cycle configuration employing an alkane as a working fluid (the most common are pentane and butane). To enhance the heat recovery from geothermal brine, more pressure levels can be employed [START_REF] Dipippo | Chapter 8 -binary cycle power plants[END_REF] and sometimes the system can include a regenerator to improve the cycle efficiency. An example of application employing a non regenerative one pressure-level cycle is depicted in Fig. 1.3.

Figure 1.3: Example of application employing a non regenerative one pressure-level cycle [START_REF] Council | Portugal: Feature on pico alto geothermal orc power plant[END_REF] Biomass represents another quite widespread ORC application: starting from the early 2000, several high-temperature ORC power plants with the size of about 1 M W e have been installed in Europe to use various types of solid biomass [START_REF] Colonna | Organic Rankine Cycle Power Systems: From the Concept to Current Technology, Applications, and an Outlook to the Future[END_REF]. Quite commonly, these plants are cogenerative, providing both electricity and heat, which is typically employed below 100 degC for district heating or for process purposes (e.g., wood drying). [START_REF] Guercio | 15 -biomass-fired organic rankine cycle combined heat and power systems[END_REF]. The majority of these systems adopts siloxanes like hexamethyldisiloxane or octamethyltrisiloxane as the working fluid in a superheated regenerative cycle, where the working fluid is warmed up by thermal oil. Usually, the rated net electrical efficiency is in the 15-20% range, while the total energy efficiency can be higher, depending on the cogeneration needs [START_REF] Bini | Operational results of the first biomass chp plant in italy based on an organic rankine cycle turbogenerator and overview of a number of plants in operation in europe since 1998, 2nd, world biomass conference; biomass for energy industry and climate protection[END_REF].

Another application with an interesting potential for all unit sizes is waste heat recovery (WHR); in fact, many are the opportunities for heat recovery from the manufacturing and process industry [START_REF] Campana | Orc waste heat recovery in european energy intensive industries: Energy and ghg savings[END_REF].

According to data plotted in Fig. 1.1 and Fig. 1.2, this should be considered an emerging field for ORCs:

on that market several solutions are available with medium and large size recovery solutions from gas turbines, internal combustion engines or industrial processes (i.e., cement factory, steel mill and glass factory). Usually in this application siloxanes, alkanes, cycloalcanes and some HFO refrigerant can be adopted in thermodynamic cycles that are generally superheated. Nowadays, this application is the one presenting the highest flexibility in the cycle configuration; a 3.2 MW plant for waste heat recovery from a steel mill is shown in Fig. 1.4.

Finally, solar applications are negligible if compared with geothermal power, biomass and WHR:

this is probably due to the fact that ORCs could be coupled with concentrating solar power plants that for the moment are more expensive than photovoltaic panels and battery systems [1]. However, because of their high reliability, availability and performance, ORCs have been identified as the optimal conversion technology in this context [START_REF] Price | Modular trough power plant cycle and systems analysis[END_REF][START_REF] Casartelli | Power block off-design control strategies for indirect solar orc cycles[END_REF].

One of the biggest strength of ORC technology, that should be considered as one of the main reasons of their dramatic success, has probably to be found in their extreme flexibility, which allows them to be conveniently adapted to customer's necessities; in fact, they are so much customizable that sometimes some tailor-made solutions proposed by the ORC suppliers can be considered as real engineering challenges in terms of design, construction, commissioning and operation. By moving from a project to another, a number of elements are likely to change; some of which are listed below.

1. The purpose of the application: pure electric power production or cogeneration of electric power and heat.

2. The type of heat source (geothermal or biomass or WHR or solar).

3. The ORC operating conditions, which may require different plant control solutions.

4. The plant components (turbine, pump, heat exchangers, valves), that from an application to Figure 1.4: Exemple of a 3.2 MW plant for waste heat recovery from a steel mill [START_REF] Council | Portugal: Feature on pico alto geothermal orc power plant[END_REF] another usually can change dimensions and configurations.

5. The plant layout, depending on availability of room for all the components.

6. The working fluid.

7. The type of environment (indoor or outdoor) and the presence of acidity in case of geothermal plants. [START_REF] Marelli | Ingredients for an innovative uncertainty quantification platform in matlab[END_REF]. Safety restrictions, depending on the country and the applications.

9. Regulations about fluids and their interactions with pollution (i.e. Montreal protocol or the European F-gas regulation).

As a matter of fact, even if the totality of ORC manufacturers tries to standardize the solutions proposed to customers, by defining series of standard product range, they still continue to supply mostly tailor-made plants since it is not possible to fully standardize ORCs. This is confirmed by inspection of commercial material (brochures, presentations or white papers) provided by the biggest ORC suppliers (for instance, see [24,[START_REF] Enertime | Corporate brochure[END_REF]26,27,[START_REF] Exergy | Efficient and tailor-made organic rankine cycle geothermal binary power plants[END_REF]).

Paradoxically, the large number of degrees of freedom left free when designing an ORC plant is at the same time a serious limitation. In fact, in order to better understand the dynamics of the ORC market, it is mandatory to identify the shareholders that usually buy ORCs to operate them: the large majority of the ORC plants that have been built until today and that are still in operation is property of privately held entities, which usually consider ORCs as an investment of a part of their resources to be analysed in terms of risk, payback time and return on investment; usually the revenues come almost exclusively from the amount of energy generated, which needs to be evaluated as accurately as possible.

Thus, despite the high level of confidence and know-how reached nowadays, usually a project concerning an ORC is still dominated by a large series of unknown variables which are related to the design, the commissioning, the operation and the decommissioning of the plant and which could reduce shareholders' enthusiasm about ORCs. It is therefore mandatory to identify innovative techniques capable to deal with the uncertainty affecting ORCs: concerning the design of ORC components and whole ORC systems, robust design optimization (RDO) can be an answer to this issue.

In the following several RDO approaches are applied to ORC systems and to ORC turbo-expanders.

The manuscript is organised as here described:

• Chapter 2 introduces the verification and validation framework for computer models of real systems and the definitions of uncertainty; afterwards, it presents the uncertainty quantification (UQ) methodology and it provides a quick outline of some techniques that are commonly used to propagate the uncertainty through the model and to quantify it.

• Chapter 3 depicts the deterministic design process of a system or of a component and it provides few elements about optimization. Afterwards, the robust design optimization is introduced and finally two RDO techniques specifically investigated in this work are introduced, namely the "two nested Bayesian Kriging" (TNBK) and the multi-fidelity (MF) approaches.

• Chapter 4 presents the application of the TNBK technique to two ORC systems. First a proof of concept of this methodology is provided, by means of a direct comparison with a "brute force" RDO based on Monte Carlo sampling: to perform this validation, a crudely simplified model of an ORC for WHR is employed. Once that this approach has been validated, it is used for the RDO of a real ORC system for geothermal application.

• Chapter 5 focuses on the application of both RDO techniques explained in Chapter 3 to the robust design optimization of turbo-expanders. First several UQ techniques are assessed for a toy test problem roughly representative of an ORC nozzle, i.e. a quasi 1D nozzle. Finally, the MF technique is applied to the RDO of a 2D cascade of stator blades for an ORC turbine.

• Chapter 6 draws conclusions and provides some perspectives. The design of complex engineering systems increasingly relies on simulations based on computational models predicting their behaviour. At a glance, these models can be considered as a black-box containing a mathematical description of the relevant phenomena. Since it is usually not possible to solve analytically the set of equations employed in it, a numerical solver is embedded in the model.

Considering that the model aims at reflecting the real world, once properly calibrated and validated with respect to experimental data, it becomes an important tool for decision-making. For instance, it can be used to explore the design space, i.e. to explore how the system response (model output) reacts to changes in the system input parameters, like geometry, operating conditions and model parameters.

VERIFICATION AND VALIDATION

Another possibility is to put an optimization loop upon the model, to find the best design according to some pre-defined criterion (cost functions). Finally, models can be employed for making predictions about future system behavior (forecasting).

The present Chapter first provides introductory concepts about model verification, validation and uncertainty quantification. A selection of methods for quantification and propagation of uncertainty of interest for the present study is also presented.

Verification and Validation

Using a model as a reflection of some real phenomena is a not-straightforward process and it should be done consciously. In fact, as depicted in Fig. 2.1, it is mandatory to pass through two different types of models: a conceptual model and a computerized one. The first is constructed by analyzing and observing the physical phenomena of interest, since it aims to contain a more or less complete description of the system, based on governing equations and mathematical modeling from data [START_REF] Oberkampf | Verification and Validation in Computational Fluid Dynamics[END_REF]. In computational physics, the conceptual model contains the partial differential equations (PDE), or sometimes the ordinary differential equations (ODE), for conservation of mass, momentum, and energy; it also includes all auxiliary equations, such as 2.1. VERIFICATION AND VALIDATION turbulence models, constitutive models for materials, and electromagnetic cross-section models, the initial conditions and the boundary conditions of the PDEs [5]. The step linking the reality with the conceptual model is referred to as model qualification and it requires the good knowledge of the real system.

To simulate the behaviour of the real system, the conceptual model is usually replaced by a numerical approximation; this provides a computerized model that is finally applied to solve engineering problems. Finally, the computerized model output is assessed against measurements of the real system (possibly affected by observation errors). Based on the results, the conceptual model may be improved, and so on.

The circular process just described includes several sources of error or uncertainty that must be taken into account, namely observational errors (which is object of metrology), modeling errors and discretization ones [6]. To control them, verification and validation have a central role in the assessment of the accuracy of the conceptual and computerized models; as a consequence, a huge bibliography exists about this topic and about methods and techniques related at them. A thorough definition of verification has been provided by AIAA, as "the process of determining that a model implementation accurately represents the developer's conceptual description of the model and the solution to the model" [START_REF] Of | Guide for the Verification and Validation of Computational Fluid Dynamics Simulations[END_REF]. Thus, the emphasis in verification is pointed to the identification, 2.1. VERIFICATION AND VALIDATION quantification, and reduction of errors in the numerical methods and algorithm used to solve the equations employed in the model, in order to provide evidence that the conceptual model is solved correctly by the discrete mathematics embodied in the computational code; to do this, verification is mainly based on benchmark and comparison with highly accurate solutions of test case that are available [START_REF] Oberkampf | Verification and Validation in Computational Fluid Dynamics[END_REF]. Fig. 2.2 shows that the verification test relies on demonstration: a comparison between the computational solution and the correct answer is carried out for the evaluation of error bounds, such as global or local error norms [5].

On the other side, validation has been defined as "the process of determining the degree to which a model is an accurate representation of the real world from the perspective of the intended uses of the model" [START_REF] Of | Guide for the Verification and Validation of Computational Fluid Dynamics Simulations[END_REF]. Here, the main focus is on the process of understanding whether the right equations have been considered in the model for describing the real world. The fundamental strategy of validation involves identification and quantification of the error and uncertainty in the conceptual and computational models, quantification of the numerical error in the computational solution, estimation of the experimental uncertainty, and finally, comparison between the computational results and the experimental data, which are considered as the most faithful reflection of reality [5]. Fig. 2.3 depicts the validation test, which is based on the comparison of the computational solution with experimental data.

As typical complex engineering systems are usually interested by multidisciplinary coupled physical phenomena occurring together, the validation test presents a remarkable amount of serious complications: the quantities of interest (QoIs) are commonly not accessible for observation and, in any case, generally only just few experimental data, characterized by low quality and subject to uncertainty, are available. Consequently, because of the unfeasibility and impracticality of conducting true validation experiments on the most of complex systems, it is usually recommended to use a building-block approach, which decomposes the system of interest into progressively simpler tiers [START_REF] Of | Guide for the Verification and Validation of Computational Fluid Dynamics Simulations[END_REF][START_REF] Cosner | CFD validation requirements for technology transition[END_REF][START_REF] Marvin | Perspective on computational fluid dynamics validation[END_REF][START_REF] Sindir | Overview of the state-of-the-practice of computational fluid dynamics in advanced propulsion system design[END_REF].

Once that all blocks the validation and verification (VV) approach for scientifically-based predictions have been here presented, the whole VV framework is represented in Fig. 2.4. For further details about this topic, the reader is addressed to [6,[START_REF]Standard for Verification and Validation in Computational Fluid Dynamics and Heat Transfer[END_REF].

ERROR AND UNCERTAINTY

Error and Uncertainty

Two key elements in the VV process just presented are the concepts of error and uncertainty, that are often used interchangeably in the colloquial language even if they differ significantly.

In a probabilistic framework, the error is "a recognizable deficiency in any phase or activity of modelling and simulation that is not due to lack of knowledge " [START_REF] Of | Guide for the Verification and Validation of Computational Fluid Dynamics Simulations[END_REF] so that its deterministic nature as a deficiency can be identified by means of examination [START_REF] Walters | Uncertainty Analysis for Fluid Mechanics with Applications[END_REF]. In computerized solutions, there are five major sources of errors:

• discretization (spatial and temporal),

• insufficient convergence of an iterative procedure in the code solver,

• computer round-off,

• computer programming errors,

• improper use of the code. Figure 2.5: Classification of errors and major sources of errors in computerized solutions [START_REF] Montomoli | Uncertainty Quantification in Computational Fluid Dynamics and Aircraft Engines[END_REF] 2.3. UQ FRAMEWORK As indicated in Fig. 2.5, these can be classified in two sub-categories: acknowledged and unacknowledged errors. If the presence of the first one can be acceptable, as long as this kind of errors can be identified and quantified, for the second one instead, it is possible to identify them only by comparing the results with codes that set the benchmark; otherwise they cannot be found and removed from the code [START_REF] Montomoli | Uncertainty Quantification in Computational Fluid Dynamics and Aircraft Engines[END_REF].

Unlike error, uncertainty is "a potential deficiency in any phase or activity of the modeling process that is due to the lack of knowledge" [START_REF] Of | Guide for the Verification and Validation of Computational Fluid Dynamics Simulations[END_REF], revealing a clear stochastic nature [START_REF] Walters | Uncertainty Analysis for Fluid Mechanics with Applications[END_REF] as long as the lack of knowledge can occur in the physical models or in the input parameters, with a potential risk to affect the reliability of the simulation.

Uncertainty can be classified in two categories: the aleatory and the epistemic uncertainty. The first one is connected to the physical variability within the system or its environment and it cannot be reduced: the only way to treat it is just to characterize it by performing additional experiments, in order to get more data modelling the variables that can be used consequently in a probabilistic approach. On the other side, the epistemic uncertainty arises from the assumptions and from the simplifications made in deriving the physical formulation: thus, it can be considered reducible by performing more experiments and using the information to improve the physical models.

As depicted in Fig. 2.6, several sources of uncertainty exist for a computerized model. With regards to aleatory uncertainty, these are the randomness of the system parameters (geometric variables, material properties, manufacturing tolerances) and the randomness of the environment (initial and boundary conditions and environmental parameters). Epistemic uncertainty is largely related to uncertainties in model selection, which typically stems from a lack of knowledge about the underlying physics or from the impossibility to use a complete and accurate model to simulate them [START_REF] Fusi | Robust shape optimization of fixed and morphing rotorcraft airfoils[END_REF].

UQ Framework

Probabilistic engineering aims to take into account the uncertainties appearing in the modeling of physical systems and to study the impact of those uncertainties onto the system response [START_REF] Sudret | Uncertainty propagation and sensitivity analysis in mechanical models -Contributions to structural reliability and stochastic spectral methods[END_REF]. To this end, a well established and universally accepted framework for the quantification of the uncertainty has been developed: it is a circular multi-step process whose guidelines are presented in Fig. 2.7 and hereafter briefly described.

UQ FRAMEWORK

Step A consists in defining the model that should be used to reflect the physical system under consideration as already described in Section 2.1: this stage gathers all the ingredients used for a classical deterministic analysis of the physical system to be analyzed. Moreover this is the step where the QoIs are selected.

Step B corresponds to the quantification of the sources of uncertainty, by identifying the input parameters that cannot be considered as well-known because they are affected by uncertainty; therefore, they are modelled in a probabilistic context, through the definition of their PDFs. Typically, stochastic inputs can be associated with the operating conditions [START_REF] Cinnella | Robust optimization of dense gas flows under uncertain operating conditions[END_REF][START_REF] Congedo | Shape optimization of an airfoil in a BZT flow with multiple-source uncertainties[END_REF][START_REF] Cinnella | Robust optimization using nested Kriging surrogates: application to supersonic ORC nozzle guide vanes[END_REF], the geometry [START_REF] Chen | A new level-set based approach to shape and topology optimization under geometric uncertainty[END_REF][START_REF] Lazarov | Topology optimization with geometric uncertainties by perturbation techniques[END_REF][START_REF] Zhang | Robust shape and topology optimization considering geometric uncertainties with stochastic level set perturbation: Robust Topology Optimization Considering Geometric Uncertainties[END_REF][START_REF] Montomoli | Geometrical uncertainty in turbomachinery: Tip gap and fillet radius[END_REF] as well as the empirical parameters involved in the physical models, e.g. in turbulence models [START_REF] Harinck | Optimization of a 3d radial turbine by means of an improved genetic algorithm[END_REF] or the equation of state (EoS) [START_REF] Cinnella | Quantification of Thermodynamic Uncertainties in Real Gas Flows[END_REF]. The final product of this step is a vector of random variables having a well established probability distribution function (PDF).

Step C is the most computationally demanding one of the process, as it consists in propagating through the model all the uncertainties in the inputs, characterizing the random response appropriately with respect to the assessment criteria defined in Step A: usually, the objective of this phase is the computation of the probability distribution function of the QoI or of its statistical moments. Monte Carlo (MC) sampling is probably the most intuitive method to carry out this task, but it requires also a lot of computational resources, because it gives an entire discrete probability distribution (histogram) of the QoI as an output. In order to reduce the computational effort, some other methods are commonly used as they are faster than MC; some of them are described in the next section.

Once that the uncertainty has been propagated through the model, it is possible to perform a sensitivity analysis to get some information about the respective impact of the random input variables on the QoIs; in fact, the more complex the system is, the larger the number of input parameters becomes, and the designer will be probably interested in understanding the effects of each of them on system response. Thus, once that the PDFs (or just the statistics) of the QoIs have been computed, they can be used to characterize the output, in order to improve the knowledge about the considered problem.

Generally, two main approaches exist for sensitivity analysis, i.e. the local and the global. The first one is probably the most intuitive method as it analyses the impact on the QoI of small perturbation of the uncertain variables around their nominal values: this is done calculating the partial derivatives of the QoI with respect to the uncertain variables. Usually these quantities are derived by means of the 2.4. OVERVIEW ABOUT UNCERTAINTY PROPAGATION One-At-Time (OAT) method (the interested reader can find thorough review of this basic technique in [START_REF] Saltelli | How to avoid a perfunctory sensitivity analysis[END_REF]). However, the local approach presents several limitations, such as the assumptions of linearity and normality and it suffers in case of dependency among the uncertain variables [START_REF] Iooss | A review on global sensitivity analysis methods[END_REF]. On the other hand, the global sensitivity approach is based on the analysis of variance (ANOVA) decomposition [START_REF] Liu | Estimating mean dimensionality of analysis of variance decompositions[END_REF], and then it estimates Sobol indexes [START_REF] Ballester-Ripoll | Sobol tensor trains for global sensitivity analysis[END_REF] associated to the uncertain inputs. In the following, this last approach is the only one used to perform sensitivity analysis. The interested reader is referred to Appendix B for more details about it.

Overview about uncertainty propagation

Various techniques for uncertainty propagation are briefly reviewed hereafter; first the focus is put on Monte Carlo and some improved sampling methodologies and then the polynomial chaos expansion techniques and some surrogate-model based approaches are presented.

Monte Carlo method and improved sampling techniques

Developed at the begin of 40s of the last century within the Manhattan project [START_REF] Metropolis | The Monte Carlo Method[END_REF][START_REF] Eckhardt | Stan Ulam, John Von Neumann, and the Monte Carlo Method[END_REF], the MC method consists of random sampling of the input variables according to their own distributions and propagation of the resulting values through the model, which is run in deterministic mode for each sample. The output of the sampling process is an entire distribution of the QoI that can be analyzed with tools from statistics [START_REF] Metropolis | The Monte Carlo Method[END_REF][START_REF] Fishman | Monte Carlo -Concepts, Algorithms, and Applications[END_REF].

As described in Fig. 2.8, Monte Carlo approach consists in drawing samples from the input PDFs, and to compute for each of them the corresponding model outputs for the QoIs. The results may be used to construct an histogram of the output QoIs or to compute means and variances. The MC method has several advantages: it is well-posed, it has good convergence properties and, since it is non-intrusive, it does not need to have access to the equations employed in the model. Unfortunately, a very large number of samples is required to converge the statistics, which makes the Monte Carlo approach inapplicable to costly models: thousands or even millions of samples (each of them involving a model run) are required in order to get a good estimation of the QoI with enough accuracy. Precisely, the MC method shows convergence to the exact stochastic solution of 1 √ N , where N is the number of samples [START_REF] Marshall | Monte Carlo methods for the solution of nonlinear partial differential equations[END_REF]. For this reason, improved sampling techniques are usually adopted. LHS [START_REF] Mckay | Comparison of three methods for selecting values of input variables in the analysis of output from a computer code[END_REF] is one of the most widely used alternative to the standard random sampling of MC method. In fact, assuming to have a d-dimensional vector ξ of uncertain input parameters, it has been just shown that MC method randomly samples N vectors ξ i as represented in Eq.2.1.

ξ i = (ξ i 1 , ξ i 2 , ..., ξ i d ) f or i = 1, 2, ..., N (2.1)
In LHS, the range of probable values for each uncertain input parameter of ξ is divided into N non-overlapping segments of equal probability. Subsequently, for each uncertain input parameter of ξ, one value is randomly selected from each interval, by following the probability distribution on the interval. In this way it is possible to obtain d sequences of N numbers and to built the ξ i vectors by pairing randomly the entries in these sequences. Obviously in order to work properly, LHS must have a memory, while this is not true for MC method. To better understand this technique an example is presented in Fig. 2.9: for a set of size N = 5 and two input random variables ξ = (ξ 1 , ξ 2 ) , where ξ 1 A more detailed description of LHS technique is available in [START_REF] Helton | Latin hypercube sampling and the propagation of uncertainty in analyses of complex systems[END_REF], where it is shown also that, since the occurrence of low probability is reduced, the convergence is faster than standard MC, providing however an optimal convergence of the parameter space.

As indicated by the name, LBS is based on lattices, which can be seen as discrete repeating arrangements of points in a grid pattern. Thus, for each variable the entire range of probable values is again divided, but in this case the discretization is made using regularly spaced points and the solutions are evaluated at those points The characteristics of the lattice are strictly related to the distribution of the input variables [START_REF] Goos | Optimal Design of Experiments: A Case Study Approach[END_REF]. An example of a lattice based sampling for a two-variable problem is shown in Fig. 2.10. Usually LBS is used for for simulations dealing with condensed matter Several more improved sampling techniques exist (i.e. Low Discrepancies Sequences [START_REF] Braaten | An improved low-discrepancy sequence for multidimensional quasi-Monte Carlo integration[END_REF][START_REF] Tuffin | On the use of low discrepancy sequences in Monte Carlo methods[END_REF]), but a detailed description of them is out of the scope of the present work; the interested reader can refer to [START_REF] Goos | Optimal Design of Experiments: A Case Study Approach[END_REF][START_REF] Myers | Response Surface Methodology: Process and Product Optimization Using Designed Experiments, ser. Wiley Series in Probability and Statistics[END_REF]. However, even if they are faster than the standard MC, they remain too much expensive for complex models and large parameter spaces. For this reason, in order to reduce computational resources, they can be used with surrogates models, as it is discussed in the next section.

Polynomial Chaos Expansion (PCE)

A large class of methods for uncertainty propagation is based on the polynomial chaos framework [START_REF] Ghanem | Stochastic Finite Elements: A Spectral Approach[END_REF]: both intrusive and non-intrusive variants exist [START_REF] Loeven | Efficient uncertainty quantification in computational fluid dynamics[END_REF]. The implementation of intrusive uncertainty propagation methods requires access to the source code of the model; this formulation has two main drawbacks: first, the availability of the source code of the solver and, second, the need to modify the code as soon as the nature or number of uncertain parameters changes, the polynomial chaos order is modified, or different modeling choices are introduced. On the other hand, non-intrusive methods use the deterministic solver as a black-box for uncertainty propagation: this allows the application 2.4. OVERVIEW ABOUT UNCERTAINTY PROPAGATION of uncertainty propagation methods in combination with commercial or in-house solver. However, the number of samples required by PCE is greatly reduced in comparison with a full Monte Carlo simulation to get the same accuracy. As both the intrusive and the non-intrusive approaches are based on polynomial chaos expansion, the PCE framework is here quickly presented.

A polynomial chaos is a polynomial of random variables instead of ordinary variables [START_REF] Loeven | Efficient uncertainty quantification in computational fluid dynamics[END_REF]; it is based on the homogeneous chaos theory of Wiener [START_REF] Wiener | The Homogeneous Chaos[END_REF], who constructed a chaos expansion using Hermite polynomials. Ghanem and Spanos [START_REF] Ghanem | Stochastic Finite Elements: A Spectral Approach[END_REF] provided the basis for the current spectral stochastic finite element methods, like the Generalized Polynomial Chaos method (gPC) [START_REF] Xiu | The Wiener-Askey Polynomial Chaos for Stochastic Differential Equations[END_REF] and the Gram-Schmidt Polynomial Chaos method [START_REF] Witteveen | Modeling Arbitrary Uncertainties Using Gram-Schmidt Polynomial Chaos[END_REF]. Since a Galerkin projection is used to obtain the polynomial chaos coefficients, these methods are referred to as Galerkin Polynomial Chaos methods: they are probably the most intuitive PCE approach, but they have also a relatively cumbersome implementation, primarily due to the fact that it requires the availability of the solver source code; furthermore the equations for the expansion coefficients are almost always coupled [START_REF] Xiu | Fast numerical methods for robust optimal design[END_REF]. In case of high complexity in the original problem, the explicit derivation of the gPC equations can be impossible.

Thus, non-intrusive implementations of PCE are generally employed: usually they rely on two different approaches, namely the pseudo-spectral form [START_REF]Efficient collocational approach for parametric uncertainty analysis[END_REF][START_REF] Ko | Sensitivity of two-dimensional spatially developing mixing layers with respect to uncertain inflow conditions[END_REF] and the stochastic collocation method [START_REF] Tatang | An efficient method for parametric uncertainty analysis of numerical geophysical models[END_REF][START_REF] Babuška | A Stochastic Collocation Method for Elliptic Partial Differential Equations with Random Input Data[END_REF][START_REF] Loeven | Efficient uncertainty quantification in computational fluid dynamics[END_REF][START_REF] Cinnella | Robust optimization of dense gas flows under uncertain operating conditions[END_REF][START_REF] Cinnella | Sensitivity analysis of dense gas flow simulations to thermodynamic uncertainties[END_REF]. However, they both face the drawback that computational cost for the estimation of the statistics increases exponentially with the number of uncertain variables: this problem, known in UQ as curse of dimensionality [START_REF] De Baar | EFFICIENT UNCERTAINTY QUANTIFICATION WITH GRADIENT-ENHANCED KRIGING: APPLICATIONS IN FSI[END_REF][START_REF] Congedo | Efficient robust optimization techniques for uncertain dense gas flows[END_REF]), limits the applicability of the method to spaces not larger than 4 or 5 dimensions. A possible solution to alleviate curse of dimensionality consists in using sparse grid methods [START_REF] Xiu | High-order collocation methods for differential equations with random inputs[END_REF][START_REF] Ganapathysubramanian | Sparse grid collocation schemes for stochastic natural convection problems[END_REF][START_REF] Garcke | Sparse Grids and Applications -Stuttgart[END_REF], simplex stochastic collocation method [START_REF] Witteveen | Simplex Stochastic Collocation with Random Sampling and Extrapolation for Nonhypercube Probability Spaces[END_REF][START_REF] Edeling | Simplex-stochastic collocation method with improved scalability[END_REF] or Smolyak-based sparse pseudospectral approximations [START_REF] Eldred | Towards Goal-Oriented Stochastic Design Employing Adaptive Collocation Methods[END_REF][START_REF] Conrad | Adaptive smolyak pseudospectral approximations[END_REF][START_REF] Winokur | Sparse pseudo spectral projection methods with directional adaptation for uncertainty quantification[END_REF].

Another way to reduce curse of dimensionality consist in using adaptive techniques: a classic approach can be found in [START_REF] Gerstner | Dimension-adaptive tensor-product quadrature[END_REF][START_REF] Nobile | An anisotropic sparse grid stochastic collocation method for partial differential equations with random input data[END_REF][START_REF] Ma | An adaptive hierarchical sparse grid collocation algorithm for the solution of stochastic differential equations[END_REF]. More advanced implementations exists; the interested reader can refer to [START_REF]An adaptive high-dimensional stochastic model representation technique for the solution of stochastic partial differential equations[END_REF][START_REF] Griebel | Dimension-wise integration of high-dimensional functions with applications to finance[END_REF][START_REF] Blatman | Adaptive sparse polynomial chaos expansion based on least angle regression[END_REF][START_REF] Resmini | Sparse grids-based stochastic approximations with applications to aerodynamics sensitivity analysis[END_REF]

Surrogate Models

The basic idea in the surrogate model approach is substantially to build a fast mathematical approximation of the real model, which is too expensive in terms of computational resources. Thus, "given these approximations, many questions can be posed and answered, many graphs can be made, code to test the ideas so generated and, if necessary, update the approximations and iterate" [START_REF] Forrester | Engineering Design via Surrogate Modelling: A Practical Guide[END_REF]. The term "surrogate model" is often used interchangeably with the terms "metamodel" or "response surface model".

There are several types of surrogate model and each one needs to be built with its proper technique, but generally the way to employ them is very similar, as shown in Fig. 2.11. The process always starts with the construction of a sampling plan, also called design of experiment (DOE): the domain space of the input parameters is covered with several points that are selected according to the distribution of each input parameter. For instance if a parameter of the problem has a uniform distribution, it will 2.4. OVERVIEW ABOUT UNCERTAINTY PROPAGATION be subject to a uniform sampling. Typically, in this part of the process, random sampling, LHS, LBS or another improved sampling techniques are used.

Afterwards, since all these techniques rely on the data coming from the DOE, the original model is run for all the selected points, as in a MC simulation: the only difference lies on the fact that a MC simulation needs a huge amount of points in order to build an output distribution, while the DOE samples are required just to train and to test the surrogate model.

Once that the cheaper surrogate model has been built, it can be used to perform visualization activities, (i.e. the exploration if the design space, sensitivity analyses or uncertainty quantification)

or optimization, (i.e. the adjustment of system parameters to maximize/minimize one or more cost functions). Independently from the application, as the surrogate is only an approximation of the full model output estimates of some techniques for evaluating the approximation error are required.

Surface response enrichment by inclusion of additional observation(s) can be used to improve the surrogate quality.

Using a surrogate allows considerable gains in computational effort: for instance, in UQ applications there is usually a difference of 2-3 orders of magnitude in computational time between surrogate-based sampling and the brute MC method to get the same accuracy in the approximation of the PDF of the QoI [START_REF] Sudret | Meta-models for structural reliability and uncertainty quantification[END_REF].

Moreover, many surrogate models are non-intrusive, since it is not necessary to have access to the equations of the original model, which can be run in black-box mode for all the points in the DOE. Finally, the surrogate model approach can also benefit of parallelization; in fact, every run is independent from the other and consequently it can be solved in parallel on a High Performance Computing Infrastructure.

With regards to quantification of uncertainty, surrogates can be divided in two main category.

The first group is composed by deterministic methods: they first determine the sensitivities of the uncertain parameters, and then they combine this information with the covariance matrix of the input parameters, to obtain an estimation of the statistical moments of the QoI PDF [START_REF] Smith | Uncertainty Quantification: Theory, Implementation, and Applications[END_REF]. The second category contains stochastic-based techniques, which relies on probabilistic approach since they first require that the uncertainty of the input parameters is characterized by a PDF (e.g. normal, beta, lognormal,etc.); afterwards random samples are generated from the PDFs of the uncertain parameters 2.4. OVERVIEW ABOUT UNCERTAINTY PROPAGATION and these samples are propagated through the numerical model [START_REF] Maitre | Spectral Methods for Uncertainty Quantification: With Applications to Computational Fluid Dynamics[END_REF]. These methods allow to construct an approximation of the QoI PDF or they can be simply employed to calculate its statistical moments.

An example of deterministic method is the so-called first-order Method of Moments (MoM) [START_REF] Hazelton | Methods of moments estimation[END_REF],

which approximate statistical moments of the QoI by Taylor series expansions. Such method can be remarkably fast if derivatives of the cost function with respect to the uncertain variables are readily available by means of a discrete or continuous adjoint solver [START_REF] Jameson | Optimum Aerodynamic Design Using the Navier-Stokes Equations[END_REF][START_REF] Giles | An introduction to the adjoint approach to design[END_REF][START_REF] Cinnella | Hybrid Adjoint-based Robust Optimization Approach for Fluid-Dynamics Problems[END_REF][START_REF] Papoutsis-Kiachagias | Robust design in aerodynamics using third-order sensitivity analysis based on discrete adjoint. Application to quasi-1d flows[END_REF][START_REF] Papadimitriou | Third-order sensitivity analysis for robust aerodynamic design using continuous adjoint[END_REF][START_REF] Papadimitriou | Aerodynamic shape optimization for minimum robust drag and lift reliability constraint[END_REF]. Nevertheless, its accuracy is limited to Gaussian or weakly non-Gaussian processes with small uncertainties, since higher-order terms become increasingly important for strongly non-Gaussian input distributions. Some improvement can be achieved by using higher-order expansions, but these require information about higher-order sensitivity derivatives, which may represent a delicate and highly intrusive task. A more complete discussion can be found in [START_REF] Padulo | Comparative analysis of uncertainty propagation methods for robust engineering design[END_REF].

On the other hand, a stochastic UQ method is the radial basis functions approach: in engineering it is usual practice to use a polynomial approximation. Polynomials are generally restricted to low order approximations, especially in high dimensional problems [START_REF] Forrester | Engineering Design via Surrogate Modelling: A Practical Guide[END_REF]. Increased flexibility in the approximation can be achieved by adding further approximations to the polynomial surface with each one centered around one of the n sample points, as discussed in [START_REF] Broomhead | Multivariable Functional Interpolation and Adaptive Networks[END_REF] or in [START_REF] Buhmann | Radial Basis Functions: Theory and Implementations[END_REF] . This approach is known as radial basis functions. However, in case of complex or parametric basis functions, the additional task of estimating any other supplementary parameters should be considered; a typical example of a situation like this is the estimation of σ term in the Gaussian basis function, that as a result is usually taken to be the same for all basis functions [START_REF] Forrester | Engineering Design via Surrogate Modelling: A Practical Guide[END_REF]. Thus, for high-dimensional problems with a large number of sample points this approach turns out to be unwieldy.

Another stochastic method for the UQ is Kriging, which has emerged just in the last years as a powerful tool for building meta-models, even if it has been used (with different names) by Wiener [START_REF] Wiener | Extrapolation, interpolation, and smoothing of stationary time series[END_REF][START_REF] Wiener | The prediction theory of multivariate stochastic processes, II[END_REF] and Kolmogorov [START_REF] Kolmogorov | Stationary sequences in Hilbert space[END_REF][START_REF] Kolmogorov | On Strong Mixing Conditions for Stationary Gaussian Processes[END_REF]. Historically Kriging opened the field of geostatistics, formalized by Matheron [START_REF] Matheron | Traité de géostatistique appliquée[END_REF][START_REF]The theory of regionalized variables and its applications[END_REF] and it was named like this in order to honor the seminal work of the South African engineer D. Krige, who initiated a statistical method for evaluating the mineral resources and reserves [START_REF] Krige | A statistical approach to some basic mine valuation problems on the Witwatersrand[END_REF]. A detailed dissertation about the theory which is behind this technique and about its employment in geostatistics is provided also in Cressie [110] and in Stein [START_REF] Stein | Interpolation of Spatial Data, ser. Springer Series in Statistics[END_REF]. Sacks et al. [START_REF] Sacks | Design and Analysis of Computer Experiments[END_REF] introduced the key idea that Kriging may also be used in the analysis of computer experiments in which:

UNCERTAINTY QUANTIFICATION METHODOLOGY

• the data is not measured but results from evaluating a computer code, i.e. a simulator such as a CFD code or a finite element code;

• the points where data is collected are not physical coordinates in a 2D or 3D space, but parameters in an abstract space of arbitrary size.

Nowadays Kriging is considerably studied and applied as a prediction tool in economics and in machine learning, fields in which it is known as Gaussian Process (GP) technique.

In contrast to polynomial chaos expansions, Kriging provides a meta-model that does not depend on the probabilistic model for the input random vector [START_REF] Sudret | Meta-models for structural reliability and uncertainty quantification[END_REF]. In fact, The basic idea of Kriging is to model some function known only at a finite number of sampling points as the realization of a Gaussian Process (GP), which is defined as a collection of random variables, any finite number of which have a joint Gaussian distribution and it is completely specified by its mean and covariance functions [START_REF] Rasmussen | Gaussian Processes for Machine Learning[END_REF].

Uncertainty quantification methodology

One should consider a set of QoIs J :

J = J (x, ξ) (2.2)
with J ∈ R m depending on a vector of deterministic design parameters x ∈ R n des and on a vector of uncertain parameters ξ ∈ R nunc . One should note that some of the design parameters may also be uncertain. The objective of the UQ problem is the accurate estimation of the expectancy and the variance of the PDF of the QoI J , namely E[J ] and var[J ], given variations in the uncertain parameters; this methodology is based on a set of N deterministic samples of the solution J * . For sake of simplicity, in the following, only the case of a single QoI (m = 1) is considered, but the approach can be extended to multiple QoIs. Moreover, since the design parameters are considered to be uncertainty-free, the dependency of J on x is for the moment neglected. The required statistics of the PDF of the QoI J are calculated by means of the deterministic MoM and the stochastic Kriging, implemented in a Bayesian framework; furthermore, the formulation of the gradient enhanced Kriging (GEK) is provided.

UNCERTAINTY QUANTIFICATION METHODOLOGY

Methods of Moments (MoM)

Among deterministic UQ methodologies, an interesting approach is the MoM, which approximates the statistical moments of the fitness function by Taylor series expansions. This method may provide fast and sufficiently accurate estimates of the QoI statistics, as long as the fitness sensitivity derivatives with respect to uncertain parameters are provided by an efficient method and complete output statistics are not required [START_REF] Hazelton | Methods of moments estimation[END_REF]. Both first-and second-order MoM formulations have been considered in the literature (see for instance [START_REF] Cinnella | Hybrid Adjoint-based Robust Optimization Approach for Fluid-Dynamics Problems[END_REF]). Here, the focus is restricted on the first-order MoM. Following [START_REF] Walters | Uncertainty Analysis for Fluid Mechanics with Applications[END_REF],

the first-order approximation for the expected value (mean) of J is simply given by:

E [J(ξ)] = J(ξ ¯) (2.3)
which is nothing but the deterministic evaluation of function J at the mean value of the input ξ ¯.

The first-order variance is:

var [J(ξ)] = N ∑ i=1 N ∑ j=1 ∂J ∂ξ i ⃓ ⃓ ⃓ ⃓ ξ ¯∂J ∂ξ j ⃓ ⃓ ⃓ ⃓ ξ ¯cov(ξ i , ξ j ) (2.4)
with cov(ξ i , ξ j ) the covariance matrix. If the variables ξ = {ξ 1 , ξ 2 , ....ξ N } are uncorrelated, the covariance is a diagonal matrix and Equation (2.4) takes the simplified form:

var [J(ξ)] = N ∑ i=1 N ∑ j=1 ∂J ∂ξ i ⃓ ⃓ ⃓ ⃓ ξ ¯∂J ∂ξ j ⃓ ⃓ ⃓ ⃓ ξ ¯σ2 i (2.5)
where σ 2 i stands for the variance of the i-th uncertain parameter.

Bayesian Kriging (BK)

In general, the unknown QoI J , depending on a set of uncertain variables ξ with dimension n unc = M , is modelled as a regression function of the form:

J(ξ) = m(ξ) + Z(ξ) (2.6)
where:

1. m(ξ) = ∑ i f i (ξ)β i is the mean of the process.
2. f i is a set of basis functions (e.g. polynomial) for the regression model of the mean of the process.

UNCERTAINTY QUANTIFICATION METHODOLOGY

3. β i are the regression coefficients to be calculated for the mean of the process.

4. Z(ξ) = GP (0, P) is a stationary zero mean Gaussian process with covariance function P, modeling the deviation between the regression function and the data, that needs to be evaluated.

This Kriging formulation can be revised in the Bayesian framework (see for instance [START_REF] Kennedy | Bayesian calibration of computer models[END_REF][START_REF] Wikle | A Bayesian tutorial for data assimilation[END_REF]116,117]), which is particularly suitable to deal with the uncertainty of the model parameters and to provide a compensation in case of availability of only few measurements [START_REF] Omre | Bayesian kriging-Merging observations and qualified guesses in kriging[END_REF][START_REF] Pilz | Why do we need and how should we implement Bayesian kriging methods[END_REF]. Therefore, considering N sampling points contained in the design of experiment (DoE) ξ * , which is a matrix with dimension N × M , and the vector of the observed data J * = J(ξ * ), with dimension N composed by the evaluation of the QoI in the N samples ξ * , this approach is founded on Bayes' theorem, that is one of the cornerstones of Bayesian inference and it is summarized by Eq. 2.7:

p(J|J * ) = p(J * |J)p(J) p(J * |ξ) (2.7)
where:

• p(J|J * ) is the posterior,

• p(J * |J) is the likelihood,

• p(J) is the prior,

• p(J * |ξ) is the evidence.

The observed data vector J * should be thought as a subset of the QoI J selected by means of the observation matrix H, whose formulation is given in Eq. 2.8.

H ij = { 1, if i = j 0, otherwise i = 1, ..., M and j = 1, ..., N (2.8)
The likelihood in Eq. 2.7 comes from the sampling distribution as it is the conditional probability to observe the N data J * given the unknowns J and, because of the assumption that it is Gaussian distributed, it is formulated as in Eq. 2.9:

p(J * |J) ∼ N (HJ * , R) (2.9)
where R is the covariance of the observation error, that in the present work is considered uniform and uncorrelated, so that R = σ 2 I, with I equal to the identity matrix and σ which is a predefined error of the observed variable values.

The prior in Eq. 2.7 expresses the prior knowledge over the vector of the unknown QoI J(ξ) and, as it is also assumed to follow a Gaussian distribution, it is considered p(J) ∼ N (0, P); moreover, the evidence in Eq. 2.7 is also Gaussian, as it is the marginalization of the likelihood over the prior [START_REF] Rasmussen | Gaussian Processes for Machine Learning[END_REF].

As the right hand side of Eq. 2.7 is Gaussian, the posterior distribution obtained by means of the Bayesian framework just presented is calculated as a Gaussian distribution whose formulation is expressed in Eq. 2.10 [START_REF] Wikle | A Bayesian tutorial for data assimilation[END_REF].

p(J|J * ) ∼ N (E [J|J * ] , var [J|J * ]) (2.10)
where:

1. the mean E [J|J * ] = PH T ( R + HPH T ) -1 J * 2. the variance var [J|J * ] = [ I -PH T ( R + HPH T ) -1 ]

P

Finally, once that the posterior is available, it can be used to predict the values J ′ = J(ξ ′ ) of the QoI on the prediction points ξ ′ : as expressed in Eq.2.11, the predictive distribution J(ξ ′ ) can be obtained averaging the output of all possible linear models with respect to the posterior distribution [START_REF] Rasmussen | Gaussian Processes for Machine Learning[END_REF].

p(J ′ |ξ ′ , J) = ∫ p(J ′ |J, ξ ′ )p(J|J * )dJ (2.11)
Generally, for the most of the models, the integral in Eq.2.11 is not easily tractable, but Kriging can be considered as an exception to this statement, because of the properties of Gaussian distribution.

In fact, the predictive distribution is again a Gaussian, with the predictive mean equal to the posterior mean and the predictive covariance equal to the sum of the posterior covariance and the covariance of the observation error R (Eq. 2.12). For a better understanding of the technique just described, the whole process of the Kriging method in Bayesian framework is depicted in Fig. 2.12. At the beginning there are no function evaluations available, and the prior is supposed to be an infinite collection of random pathlines; to fix ideas, in Fig. 2.12a just eight of these pathlines have been depicted. The evaluation of the function J(ξ) in some points allows to improve the knowledge about J and to compute the posterior by means of Eq. 2.10; in Fig. 2.12b only eight posterior realization are shown. Finally, once that the information about the posterior is available, it can be used to perform the prediction considering also the uncertainty about the model; in Fig. 2.12c this uncertainty is represented by a 6σ confidence interval.

p(J ′ |ξ ′ , J) ∼ N (E(J|J * ), var(J|J * ) + R) (2.12)
To complete the description of the BK method, some information about the covariance function P must be provided. This must be semidefinite positive and it is sometimes also referred to as Kernel.

The most widely-used kernel within the kernel machines field is probably the squared exponential (SE), whose formulation is provided in Eq.2.13.

cov SE (ξ i , ξ j , θ) = exp ( - (ξ i -ξ j ) 2 2θ 2 ) (2.13)
where:

• ξ is the generic coordinate of the prediction y.

• ξ i -ξ j is the correlation range.

• θ is a parameter defining the characteristic length-scale (positive).

θ is a hyperparameter, that needs to be estimated using the data available from the prior distribution. The shape of the SE kernel as a function of hyperparameter θ is represented in Fig. 2.13a, while The SE covariance function is infinitely differentiable, which means that the GP with this covariance function has mean square derivatives of all orders, and is thus very smooth [START_REF] Rasmussen | Gaussian Processes for Machine Learning[END_REF]. Stein [START_REF] Stein | Interpolation of Spatial Data, ser. Springer Series in Statistics[END_REF] argues that such strong smoothness assumptions are unrealistic for modeling many physical processes, and he recommends therefore the Matérn class of covariance functions (Eq.2.14).

cov M atern (ξ i , ξ j , l, ν) = 2 1-ν Γ(ν) ( √ 2ν(ξ i -ξ j ) l ) ν K ν ( √ 2ν(ξ i -ξ j ) l )
(2.14) • ν and l are 2 hyperparameters (both positive).

• Γ is the Gamma function.

• K ν is the modified Bessel function of the second kind. A multitude of other possible families of covariance functions exists. The interested reader can find more information about this topic in [START_REF] Rasmussen | Gaussian Processes for Machine Learning[END_REF]. To estimate hyperparameters, the Maximum Likelihood )} (2.15) where A = (R+HPH T ). Actually, problem at Eq. 2.15 represents the bottleneck for the technique, since it is performed over a multidimensional space whose cardinality is given by the number of uncertain variables and the cost varies as o(n s × N 3 ), where n s is the number of optimization steps and N is the number of observed samples. As an alternative to MLE, a full Bayesian approach for determining the hyperparameters is also possible [START_REF] Rasmussen | Gaussian Processes for Machine Learning[END_REF], even if it is usually less common because it requires more resources. More efficient methodologies for estimating the hyperparameters could be implemented [START_REF] Bufi | Robust optimization of ORC turbine expanders[END_REF][START_REF] De Baar | Speeding up Kriging through fast estimation of the hyperparameters in the frequency-domain[END_REF]. Once that the hyperparameter is estimated, the Kriging response surface is used to approximate the PDF of the QoI, allowing hence the estimation of the required statistical moments (see for instance [START_REF] Cinnella | Robust optimization using nested Kriging surrogates: application to supersonic ORC nozzle guide vanes[END_REF]).

Gradient Enhanced Kriging (GEK)

Bayesian kriging can be accelerated by adding gradient information. The formulation described hereafter is the form of so-called gradient-enhanced co-kriging, denoted GEK for simplicity, where the gradient data are added as covariables. This approach has been shown to be simpler and more robust for indirect gradient-enhanced kriging, where the gradient information is included in the surrogate by adding finite difference-based values at small distances from the sample locations, and then, a standard kriging is performed on the resulting augmented sample [START_REF] Chung | Using gradients to construct cokriging approximation models for high-dimensional design optimization problems[END_REF][START_REF] Laurenceau | Building efficient response surfaces of aerodynamic functions with kriging and cokriging[END_REF][START_REF] Dwight | Efficient uncertainty quantification using gradient-enhanced kriging[END_REF].

For the sake of clarity, first GEK is described in the one-dimensional case: J = J(ξ). Assuming that the derivatives of dJ/dξ are known at the M sampling points, the observation vector is redefined by concatenating the gradient information to the function values:

J * = ( J 1 * , J 2 * , . . . , J M * , dJ dξ * ⃓ ⃓ ⃓ ⃓ 1 * , dJ dξ * ⃓ ⃓ ⃓ ⃓ 2 * , . . . , dJ dξ * ⃓ ⃓ ⃓ ⃓ M * ) T (2.16)
The covariance matrix is subsequently modified to account for the observational errors on the derivatives, also assumed to be uniform and uncorrelated. This results in a diagonal matrix of dimension

P 11 = ∂ 2 cov(h ij , θ) ∂h 2 ij .
In the above, P 00 is the covariance of the function values, P 11 is the covariance of the derivatives, and P 01 and P 10 are the cross-covariances.

The GEK formulation is extended to M -dimensional problems in a straightforward manner: the observation vector is expanded by adding N M -dimensional gradients, which results in the N ×(M +1)

compiled observations [START_REF] De Baar | Stochastic Surrogates for Measurements and Computer Models of Fluids[END_REF]. The prior covariance P becomes:

P = ⎡ ⎢ ⎢ ⎢ ⎢ ⎣
P 00 P 01

P M 0 P M 1 • • • P M M ⎤ ⎥ ⎥ ⎥ ⎥ ⎦
(2. [START_REF] Guercio | 15 -biomass-fired organic rankine cycle combined heat and power systems[END_REF] where:

P kl ij = ∂ 2 cov(ξ k,j , ξ l,i ) ∂ξ k ∂ξ l (2.20)
Due to the use of additional gradient information, GEK achieves a given level of accuracy with a number of samples that is consistently lower than BK. However, the solution quality may be highly dependent on the accuracy of the computed gradients [START_REF] Dwight | Efficient uncertainty quantification using gradient-enhanced kriging[END_REF], which must be properly accounted for when building the surrogate.

Preliminary validations

In the following both BK and GEK have been tested: first their performance have been compared

in the prediction of some simple 1D and 2D test functions for which analytical solution is available.

Afterwards, they are both used for the UQ of a simple ORC; in this last case, they have been compared with a MC, which has been used as a reference.

Simple analytic test functions

The simple 1D test function J = ξ 2 is considered, with variable x varying with uniform distribution in a range [-10, 10]. A Kriging approximation with just 3 samples of the test function is represented in Fig. 2.16. The red line is the Kriging response surface built by applying Eq. 2.11 and it presents important deviation from the blue line, which is the real solution; the red area is the 6σ confidence interval. In order to increase the accuracy, a larger number of samples should be used. Fig. 2.17

shows how Kriging approximation improves with 4, 5 and 6 samples. • the second-order 2D test function

J = ξ 2 1 + ξ 2 2 ,
• the third-order 2D test function

J = ξ 3 1 + ξ 3 2 ,
• the fourth-order 2D test function

J = ξ 4 1 + ξ 4 2 , • the Rosenbrock function [129] J = (a -ξ 1 ) 2 + b(ξ 2 -ξ 2 1 ) 2 , with a = 1, b = 100.
For all test functions are considered dependent from two variables ξ 1 and ξ 2 both varying with uniform distribution in a range [-10, 10]. The results of these validations are presented respectively in Each case confirm that GEK is very interesting as it can outperform BK in terms of observations;

in fact, with only 20 samples it reaches an acceptable accuracy, while BK requires always at least 30 2.6. PRELIMINARY VALIDATIONS samples to provide similar results. In Fig. 2.23d some noise appears in the GEK reconstruction: this is a possible issue for such a technique as a consequence of the increased size of the covariance matrix, leading to a higher conditioning number w.r.t. Bayesian Kriging. A little nugget can be applied to mitigate the magnitude of this side effect [START_REF] Laurent | An overview of gradient-enhanced metamodels with applications[END_REF]. However in the present calculation no nugget has been used. Finally, even if the employed methodologies deal with sampling techniques, a randomness analysis has not been considered in the present work. Once that BK and GEK have been tested with the fitting of some simple 1D and 2D analytical functions, the two methods have been compared on a more realistic test case: the quantification of the uncertainty in an ORC. However, several important simplifications have been employed in the ORC model for this phase. Such a choice has been done mainly for two reasons. First, the objective of this work is mainly to explore the capabilities of BK and GEK, comparing their predictions with the true values calculated directly with the code simulating the system; therefore a costly ORC cannot be used in this step. Moreover, GEK requires the gradient of the QoI with respect to the uncertain parameters:

in order to avoid the problems about gradient calculation described in Appendix A, in this step it has been decided to compute this quantity analytically. Consequently, an extremely simplified ORC model has been mandatory in this phase.

Therefore, it has been considered a 1 MW ORC for waste heat recovery employing R245fa For the ORC subject of the present work some assumptions have been made:

• Absence of pressure losses at the condenser; thus, T 1 = T 7 = Temperature at the Condenser (in the followings also called T con ).

• Absence of pressure losses in the heat exchangers; thus, P 2 = P 3 = P 4 = P 5 = Pressure at the Evaporator (in the followings also called P eva ).

• The work of the pump and off all other self-consumptions are not considered.

• Mechanical losses at the shaft of both turbine and pump are neglected.

• Disk windage losses in the turbine are not considered.

• Both mechanical and electrical losses at the generator are neglected.

After having defined the following quantities:

• Turbine isentropic efficiency η t = T urbine W ork T urbine Isentropic W ork • Pump isentropic efficiency η p = P ump Isentropic W ork P ump W ork

• Superheating ∆T sh = T 5 -T 4
A UQ study has been carried out in order to understand the effect of 5 uncertain parameters on the overall cycle efficiency η ORC , that has been chosen as QoI. The 5 uncertain parameters are:

• Temperature at the condenser T con varying with uniform distribution between 288 K and 315 K.

• Pressure at the evaporator P eva varying with uniform distribution between 18.0 bar and 23 bar.

• Turbine isentropic efficiency η t varying with uniform distribution between 75.0 % and 90.0 %.

• Superheating ∆T sh varying with uniform distribution between 1 K and 16 K.

• Pump isentropic efficiency η p varying with uniform distribution between 70.0 % and 85.0 %.

PRELIMINARY VALIDATIONS

The definition of the QoI is provided in Eq. 2.21. By mean of thermodynamic relations [START_REF] Moran | Fundamentals of Engineering Thermodynamics[END_REF] is it possible to express the QoI as in Eq. 2.22.

η
η ORC = m ̇(h 5 (P 5 , T 5 ) -h 6 (P 6 , s 6 ) -(h 2 (P 2 , T 2 ) -h 1 (P 1 , T 1 ))) m ̇(h 5 (P 5 , T 5 ) -h 2 (P 2 , s 2 )) (2.22)
where h is enthalpy (specific to mass), s is entropy (specific to mass), P is pressure and T is temperature. Since gradient of the QoI with respect to the five uncertain parameters has been computed

analytically, an equivalent formulation of the QoI has been considered (Eq. 2.23).

η ORC = ( h 5 (P 5 , T 5 ) -h 6 isos (P 6 , s 5 ) ) η t - P 2 -P 1 ρ 1 η p h 5 (P 5 , T 5 ) - ( h 1 (P 1 , T 1 ) + P 2 -P 1 ρ 1 η p ) (2.23)
where h 6 isos is the isentropic enthalpy at the discharge of the turbine, computed at turbine outlet pressure and at turbine inlet entropy.

The analytical formulation of gradient is provided by Eq. 2.24, Eq. 2.25, Eq. 2.26 and Eq. 2.27, where the thermodynamic derivatives can be computed both analytically [134] or with the thermodynamic library.

∂η ORC ∂P eva = ( η t (( ∂h 5 ∂P 5 ⃓ ⃓ ⃓ ⃓ T ∂P 5 ∂P eva + ∂h 5 ∂T 5 ⃓ ⃓ ⃓ ⃓ T ∂T 5 ∂P eva ) - ∂h 6 isos ∂s 5 ⃓ ⃓ ⃓ ⃓ P ∂s 5 ∂P eva - 1 ρ 1 η p ) × ( h 5 -h 1 - P 2 -P 1 ρ 1 η p ) - ( (h 5 -h 6 isos )η t ( ∂h 5 ∂P 5 ⃓ ⃓ ⃓ ⃓ T ∂P 5 ∂P eva + ∂h 5 ∂T 5 ⃓ ⃓ ⃓ ⃓ T ∂T 5 ∂P eva - 1 ρ 1 η p )))( h 5 -h 1 - P 2 -P 1 ρ 1 η p ) -2
(2.24)

∂η ORC ∂T con = (( - ∂h 6 isos ∂P 6 ⃓ ⃓ ⃓ ⃓ s ∂P 6 ∂T con η t - ∂( 1 ρ 1 ) ∂T 1 ⃓ ⃓ ⃓ ⃓ P ∂T 1 ∂T con P 2 -P 1 η p - ∂P 1 ∂T 1 ⃓ ⃓ ⃓ ⃓ P ∂T 1 ∂T con 1 ρ 1 η p ) × ( h 5 -h 1 - P 2 -P 1 ρ 1 η p ) + ( (h 5 -h 6 isos )η t ( ∂h 1 ∂P 1 ⃓ ⃓ ⃓ ⃓ T ∂P 1 ∂T con + ∂h 1 ∂T 1 ⃓ ⃓ ⃓ ⃓ P ∂T 1 ∂T con + ∂( 1 ρ 1 ) ∂T 1 ⃓ ⃓ ⃓ ⃓ P ∂T 1 ∂T con P 2 -P 1 η p - ∂P 1 ∂T 1 ⃓ ⃓ ⃓ ⃓ P ∂T 1 ∂T con 1 ρ 1 η p ))( h 5 -h 1 - P 2 -P 1 ρ 1 η p ) -2
(2.25) 2.6. PRELIMINARY VALIDATIONS

∂η ORC ∂η p = (( P 2 -P 1 ρ 1 η 2 p ) - (( (h 5 -h 6 isos )η t - P 2 -P 1 ρ 1 η p )( P 2 -P 1 ρ 1 η 2 p ))) × ( h 5 -h 1 - P 2 -P 1 ρ 1 η p ) -2
(2.26)

∂η ORC ∂∆T sh = ( η t ( ∂h 5 ∂T 5 ⃓ ⃓ ⃓ ⃓ P ∂T 5 ∂δT sh - ∂h 6 isos ∂s 5 ⃓ ⃓ ⃓ ⃓ P ∂s 5 ∂T 5 ⃓ ⃓ ⃓ ⃓ P ∂T 5 ∂δT sh ) × ( h 5 -h 1 - P 2 -P 1 ρ 1 η p ) -(h 5 -h 6 isos )η t ( ∂h 5 ∂T 5 ⃓ ⃓ ⃓ ⃓ P ∂T 5 ∂δT sh )) × ( h 5 -h 1 - P 2 -P 1 ρ 1 η p ) -2 (2.27)
Once that the problem has been formulated as described so far, the response surfaces of BK and GEK are constructed by using several experimental designs with different numbers of samples of the 5 uncertain parameters (Latin Hypercube). Afterwards, the surfaces are evaluated on a full factorial grid of 5 5 points in order to get the predicted values, which, because of the low complexity of the present test case, can be compared with the exact results obtained through direct simulation of the model by mean of Eq. 2.21. This comparison is depicted in Fig. 2. [START_REF] Enertime | Corporate brochure[END_REF], where values predicted by surrogate models (yellow points) should lie on the quadrant bisector (blue dotted line) when they correspond to the exact solution.

For both surrogate models the convergence is verified by monitoring the trend of the average value of the coefficient of variation (defined as the ratio of the standard deviation σ to the mean value µ)

on the response surface with respect to the number of samples in the experimental design; this trend is shown in Fig. 2 The design of a component is the process by which engineers select the technical specification of a system or component allowing to meet a set of requirements in terms of performance, cost, durability.

The search of the best design parameters is to all intents and purposes an optimization [START_REF] Beyer | Robust optimization -A comprehensive survey[END_REF]. The goal of the present section is to present the robust design optimization strategies employed in the next chapters. For this purpose, a short general overview of optimization techniques is first provided in Section 3.1 and the deterministic design approach is briefly introduced in Section 3.2. Second, optimization under uncertainty is introduced in Section 3.3. Section 3.4 introduces the notion of multifidelity (MF) surrogate model. Finally, the robust optimization selected for the present study are presented in Section 3.5, including a MF strategy combining two UQ methods of different fidelities introduced in the preceding chapter.

INTRODUCTORY NOTIONS IN OPTIMIZATION

Introductory notions in optimization

A general optimization problem can be formulated as in Eqs. 3.1 -3.4.

minimize : J (x) (3.1)
subject to :

g(x) ≤ 0 (3.2) h(x) = 0 (3.3)
where :

x = ⎛ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎝ x 1 : x L 1 ≤ x 1 ≤ x U 1 x 2 : x L 2 ≤ x 2 ≤ x U 2 .
. .

x i : x L i ≤ x i ≤ x U i x N : x L N ≤ x N ≤ x U N ⎞ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎠ (3.4)
The target is to minimize the function J , which is usually called cost function or objective function [START_REF] Kochenderfer | Algorithms for Optimization[END_REF]. J can be a quantity like a cost, a mass, a stress, a deformation, an efficiency or even a combination of them [137]. If J is a scalar the problem is referred to as single-objective optimization, while if it is a vector containing multiple cost functions, the problem is classified as a multi-objective optimization.

J depends on N design variables contained in the N-dimensional vector x, each one defined in its own validity range, as indicated in Eq. 3.4. Typically, in engineering it is possible to find optimization problems counting from one or few couples of design variables to several tens: they can be seen as the degrees of freedom of the problem. For the optimization formalized in Eqs. 3.1 -3.4, the idea is to find the best value of the design variables in x in order to minimize the cost function J . In case of a maximization problem, usually this is converted in an equivalent one, aiming to minimize the inverse problem -J (x); otherwise, another approach consists in the minimization of the reciprocal problem 1/J (x). However, this second strategy results to be less stable of the first one, as it can diverge whether J (x) = 0 and it is therefore less common.

INTRODUCTORY NOTIONS IN OPTIMIZATION

Optimization problems may be subject to a set of constrains under the form of an inequality (Eq. 3.2) or an equality (Eq. 3.3). Because of computers floating point arithmetic [START_REF] Engineers | Ieee standard for floating-point arithmetic[END_REF], equality constraints are much harder to satisfy than inequality constraints; therefore, they are usually converted in inequality constraints, by introducing a tolerance value. In any case, the constrained optimization problem as defined in Eqs. 3.1 -3.4 is usually converted in an unconstrained equivalent one, where some penalties are applied to the cost function whether the selected point falls outside the feasible region.

In case of multi-objective optimization, the problem formulation provided in Eqs. 3.1 -3.4 is still valid: the only difference is that J is a vector of (possibly conflicting) objectives [137]. Therefore, finding a single optimal solution optimizing all objectives simultaneously is usually an impossible operation, and the optimization process aims instead at identifying a set of efficient solutions whose objectives cannot be improved in one dimension without being worsened in the others [START_REF] Legriel | Approximating the pareto front of multi-criteria optimization problems[END_REF]: these points are called Pareto solutions and together they form the Pareto front, that is the subset of the solutions, for which all objective functions are equal or better than for all other solutions. A typical Pareto front for a two-objective optimization (minimization) problem is depicted in red in fig. 3.1:

the solution C is dominated both by solutions A and B, because

J 1 (A) < J 1 (C), J 2 (A) < J 2 (C), J 1 (B) < J 1 (C) and J 2 (B) < J 2 (C).
On the other hand, solution A and solution B do not dominate each other as they both are Pareto solutions.

In such a situation, expert judgement is required to select a single "best" solution from the Pareto front; this is not a simple task, since the cost functions are not dimensionally homogeneous quantities and comparisons among them are therefore impossible.

Several techniques exist to solve the optimization problem defined in Eq. 3.1 -3.4. Usually they are divided in two broad categories: the group of derivative-free methods and the one of derivative-based methods. For a review, the interested reader can refer to [START_REF] Kochenderfer | Algorithms for Optimization[END_REF]. Methods employing derivatives are generally supposed to converge to local optima, while derivative-free methods are usually reputed to reach global optima, as they explore the design space [START_REF] Conn | Use of quadratic models with mesh-adaptive direct search for constrained black box optimization[END_REF][START_REF] Wild | Orbit: Optimization by radial basis function interpolation in trust-regions[END_REF]. Moreover, if the problem is modelled as a black-box function, it could be complicated to calculate derivatives and usually the only available choice is the expensive finite differences approach (for an overview about the methods to calculate gradients, the reader is addressed to Appendix A). On the other hand, derivative-free optimization methods always require a larger number of evaluations in comparison with derivative-based methods (usually, depending on the problem from 10 times to 100 times more [START_REF] Larson | Derivative-free optimization methods[END_REF][START_REF] Zingg | A comparative evaluation of genetic and gradientbased algorithms applied to aerodynamic optimization[END_REF]) and they tend to suffer from curse of dimensionality problem, as with the increasing number of parameters they present a more difficult convergence towards the optimal solutions [START_REF] Bellman | Adaptive Control Processes: A Guided Tour, ser[END_REF]. Therefore, the selection of an optimization method is always a compromise among these issues.

Design and optimization

As already mentioned in the introduction of the present chapter, designing a component means to perform a constrained optimization of some (one or more) objective functions, that are a mathematical translation of all desired goals and requirements; a schematic representation of this iterative process Once that the numerical model is available, it should be fed with the environmental variables, which are the boundary conditions of the problem, and with the design parameters, defining the design space and they are the degree of freedom of the problem: the aim of design is to find the best value for all design parameters in order to optimize the value of the cost functions. With these ingredients is possible to perform a simulation, whose results can be post-processed to evaluate the cost functions.

When this point of the process is reached, one should perturb the design point and re-run the whole process. This operation is generally carried out by the optimization algorithm.

Usually, this is a deterministic process: for instance, in a case of single-objective optimization,

DESIGN AND OPTIMIZATION

given a set of constant deterministic input parameters, the outcome of the optimization is a single result, that is a vector containing the deterministic optimum of the design parameters.

However, as already seen in Chapter 2, the whole design process depicted in Fig. 3.2 can be affected by several sources of uncertainty, acting on:

• the environmental variables of the problem,

• the design parameters,

• the approximate numerical model employed to simulate the real system.

As a consequence, the deterministic design process just presented is no longer applicable. To overcome this difficulty, a classical technique is to consider the average values for the boundary conditions and then to apply during the design some safety margins; alternatively, a more modern approach requires to switch from the deterministic viewpoint to the stochastic one and to perform the quantification of uncertainty.

Moreover, the standard methodology of the deterministic optimization presents other limits that should be taken into account, as systems optimized in the classical sense can usually be very sensitive to small changes [START_REF] Bufi | Robust optimization of ORC turbine expanders[END_REF]. Therefore, because of the effect of modeling errors, it can happen that the found optimal solution differs from the true one. Finally, even in the case of a perfect model devoid of errors, the realization of the predicted optimum design could be affected by all sources of uncertainties, both aleatory (i.e. environmental conditions and geometric tolerances) and epistemic (e.g. imperfect knowledge of the input distribution), leading to an optimized configuration where actual performance can be far from the predicted one.

In such a context, non-deterministic methods should be considered as a strategy. Without the aim to present an exhaustive review about these techniques for optimizing under uncertainties, one can classify them into two main broad approaches, namely, reliability-based design optimization (RBDO)

and robust design optimization (RDO) [START_REF] Schuëller | Computational methods in optimization considering uncertainties -an overview[END_REF]. The first one concentrates on finding an optimum with a given maximum probability of failure that should be met by the design and it is usually adopted in structural problems [START_REF] Dersjö | Efficient design of experiments for structural optimization using significance screening[END_REF], while the second one aims at identify a design minimizing the sensitivity of objective and constraints to to uncertainties. Several formulations exist for both strategies (see for instance [START_REF] Deng | Methods for robust and reliability-based design optimization of electromagnetic devices[END_REF][START_REF] Kang | Robust Design Optimization of Structures under Uncertainties[END_REF]). In the works presented in the following, only the RDO approach has been 3.3. ROBUST DESIGN OPTIMIZATION considered.

Robust Design Optimization

In recent years, robust design optimization (RDO) [START_REF] Beyer | Robust optimization -A comprehensive survey[END_REF] has received increasing interest in engineering applications, due to its ability to provide efficient designs with a stable behavior under uncertainties of a diverse nature, such as randomly fluctuating operating conditions, geometric tolerances, and model uncertainties. Taguchi's method [START_REF] Taguchi | System of experimental design: engineering methods to optimize quality and minimize costs[END_REF], relying on the simultaneous optimization of the average and variance of the stochastic cost functions, is by far the most popular RDO method, although approaches allowing accounting for rare events, such as the low-quantile [START_REF] Maliki | Quantile-based optimization under uncertainties using adaptive Kriging surrogate models[END_REF][START_REF] Razaaly | Quantile-based robust optimization of a supersonic nozzle for organic rankine cycle turbines[END_REF] or the "horsetail matching" [START_REF] Cook | Horsetail matching: a flexible approach to optimization under uncertainty[END_REF] methods, have been paid significant interest recently.

The main ingredient for RDO is an uncertainty quantification method, allowing characterizing the probability distribution functions or, at least, the lower order statistics of the cost functionals for each proposed design, in order to select those that guarantee the best possible average performance while avoiding critical deviations when nominal design conditions are not matched. According to the RDO method in use, a single objective deterministic design problem is generally converted into a multiobjective (Pareto front) one, with the aim to optimize the average performance while avoiding critical performance loss at off-design conditions. For this reason, RDO often combines an UQ solver with evolutionary algorithms (typically, multi-objective genetic algorithms (MOGA) [START_REF] Deb | Optimization for Engineering Design -Algorithms and Examples[END_REF][START_REF] Kochenderfer | Algorithms for Optimization[END_REF]), which are naturally suited for providing a full set of compromise solutions among the multiple objectives. On the other hand, evolutionary optimizers are generally very demanding in terms of cost function evaluations, which may require in the end a prohibitive computational effort for problems described by costly computer models, such as those encountered in computational fluid dynamics (CFD), despite the use of massive parallelization [START_REF] Congedo | Efficient robust optimization techniques for uncertain dense gas flows[END_REF]117,[START_REF] Cinnella | Robust optimization using nested Kriging surrogates: application to supersonic ORC nozzle guide vanes[END_REF]. In order to reduce the number of costly function calls, it is crucial to select parsimonious UQ methods and optimizers, the overall cost of RDO being typically the product of the cost of the two approaches [START_REF] Congedo | Shape optimization of an airfoil in a BZT flow with multiple-source uncertainties[END_REF]. Past examples of RDO on expensive CFD simulations include various forms of UQ solvers based on non-intrusive polynomial chaos expansion [START_REF] Cinnella | Robust optimization of dense gas flows under uncertain operating conditions[END_REF][START_REF] Congedo | Shape optimization of an airfoil in a BZT flow with multiple-source uncertainties[END_REF] or surrogate models such as simplex stochastic collocation [START_REF] Congedo | A simplex-based numerical framework for simple and efficient robust design optimization[END_REF] or kriging [117]. All of them require a number of CFD solves quickly increasing with the number of uncertain parameters, and their direct coupling with MOGA optimizers is not computationally affordable for industrial applications, especially if massively parallel computers are not available.

MULTI-FIDELITY SURROGATE MODELS

An interesting option for reducing the cost of UQ solves is to use gradient information. It has been already seen that a simple method for approximating statistical moments of the cost function by Taylor series expansions is the first order MoM. An alternative to MoM, better suited for high uncertainty ranges and generic pdf, consists of leveraging gradient information to construct a highquality surrogate from a reduced number of samples. Such an approach is used for instance in GEK surrogates.

Massive parallelization is of great help for speeding up the RDO process [START_REF] Congedo | Shape optimization of an airfoil in a BZT flow with multiple-source uncertainties[END_REF][START_REF] Hercus | Robust shape optimization of uncertain dense gas flows through a plane turbine cascade[END_REF][START_REF] Cinnella | Robust optimization using nested Kriging surrogates: application to supersonic ORC nozzle guide vanes[END_REF], but it is not promptly applicable for routine industrial use. A way of drastically reducing the required number of function calls consists of replacing the costly CFD or UQ solvers with surrogate models, such as radial basis functions [START_REF] Cinnella | Optimal airfoil shapes for viscous transonic flows of bethe-zel'dovich-thompson fluids[END_REF], artificial neural networks [158], and kriging [117,[START_REF] Cinnella | Robust optimization using nested Kriging surrogates: application to supersonic ORC nozzle guide vanes[END_REF], approximating variations of the cost functions through the design space. Such an approach is called a surrogate-based multi-objective genetic algorithm (SMOGA).

Finally, further reductions of computational time can be achieved by combining models with various levels of fidelity during the optimization.

Multi-fidelity surrogate models

In the last decade multi-fidelity (MF) surrogate models have been paid significant interest for optimization, UQ and statistical inference as they can allow to achieve the desired accuracy at lower cost, if compared with the traditional techniques [START_REF] Fernández-Godino | Issues in deciding whether to use multifidelity surrogates[END_REF]. To come to this end, MF methods combine information from multiple models, with different cost and accuracy (fidelity).

In fact, despite increased availability of computational power, using exclusively (costly) high-fidelity for UQ or optimization may exceed the available resources. Thus, MF methods aims to provide good results and to reduce the cost by leveraging the use of low-fidelity (LF) models, which are cheap but they can also be inaccurate (as depicted in Fig. 3.3) [START_REF] Peherstorfer | Survey of multifidelity methods in uncertainty propagation, inference, and optimization[END_REF]; the recourse to the HF model is tough to preserve the accuracy.

LF models can be for instance coarse-grid approximations (see for instance [START_REF] Choi | Multi-Fidelity Design Optimization Studies for Supersonic Jets Using Surrogate Management Frame Method[END_REF]160,[START_REF] Kennedy | Bayesian calibration of computer models[END_REF]), data-fit interpolation and regression models (for example in [161]), projection-based reduced models [START_REF] Rozza | Reduced Basis Approximation and a Posteriori Error Estimation for Affinely Parametrized Elliptic Coercive Partial Differential Equations[END_REF][START_REF] Antoulas | Interpolatory Model Reduction of Large-Scale Dynamical Systems[END_REF], machine-learning-based models [START_REF] Vapnik | Statistical Learning Theory[END_REF] or simplified models relying on approximations of the underlying physics [START_REF] Ng | Multifidelity approaches for optimization under uncertainty[END_REF]. Concerning point 3, several choices are possible. The simpler approach is to link the HF and the LF models by using an additive correlation, as proposed in [START_REF] Lewis | A multigrid approach to the optimization of systems governed by differential equations[END_REF]: given a LF model f LF (ξ) and a HF model f HF (ξ), this is presented in Eq. 3.5.

f HF (ξ) = f LF (ξ) + δ(ξ). (3.5)
δ(ξ) is a function defined in R that can be used to capture the difference between the two levels of fidelity. This approach is accurate enough when HF and LF models are expected to have the same scale and a good correlation, as for instance in case of a coarse-grid approximations [START_REF] Meliani | Multi-fidelity efficient global optimization: Methodology and application to airfoil shape design[END_REF].

MULTI-FIDELITY SURROGATE MODELS

An alternative is provided by the multiplicative approach, used for example in [START_REF] Forrester | Optimization using surrogate models and partially converged computational fluid dynamics simulations[END_REF][START_REF] Han | Improving variable-fidelity surrogate modeling via gradient-enhanced kriging and a generalized hybrid bridge function[END_REF][START_REF] Zhang | Multifidelity surrogate based on single linear regression[END_REF],

which is described in Eq. 3.6.

f HF (ξ) = ρ(ξ)f LF (ξ) (3.6)
Here ρ(ξ) is in most cases a constant scalar multiplier to the LF function instead of a real function depending on ξ.

A comprehensive formulation combining the approach presented in Eq. 3.5 and Eq. 3.6 has been proposed in [START_REF] Kennedy | Bayesian calibration of computer models[END_REF] and it is shown in Eq. 3.7.

f HF (ξ) = ρf LF (ξ) + δ(ξ) (3.7)
with f LF (ξ) and δ(ξ) independent. This approach is considerably more robust of the previous two and for this reason this method has been extensively adopted. In particular, an interesting application is the combination of GP with a Bayesian identification of ρ and δ(ξ) by using the maximum likelihood estimation (see for instance [START_REF] Forrester | Multi-fidelity optimization via surrogate modelling[END_REF][START_REF] March | Gradient-based multifidelity optimisation for aircraft design using bayesian model calibration[END_REF][START_REF] Gratiet | Cokriging-based sequential design strategies using fast crossvalidation techniques for multi-fidelity computer codes[END_REF][START_REF] Park | Remarks on multi-fidelity surrogates[END_REF]).

A more complex correlation can be found in [START_REF] Gano | Hybrid variable fidelity optimization by using a kriging-based scaling function[END_REF] and it has the expression reported in Eq. 3.8.

f HF (ξ) = ω(ξ)ρ(ξ)f LF (ξ) + (1 -ω(ξ))(f LF (ξ) + δ(ξ)) (3.8)
where ω(ξ) is a weighting function. For the most of the applications, this method is usually too complex and therefore the model in Eq. 3.7 is generally preferred; however the employment of model Eq. 3.8 can be found in [START_REF] Zheng | A hybrid variable-fidelity global approximation modelling method combining tuned radial basis function base and kriging correction[END_REF][START_REF] Fischer | Bayesian Low-Fidelity Correction Approach to Multi-Fidelity Aerospace Design[END_REF].

Concerning point 4, i.e. the strategy for generating the LF and HF DoEs, an extremely broad range of alternative options exist. Without the aim to provide a full review, the three simplest strategies are here quickly outlined. For a thorough presentation of more complex strategies the reader is addressed to [START_REF] Fernández-Godino | Issues in deciding whether to use multifidelity surrogates[END_REF]. The first, and probably most common one, is to generate independently the LF and HF sampling points. Alternatively, is it possible to generate HF data points as a subset of LF data points:

this method is referred to as Nested sampling design and it is depicted in Fig. 3.4, where LF points (in blue) are placed first and then some HF points (in orange bubbles) are selected among them. Usually this operation is carried out by means of the D-optimal design technique [START_REF] De Aguiar | D-optimal designs[END_REF]. Finally, a third possible technique is given by the nearest neighbor sampling: as presented in Fig.

Robust Design Optimization Strategies

Following Taguchi's RDO method, the objective is to look for a methodology that allows to optimize a set of QoIs, J = J (x, ξ), J ∈ R m depending on a vector of deterministic design parameters x ∈ R n des and on a vector of uncertain parameters ξ ∈ R nunc . Note that some of the design parameters may also be uncertain. The RDO problem is therefore formulated by using the expectancy and the variance of J as measures of robustness, which leads to the solution of the two objective deterministic optimization problem in Equation (3.9).

min x E[J ], min x var[J ] (3.9)
To solve the preceding optimization problem, one should take into account that:

• evaluating J can be very costly,

• derivatives of J with respect to x and ξ can be not easily accessible,

ROBUST DESIGN OPTIMIZATION STRATEGIES

• J can be a non-convex function.

Since the RDO problem defined in Eq. 3.9 is an intrinsically multi-objective problem and considering that the aim is the pursuit of global optima, a MOGA is selected for searching the design space.

More precisely, following some previous studies [START_REF] Cinnella | Robust optimization of dense gas flows under uncertain operating conditions[END_REF][START_REF] Cinnella | Robust optimization using nested Kriging surrogates: application to supersonic ORC nozzle guide vanes[END_REF], the non-dominated sorting genetic algorithm (NSGA-II) of Deb et al. [START_REF] Deb | A fast and elitist multiobjective genetic algorithm: NSGA-II[END_REF] is used: this provides an approximated Pareto front of optimal solutions corresponding to different trade-offs between average performance and robustness for the various QoIs at hand. For simplicity, in the following, only the case of a single QoI, m = 1, is considered, but the approach can be extended to multiple QoIs. The required statistics of the QoIs are calculated by means of some non-intrusive UQ methods, already presented in Chapter 2. Thus, the first ingredient of the RDO process is an efficient UQ approach, which provides accurate approximations of E[J ] and var[J ] based on a set of N deterministic samples of the solution.

Direct coupling of the MOGA with the UQ solver is overly costly for complex expensive simulations code (i.e. CFD), due to the high number of function evaluations. For instance, running the MOGA with a population of n pop individuals over n gen generations and using N samples for UQ lead to an overall number of the QoI evaluations of about N ×n pop ×n gen . The computational cost can be greatly alleviated by running N × n pop deterministic runs in parallel at each NSGA generation [START_REF] Hercus | Robust shape optimization of uncertain dense gas flows through a plane turbine cascade[END_REF][START_REF] Congedo | Efficient robust optimization techniques for uncertain dense gas flows[END_REF], but:

(i) the required number of computational cores may exceed the computational resources available, and

(ii) even with a perfect parallel scaling at each generation, the turn-around time of the RDO equals at least the average cost of a simulation run multiplied by n gen . To reduce the computational cost, a second (external) surrogate model is introduced to predict the response of the cost functions to the design parameters.

In the following, two RDO methodology are presented: both are built on a SMOGA-based RDO approach introduced in [START_REF] Cinnella | Robust optimization using nested Kriging surrogates: application to supersonic ORC nozzle guide vanes[END_REF], relying on the coupling of two nested Bayesian kriging (BK) surrogates:

the first one is used to compute the required statistics of the objective functions in the uncertain parameter space, while the second one is used to model the response of these statistics to the design variables. Such an approach is sometimes also called "combined kriging" [START_REF] Keane | Robust design optimization using surrogate models[END_REF]. An infill criterion is used to update the second kriging surrogate during convergence towards the optimum. This technique has been successfully applied to the design of turbine blades for organic Rankine cycles [117] and to the RDO of the thermodynamic cycle [START_REF] Serafino | Assessment of an Innovative Technique for the Robust Optimization of Organic Rankine Cycles[END_REF]. Assuming that each kriging surrogate requires a number of available, GEK surrogates can be used to reduce the number of samples for the UQ step; however a GEK-based MOGA is not straightforward in the context of RDO problems, since it requires also the gradient of the statistical moments of the QoI's pdf with respect to the design variables. Obtaining such a piece of information by using efficient adjoint methods is not a trivial task; on the other hand, finite difference approximations are easily applicable, but at the price of a considerable computational expense for high-dimensional design spaces. This is why a new multi-fidelity strategy for RDO is here also proposed: this can drastically reduce the required number of function calls by leveraging an inexpensive (but low-accuracy) first-order MoM with an higher fidelity BK or GEK. Therefore, two different fidelity models are fused together by using a methodology similar to [START_REF] Gratiet | Multi-fidelity Gaussian process regression for computer experiments[END_REF] to generate a surrogate model for the MOGA optimization.

Details of the two strategies are given in the following.

BK-Based Robust Design Optimization

In [117,[START_REF] Cinnella | Robust optimization using nested Kriging surrogates: application to supersonic ORC nozzle guide vanes[END_REF], an external BK surrogate was constructed to describe variations of the mean and variance of the QoIs in the design space, by using n init samples chosen according to a preliminary design of experiments (DOE). In those references, a first-level UQ BK using N samples was used to evaluate the cost functions. Therefore, as in this configuration the method is based on two nested Bayesian Kriging surrogates, the method has been named "TNBK" in [START_REF] Serafino | Assessment of an Innovative Technique for the Robust Optimization of Organic Rankine Cycles[END_REF].

To further reduce the computational cost, a similar approach can be used, whereby a GEK surrogate or the MoM may replace the first-level BK for the UQ step; a flowchart of the resulting RDO process is provided in Figure 3.6. BK is constructed and coupled with the MOGA. In order to control the accuracy of the approximated cost functions, an adaptive infill strategy is adopted to enrich the external BK surrogate during the evolution. For this purpose, at the end of the MOGA iteration a new sample is selected, evaluated and added to the initial DOE: afterwards, the BK model is retrained.

Several adaptive infill criteria exist, like for instance:

• random search,

• prediction minimization,

• maximum entropy,

• lower confidence bound,

• probability of improvement,

• expected improvement,

• cyclic infill.

For an overview about all these techniques the reader is addressed to [184,[START_REF] Rojas-Gonzalez | A survey on kriging-based infill algorithms for multiobjective simulation optimization[END_REF]. Among them, the expected improvement (EI) criterion [START_REF] Močkus | On bayesian methods for seeking the extremum[END_REF] has been selected because it provides a good trade-off between exploitation and exploration.A brief description of the EI infill technique is given below.

To fix ideas, one should consider the minimization of a single-objective function as depicted in Fig. 3.7: initially six samples of the function are available and they are used to build an interpolating response surface and the 95% confidence interval. In this step the minimum of the function is the point marked as y min . To identify which is the most convenient point to be further sampled, a strategy can be to quantify the probability of improvement P (I) as the probability that the new sample provides a lower function response with respect to the current best y min . This P (I) is the shaded area that can be identified in the y -p chart in Fig. 3.7 and it can be considered as an indicator of the probability that the new point can be lower than y min . Therefore, it can be possible to find the point with the global maximization of this area in the design space. However, this approach has the drawback that P (I) is an area and it has the properties of the probability; thus, it can be difficult to evaluate and it can lead to wrong conclusions. As a consequence of this, instead of maximizing a probability it can be more convenient to evaluate its expectancy: this is the EI, whose formulation is provided in Eq.

[187]

.

EI = ∫ ∞ 0 I • P (I)dI (3.10)
Since RDO is intrinsically a multi-objective optimization process, the EI function is a surface of a hyper-dimensional parameter space and a more complex formulation is developed, referred-to as the multi objective expected improvement (MOEI) approach (see [START_REF] Keane | Statistical improvement criteria for use in multiobjective design optimization[END_REF]). The accuracy of the surrogate is rapidly improved with few MOEI infills, so that n MOEI << n gen . A MOEI update implies running the UQ algorithm for the additional sample to be integrated to the second-level BK. The N deterministic direct simulation and adjoint (if needed) runs required for the BK and GEK UQ solvers can be carried out in parallel. Finally, the cost of the RDO based on the adaptive infill strategy in terms of cost-function evaluation is given by N × (n init + n M OEI ) direct CFD calculations for BK, N × (n init + n M OEI ) direct and adjoint CFD calculations for GEK, and (n init + n M OEI ) direct and adjoint calculations for MoM. The turn-around time, in the case of a perfectly-scaling parallel implementation, is approximately equal to n M OEI + 1 runs. Based on numerical experiments [START_REF] Bufi | Robust optimization of ORC turbine expanders[END_REF], a MOEI adaptation every three to five MOGA generations is generally sufficient to achieve an accurate approximation of the optimum, as shown in the next two Chapters.

The full BK-based RDO loop is described in the pseudo-code provided hereafter in Algorithm 1.

Algorithm 1 BK RDO loop with the MOEI adaptation.

1. Initialization: LHS DOE with n init samples S DOE = { s (1) , s (2) , ...,

s (n init ) } While i < n init Select design s (i) in S DOE
Run UQ solver to compute the statistics µ s (i) ,σ 2

s (i)

Multi-Fidelity RDO

Although the BK-based SMOGA allows a considerable reduction of function calls during the RDO loop, the overall computational cost remains significant for high-dimensional design and uncertain spaces. While the computational cost is higher for BK or GEK solvers than for MoM, the former may provide an accurate estimate of the QoI statistics, which is not the case for the first-order MoM.

With the aim of achieving a compromise between the accuracy of kriging-surrogate UQ solvers and the computational efficiency of the MoM, in this section, an advanced SMOGA is introduced, based on a multi-fidelity surrogate model of the design space replacing the preceding external BK surrogate (MF-SMOGA).

In the MF surrogate, the MoM UQ solver is used as the low-fidelity (LF) model, and the GEK (or the BK, if gradients of the cost functions with respect to uncertain variables are not available) UQ solver is the high-fidelity (HF) one. In the present implementation, an initial DoE is run at the first 3.5. ROBUST DESIGN OPTIMIZATION STRATEGIES iteration of the SMOGA, where the HF and LF sampling points are generated independently. The auto-regressive correlation in Eq. 3.7 is used for correcting the discrepancy between the LF and HF models. More specifically, the implementation proposed in [START_REF] Gratiet | Multi-fidelity Gaussian process regression for computer experiments[END_REF] is here adopted, and the LF model is used as a basis function for the regression term expression of a universal kriging model; therefore, the term m(ξ) in Eq. 2.6 becomes:

m(ξ) = ∑ i f i (ξ)β i + f LF β ρ (3.11)
where β ρ is an estimate of the coefficient ρ of Equation (3.7) by means of a GP regression. Assuming that the LF and HF models are independent, the mean µ HF and variance σ 2 HF of the high-fidelity model are given by:

{ µ HF = ρµ LF + µ δ σ 2 HF = ρ 2 σ 2 LF + σ 2 δ (3.12)
with µ LF , σ 2 LF , µ δ , and σ 2 δ the means and variances of the LF model and of the discrepancy function δ, respectively. In order to improve the MF surrogate accuracy during SMOGA convergence, adaptive infill based on the MOEI criterion is used. In this case, however, either the LF or the HF model can be used for the infill. In the following calculations, the strategy of [START_REF] Meliani | Multi-fidelity efficient global optimization: Methodology and application to airfoil shape design[END_REF] is adopted: for the infill, priority is given to the less expensive LF model, and the HF one is used only when improvement achieved with the lower level of fidelity is below a given tolerance, tol = 10 -4 in the present calculations. Either way, re-sampling at the same location is avoided.

Pseudo-code for the MF model is provided below in Algorithm 2; for more information about it, the reader is addressed to [START_REF] Gratiet | Multi-fidelity Gaussian process regression for computer experiments[END_REF].

ROBUST DESIGN OPTIMIZATION STRATEGIES

HF , s (2)
HF , ..., s

(n HF,init ) HF } While i < n LF,init
Select design s

(i) LF in S DOE,LF
Run LF solver to compute the statistics µ(s

(i) LF ), σ 2 (s (i) LF ) Train LF surrogates S LF,µ (S DOE,LF ), S LF,σ 2 (S DOE,LF ) While i < n HF,init
Select design s

Summary of the chapter

In this chapter, after some introductory notions in optimization, two strategies are presented to reduce the RDO computational burden. Both are built on a two nested Bayesian Kriging (TNBK) approach, which is a SMOGA-based RDO approach, relying on the coupling of two nested BK surrogates: the first BK surrogate is used to compute the required statistics of the QoIs in the uncertain parameter space, while the second one is coupled with the optimizer, to model the response of these statistics to the design variables. As an improvement of the TNBK approach, the first proposed RDO technique is the BK-based RDO and it can be considered as a natural extension of the TNBK approach: to reduce the computational cost, the GEK or the MoM can be used to replace the first-level BK for the UQ step, while the second-level surrogates is always a BK directed coupled with the MOGA whose knowledge about the prediction is step by step increased by means of the MOEI criterion. The second proposed RDO technique is a MF SMOGA approach: the MoM UQ solver is used as the low-fidelity model, and the GEK (or the BK, if gradients of the QoIs with respect to uncertain variables are not available) UQ solver is the high-fidelity one; a second-level BK is coupled with the MOGA and its accuracy is improved during the SMOGA convergence with an adaptive infill based on the MOEI criterion. Both RDO strategies are used in the following for the RDO of ORC systems and turbo-expanders. ORC design typically relies on a mathematical model of the cycle, allowing to evaluate a set of performance parameters (or cost functions) given a set of design variables. The model is supposed to provide an accurate enough description of the ORC system over a range of operating conditions. Several examples of thermodynamic and techno-economic optimization can be found in [START_REF] Macchi | Organic Rankine Cycle (ORC) Power Systems[END_REF]. Historically, a black-box strategy has been typically applied, whereby the cycle performance is computed with a simulation code, while the design parameters are optimized with an evolutionary algorithm (e.g. simulated annealing, particle swarm, artificial bees colony or genetic algorithms); some examples can be found in [START_REF] Dai | Parametric optimization and comparative study of organic Rankine cycle (ORC) for low grade waste heat recovery[END_REF][START_REF] Wang | Fluid selection and parametric optimization of organic Rankine cycle using low temperature waste heat[END_REF][START_REF] Wang | Thermodynamic analysis and optimization of an (organic Rankine cycle) ORC using low grade heat source[END_REF][START_REF]Multi-objective optimization of an organic Rankine cycle (ORC) for low grade waste heat recovery using evolutionary algorithm[END_REF][START_REF] Xi | Parametric optimization of regenerative organic Rankine cycle (ORC) for low grade waste heat recovery using genetic algorithm[END_REF][START_REF] Pierobon | Multi-objective optimization of organic Rankine cycles for waste heat recovery: Application in an offshore platform[END_REF][START_REF] Andreasen | Selection and optimization of pure and mixed working fluids for low grade heat utilization using organic Rankine cycles[END_REF]. For several years, the analysis has been carried out solely at a fixed design point, corresponding to the nominal working conditions of the ORC plant. Only recently some thermodynamic and techno-economic optimizations have been carried out considering also part-load performance: a first contribution has been given by [START_REF] Lecompte | Part load based thermo-economic optimization of the Organic Rankine Cycle (ORC) applied to a combined heat and power (CHP) system[END_REF]. Other interesting applications 4.1. PRELIMINARY ASSESSMENT OF THE TNBK RDO STRATEGY FOR A SIMPLIFIED ORC can be found in [START_REF] Manente | An Organic Rankine Cycle off-design model for the search of the optimal control strategy[END_REF][START_REF] Pierobon | Design methodology for flexible energy conversion systems accounting for dynamic performance[END_REF][START_REF] Maraver | Systematic optimization of subcritical and transcritical organic Rankine cycles (ORCs) constrained by technical parameters in multiple applications[END_REF][START_REF] Walraven | Optimum configuration of shell-and-tube heat exchangers for the use in low-temperature organic Rankine cycles[END_REF][START_REF] Larsen | A comparison of advanced heat recovery power cycles in a combined cycle for large ships[END_REF][START_REF] Walraven | Economic system optimization of air-cooled organic Rankine cycles powered by low-temperature geothermal heat sources[END_REF][START_REF] Martelli | Numerical optimization of Combined Heat and Power Organic Rankine Cycles -Part A: Design optimization[END_REF][START_REF] Capra | Numerical optimization of combined heat and power Organic Rankine Cycles -Part B: Simultaneous design & part-load optimization[END_REF]. Part-load performance is evaluated by means of dynamic simulations [START_REF] Macchi | Organic Rankine Cycle (ORC) Power Systems[END_REF] or, more often, through steady state calculations at off design conditions [START_REF] Martelli | Numerical optimization of Combined Heat and Power Organic Rankine Cycles -Part A: Design optimization[END_REF]. In all cases, the optimization process is deterministic in the sense that the design variables, the operating conditions and any other input required by the model, as well as the model itself, are supposed as perfectly known.

Because of its high cost, this approach is as far as possible avoided and replaced by some assumptions and approximations. For instance, in case of geothermal energy, WHR from continuous industrial processes, or a concentrating solar plant with a large solar multiple and a large storage, the ORC is usually thought to work close to nominal conditions and it is common practice to use average values over a period and to account for off-design penalization by reducing the number of equivalent working hours [START_REF] Macchi | Organic Rankine Cycle (ORC) Power Systems[END_REF]. However, as already seen, several can be the sources of uncertainty affecting ORC systems.

In the present chapter, a preliminary applications of RDO to two ORCs is presented. The first one is described using a very simple (and computationally cheap) model, which allows direct comparisons with a "brute force" method relying on pure Monte Carlo sampling. The second one is an ORC of practical interest for geothermal applications. It is described by a more complete (and costly) model.

In both cases, only the TNBK RDO approach is used, since the cost of a model run remains relatively low compared, e.g. to CFD simulations. The results presented in the following have been published in [START_REF] Serafino | Assessment of an Innovative Technique for the Robust Optimization of Organic Rankine Cycles[END_REF] and [START_REF] Serafino | Robust optimization of an organic rankine cycle for geothermal application[END_REF].

Preliminary assessment of the TNBK RDO strategy for a simplified ORC

Since for the moment the aim is the proof of concept assessing and validating the TNBK RDO approach, this methodology has been here compared with MC sampling: this is considered as a reference technique but it has the huge drawback of being seriously costly, therefore a simplified ORC model has been used as a test case in this step. Q in = m gas c p gas (T in gas -T out gas ) (4.1)

1. the mass flow rate of the flue gas m gas , 2. its heat capacity c p gas , 3. the temperature difference (T in gas -T out gas ), where T in gas is the flue gas temperature at the inlet of the thermal loop (see Fig. 4.1) which is here considered uncertain, while T out gas is the flue gas temperature at the outlet of the thermal loop, before being discharged to the chimney, which is here considered deterministically constrained at 430 K.

The mass flow rate of the flue gas m gas and the flue gas inlet temperature T in gas are here considered sources of aleatory uncertainty, while the flue gas heat capacity c p gas is a source of epistemic uncertainty: they all affect the heat exchange between flue gas and thermal oil, occurring in the thermal loop as shown in Fig. 4.3.

To complete the model of this component few further assumptions are taken into account: the heat exchanger is supposed to operate with two imposed deterministic ∆T between the fluids on both its extremities, which are respectively 110 K on the side where hot flue gas enters and 40 K on the other side; moreover, all thermodynamic properties of thermal oil are here considered deterministically known.

Since an Air Cooled Condenser (ACC) has been here considered, another source of aleatory uncertainty that must be taken into account is the environmental temperature T env which has a large effect on the ORC condensing temperature T cond , modelled in Eq. 4.2 as the sum of three elements:

T cond = T air in + ∆T air + ∆T pp cond (4.2)
1. T air in is the air inlet temperature in the ACC and it is considered equal to the environmental temperature T env and hence uncertain, 2. ∆T air is the air ∆T in the ACC, SIMPLIFIED ORC For sake of simplicity, the ACC is supposed to be always ON and the driver of the fan does not allow to regulate its speed: in this case, the lower is the air density, the lower is the air mass flow rate available at the condenser. Thus, ∆T air is not given by a complex multi-parameter function as usual, but it is here approximated with a deterministic linear function depending only on T env whose formulation is provided in Eq.4.3.

∆T air = 15.0 + 1 3 (T env [K] -275) (4.3)
Finally, another source of epistemic uncertainty here taken into account is the pump isentropic efficiency η pump .

As each uncertain quantity can be modelled by means of an elementary pdf or a combination of some of them [START_REF] Montgomery | Applied statistics and probability for engineers[END_REF], it is assumed that the uncertain nature of m gas , T in gas , T env and η pump is totally

Gaussian ∼ N(50 kg s -1 , 25 kg 2 s -2 ) 2 T in gas Gaussian ∼ N(640 K, 25 K 2 ) 3 c p gas Uniform in [1100 -1400] J kg -1 K -1 4 T env Gaussian ∼ N(290 K, 25 K 2 ) 5 η pump Gaussian ∼ N(0.7, 0.0016)
The TNBK RDO approach is general and more sources of uncertainty could be considered, like for instance the thermodynamic properties of the working or turbine performance; however, since the aim of the present analysis is to validate the methodology for the RDO of ORCs, only the five uncertainties listed in Tab. 4.1 have been here considered. Other employed assumptions are listed below:

1. at the outlet of the recuperator, temperature T 8 * is fixed 3 K higher than temperature T 2 , 2. for each heat exchanger, on the working fluid side a fix pressure drop (1% of the inlet pressure) is imposed, 3. for each heat exchanger, no thermal losses are considered, 4. condenser auxiliaries are neglected, as the electrical work required by the ACC is always constant and independent from T env .

5. the turbine isentropic efficiency η turb is assumed to be always 85% at the design point, while for the off-design simulation it is assumed to be just a deterministic function of the pressure ratio

η ORC = [h 6 (P 6 , T 6 ) -h 7 (P 7 , T 7 )] -[h 2 (P 2 , T 2 ) -h 1 (P 1 , T 1 )] h 6 (P 6 , T 6 ) -h 3 (P 3 , T 3 ) (4.4)
An in-house code written in Python has been utilized to simulate the ORC as just described; the thermodynamic properties of the working fluid are calculated with an equation of state written in terms of the Helmholtz energy [START_REF] Mondéjar | Thermodynamic properties of trans-1-chloro-3,3,3-trifluoropropene (r1233zd(e)): Vapor pressure, (p, rho, t) behavior, and speed of sound measurements, and equation of state[END_REF] by means of Coolprop 6.1.0 [132] and they are all assumed devoid of every kind of uncertainty.

Methodology

The only QoI of the rdo for the ORC just presented is the cycle efficiency η ORC : given all deterministic and uncertain inputs listed before, the goal of the optimization is to find the values of the pressure at the evaporator P ev and the superheating ∆T sh = T 6 -T 5 that maximize the mean value of the pdf of η ORC minimizing at the same time its variance.

The "brute-force" approach, described and applied for instance in [START_REF] Bufi | Robust optimization of an Organic Rankine Cycle for heavy duty engine waste heat recovery[END_REF], is generally considered the most reliable strategy to perform a rdo; unfortunately, it usually proves to be also computationally

PRELIMINARY ASSESSMENT OF THE TNBK RDO STRATEGY FOR A SIMPLIFIED ORC

too expensive to be used. In fact, it deals with all the uncertainties by means of a MC method which needs, considering the high-dimensionality of the problem, a Design of Experiment (DOE) composed at least by hundreds thousands samples to build the distribution of the QoI that will be afterwards used to compute the statistical moments (mean and variance) of the ORC efficiency. As a consequence of this, assuming for example that for each simulation a calculation time of 0.5 seconds is required on a standard PC of a small/middle size company, just one MC takes an order of magnitude of 10 5

seconds to generate just one PDF of the QoI. Moreover, as already said, for this kind of problems a MOGA is usually used as optimizer, which requires a large number of individuals (depending on the dimensionality of the problem) for several generations, where each individual is a MC simulation. In order to speed-up the robust design methodology, the TNBK approach is here considered.

Moreover, since the NSGAII used as a MOGA is an evolutive optimization algorithm, design parameters have been bounded as follows: between 15.0 bar and 28.0 bar for P eva and between 1.0 K and 50.0 K for ∆T sh .

To validate the TNBK methodology, which has never been employed for ORCs, its results are here compared with the ones obtained by means of the "brute-force " approach; in order to make this comparison possible, an inexpensive simplified ORC model must be used, justifying all the assumptions already presented in the previous section.

Results

First the BK used for the uncertainty quantification step has to be validated. Therefore, several Kriging response surfaces are built for various combinations of couples of the design parameters: for each design, the BK surrogate is constructed from an initial LHS DOE containing different numbers of samples, in order to investigate the influence of this parameter on the approximation error. The error is therefore calculated using a full factorial grid of 10 5 points to evaluate the prediction on the BK response surface of the QoI, which is the ORC efficiency: these predicted values are compared with the exact ones, that are directly obtained running the real ORC model on the same grid.

This validation results to be especially satisfactory for DOEs with 120 samples allowing to obtain comparisons like the one depicted in Fig. 4.5, where values predicted by surrogate models (red dots) should lie on the quadrant bisector (blue dotted line) when they correspond to the exact solution. As depicted in Fig. 4.8, both Pareto fronts are very similar and consequently they agree providing the same result: for the problem object of the present work is it possible to find an high-efficiency design which is also robust. Actually, the shape of the Pareto fronts suggests a further observation:

for the test case here studied it is likely that a deterministic optimization of the ORC can converge to the same results obtained with the RDO; this outcome is probably due to the strong assumptions used to model the ORC. However the main goal of the present analysis is the validation of the TNBK approach for the RDO of ORC systems; therefore, no deterministic optimization has been carried out.

Finally, in Fig. 4.9 and Fig. 4.10 the optimal solutions from both the "brute-force" method and 4.2, the TNBK approach requires a DOE for the second-level BK composed by only 15 samples, where each of them use a 120 samples DOE to train the first-level BK for the UQ process: with just 1800 evaluations of the ORC model this method can provide a satisfactory prediction completely in agreement with the "brute force"approach, which on the other side needs for the optimization 30 individuals evolving for 60 generations and each one requires a 10000 samples DOE to quantify the uncertainty with the MC method. Consequently, in this case the potentiality of the TNBK technique with respect to the "brute force" method lies in its capability to speed-up the RDO process theoretically up to 10000 times.

ROBUST DESIGN OPTIMIZATION OF AN ORC FOR GEOTHERMAL APPLICATIONS

Robust design optimization of an ORC for geothermal applications

The development of ORCs in the last fifteen years has contributed significantly to the enlargement of the commercial exploitation of geothermal power. Specifically, the ORC technology has fostered the use of low-temperature geothermal resources. Data collected by Bertani [START_REF] Bertani | Geothermal power generation in the world 2005-2010 update report[END_REF][START_REF]Geothermal power generation in the world 2010-2014 update report[END_REF] show a noticeable growth in the worldwide installed geothermal capacity of binary plants (almost exclusively based on ORC technology), which in just 5 years has almost doubled passing from 11% to 14.2% of the overall geothermal applications. Among the 780 MW of new geothermal capacity installed in 2016 [START_REF]Geothermal power: Technology brief[END_REF], almost 30% consists in ORC systems [1].

However, even if nowadays geothermal power is a mature, commercially available and well known technology providing low-cost base load capacity, several sources of uncertainty are hidden in the whole process adopted for its exploitation. Some major sources of epistemic uncertainties in geothermal power generation [START_REF] Barr | Coping With Uncertainty in Geothermal Field Development[END_REF] are listed below:

• temperature of the geothermal reservoir;

• field extent for the definition of the drilling area;

• soil permeability;

• average well production;

• re-injection cost;

• phenomena like quenching, chemical clogging and corrosion.

Among them, just a few can be reduced, investing a considerable amount of time and capital resources in preliminary discovery and exploration activities, while some others can be reduced only

ROBUST DESIGN OPTIMIZATION OF AN ORC FOR GEOTHERMAL APPLICATIONS

through long term operation of the field (about a decade, as an order of magnitude). All of them are due to a lack of knowledge, and for that reason they are qualified as "epistemic" uncertainties [START_REF] Of | Guide for the Verification and Validation of Computational Fluid Dynamics Simulations[END_REF][START_REF]Standard for Verification and Validation in Computational Fluid Dynamics and Heat Transfer[END_REF].

A large amount of literature has been written since the 70s to identify the sources of uncertainty affecting geothermal power and to quantify their effects. A detailed overview can be found for instance in [START_REF] Robertson-Tait | Managing Geothermal Resource Risk -Experience From the United States[END_REF]. These sources of uncertainty can lead to the scenario presented in Fig. 4.11, showing that from an investor viewpoint a geothermal project can be a risky and capital-intensive investment. The risk derives from the fact that an important part of the capital is required for preliminary activities like pre-survey, exploration and test drilling without any certainty about the presence of an exploitable geothermal resource. For instance, the test drilling can account alone for up to 15% of the overall capital cost [START_REF]The World Bank -Energy Sector Management Assistance Program[END_REF], before that the project profitability can be determined. Moreover, these sources of uncertainty can result in the oversizing of the geothermal plant, with a significant reduction of profits and a possible failure of the whole geothermal project. Considering that, historically, the majority of the worldwide geothermal installed capacity has been funded mainly through private financing [START_REF] Sanyal | Resource risk and its mitigation for the financing of geothermal projects[END_REF],

such a scenario can deter investments in this technology.

ROBUST DESIGN OPTIMIZATION OF AN ORC FOR GEOTHERMAL APPLICATIONS

Even if there is no reason to believe that the uncertainty and the risk associated with geothermal fields are any greater than those for other forms of electricity generation [START_REF] Barr | Coping With Uncertainty in Geothermal Field Development[END_REF], the accurate quantification of geothermal resource risk is of a paramount importance in the financing of geothermal projects [START_REF] Sanyal | Quantification of geothermal resource risk -A practical perspective[END_REF]. The objective of the present work is to employ a promising methodology aiming to perform the RDO of an ORC for a geothermal application, which can take into account all these uncertainties in order to reduce the associated risk.

In the context of ORC design under uncertain operating conditions, RDO represents a promising alternative to standard design methods [START_REF] Bufi | Robust optimization of an Organic Rankine Cycle for heavy duty engine waste heat recovery[END_REF], allowing to ensure a more stable performance over a range of randomly varying operating conditions. In the following the ORC model for geothermal applications is presented in Section 4.

ORC model

The object of the RDO carried out in this study is an ORC for geothermal application, whose plant layout is depicted in Fig. 4.12.

Figure 4.12: Layout of the geothermal plant

The plant exploits geothermal brine in a single-pressure ORC employing iC4 as the working fluid.

ROBUST DESIGN OPTIMIZATION OF AN ORC FOR GEOTHERMAL APPLICATIONS

As depicted in Fig. 4.12, the brine is cooled-down before the re-injection by a set of heat exchangers (HEXs). These include a pre-heater (PRE), an evaporator (EVA) and, eventually, a super-heater (SH)

where heat is transferred to the ORC working fluid. The PRE is a shell and tube heat exchanger, while the EVA is a kettle reboiler. The vapour of the working fluid, which can be saturated or superheated (depending on the outcomes of the RDO), expands in an axial turbine directly connected with the generator by means of a coupling disc. The condenser (CD) is a shell and tube heat exchanger using water, which is cooled in a wet cooling tower to condensate the vapour coming from the discharge of the turbine. Finally, the pump is of the centrifugal type and it is connected to an electrical driver equipped with a variable-frequency drive.

An in-house deterministic model for ORC systems written in Python is used for the design and the off-design simulation of the cycle, described in [START_REF] Fontaine | Design and Simulation Model of Medium Scale Organic Rankine Cycles -Validation on Waste Heat Recovery Plant and Case Studies[END_REF]. A brief overview of both design and simulation algorithms is provided hereafter for completeness.

The design algorithm, which is sketched in Fig. 4.13, defines the proper size of all ORC components, in order to evaluate the performance and the cost of the whole system. Therefore, first the nominal operating conditions of the cycle should be assigned; these include the characteristics of both hot and cold energy sources and, more specifically, the nominal temperature, mass flow rate and specific heat capacity of the geothermal brine, the brine re-injection temperature (which is fixed at 353.15 K, to avoid corrosion and deposition issues), the inlet temperature and the mass flow rate of the cooling water. Additional inputs required by the model are:

• the subcooling,

• the superheating (∆T SH ),

• the isentropic efficiency of the turbine,

• the efficiency of the pump,

• the pressure drops across the HEXs,

• the pinch points in the evaporator and in the condenser.

• the working fluid. The turbine is designed to perform 85% isentropic efficiency at the design point, while the pump is chosen with a 70% overall efficiency. The PRE, the EVA and the SH (considered all together

ROBUST DESIGN OPTIMIZATION OF AN ORC FOR GEOTHERMAL APPLICATIONS

as a single block of HEXs) and the condenser are modelled with the generalized moving-boundary algorithm proposed in [START_REF] Bell | A generalized moving-boundary algorithm to predict the heat transfer rate of counterflow heat exchangers for any phase configuration[END_REF]; the ϵ-NTU method [START_REF] Incropera | Fundamentals of Heat and Mass Transfer[END_REF] is then used to estimate the exchange area at given by Enertime are used to calculate the pump and turbine off-design isentropic efficiencies. The turbine simulation allows to compute the ORC mass flow rate; with this last piece of information, one can simulate the pipelines between the turbine and the condenser and the one between the pump and the HEXs, obtaining an estimation of pressure drops in each line knowing its own diameter.

Furthermore, the HEXs and the condenser can be also simulated: the same generalized movingboundary algorithm already employed in the design mode, is now used to define the thermodynamic state at the boundaries and inside the HEXs block and the condenser, while the ϵ-NTU method is still utilized to compute the exchange area of both these two items. Finally, the areas of the HEXs and of the condenser are used to check convergence: if the values calculated in the simulation mode differ from the design values at nominal conditions by more than 10 Diameter of the pipe between CD and pump Diameter of the pipe between turbine and CD Diameter of the pipe between pump and EVA COBYLA optimization algorithm [START_REF] Powell | A view of algorithms for optimization without derivatives[END_REF] available in the Python library Scipy is used.

In the following, the deterministic ORC model is coupled with an uncertainty quantification method, described in the next section, to account for uncertainties in the nominal conditions.

ORC Uncertainty quantification

For the ORC model here considered, several can be the sources of uncertainty and the uncertain output QoIs. Here below the selection of performance parameters that will serve as cost functions for the RDO is first presented; afterwards the sources of uncertainty are identified. This information is mandatory before propagating the uncertainty through the model with the BK technique. the auxiliaries like the cooling tower system, the oil system, the air compressors and the control room system. W ORC,net is here preferred to the global cycle efficiency, defined as the ratio of the ORC net power to the total heat content of the heat source: in fact, for a randomly varying thermal input, global efficiency does not necessarily correspond to the maximum power output operation of the ORC system. This happens because during the optimization, in case of an increase in the evaporation pressure, an increased cycle efficiency can be obtained at the cost of a decrease in the amount of heat extracted from the heat source. Therefore, instead of considering a possible trade-off between these two quantities, here it is chosen to work only with W ORC,net , as it is an extensive quantity which is easier to be used in the optimization with respect to any possible constraint about the size of the plant.

To evaluate the influence of uncertain inputs on the selected QoI, the ORC model is first run in

ROBUST DESIGN OPTIMIZATION OF AN ORC FOR GEOTHERMAL APPLICATIONS

the design mode, with nominal conditions fixed to the average values of the input parameters. This allows to determine the geometrical characteristics of the ORC. Afterwards, the ORC model is run several times in the simulation mode for a set of conditions sampled from the input PDFs, and the corresponding values of the ouput QoI are used to build an approximation of its PDF and of its main statistical moments (mean and variance).

Identification of sources of uncertainty

The ORC system under investigation is affected by two major uncertainties: the first is of epistemic nature, as it is due to the lack of complete information about the characteristics of the geothermal heat source; the second one is aleatory and associated with the random variations of the temperature of the cold medium at the condenser.

About the first uncertainty, Sanyal and Morrow [START_REF] Sanyal | Quantification of geothermal resource risk -A practical perspective[END_REF] conducted a local sensitivity study on the internal rate of return for a typical geothermal project, and illustrated that the most sensitive variables are resource parameters. Hereafter this uncertainty is considered as the outcome of the interaction of three sources of epistemic uncertainty, namely:

1. mass flow rate of the geothermal brine m ̇geo , 2. temperature of the geothermal brine T geo , 3. heat capacity of the geothermal brine c geo .

Based on information available in the technical literature (see for instance [START_REF] Vogt | Quantifying Uncertainty in Geothermal Reservoir Modeling[END_REF][START_REF] Foerster | Analysis of borehole temperature data in the Northeast German Basin: Continuous logs versus bottom-hole temperatures[END_REF][START_REF] Lee | Heat flow and heat production in the Arkoma Basin and Oklahoma Platform, southeastern Oklahoma[END_REF][START_REF] Olsen | Geothermal reservoir assessment based on slim hole drilling[END_REF]), the preceding quantities are modelled as independent Gaussian random variables, whose parameters are listed in Tab. 4.4.

Regarding the temperature of the cooling water T cw in at the inlet of the condenser, the observation of historical data could be used to infer a probability distribution modelling for this source of aleatory uncertainty. However, for the present work, no data were available; therefore, some information found in literature [START_REF] Olsen | Geothermal reservoir assessment based on slim hole drilling[END_REF][START_REF] Ashrae | HVAC Systems and Equipment: SI Edition, ser. ASHRAE Handbook of Heating, Ventilating and Air-Conditioning Systems and Equipment SI[END_REF] was used to model this quantity as a Gaussian random variable, with parameters also reported in Tab. 4.4.

ROBUST DESIGN OPTIMIZATION OF AN ORC FOR GEOTHERMAL APPLICATIONS

ORC optimization strategy

In order to design an ORC with a stable performance under the uncertain inputs, a RDO strategy is considered: the ORC model is coupled with an UQ method and an optimization algorithm to select design parameters that maximize (or depending on the problem, minimize) the expected (average) value of a set of design criteria while optimizing some of their statistical properties, which are then added to the set of the cost functions.

In the following, Taguchi's RDO criterion is adopted an it is applied by means of the TNBK approach. For the optimization, the following cost functions (CFs) are taken in account: The RDO algorithm, whose flowchart is depicted in Fig. 4.17, proceeds as follows:

1. A nominal design point is randomly sorted from the input parameter distribution.

2. The ORC model is run in the design mode with deterministic inputs, determined at Point 1.

3. As an outcome of Point 2, the ORC characteristics are defined, including A ORC .

4. The operating conditions are then randomly perturbed, according to their PDFs.

5. The random samples are propagated through the ORC model in simulation mode.

ROBUST DESIGN OPTIMIZATION OF AN ORC FOR GEOTHERMAL APPLICATIONS

6. The cost functions CF1 and CF2 are evaluated as statistical moments of the QoI W ORC,net .

7. Based on the values of the cost functions, a new set of design variables is selected and the process is repeated from Point 2.

The NSGAII algorithm has been chosen as the MOGA; therefore the design parameters space is defined by the support of the PDFs for the uncertain parameters and the interval [0, 4] K for ∆T SH .

In the genetic algorithm, an initial population of 40 individuals (alternative nominal designs) is let to evolve over 80 generations. The evaluation of the cost functions for each individual corresponds to an UQ calculation, as described in Section 4.2.2. To alleviate the computational cost of the optimization process, the RDO is conducted by using the TNBK approach, as in [117,[START_REF] Serafino | Assessment of an Innovative Technique for the Robust Optimization of Organic Rankine Cycles[END_REF]. The MOEI technique has been applied each 5 generations to increase the accuracy of the external BK surrogate.

UQ results and global sensitivity analysis

The sensitivity of the selected QoI, i.e. the ORC net power W ORC,net , to random variations of the uncertain parameters described in Section 4. The UQ results are also used to carry out a global sensitivity analysis of the QoI to the uncertain parameters. For that purpose, the analysis of variance (ANOVA) decomposition [START_REF] Liu | Estimating mean dimensionality of analysis of variance decompositions[END_REF] is applied and the Monte Carlo sampling is used to estimate Sobol indexes associated to the uncertain inputs. More specifically, for each input parameter the first-order and total Sobol indexes are extracted. carried out in three design points corresponding respectively to average, minimum and maximum conditions for the geothermal source, according to the PDFs in Tab.4.4, since it allows to maximize the plant profitability avoiding a dramatic increase in the variance of the QoI PDF.

ROBUST DESIGN OPTIMIZATION OF AN ORC FOR GEOTHERMAL APPLICATIONS

SUMMARY OF THE CHAPTER

Summary of the chapter

In this chapter the TNBK approach has been validated by means of a comparison with a "bruteforce" approach for a simplified ORC cycle. More specifically, it has been shown that the TNBK technique provides the same results of the "brute-force" MC approach with a potential speed-up of the optimization process by a factor 10000. Secondly, the TNBK strategy has been successfully applied to an ORC for a geothermal application, affected both by epistemic and aleatory uncertainty: a preliminary global sensitivity analysis shows that the most influential parameters are T geo , m ̇geo and T cd , while c geo can be neglected. These results are used to reduce the number of optimization parameters that are considered in the RDO of the ORC. Afterwards, TNBK RDO of the cycle is carried out to generate a Pareto front of solutions, and the one offering the best trade-off among all optimization targets is selected. This is then compared with the optimal solutions from the deterministic optimization, which are always outperformed. In particular, the best RDO solution wins against the best deterministic optimum, increasing the mean value of the PDF of the QoI (the power production) by 1.5%, while its standard deviation is reduced by 8.5% and the surface of the heat exchangers by 34%. In the RDO of ORCs performed in the previous Chapter, the TNBK approach has allowed a dramatic reduction in the required computational resources in comparison to the "brute-force" MC approach; however, this strategy continues to be too expensive for extremely expensive simulation codes, like for instance CFD applications. In this case, a first contribution can be the adoption of UQ methods using the gradient information, i.e. GEK or MoM instead of BK. In addition, more advanced RDO techniques should be adopted. Therefore, an assessment of the MF-based RDO methodology (see Chapter 3) is here provided. The results have been published in [START_REF] Serafino | Multi-fidelity gradient-based strategy for robust optimization in computational fluid dynamics[END_REF]. Finally, the shape optimization of an ORC turbine blade is carried out by means of the validated MF-based RDO strategy.

ASSESSMENT OF GRADIENT-BASED RDO STRATEGIES FOR AN IDEALIZED TURBINE NOZZLE

Assessment of gradient-based RDO strategies for an idealized turbine nozzle

Several RDO techniques already presented in Chapter 3 are assessed against an inexpensive test problem (also studied in [START_REF] Cinnella | Hybrid Adjoint-based Robust Optimization Approach for Fluid-Dynamics Problems[END_REF][START_REF] Cinnella | Robust optimization using nested Kriging surrogates: application to supersonic ORC nozzle guide vanes[END_REF]), namely, the inverse design of a supersonic quasi-1D diverging nozzle.

This work is organized as follows: in section 5.1.1 the test-case is presented; some information about the UQ methods considered and the SMOGA strategies are provided in Section 5.1.2. In section 5.1.3, first various UQ methods are applied to the test configuration and compared in terms of accuracy and computational costs; afterwards, such methods are combined with a SMOGA or MF-SMOGA and assessed their efficiency in solving the RDO problem.

Problem definition

The nozzle geometry is assigned through the area distribution S(x) along the longitudinal axis x.

This is chosen to be of the form:

S(x) = a + b tanh(cx -d) ( 5
J = 1 2 ∫ L 0 (P -P des ) 2 dx (5.2)
The optimization goal is to determine the design parameters a, b, c, d in Eq. 5.1 providing the best fit to the target pressure distribution under multiple uncertainties, in the sense of Eq. (3.9).

The flow is assumed to be governed by the Euler equations for quasi-1D flows (Eq. 5.3):

R(w, x) = ∂(f S) ∂x + K dS dx = 0 (5.3)
where w = [ρ, ρv, ρe t ] T and f (w) = [ρv, ρv 2 + P, ρvh t ] T are, respectively, the conservative variable and the physical flux vectors [START_REF] Hirsch | Numerical Computation of Internal and External Flows[END_REF]. In the preceding equations, ρ is the fluid density, v is the velocity along the nozzle axis, e t and h t are the total specific energy and enthalpy, and K = [0, -P, 0] T . The system of equations is supplemented by the equation of state for thermally and calorically perfect gases,

P = (γ -1)ρ(e T -v 2 /2).
The governing equations are discretized by a cell-centered finite volume formulation, using Rusanov's first-order upwind scheme for space integration and four-stage explicit Runge-Kutta time-stepping [START_REF] Hirsch | Numerical Computation of Internal and External Flows[END_REF]. Characteristic boundary conditions based on Riemann invariants are imposed at nozzle inlet and outlet. Sonic flow conditions are prescribed at the inlet, so that all Riemann invariants enter the domain. The range of variation of the total pressure is such that a shock is always created in the divergent. As a consequence, outlet flow conditions are always subsonic. In this case, we impose the outlet static pressure, which is treated as deterministic, and fixed to P out = 0.6

bar. Based on a preliminary mesh study, a computational grid of 300 uniformly spaced cells is used in all of the following calculations. Despite the presence of the shock, a first order model has been preferred to more accurate high-resolution schemes, to unsure a satisfactory robustness in the solution of both direct and adjoint problems.

The system is assumed to be subject to uncertainties of various nature, specifically:

• geometric tolerances on the nozzle shape, modelled by treating the shape parameters a, b, c, d as normally distributed random variables, with mean µ and coefficient of variation CoV = σ/µ,

ASSESSMENT OF GRADIENT-BASED RDO STRATEGIES FOR AN IDEALIZED TURBINE NOZZLE

this issue and to facilitate the use of gradient-based UQ, including MoM, it is chosen to construct instead a separate surrogate model mapping a variable in the design space to the cost function space.

Therefore, two RDO strategies are first performed: these are the GEK+BK and the MoM+BK, which are both compared with the already validated TNBK approach. Finally, a MF based RDO method is built, using the MoM UQ solver as the low-fidelity (LF) model and the GEK UQ solver is the high-fidelity (HF) one.

The gradient information required for the GEK and MoM approaches is calculated by means of the adjoint method. Both the continuous and the discrete adjoint formulations are used in the following.

Their derivation for the present quasi-1D nozzle problem is described in Appendix A.

Results

The quasi-1D supersonic nozzle test problem presented in Section 5.1.1 is first used to assess the BK, GEK and MoM UQ methods against reference MC sampling. Afterwards, the methods are applied to a sample of nozzle geometries and used to build single or multi-fidelity surrogates used in the SMOGA RDO loop.

Preliminary assessment of UQ methods

One of the nozzle geometries is selected in the design space by assigning the geometric coefficients fixed normal pdfs with standard deviation equal to 1% of the mean. The pdf parameters for the geometric variables are provided in Tab. 5.2. The operation and thermodynamic parameters are assigned the same pdf as in Tab. 5.1. The pdfs are sampled in order to build BK and GEK surrogates, while the MoM is applied by computing the QoI and its gradient at the expected value of the input parameters.

ASSESSMENT OF GRADIENT-BASED RDO STRATEGIES FOR AN IDEALIZED TURBINE NOZZLE

A summary of the UQ results is given in Tab. 5.3, where the approximated mean and variance of the QoI J according to the various UQ methods are reported. A reference calculation based on MC integration over 10 5 samples is carried out to provide a reference solution. For the present cheap test problem, the CPU time required for MC sampling is approximately 3 × 10 6 seconds on a personal computer having an Intel(R) Xeon(R) CPU E5-1620 v3 at 3.50 GHz. Results are also reported for MC sampling on the BK surrogate. The latter uses 60 function evaluations at points selected according to a latin hypercube sampling (LHS) of the six-dimensional uncertain space. The sample size corresponds to the empirical rule N = 10n unc . This is already sufficient to achieve extremely low errors with respect to the MC mean and variance, while reducing the overall computational cost of the UQ by three orders of magnitude. The GEK UQ based on the discrete adjoint solver provides an accuracy similar to BK by using only 15 samples. Despite the additional adjoint solves for gradient computations, the computational cost is reduced by more than 1/3 with respect to BK. The accuracy of GEK is confirmed by inspection of figure 5.2, showing the full empirical pdfs of J computed with PC, BK and GEK. A very good agreement between MC, BK and GEK is observed.

The accuracy is less satisfactory for the continuous adjoint approach, due to the less accurate computation of gradients. This reflect numerical errors introduced by the discretization of the adjoint equations and the treatment of boundary conditions, which are not completely consistent with the discretization errors introduced by the direct CFD solver. On the other hand, the continuous adjoint solver is developed independently on the direct solver and, in this sense, it is non intrusive. The computational cost of GEK samples using direct or continuous adjoint is approximately the same.

The MoM method is obviously less accurate than the other UQ solvers. Nevertheless, the errors on both the computed mean and variance remain very reasonable, despite the presence of a shock in the divergent (where the Taylor-series expansion is not defined) and the relatively large uncertainty range on the reservoir pressure (approximately 20%). The shock is in practice regularized by numerical diffusion both in the direct and the adjoint solvers, while the uncertainty ranges of all other parameters are reasonably small. While the discrete and continuous-adjoint MoM calculations provide strictly the same results for E[J] (which does not use gradient information), they predict different results for the standard deviation. The slightly lower errors obtained for the continuous adjoint method is the effect of a compensation of errors for the case at hand. Since MoM requires only one direct and one adjoint CFD calculation, its computational cost is essentially 1/15 of the GEK sampling. Overall, the first-

ASSESSMENT OF GRADIENT-BASED RDO STRATEGIES FOR AN IDEALIZED TURBINE NOZZLE

order MoM provides a reasonably accurate estimate of the lower order statistics of the QoI and a very good tradeoff between cost and accuracy for the present problem. Its accuracy is however expected to decrease for more severe uncertainty ranges. For this reason, MoM is categorized as a lower fidelity model than BK or GEK. The UQ results for MC, BK and GEK can be used to carry out a global sensitivity analysis and to identify the random parameters contributing the most to the variance of the QoI by means of an ANOVA decomposition. Specifically, exact or surrogate-based MC samples are employed to calculate the Sobol Indexes in the parameter space defined by all uncertain inputs. For this purpose the SALib Python library [START_REF] Herman | SALib: An open-source python library for sensitivity analysis[END_REF] has been adopted. Only the discrete GEK results (reported in Fig. 5.3) are considered here, the MC, BK, and continuous GEK methods leading to similar conclusions. The figure reports the first order Sobol indexes with respect to each uncertain parameter, and the sum of the higher-order indexes, corresponding to interaction between parameters when these are changed

ASSESSMENT OF GRADIENT-BASED RDO STRATEGIES FOR AN IDEALIZED TURBINE NOZZLE

simultaneously. The inlet total pressure P T,in appears to be the most influential parameter when taken alone. This is consistent with the much larger range of its pdf. Nevertheless, the interactions terms are very significant. This is due to the highly nonlinear nature of the compressible CFD problem. Due to such strong inter-dependency of the QoI on the input parameters, it is not possible to neglect any of them. As a consequence, the RDO calculations of the following section treat all of the six input parameters as random variables.

ASSESSMENT OF GRADIENT-BASED RDO STRATEGIES FOR AN IDEALIZED TURBINE NOZZLE

displays the values of the objective functions computed on the surrogate models. BK-BK, GEK-BK and MF surrogates lead to close estimates of the objective functions, with discrepancies of less than 1% on the mean and less than 2% for the variance. The CPU cost of the different RDO algorithms on a single processor personal computer are reported in the same table. Due to the reduced number of samples used in the UQ runs, the GEK algorithm allows to gain a factor 2 to 2.5 over the BK.

A further reduction of a factor 15 or larger is obtained using MoM, but caution must be taken in systematically preferring such a method for RDO problems, because of its lower accuracy. Finally, the MF SMOGA RDO has a computational time 8 times smaller than BK and 4 times smaller than GEK while ensuring similar accuracy to the BK SMOGA.

Table 5.4: Optimal solutions of the RDO according to the various strategies (DA = Discrete adjoint, CA = Continous adjoint).

UQ method E[a] E[b] E[c] E[d] E[J] var[J]

Optimization As a further verification of the RDO results, the optimal design is recomputed using MC sampling, resulting in E[J] = 0.32325 and var[J] = 0.09487, which is in rather good agreement with the SMOGA estimates, appearing to be slightly over-optimistic in predicting the objective functions. Indeed, the full pdf of the QoI computed by propagating the uncertain parameters through the CFD solution for the optimal nozzle geometry by using the BK and GEK methods (shown in Fig. 5.4 alongside the MC distribution) exhibits moderate deviations from the MC one. Globally, the optimal distributions are much closer to 0 in average and they exhibit a smaller variance than the baseline geometry investigated in the preceding Section, showing that the RDO effectively improves both criteria.

Finally, the optimized geometry is depicted for all of the employed methodologies in Fig. 5.

RDO OF AN ORC TURBINE BLADE

RDO of an ORC turbine blade

Internal flows in turbomachinery components represent a challenge for numerical models. Such flows are typically highly three-dimensional, unsteady and characterized by strong interactions between vortexes and boundary layer. Their complexity increases at higher Mach numbers, due to the appearance of shock waves and shock-induced separations. A simplified approach for analyzing turbomachinery flows is described in [START_REF] Wu | NACA-TN-2604: A General Theory of Three Dimensional Flow in Subsonic and Supersonic Turbomachines of Axial-, Radial-and Mixed-Flow Types[END_REF]; it relies on the decomposition of the three dimensional S1 is the blade-to-blade plane; considering a cylindrical coordinate system with radial distance r, azimuthal angle θ and axial direction x, this surface has coordinates (x, θ) and it is parametric in r. S2

is the meridional plane, with coordinates (x, r) and parametric in θ in a cylindrical coordinate system.

This surface is the place where the radial equilibrium [START_REF] Lakshminarayana | Fluid Dynamics and Heat Transfer of Turbomachinery[END_REF] is studied and its flow can be simulated by means of a throughflow models; at a first approximation it is parallel to the machine axis of rotation.

Finally, S3 is called secondary plane and it is orthogonal to the other two.

RDO OF AN ORC TURBINE BLADE

The framework just presented has been widely used for understanding and analyzing loss mechanisms in turbomachinery. The classical approach consists in separating the vortex structures and to define their behavior independently from each other. Therefore, the historical breakdown of losses into "profile losses", "end-wall losses" and "leakage losses" is still quite common, even if nowadays it has been recognized that each loss mechanism is seldom really independent [START_REF] Denton | The 1993 IGTI Scholar Lecture: Loss Mechanisms in Turbomachines[END_REF].

Profile losses are considered to be generated in the blade-to-blade plane [START_REF] Schobeiri | Turbomachinery Flow Physics and Dynamic Performance[END_REF]. These are the losses produced by the energy degradation due to the boundary layer developing on the blade surfaces as an effect of fluid viscosity; generally, with equal number of blades, the greater the flow deflections, the greater the losses. At the same time, with the same deflection, the greater the number of blades and the greater the losses [START_REF] Carter | Three-dimensional-flow theories for axial compressors and turbines[END_REF]. The extra losses generated at the blade trailing edge are usually considered as a part of profile losses. End-wall losses are also called "secondary" since they are generated by phenomena like the passage vortex, the horseshoe vortex, the corner vortex, the shed vortex and the scraping vortex, which are mainly led-back to the flow on the secondary surface S3. However, end-wall losses are due to a combination of many factors [START_REF] Yamamoto | Production and Development of Secondary Flows and Losses in Two Types of Straight Turbine Cascades: Part 1-A Stator Case[END_REF]. Tip leakage losses are due to the presence of a clearance between fixed and rotating parts of the turbomachine; the resulting tip leakage vortex interacts with end-wall loss and the energy degradation can be particularly strong in case of unshrouded blades. The relative magnitude of these three categories of losses depends several details about the design of the turbomachine. However, generally for big-sized turbomachinery these three can be considered comparable, each accounting for about 1/3 of the total loss [START_REF] Denton | The 1993 IGTI Scholar Lecture: Loss Mechanisms in Turbomachines[END_REF]. An extensive discussion about the flow in turbomachinery and the mechanisms generating losses is out of the scope of the present work; the interested reader can find relevant information about these topics in [START_REF] Lakshminarayana | Fluid Dynamics and Heat Transfer of Turbomachinery[END_REF][START_REF] Schobeiri | Turbomachinery Flow Physics and Dynamic Performance[END_REF][START_REF] Baskharone | Principles of Turbomachinery in Air-Breathing Engines[END_REF][START_REF] Dixon | Fluid Mechanics and Thermodynamics of Turbomachinery (Seventh Edition)[END_REF].

Performing the fluid-dynamic optimization of a cascade of a turbomachine is an articulate engineering process aiming to improve the performance of the component by reducing the losses just listed above. Generally this problem depends on a great amount of design parameters, which are:

• the characteristics of the working fluid, i.e. density, viscosity and speed of sound,

• operating conditions, i.e. the mass flow rate, the speed of the machine and the isentropic specific work,

• the geometry of the cascade, given by the diameter, the blade profile and its height.

RDO OF AN ORC TURBINE BLADE

Whether all these degrees of freedom are constrained, with the exception of the blade profile, the optimization of the cascade coincides with a shape optimization problem [137]. With regards to shape optimization of ORC turbomachinery the SMOGA approach has been widely used with the aim of reducing computational efforts: Harinck [START_REF] Harinck | Super-and transcritical fluid expansions for next-generation energy conversion systems[END_REF] performed a deterministic global optimization on a BZT nozzle employing D6 [START_REF] Colonna | Multiparameter equations of state for siloxanes: [(ch3)3-si-o1/2]2-[o-si-(ch3)2]i=1,. . . ,3, and [o-si-(ch3)2]6[END_REF] coupling a MOGA with an artificial neural network used as an inexpensive predictive method for the optimization. This approach has been later adopted for the deterministic shape optimization of a 3D ORC radial inflow turbine using Toluene as working fluid [START_REF] Harinck | Performance improvement of a radial organic rankine cycle turbine by means of automated computational fluid dynamic design[END_REF]. A comparison between Kriging and an artificial neural network as surrogates for the SMOGA is provided in [START_REF] Pasquale | Shape Optimization of an Organic Rankine Cycle Radial Turbine Nozzle[END_REF], where they are both used for the deterministic shape optimization of a 2D stage of an ORC radial inflow turbine.

To obtain more efficient and improved strategies to get more accurate surrogate predictions in the SMOGA approach, some authors have proposed adaptive surrogate training methodologies: Rodriguez-Fernandez and Persico [START_REF] Rodríguez-Fernández | Automatic design of orc turbine profiles using evolutionary algorithms[END_REF] developed an automatic design tool for ORC supersonic turbines by means of a Kriging-based SMOGA with a surrogate improvement during iterations. This approach has been successfully applied for the deterministic shape optimization of a nozzle for an ORC radial outflow turbine [START_REF] Persico | Evolutionary optimization of centrifugal nozzles for organic vapours[END_REF] and for an ORC centrifugal turbine stage [START_REF]Optimal aerodynamic design of a transonic centrifugal turbine stage for organic rankine cycle applications[END_REF].

In opposition to SMOGA, a different strategy, allowing to reduce the computational burden of a shape optimization problem for an ORC turbine, is the gradient-based optimization proposed in [START_REF] Pini | Adjoint Method for Shape Optimization in Real-Gas Flow Applications[END_REF] for the search of the nearest optimum: discrete adjoint is used to compute the sensitivities of the objective function with respect to design parameters. The drawback of this technique is that it is intrusive, as to develop the adjoint code, the CFD code must be differentiated with automatic differentiation techniques.

Concerning RDO for ORC applications, a surrogate-based approach is mandatory as a consequence of the high computational cost of the CFD simulations, especially when they employ accurate EOS to simulate real gas effects, and because of the high number of possble sources of uncertainty. Therefore, in this scenario, Congedo et al. [START_REF] Congedo | Tsi metamodels-based multi-objective robust optimization[END_REF] developed the robust optimization of an ORC cascade using the VKI LS-59 as a baseline and reconstructing response surfaces for sensitivity indexes in the design variables plan. Bufi and Cinnella [117] developed a two nested Bayesian Kriging with adaptive infill sampling and they used it to perform the RDO of an ORC turbine vane following Taguchi's criteria.

A similar approach (Kriging used both for UQ and for the SMOGA) has been proposed by Razaali

RDO OF AN ORC TURBINE BLADE

and al. [START_REF] Razaaly | Quantile-based robust optimization of a supersonic nozzle for organic rankine cycle turbines[END_REF] to perform a quantile-based RDO of an ORC turbine nozzle.

To further reduce the computational cost of the shape optimization, in the following a MF approach is adopted to carry out the RDO of a blade profile for an ORC axial turbine: this methodology, already described in Chapter 3 combines the information from a low-fidelity UQ model and an high-fidelity one to predict the statistics of the QoIs and it uses a BK-based MOGA to speed-up the convergence towards the global optimal solutions.

Problem definition

The objective of the RDO here carried out is the shape optimization of the blade profile for the first statoric cascade of an ORC turbine to be used for WHR applications at high temperature. For the present test case, the working fluid is Cyclopentane [START_REF] Gedanitz | Speed of sound measurements and a fundamental equation of state for cyclopentane[END_REF] and the cascade is supposed to work under uncertain operating conditions: the inlet total pressure P T,in , the inlet total temperature T T,in and the outlet static pressure P out are sources of aleatory uncertainty and they are modelled as Gaussian PDFs. Tab. 5.5 reports the mean value of these three PDFs and their coefficient of variation CoV , defined as the ratio between the mean value and the standard deviation. These uncertainties are also depicted in the h-s chart in Fig. 5.9: the lack of knowledge about P T,in and T T,in determines an uncertainty at the inlet of the cascade that is identified with the blue hatched area. In the same chart, the uncertainty on P out is represented by the red hatched area. On the lower left margin of the chart, it is possible to find the saturation curve, meaning that at the cascade inlet the working fluid is highly superheated. The cascade is composed by 50 blades with an axial chord of 20 mm and a stagger angle of 40 deg; the hub diameter is 372 mm. Since the ratio between the height of the blade and blade mean diameter is below the 2%, the cascade has been simulated as a 2D blade.

Concerning the blade profile, the VKI LS89 [13] has been chosen as a baseline geometry, to be this decision yields the optimization parameters, i.e. the quantities that must be optimized in order to get a performing design. Hence, the parameterization should be efficient, since it has to require a low number of design parameters, and flexible, because it has to represent different shapes [START_REF] Harinck | Performance improvement of a radial organic rankine cycle turbine by means of automated computational fluid dynamic design[END_REF];

moreover the geometry should be smooth, without useless edges. Several works exist proposing several parametrization techniques (see for instance [START_REF] Abbott | Theory of Wing Sections, Including a Summary of Airfoil Data[END_REF][START_REF]A Parametric Method of Turbine Blade Profile Design[END_REF][START_REF] Braibant | Shape optimal design using b-splines[END_REF][START_REF] Hicks | Wing design by numerical optimization[END_REF][START_REF] Sobieczky | Parametric Airfoils and Wings[END_REF][START_REF] Anderson | Geometry parameterization method for multidisciplinary applications[END_REF]). For a review the reader can also refer to [START_REF] Masters | Review of Aerofoil Parameterisation Methods for Aerodynamic Shape Optimisation[END_REF].

For the present work the free-form deformation (FFD) [255] method has been employed as a geometry parametrization strategy: this method allows to keep a low number of design parameters ensuring at the same time the successful modelling of complex geometries. As a drawback, it does not provide the possibility to control some engineering quantities, like for instance the blade thickness, that are usually treated as optimization constraints. Despite this point, it is widely adopted in CFD shape optimization [START_REF] Vitale | Multistage turbomachinery design using the discrete adjoint method within the open-source software su2[END_REF], because of its flexibility and ease of use. This parametrization leads therefore to a shape optimization problem with eight design points, which are the vertical displacements of the eight FFD control points free to move deforming the baseline profile to generate an optimal blade with improved performance. Since the optimization is performed on a 2D plane that coincides with the blade-to-blade surface, an approximation is expected in the evaluation of the performance of the system: secondary flows and end-wall losses are almost

RDO OF AN ORC TURBINE BLADE

RDO OF AN ORC TURBINE BLADE

totally neglected in the simulation of the flow flied within the cascade. Despite this limitation, this simplification avoids an increase in the number of design parameters and it can be reliably used to understand many of the phenomena involved in turbomachinery flow [START_REF] Denton | The exploitation of three-dimensional flow in turbomachinery design[END_REF].

To study the domain, three unstructured computational grids have been generated with the open source code GMSH [START_REF] Geuzaine | Gmsh: A 3-d finite element mesh generator with built-in preand post-processing facilities[END_REF]: they have been built with the same strategy, which is depicted in Fig. 5.11, and the only difference among them is the refinement. Thus, the coarse, the medium and the fine mesh present an inlet surface (identified as A in Fig. 5.11) an outlet surface (indicated as B) and two periodic surfaces (called C and D); in the middle of the region bounded by these four surfaces, it is possible to find the blade surface, which in Fig. 5.11 is referred to as E.

In each mesh, to resolve viscous layers close to the blade walls, a boundary layer infill has been adopted; a close-up view of it can be observed in Fig. 5.12, which is taken from the medium grid.

Moreover, as the trailing edge effects are important in the account of total losses, this region has been particularly refined in each computational grid adopted. In Fig. 5.13 a detail of the refinement in the trailing edge area is shown for the fine mesh. The quality of the meshes has been evaluated by means of the GMSH gamma quality parameter, defined as the ratio of the element inscribed radius and the element circumscribed radius [259]: all meshes show a similar distribution of this quality estimator. As an example, in Fig. 5.14 the element quality distribution for the medium mesh is depicted; moreover, in Tab.5.6 the size and the mean value of the element quality distribution are provided for every grid.

The maximum Y + value is always lower then 1.2, as further discussed at Par. 5.2.4.

∂ω ∂t + ∇ • F c -∇ • (µ vk F vk ) = 0 (5.4)
where:

• ω is the state vector containing the conservative variables:

ω = ⎡ ⎢ ⎣ ρ ρu ρe ⎤ ⎥ ⎦
(5.5)

• F c is the convective (inviscid) fluxes vector:

F c = ⎡ ⎢ ⎣ ρu ρu ⊗ u + IP ρuh ⎤ ⎥ ⎦
(5.6)

• F v1 and F v2 are the viscous fluxes vectors:

F v1 = ⎡ ⎢ ⎣ 0 τ τ • u ⎤ ⎥ ⎦
(5.7)

F v2 = ⎡ ⎢ ⎣ 0 0 c p ∇T ⎤ ⎥ ⎦ (5.8)
• ρ is the density.

• u = u 1 , u 2 ∈ R 2 is the flow speed in the 2D Cartesian system of reference.

• e is the total energy (specific to mass).

• P is the pressure, depending on ρ and e by means of an equation of state P (ρ, e).

• h is the total enthalpy (specific to mass), that can be calculated from an equation of state h(ρ, e).

• T is the temperature, that can be obtained with an equation of state T (ρ, e).

RDO OF AN ORC TURBINE BLADE

In the present work, the focus is put on entropy, because it is a convenient quantity in the evaluation of the cascade performance. In fact, from thermodynamics, any irreversibility creates inevitable a welldefined amount of entropy, reducing the process efficiency. Since the flow in general turbomachinery is closely adiabatic, entropy cannot decrease, but it can only increment whether some dissipation is present. Moreover, unlike other quantities like total enthalpy, total pressure or kinetic energy, entropy does not change if the reference frame is stationary or rotating. Because of these reasons, in [START_REF] Denton | The 1993 IGTI Scholar Lecture: Loss Mechanisms in Turbomachines[END_REF] the author concludes that "the only rational measure of loss in an adiabatic machine is entropy creation" and he suggests to consider entropy like "smoke" that is created within the flow whenever something deleterious to efficiency is taking place, like for instance in blade boundary layers or in shock waves, and whose concentration at the exit of the machine is proportional to the loss encountered inside it. Therefore, here the difference between the average entropy at the outlet of the domain and the entropy average value at the inlet is considered as the QoI J for the UQ carried out within the RDO.

The formal definition of J is provided in eq. 5.10. J = s out -s in (5.10) Corncerning the RDO, the objective is to identify the best values of vertical displacements for the FFD control points parametrizing the blade profile which minimize both the mean value and the variance of the QoI J. As the RDO exploits a MOGA optimizer, the maximum deformation allowed for the FFD control points must be specified: in the present optimization this has been set to 20% in the lattice reference system.

Validation of the CFD model

Before performing the RDO, a two step validation of the CFD model has been carried out. First, the experimental data from the reliable turbomachinery test cases [13] for dry air at ideal-gas conditions have been compared with the results from a simulation with the SU2 code: the selected test case is the MUR-43 subsonic flow scenario, where the isentropic outlet Mach number is 0.84, which has been modelled with the same modelling strategy presented in Section 5.2.1. The adopted mesh is depicted in Fig. 5.15.

The results of this first comparison are presented in Fig. 5.16: the static pressure values along the weakly supersonic at the pressure side; when they meet downstream of the trailing edge, they must turn to a common flow direction, generating a classical "fish-tail shock system" which propagates in the flow with large effect at the resulting flow field. The shock generated on the pressure side reaches the suction side of the adjacent blade, and it is reflected, increasing its strength merging with the shock wave generated on the suction side side, that interacts also with the blade wake coming from the other blades.

As a consequence of this extremely dissipative phenomenon, a value of 6.582 J kg -1 K -1 has been calculated for the QoI J, which is the entropy increase ∆s between the inlet at the outlet sections of the domain; this leads to a total pressure drop ∆P t = 173029 P a, corresponding to a pressure loss factor γ P = ∆P t /(P t,in -P out ) = 9.38%. All these performance estimation are reported in Tab 5.8. The present analysis has been performed on the coarse, the medium and the fine mesh, to check the grid convergence. All simulations agree with the analysis just presented. Fig. 5.21 shows that the estimation of the QoI J on the fine grid is not completely converged: this is probably due to the high sensibility of entropy with regards to variations in pressure and temperature. Hence, a further check of grid convergence is performed on the average outlet Mach number M out , which is depicted in Fig. 5.22. Even if a more refined mesh would have been required to obtain a more grid-independent solution, the medium mesh has been chosen for the RDO to limit the computational cost.

RDO OF AN ORC TURBINE BLADE

Moreover, a common approach has been employed for the wall treatment in the k-omega turbulence

Robust design optimization of the cascade

To propagate the uncertainty through the model both BK and the MoM are used. The DOE used by the first one is composed by 32 points, randomly sampled with a LHS strategy from the PDFs defined in Tab. 5.5. Therefore, the computational cost of a UQ BK is 32 CFD simulations.

Even if BK needs a larger DOE, it has been here preferred to GEK in the present case, since the adjoint formulation implemented in SU2 does not allow to compute the gradient of the QoI J with respect to the uncertain variables reported in Tab. 5.5: the computation of this quantity with a finite-difference approach on a DOE composed by 16 samples (half of a BK DOE) would have been required the computation of 96 further CFD simulations, which is extremely expensive in terms of computational resources and time. On the other hand, this extra-cost has been considered acceptable for the MoM; in fact, if the gradient of the QoI was available from the adjoint code, only one CFD simulation and one adjoint simulation would had been enough to compute the statistic moments of the PDF of the QoI J. However, since in the present case the adjoint strategy cannot be used to this end, the finite-difference approach can be adopted, with a cost of 7 CFD similation. The optimal solution found has been therefore compared with the baseline geometry: Fig. 5.31

shows both the profiles. It appears that the RDO has led to a blade that is more front-loaded than the baseline, with an higher load in the initial part of the profile and a lower one at the end. The blade is severely transformed by the optimizer, while retaining a regular and smooth shape. However, the change in curvature on the LE and on the initial part of the suction side seems to be severe, with the risk of local separations; moreover, the high curvature at the LE on the pressure side is likely to Both the statistical moments of these two PDF have been calculated and they are reported in Tab.

5.10: this shows that the RDO effectively improves both the quantities, reducing the expected value of the QoI PDF E[J] by 24.3% and its standard deviation √ var[J] by 18.6%. Moreover, the reader can here observe that the values of E[J] and var[J] for the RDO optimal solution are in rather good agreement with the SMOGA estimations, even if the latter appear to be slightly over-optimistic in their prediction. Tab. 5.11 reports a summary of performance estimation for the deterministic simulation carried out on the optimal geometry: the entropy increase ∆s between the inlet at the outlet sections of the domain is equal to 4.618 J kg -1 K -1 , which is 29.8% lower than in the initial configuration, while total pressure drop ∆P t = 118300 P a, decreased by 31.6% than the drop calculated for the baseline geometry; finally, the pressure loss factor γ P has been reduced to 6.41%.

Summary of the chapter

In this chapter first RDO techniques employing gradient-based UQ methods, i.e. GEK or MoM, have been successfully validated on the quasi 1D test case: both approaches allow to obtain similar accuracy than the TNBK, with a dramatic reduction in the required computational cost. On the same test case, the MF RDO methodology has been assesses: the MoM is used as low fidelity UQ model, while the GEK as high fidelity one. This approach is able to provide results which are as accurate as in a BK or in a GEK with a computational cost that is close to the MoM one. After this validation, the MF-based RDO techniques is adopted for the shape optimization of an ORC turbine 2D stator blade profile: the MoM is chosen as LF UQ model while BK is the HF one. The BK has been preferred to the GEK because the code employed doe not allows the computation of the gradients of the QoI with respect to uncertain variable using the adjoint method; hence, such a calculation would have been too expensive with the finite-difference approach. The whole RDO process requires 1271 CFD simulations, solved in 67 hours, and it effectively improves both the quantities, reducing the expected value of the QoI PDF E[J] by 24.3% and its standard deviation √ var[J] by 18.6%.

Chapter 6 Conclusion

In this thesis some robust design optimization (RDO) methods have been developed and employed for organic Rankine cycles (ORC) applications. In fact, despite the high level of confidence and know-how reached about this technology, ORC projects are still dominated by a number of uncertain variables affecting the design, the commissioning, the operation and the decommissioning of the plant. It is therefore mandatory to identify innovative techniques capable to deal with the uncertainty influencing these systems: RDO can be an answer to this issue.

In the first part of this work, an overview about uncertainty quantification, optimization and RDO methods is provided: RDO approaches studied in the literature most often combine evolutionary optimization algorithms such as genetic algorithms, and techniques for uncertainty quantification (UQ). Unfortunately, UQ for complex systems described by expensive numerical models, which is the case in Computational Fluid Dynamics (CFD), comes up against the high computational cost of the model for configurations of practical interest. With the aim of making progress toward the application of RDO to configurations of industrial interest, this research has investigated RDO techniques wellsuited for costly applications. The main requirements are: 1) greatest possible parsimony in terms of costly CFD simulations (for a given accuracy level) and 2) non intrusiveness (as in industry the use of black-box codes can be usual).

In this respect, in Chapter 3 a recently developed non-intrusive RDO technique has been introduced, based on the coupling of two nested Bayesian Kriging surrogates (TNBK): the first one is used to compute the required statistics of the objective functions in the uncertain parameter space, while the second one is used to model the response of these statistics to the design variables. An expected improvement criterion is used to update the second BK surrogate during convergence towards the optimum. The approach was shown to require O(10 × n × d) function evaluations (with n the number of uncertain parameters and d the number of design variables) at the first generation of the multiobjective genetic algorithm (MOGA), in order to generate the BK surrogates, which can be efficiently parallelized on a high-performance computer, while additional O(10 × n) evaluations are required for each update of the external BK surrogate. Despite considerable efficiency gains with respect to other surrogate models like Polynomial Chaos Expansion or Probabilistic Collocation Method, the TNBK approach remains too expensive for industrial problems, especially if a large number of uncertain parameters has to be taken into account and if massively parallel computers are not available. Thus, TNBK RDO is limited to uncertain spaces of low to moderate dimensionality (up to about 8 uncertain parameters). With the aim of enabling the application of RDO to higher dimensional parameter spaces and/or to more realistic and complex fluid flow problems, a further reduction of the computational cost is mandatory. In order to achieve this goal, the TNBK approach can be modified by replacing the inner BK surrogate used for the UQ step with gradient based approaches. More specifically, a Gradient-Enhanced Kriging (GEK) surrogate or a first-order Method of Moments (MoM) are considered. An external BK is coupled to the MOGA in all cases, to speed up the exploration of the design space. The UQ GEK approach provides as accurate results as BK with fewer function calls thanks to the use of gradient information, but the overall number of CFD runs necessary for a full RDO loop is still high. UQ MoM is inexpensive, since it requires only a function call and the calculation of the derivatives, but it provides a rather rough approximation of the moments of the probability distribution function for the quantities of interest. In both GEK and MoM, the availability of accurate and cheap estimates of the gradients is crucial. These can be computed efficiently using continuous or discrete adjoint solvers, when available, or by means of finite difference approximations. Finally, to combine the advantages of the both gradient based methods, a multi-fidelity RDO strategy is proposed, where the GEK and MoM UQ solvers are fused together to build a multi-fidelity surrogate with adaptive infill enrichment for the SMOGA optimizer: this method can ensure a good tradeoff between cost and accuracy, thus representing an efficient approach for complex RDO problems.

In Chapter 4, the TNBK approach has been assessed by means of a direct comparison with the "brute-force" MC approach for the RDO under aleatory and epistemic uncertainty of a (simplified) ORC for WHR: the technique has shown all its potential speeding the optimization process up by a factor 10000. Afterwards, the TNBK method is successfully applied to an ORC for geothermal applications, affected both by epistemic and aleatory uncertainty. The first comes from the lack of information about the geothermal source and it is considered as the outcome of the interaction of the mass flow rate of the geothermal brine m ̇geo , the temperature of the geothermal brine T geo , and the heat capacity of the geothermal brine c geo ; the aleatory uncertainty is due to the yearly variation of the cooling water temperature at the condenser T cd . Both these uncertainties can lead to an economic risk involving the whole geothermal project. A preliminary global sensitivity analysis shows that the the most influential parameters are T geo , m ̇geo and T cd , while c geo can be neglected. These results are used to reduce the number of optimization parameters that are considered in the RDO of the ORC. Afterwards, a TNBK RDO is carried out to generate a Pareto front of alternative optimal solutions, and a single design, offering the best trade-off among all optimization targets is selected. This solution is compared with the optimal solutions from the deterministic optimization, which are always outperformed. In particular, the best RDO solution wins against the best deterministic optimum, increasing the mean value of the PDF of the QoI (the power production) by 1.5%, while its standard deviation is reduced by 8.5% and the surface of the heat exchangers by 34%.

In Chapter 5, more advanced RDO strategies using gradient-based UQ techniques are first evaluated for an inexpensive test problem, representative of a supersonic divergent nozzle, for which the results can be compared with well-converged Monte Carlo sampling. The results show that GEK allows computational gains of a factor 2 or more with respect to a Bayesian Kriging surrogate not using gradient information, when the gradients are efficiently computed using an adjoint method. For the present test problem, both a discrete and a continuous adjoint method where used for building GEK surrogates. The first-order moment method, based on either discrete or continuous adjoint gradients, is less accurate than GEK but it still provides satisfactory estimates of the QoI interest for the present shocked flow problem, due to the relatively small variation ranges of the uncertain parameters. The UQ methods are then combined with a genetic algorithm for solving the RDO problem. The computations are sped up by constructing a Bayesian Kriging surrogate model of the design space. The surrogate is enriched during GA iterations by means of a multi-objective expected improvement (MOEI) infill criterion. For the test problem at hand, RDO results are found to be similar for BK, GEK and MoM methods, the latter being much cheaper in terms of CPU time but slightly less accurate than the former ones. In order to benefit of the computational efficiency of the MoM and the accuracy of the GEK UQ solvers at the same time, a multi-fidelity surrogate model is build by fusing together the low-fidelity MoM and the high-fidelity GEK. A MOEI infill is used again to enrich the surrogate during convergence, with preference for low-fidelity infills. The multi-fidelity approach successfully identifies the RDO optimum, while dividing by a factor 3 ÷ 4 the computational cost with respect to GEK. Such an approach is then identified as a promising candidate for more complex RDO problems using CFD models. The multi-fidelity RDO approach is finally applied to the shape optimization of an ORC turbine 2D stator blade profile, parametrized with the free-form deformation technique. The average entropy increase across the cascade is considered as the QoI for the RDO. The MoM is chosen as LF UQ model while BK is the HF one. The BK has been preferred to the GEK because the CFD solver in use does not allow the computation of the gradients to the QoI with respect to uncertain variable with the adjoint method; for the MoM, finite difference approximations of the gradients have been used. The whole RDO process requires 1271 CFD simulations, solved in 67 hours, and it effectively improves both the quantities, reducing the expected value of the QoI PDF E[J] by 24.3% and its standard deviation √ var[J] by 18.6%. This result is possible since the geometry found with the RDO presents a dramatic reduction of the main shock strength and lower Mach numbers. Therefore, even if the shock pattern is similar to the one observed in the baseline geometry, the severe pressure gradient observed in the previous analysis is consistently reduced and the interactions between shocks and the blade wake are weakened. As a result, the overall entropy generation is significantly reduced, if compared with the baseline geometry.

As future perspectives, the development of a gradient-adjoint multi-fidelity RDO approach could be a mandatory step to further reduce the computational cost of the RDO: a multi-fidelity strategy employing the MoM as low-fidelity UQ code and the GEK as high-fidelity one can potentially further reduce the RDO computational cost by more than 50%, if the gradients of the QoI are calculated with the adjoint method. This would probably allow the RDO of a 3D turbine blade or of a complex multi-cascade system, i.e. a stator-rotor configuration with a mixing plane interface.

Another improvement could be the replacement of the external BK coupled with the MOGA with the GEK technique. However, this strategy would require the gradient of the statistical moments of the PDF of the QoI with respect to the design variables; such a result can be obtained by means of the derivation of the underlying UQ method used to compute the statistics, MoM, BK or GEK. This operation is possible in theory, but not straightforward using adjoint; on the other hand, finite differ-ence approximations are easily applicable, but at the price of a considerable computational expense for highly dimensional design spaces.

Moreover, with regard to MF-based RDO for CFD applications, the combination of Euler and RANS equations could be thought as a way to lower the computational burden; another option could be to combine RANS solutions on meshes with different resolution (multi-level approach).

Finally, another option to obtain this result could be the implementation of a three-levels multifidelity based RDO, where BK is used as high-fidelity UQ code only once or twice in all the RDO, while the lighter low-fidelity and the medium-fidelity codes are calculated more often; however, in this case more complex adaptive infill strategies should be considered.

Résumé étendu

Le cycle organique de Rankine (ORC) est une technologie viable pour l'exploitation des énergies renouvelables comme l'énergie solaire concentrée, la géothermie, la biomasse ou la récupération de la chaleur résiduelle. Dans ces applications, il surpasse généralement les cycles à vapeur classiques par sa simplicité, ses coûts d'exploitation plus faibles et son efficacité thermodynamique plus élevée [START_REF] Macchi | Organic Rankine Cycle (ORC) Power Systems[END_REF].

Les ORC sont des cycles de Rankine utilisant comme fluide de travail des composés organiques complexes (hydrocarbures, huiles de silicone ou réfrigérants), au lieu de la vapeur : ce sont des cycles fermés comprenant au moins une pompe, qui comprime le fluide de travail, un groupe d'échangeurs de chaleur chauds, généralement composé d'un ou plusieurs préchauffeurs, un évaporateur et parfois aussi un surchauffeur, où le fluide de travail est chauffé par une source de chaleur externe pour devenir de la vapeur ; ensuite, il y a une turbine, qui convertit la puissance thermodynamique du fluide en puissance mécanique, et un condenseur, où la chaleur résiduelle est libérée dans l'environnement permettant au fluide de revenir à l'état liquide. La puissance mécanique à l'arbre de la turbine est convertie en électricité par un générateur.

Au cours des dix dernières années, les cycles organiques de Rankine (ORC) sont devenus une solution technique compétitive pour l'exploitation de sources de chaleur à basse et moyenne température de capacité limitée [START_REF] Macchi | Organic Rankine Cycle (ORC) Power Systems[END_REF], entraînant une croissance extraordinaire de leur marché, en particulier pour les applications géothermiques, de biomasse et de récupération de la chaleur fatale (WHR). Les données de la Fig. La biomasse représente une autre application assez répandue de l'ORC : depuis le début de l'année 2000, plusieurs centrales électriques ORC à haute température d'une taille d'environ 1 M W e ont été installées en Europe pour utiliser différents types de biomasse solide [START_REF] Colonna | Organic Rankine Cycle Power Systems: From the Concept to Current Technology, Applications, and an Outlook to the Future[END_REF]. Très souvent, ces centrales sont cogénératives, fournissant à la fois de l'électricité et de la chaleur, qui est généralement R ÉSUM É ÉTENDU utilisée en dessous de 100 degrés Celsius pour le chauffage urbain ou à des fins de traitement (par exemple, le séchage du bois). [START_REF] Guercio | 15 -biomass-fired organic rankine cycle combined heat and power systems[END_REF]. La majorité de ces systèmes adoptent des siloxanes comme l'hexaméthyldisiloxane ou l'octaméthyltrisiloxane comme fluide de travail dans un cycle régénératif surchauffé, où le fluide de travail est réchauffé par de l'huile thermique. Habituellement, le rendement électrique net nominal se situe entre 15 et 20 %, tandis que le rendement énergétique total peut être plus élevé, en fonction des besoins de la cogénération [START_REF] Bini | Operational results of the first biomass chp plant in italy based on an organic rankine cycle turbogenerator and overview of a number of plants in operation in europe since 1998, 2nd, world biomass conference; biomass for energy industry and climate protection[END_REF].

Une autre application présentant un potentiel intéressant pour toutes les tailles d'unités est la récupération de la chaleur fatale (WHR) ; en fait, les possibilités de récupération de la chaleur dans l'industrie manufacturière et de transformation sont nombreuses [START_REF] Campana | Orc waste heat recovery in european energy intensive industries: Energy and ghg savings[END_REF]. Selon les données tracées dans la Fig. 6.1 et la Fig. 6.2, cela devrait être considéré comme un domaine émergent pour les ORC : sur ce marché, plusieurs solutions sont disponibles avec des solutions de récupération de moyenne et grande taille à partir de turbines à gaz, de moteurs à combustion interne ou de processus industriels (par exemple, cimenterie, aciérie et verrerie). Habituellement, dans cette application, les siloxanes, les alcanes, les cycloalcanes et certains réfrigérants HFO peuvent être adoptés dans des cycles thermodynamiques qui sont généralement surchauffés. Aujourd'hui, cette application est celle qui présente la plus grande flexibilité dans la configuration du cycle ; une installation de 3,2 MW pour la récupération de la chaleur résiduelle d'une aciérie est présentée dans la figure 6.4.

Enfin, les applications solaires sont négligeables si on les compare à l'énergie géothermique, à la biomasse et à la WHR : cela est probablement dû au fait que les ORC pourraient être couplées à des centrales solaires à concentration qui, pour le moment, sont plus coûteuses que les panneaux photovoltaïques et les systèmes de batteries [1]. Toutefois, en raison de leur grande fiabilité, disponibilité et performance, les ORC ont été identifiées comme la technologie de conversion optimale dans ce contexte [START_REF] Price | Modular trough power plant cycle and systems analysis[END_REF][START_REF] Casartelli | Power block off-design control strategies for indirect solar orc cycles[END_REF].

L'une des plus grandes forces de la technologie ORC, qui doit être considérée comme l'une des principales raisons de leur succès spectaculaire, réside probablement dans leur extrême flexibilité, qui leur permet de s'adapter facilement aux besoins des clients ; en fait, elles sont tellement personnalisables que parfois certaines solutions sur mesure proposées par les fournisseurs d'ORC peuvent être considérées comme de véritables défis d'ingénierie en termes de conception, de construction, de mise en service et d'exploitation. En passant d'un projet à l'autre, un certain nombre d'éléments sont susceptibles de changer ; certains sont énumérés ci-dessous. 2. Le type de source de chaleur (géothermique ou biomasse ou WHR ou solaire).

3. Les conditions de fonctionnement de la ORC, qui peuvent nécessiter différentes solutions de contrôle de la centrale.

4. Les composants de l'installation (turbine, pompe, échangeurs de chaleur, vannes), qui d'une application à l'autre peuvent généralement changer de dimensions et de configurations. En fait, même si l'ensemble des fabricants de ORC s'efforce de normaliser les solutions proposées aux clients, en définissant des séries de produits standard, ils continuent à fournir principalement des installations sur mesure, car il n'est pas possible de normaliser entièrement les ORC. Ceci est confirmé par l'inspection du matériel commercial (brochures, présentations ou livres blancs) fourni par les plus grands fournisseurs de ORC (par exemple, voir [24, [START_REF] Enertime | Corporate brochure[END_REF]26,27,[START_REF] Exergy | Efficient and tailor-made organic rankine cycle geothermal binary power plants[END_REF]).

Paradoxalement, le grand nombre de degrés de liberté laissés libres lors de la conception d'une usine ORC constitue en même temps une sérieuse limitation. En fait, pour mieux comprendre la dynamique du marché des ORC, il est obligatoire d'identifier les actionnaires qui achètent habituellement des ORC pour les exploiter : la grande majorité des centrales ORC qui ont été construites jusqu'à aujourd'hui et qui sont toujours en exploitation sont la propriété d'entités privées, qui considèrent généralement les ORC comme un investissement d'une partie de leurs ressources à analyser en termes de risque, de temps de retour sur investissement et de rentabilité ; en général, les revenus proviennent presque exclusivement de la quantité d'énergie produite, qui doit être évaluée le plus précisément possible.

Ainsi, malgré le niveau élevé de confiance et de savoir-faire atteint de nos jours, un projet concernant une ORC est généralement encore dominé par une grande série de variables inconnues qui sont liées à la conception, à la mise en service, à l'exploitation et au déclassement de la centrale et qui pourraient réduire l'enthousiasme des actionnaires pour les ORC. Il est donc obligatoire d'identifier des techniques innovantes capables de faire face à l'incertitude qui affecte les ORC : en ce qui concerne la conception des composants et des systèmes complets de ORC, l'optimisation robuste de la conception (RDO) peut être une réponse à cette problématique. L'incertitude est une déficience potentielle dans toute phase ou activité du processus de modélisation qui est due au manque de connaissances, révélant une nature stochastique claire [START_REF] Walters | Uncertainty Analysis for Fluid Mechanics with Applications[END_REF] En ce qui concerne l'incertitude aléatoire, il s'agit du caractère aléatoire des paramètres du système (variables géométriques, propriétés des matériaux, tolérances de fabrication) et du caractère aléatoire de l'environnement (conditions initiales et limites et paramètres environnementaux). L'incertitude épistémique est largement liée aux incertitudes dans le choix du modèle, qui découlent généralement R ÉSUM É ÉTENDU d'un manque de connaissances sur la physique sous-jacente ou de l'impossibilité d'utiliser un modèle complet et précis pour les simuler [START_REF] Fusi | Robust shape optimization of fixed and morphing rotorcraft airfoils[END_REF].

L'ingénierie probabiliste vise à prendre en compte les incertitudes apparaissant dans la modélisation des systèmes physiques et à étudier l'impact de ces incertitudes sur la réponse du système [START_REF] Sudret | Uncertainty propagation and sensitivity analysis in mechanical models -Contributions to structural reliability and stochastic spectral methods[END_REF]. À cette fin, un cadre bien établi et universellement accepté pour la quantification de l'incertitude a été élaboré : il s'agit d'un processus circulaire à plusieurs étapes dont les lignes directrices sont présentées dans la figure 6.6 et brièvement décrites ci-après. L'étape B correspond à la quantification des sources d'incertitude, en identifiant les paramètres d'entrée qui ne peuvent être considérés comme connus parce qu'ils sont affectés par l'incertitude ; ils sont donc modélisés dans un contexte probabiliste, par la définition de leurs PDF. En général, les entrées stochastiques peuvent être associées aux conditions de fonctionnement [START_REF] Cinnella | Robust optimization of dense gas flows under uncertain operating conditions[END_REF][START_REF] Congedo | Shape optimization of an airfoil in a BZT flow with multiple-source uncertainties[END_REF][START_REF] Cinnella | Robust optimization using nested Kriging surrogates: application to supersonic ORC nozzle guide vanes[END_REF], à la géométrie [START_REF] Chen | A new level-set based approach to shape and topology optimization under geometric uncertainty[END_REF][START_REF] Lazarov | Topology optimization with geometric uncertainties by perturbation techniques[END_REF][START_REF] Zhang | Robust shape and topology optimization considering geometric uncertainties with stochastic level set perturbation: Robust Topology Optimization Considering Geometric Uncertainties[END_REF][START_REF] Montomoli | Geometrical uncertainty in turbomachinery: Tip gap and fillet radius[END_REF] ainsi qu'aux paramètres empiriques impliqués dans les modèles physiques, R ÉSUM É ÉTENDU e. g. dans les modèles de turbulence [START_REF] Harinck | Optimization of a 3d radial turbine by means of an improved genetic algorithm[END_REF] ou l'équation d'état (EoS) [START_REF] Cinnella | Quantification of Thermodynamic Uncertainties in Real Gas Flows[END_REF]. Le produit final de cette étape est un vecteur de variables aléatoires ayant une fonction de distribution de probabilité (PDF)

bien établie.

L'étape C est la plus exigeante sur le plan des calculs, car elle consiste à propager dans le modèle toutes les incertitudes des entrées, en caractérisant la réponse aléatoire de manière appropriée par rapport aux critères d'évaluation définis à l'étape A : généralement, l'objectif de cette phase est le calcul de la fonction de distribution de probabilité de la QoI ou de ses moments statistiques. L'échantillonnage Monte Carlo (MC) est probablement la méthode la plus intuitive pour mener à bien cette tâche, mais elle nécessite également beaucoup de ressources de calcul, car elle donne en sortie une distribution de probabilité discrète (histogramme) entière de la QoI. Afin de réduire l'effort de calcul, d'autres méthodes sont couramment utilisées car elles sont plus rapides que la MC ; certaines d'entre elles sont décrites dans la section suivante.

Une fois que l'incertitude s'est propagée dans le modèle, il est possible d'effectuer une analyse de sensibilité pour obtenir des informations sur l'impact respectif des variables d'entrée aléatoires sur les QoI ; en fait, plus le système est complexe, plus le nombre de paramètres d'entrée augmente, et le concepteur sera probablement intéressé à comprendre les effets de chacun d'entre eux sur la réponse du système. Ainsi, une fois que les PDF (ou simplement les statistiques) des QoI ont été calculées, elles peuvent être utilisées pour caractériser la sortie, afin d'améliorer la connaissance du problème considéré.

En général, il existe deux approches principales pour l'analyse de sensibilité, à savoir l'approche locale et l'approche globale. La première est probablement la méthode la plus intuitive car elle analyse l'impact sur le QoI d'une petite perturbation des variables incertaines autour de leurs valeurs nominales : elle consiste à calculer les dérivées partielles du QoI par rapport aux variables incertaines. Habituellement, ces quantités sont dérivées au moyen de la méthode "One-At-Time" (OAT) (le lecteur intéressé peut trouver un examen approfondi de cette technique de base dans [START_REF] Saltelli | How to avoid a perfunctory sensitivity analysis[END_REF]). Cependant, l'approche locale présente plusieurs limites, telles que les hypothèses de linéarité et de normalité et elle souffre en cas de dépendance parmi les variables incertaines [START_REF] Iooss | A review on global sensitivity analysis methods[END_REF]. En revanche, l'approche de sensibilité globale est basée sur l'analyse de la décomposition de la variance (ANOVA) [START_REF] Liu | Estimating mean dimensionality of analysis of variance decompositions[END_REF], puis elle estime les indices de Sobol [START_REF] Ballester-Ripoll | Sobol tensor trains for global sensitivity analysis[END_REF] associés aux entrées incertaines. Dans ce qui suit, cette dernière approche est la seule utilisée pour effectuer une analyse de sensibilité. Le lecteur intéressé est invité à se reporter à l'annexe B pour plus de détails à ce sujet.

Dans le Ch.2, diverses techniques de propagation de l'incertitude sont brièvement passées en revue.

Pour le RDO effectué dans le cadre du présent travail, certains modèles substitut ont été utilisés pour la quantification de l'incertitude : il s'agit du krigeage bayésien (BK), du krigeage à gradient amélioré et de la méthode des moments du premier ordre.

Concevoir un composant signifie effectuer une optimisation sous contrainte de certaines (une ou plusieurs) fonctions objectives, qui sont une traduction mathématique de tous les buts et exigences souhaités ; une représentation schématique de ce processus itératif est fournie dans la Fig. 6.7. Dans un tel contexte, les méthodes non déterministes doivent être considérées comme une stratégie.

Sans vouloir présenter un examen exhaustif de ces techniques d'optimisation dans des conditions d'incertitude, on peut les classer en deux grandes approches principales, à savoir l'optimisation de la conception basée sur la fiabilité (RBDO) et l'optimisation de la conception robuste (RDO) [START_REF] Schuëller | Computational methods in optimization considering uncertainties -an overview[END_REF]. La première se concentre sur la recherche d'un optimum avec une probabilité maximale de défaillance donnée qui devrait être atteinte par la conception et elle est généralement adoptée dans les problèmes structurels [START_REF] Dersjö | Efficient design of experiments for structural optimization using significance screening[END_REF], tandis que la seconde vise à identifier une conception minimisant la sensibilité de l'objectif et des contraintes aux incertitudes. Plusieurs formulations existent pour les deux stratégies (voir par exemple [START_REF] Deng | Methods for robust and reliability-based design optimization of electromagnetic devices[END_REF][START_REF] Kang | Robust Design Optimization of Structures under Uncertainties[END_REF]). Dans les travaux présentés ci-après, seule l'approche RDO a été prise en compte.

Ces dernières années, l'optimisation de la conception robuste (RDO) [START_REF] Beyer | Robust optimization -A comprehensive survey[END_REF] a suscité un intérêt croissant pour les applications d'ingénierie, en raison de sa capacité à fournir des conceptions efficaces avec un comportement stable dans des incertitudes de nature diverse, telles que des conditions de fonctionnement fluctuant de manière aléatoire, des tolérances géométriques, et des incertitudes de modèle. La méthode de Taguchi [START_REF] Taguchi | System of experimental design: engineering methods to optimize quality and minimize costs[END_REF], qui repose sur l'optimisation simultanée de la moyenne et de la variance des fonctions de coût stochastiques, est de loin la méthode RDO la plus populaire, bien que les approches permettant de tenir compte des événements rares, telles que les méthodes low-quantile [START_REF] Maliki | Quantile-based optimization under uncertainties using adaptive Kriging surrogate models[END_REF][START_REF] Razaaly | Quantile-based robust optimization of a supersonic nozzle for organic rankine cycle turbines[END_REF] ou la méthode "horsetail matching" [START_REF] Cook | Horsetail matching: a flexible approach to optimization under uncertainty[END_REF], ont récemment reçu un intérêt important.

Le principal élément de la RDO est une méthode de quantification de l'incertitude, permettant de caractériser les fonctions de distribution de probabilité ou, au moins, les statistiques d'ordre inférieur des fonctions de coût pour chaque conception proposée, afin de sélectionner celles qui garantissent la meilleure performance moyenne possible tout en évitant les écarts critiques lorsque les conditions nominales de conception ne sont pas respectées. Selon la méthode RDO utilisée, un problème de conception déterministe à objectif unique est généralement converti en un problème de conception multi-objectif (front de Pareto), dans le but d'optimiser la performance moyenne tout en évitant la perte de performance critique dans des conditions hors conception. Pour cette raison, RDO combine souvent un solveur UQ avec des algorithmes évolutionnaires (typiquement, des algorithmes génétiques multi-objectifs (MOGA) [START_REF] Deb | Optimization for Engineering Design -Algorithms and Examples[END_REF][START_REF] Kochenderfer | Algorithms for Optimization[END_REF]), qui sont naturellement adaptés pour fournir un ensemble complet R ÉSUM É ÉTENDU de solutions de compromis parmi les multiples objectifs. D'autre part, les optimiseurs évolutifs sont généralement très exigeants en termes d'évaluation des fonctions de coût, ce qui peut nécessiter en fin de compte un effort de calcul prohibitif pour les problèmes décrits par des modèles informatiques coûteux, tels que ceux rencontrés en dynamique des fluides numérique (CFD), malgré l'utilisation de la parallélisation massive [START_REF] Congedo | Efficient robust optimization techniques for uncertain dense gas flows[END_REF]117,[START_REF] Cinnella | Robust optimization using nested Kriging surrogates: application to supersonic ORC nozzle guide vanes[END_REF] Pour résoudre le problème d'optimisation précédent, il faut tenir compte de certains éléments :

• évaluer J peut être très coûteux,

• les dérivés de J par rapport à x et ξ peuvent ne pas être facilement accessibles,

• J peut être une fonction non convexe.

Comme le problème de RDO défini dans l'équation 6.1 est un problème intrinsèquement multiobjectifs et considérant que le but est la poursuite d'une optima globale, un MOGA est sélectionné pour la recherche de l'espace de conception. Plus précisément, à la suite de certaines études antérieures [START_REF] Cinnella | Robust optimization of dense gas flows under uncertain operating conditions[END_REF][START_REF] Cinnella | Robust optimization using nested Kriging surrogates: application to supersonic ORC nozzle guide vanes[END_REF], l'algorithme génétique de tri non dominé (NSGA-II) de Deb et al. [START_REF] Deb | A fast and elitist multiobjective genetic algorithm: NSGA-II[END_REF] Il existe plusieurs critères de remplissage adaptatifs, comme par exemple :

• random search,

• prediction minimization,

• maximum entropy,

• lower confidence bound,

• probability of improvement,

• expected improvement,

• cyclic infill. Le pseudo-code pour le modèle MF est fourni ci-dessous dans l'algorithme 2 ; pour plus d'informations à ce sujet, le lecteur est prié de s'adresser à [START_REF] Gratiet | Multi-fidelity Gaussian process regression for computer experiments[END_REF].

Dans le présent travail, la méthode TNBK a été utilisée pour les RDO de deux ORC. La première application est décrite à l'aide d'un modèle très simple (et peu coûteux sur le plan du calcul), qui permet des comparaisons directes avec une méthode de "force brute" reposant sur un échantillonnage Monte Carlo pur. La seconde est un ORC présentant un intérêt pratique pour les applications géothermiques.

Il est décrit par un modèle plus complet (et coûteux). Les résultats présentés ci-après ont été publiés dans [START_REF] Serafino | Assessment of an Innovative Technique for the Robust Optimization of Organic Rankine Cycles[END_REF] et [START_REF] Serafino | Robust optimization of an organic rankine cycle for geothermal application[END_REF]. Une discussion détaillée de ces deux cas est fournie dans le Ch. 4. L'approche TNBK permet de réduire considérablement les ressources de calcul nécessaires par rapport à l'approche "brute-force" MC ; cependant, cette stratégie reste trop coûteuse pour des codes de simulation extrêmement coûteux, comme par exemple les applications CFD. Dans ce cas, une première contribution peut être l'adoption de méthodes UQ utilisant les informations de gradient, c'est-à-dire GEK ou MoM au lieu de BK. En outre, des techniques RDO plus avancées devraient être adoptées. Par conséquent, une évaluation de la méthodologie RDO basée sur la MF (voir chapitre 3) est ici fournie. Les résultats ont été publiés dans [START_REF] Serafino | Multi-fidelity gradient-based strategy for robust optimization in computational fluid dynamics[END_REF].

La géométrie de la tuyère est attribuée par la distribution de la surface S(x) le long de l'axe longitudinal x. Celle-ci est choisie pour être de la forme : Malgré la présence du choc, un modèle du premier ordre a été préféré à des schémas plus précis à haute résolution, pour ne pas garantir une robustesse satisfaisante dans la solution des problèmes directs et adjoints.

S(x) = a + b tanh(cx -d) (6.4) R ÉSUM É ÉTENDU K = [0, -P, 0] T .
Le système est supposé être soumis à des incertitudes de nature diverse, notamment :

• les tolérances géométriques des éléments sur la forme de la tuyère, modélisées en traitant les paramètres de forme a, b, c, d comme des variables aléatoires normalement distribuées, avec la moyenne µ et le coefficient de variation CoV = σ/µ, avec σ comme écart type ;

• l'incertitudes sur la pression totale à l'entrée P T,in décrit comme une variable aléatoire uniformément distribuée avec des limites inférieure et supérieure imposées ;

• les incertitudes sur les propriétés du gaz, représentées ici par le rapport de chaleur spécifique γ, qui est également supposé être uniformément réparti.

Les caractéristiques des paramètres aléatoires sont énumérées dans le tableau 6.1. Dans le processus de conception inverse, les paramètres géométriques incertains sont également des variables de conception (incertaines) : pour cette raison, leur moyenne n'est pas fixe, mais varie dans des plages correspondant aux limites de l'espace de conception. Cela signifie que, même pour les conceptions correspondant aux limites supérieures/inférieures, une réalisation de la géométrie de la tuyère peut se situer en dehors des limites prescrites, en raison des tolérances géométriques.

∂ω ∂t + ∇ • F c -∇ • (µ vk F vk ) = 0 (6.7)
où:

• ω est le vecteur d'état contenant les variables conservatrices :

ω = ⎡ ⎢ ⎣ ρ ρu ρe ⎤ ⎥ ⎦ (6.8)
• F c est le vecteur des flux convectifs ( non visqueux ) :

F c = ⎡ ⎢ ⎣ ρu ρu ⊗ u + IP ρuh ⎤ ⎥ ⎦ (6.9)
• F v1 et F v2 sont les vecteurs de flux visqueux :

F v1 = ⎡ ⎢ ⎣ 0 τ τ • u ⎤ ⎥ ⎦ (6.10) F v2 = ⎡ ⎢ ⎣ 0 0 c p ∇T ⎤ ⎥ ⎦ (6.11)
• ρ est la densité.

• u = u 1 , u 2 ∈ R 2 est la vitesse d'écoulement dans le système de référence cartésien 2D.

• e est l'énergie totale (spécifique à la masse).

• P est la pression, en fonction de ρ et e au moyen d'une équation d'état P (ρ, e).

• h est l'enthalpie totale (spécifique à la masse), qui peut être calculée à partir d'une équation d'état h(ρ, e).

• T est la température, qui peut être obtenue avec une équation d'état T (ρ, e).

• τ est le tenseur des contraintes visqueuses qui peut être écrit en notation vectorielle comme :

τ = ∇u + ∇u T - 2 3 I(∇ • u) (6.12)
• c p est la capacité thermique à pression constante, en fonction de ρ et e au moyen d'une équation d'ét c p (ρ, e).

• µ v1 = µ dyn + µ turb .

• µ v2 = µ dyn P r d + µ turb P r t .

• µ dyn est la viscosité dynamique du fluide.

• µ turb est la viscosité turbulente, calculée avec un modèle de turbulence RANS.

• P r d and P r t sont le nombre de Prandtl dynamique et turbulent, qui est ici considéré comme égal à 0.9.

Les équations directrices sont résolues en utilisant le code SU2 [START_REF] Economon | SU2: An Open-Source Suite for Multiphysics Simulation and Design[END_REF] : pour les analyses présentées dans ce qui suit, le schéma central de Les résultats de cette première comparaison sont présentés dans la Fig. 6.23 : les valeurs de pression statique le long de la corde axiale de la pale calculées avec SU2 sont représentées par une Après cette première étape de la validation, la stratégie de modélisation CFD a été testée sur une étude de cas plus proche de l'objet problématique du RDO qui est ici réalisé : le flux inviscide de propane gazeux à travers la grille LS89. Pour ce cas, les résultats du RANS SU2 peuvent être comparés aux simulations d'écoulement non visqueux disponibles sur [START_REF] Harinck | The influence of molecular complexity on expanding flows of ideal and dense gases[END_REF]. Pour le calcul de SU2, la grille de calcul est le maillage de la Fig. générant un système de choc classique "fish-tail shock system" qui se propage dans l'écoulement avec un grand effet au niveau du champ d'écoulement résultant. Le choc généré du côté pression atteint le côté aspiration de l'aube adjacente, et il est réfléchi, augmentant sa force en fusionnant avec l'onde de choc générée du côté aspiration, qui interagit également avec le sillage de l'aube provenant des autres aubes.

En conséquence de ce phénomène extrêmement dissipatif, une valeur de 6.582 J kg -1 K -1 a été calculé pour le QoI J, qui est l'augmentation de l'entropie ∆s entre l'entrée et les sections de sortie du domaine ; ceci conduit à une perte de pression totale ∆P t = 173029 P a, correspondant à un facteur de perte de pression γ P = ∆P t /(P t,in -P out ) = 9, 38%. Toutes ces estimations de performance sont reportées dans la Table 6.8. R ÉSUM É ÉTENDU Figure 6.27: Distribution de la pression dans la géométrie de base Table 6.8: Paramètres de performance de la grille Performance parameter ∆s ∆P t γ P 6.582 Jkg -1 K -1 173029 P a 9.38% La présente analyse a été effectuée sur les mailles grossières, moyennes et fines, afin de vérifier la convergence des grilles. Toutes les simulations concordent avec l'analyse présentée ci-dessus. La figure 6.28 montre que l'estimation du QoI J sur la grille fine n'est pas complètement convergente : cela est probablement dû à la grande sensibilité de l'entropie par rapport aux variations de pression et de température. Par conséquent, une vérification supplémentaire de la convergence de la grille est effectuée sur le nombre moyen de Mach de sortie M out , qui est représenté sur la figure 6.29. Même si un maillage plus fin aurait été nécessaire pour obtenir une solution plus indépendante de la grille, le maillage moyen a été choisi pour le RDO afin de limiter le coût de calcul. Table 6.11: Paramètres de performance de la grille optimisée Performance parameter s P t γ P 4.618 Jkg -1 K -1 118300 P a 6.28% Tab. 6.11 présente un résumé de l'estimation des performances pour la simulation déterministe effectuée sur la géométrie optimale : l'augmentation de l'entropie ∆s entre l'entrée aux sections de sortie du domaine est égale à 4.618 J kg -1 K -1 , soit 29.8% de moins que dans la configuration initiale, tandis que la perte de pression totale ∆P t = 118300 P a, a diminué de 31.6% par rapport à la perte une stratégie multifidélité employant la MoM comme code UQ basse fidélité et la GEK comme haute fidélité peut potentiellement réduire davantage le coût de calcul de la RDO de plus de 50%, si les gradients de la QoI sont calculés avec la méthode adjointe. Cela permettrait probablement de calculer la RDO d'une aube de turbine en 3D ou d'un système complexe en plusieurs cascades, c'est-à-dire une configuration stator-rotor avec une interface de plan de mélange.

Une autre amélioration pourrait être le remplacement du BK externe couplé au MOGA par la technique GEK. Toutefois, cette stratégie nécessiterait le gradient des moments statistiques de la PDF de la QoI par rapport aux variables de conception ; un tel résultat peut être obtenu au moyen de la dérivation de la méthode UQ sous-jacente utilisée pour calculer les statistiques, MoM, BK ou GEK.

Cette opération est possible en théorie, mais elle n'est pas simple en utilisant l'adjectif "adjoint" ; en

A.1 The Gradient Evaluation Problem

In general, the gradient evaluation problem may be stated as follows: let us consider J as a vector of objective functions and constraints of dimension N j (Eq.A.1), that must comply to a set R of governing equations (Eq.A.2). J = J(ω(α), α) (A.1)

R(ω(α), α) = 0 (A.2)
where:

• α is a vector of control/design variables of dimension N α which parametrize the problem.

• ω is a vector of state variables depending on α with dimension N ω Usually, in CFD applications the governing equations R stated in Eq.A.2 correspond to Navier-Stokes equations, ω are the flow field variables and α the design parameters parametrizing the geometry. The final aim of the problem is to evaluate the gradient of the objective functions and constraints expressed in J with respect to the design variables α, as defined in Eq.A.3.

dJ dα = [ ∂J T ∂ω ] I dω dα + [ ∂J T ∂α ] II (A.3)
where the subscripts I and II are conventionally used to distinguish the contributions due to the variations of ω (δω) and of α (δα); this notation will be useful to ease the comprehension of the techniques that will be explained afterwards.

Except for the finite differences method, other numerical methods that can be applied to the gradient evaluation problem may be classified in Direct Methods and Adjoint Methods (AM), whose efficiency depends on the values of N j and N α [START_REF] Peter | Numerical sensitivity analysis for aerodynamic optimization: A survey of approaches[END_REF]. The first approach consists in the direct derivation A.4. ADJOINT METHODS

The adjoint method has been introduced in fluid mechanics by Lions [START_REF] Lions | Optimal Control of Systems Governed by Partial Differential Equations, ser. Grundlehren der mathematischen Wissenschaften[END_REF] and Pironneau [298],

who applied it to flows governed by elliptic PDEs [299]. The first application of the continuous adjoint comes from Jameson, who used it for the optimal design of aerodynamic shapes in transonic flows:

he started with the inverse design of airfoils and wings in inviscid flows ( [START_REF] Jameson | Aerodynamic design via control theory[END_REF], [START_REF] Jameson | Control theory based airfoil design using the Euler equations[END_REF], [START_REF] Jameson | Optimum aerodynamic design using CFD and control theory[END_REF]) and then he extended this formualtion to viscous flows [START_REF] Jameson | Optimum Aerodynamic Design Using the Navier-Stokes Equations[END_REF], by using the Navier-Stokes equations as governing state equations.

Since 1988, Jameson applied the adjoint approach, preferentially the continuous formulation, in a widely manifold collection of works including the inverse design of wings in subsonic, transonic and supersonic flows using multiblock structured grids ( [START_REF] Jameson | Automatic aerodynamic optimization on distributed memory architectures[END_REF], [START_REF] Reuther | Aerodynamic shape optimization of supersonic aircraft configurations via an adjoint formulation on distributed memory parallel computers[END_REF]), the minimization of drag and/or maximization of lift in inviscid [305] and viscous flows [306], sonic boom reduction for supersonic flows ( [START_REF] Alonso | Advanced algorithms for design and optimization of Quiet Supersonic Platforms[END_REF], [START_REF] Jameson | An adjoint method for the calculation of remote sensitivities in supersonic flow[END_REF]), shock wave reduction in external aerodynamics [START_REF] Harbeck | Exploring the Limits of Shock-free Transonic Airfoil Design[END_REF], unsteady aerodynamic design of isolated airfoils and wings ( [START_REF] Nadarajah | Optimum Shape Design for Unsteady Three-Dimensional Viscous Flows Using a Nonlinear Frequency-Domain Method[END_REF], [START_REF] Nadarajah | Non-Linear Frequency Domain Based Optimum Shape Design for Unsteady Three-Dimensional Flow[END_REF]), aerodynamic design of full aircraft configurations [START_REF] Jameson | Aerodynamic Shape Optimization of Complete Aircraft Configurations using Unstructured Grids[END_REF], multipoint drag minimization [START_REF] Leoviriyakit | Multipoint Wing Planform Optimization via Control Theory[END_REF] and multidisciplinary optimization ([314], [START_REF] Leoviriyakit | Aero-Structural Wing Planform Optimization Using the Navier-Stokes Equations[END_REF]).

In the same period, Giles made a considerable improvement in the development of the discrete adjoint approach [15]: he firstly analized the properties of solutions for the adjoint equations ( [START_REF] Giles | Adjoint equations in CFD -Duality, boundary conditions and solution behaviour[END_REF], [START_REF] Giles | On the Properties of Solutions of the Adjoint Euler Equations[END_REF]) and then he presented the analytical solutions of the adjoint equations using Green's functions [START_REF]Analytic adjoint solutions for the quasi-one-dimensional Euler equations[END_REF]; moreover, he provided an exact approach for the solution of the adjoint equations in [319]; finally he used the DAM to study a shock wave reduction method in external aerodynamics ( [START_REF] Giles | Discrete Adjoint Approximations with Shocks[END_REF]) and the harmonic approach to turbomachinery steady and unsteady designs [START_REF] Campobasso | Adjoint Calculation of Sensitivities of Turbomachinery Objective Functions[END_REF].

Since then, several improvements for both the continuous and the discrete approaches have been implemented ( [START_REF] Burgreen | Three-dimensional aerodynamic shape optimization using discrete sensitivity analysis[END_REF], [START_REF] Anderson | Aerodynamic design optimization on unstructured grids with a continuous adjoint formulation[END_REF], [START_REF] Hazra | Aerodynamic shape optimization using simultaneous pseudo-timestepping[END_REF], [325]), all related to the computation of first order sensitivity derivatives and their application to optimization in CFD; as a consequence of this, in the past decade the adjoint methods have become an active topic widely used in the gradient-based design optimization of turbomachinery blades in steady ( [START_REF] Luo | Three-Dimensional Aerodynamic Design Optimization of a Turbine Blade by Using an Adjoint Method[END_REF], [START_REF] Luo | Turbine Blade Row Optimization Through Endwall Contouring by an Adjoint Method[END_REF], [START_REF] Pini | Turbomachinery design optimization using adjoint method and accurate equations of state[END_REF], [START_REF] Luo | Multi-Objective Design Optimization of a Transonic Compressor Rotor Using an Adjoint Equation Method[END_REF], [329]) and unsteady applications ( [START_REF] Kavvadias | Optimal Flow Control and Topology Optimization Using the Continuous Adjoint Method in Unsteady Flows[END_REF], [START_REF] Carnarius | Optimal Control of Unsteady Flows Using a Discrete and a Continuous Adjoint Approach[END_REF]).

Recently, the adjoint method started to be applied also to other disciplines than aerodynamics, like fluid-structure interaction [START_REF] Bazilevs | Adjoint-based Control of Fluid-Structure Interaction for Computational Steering Applications[END_REF] or heat transfer [START_REF] Gkaragkounis | The continuous adjoint method for shape optimization in Conjugate Heat Transfer problems with turbulent incompressible flows[END_REF]; in these kinds of problems higher order derivatives are required. Second-order sensitivity analysis for CFD optimization using the discrete adjoint has been based on automatic differentiation [START_REF] Sherman | First-and Second-Order Aerodynamic Sensitivity Derivatives via Automatic Differentiation with Incremental Iterative Methods[END_REF] or hand differentiation ( [START_REF] Papadimitriou | Direct, adjoint and mixed approaches for the computation of Hessian in airfoil design problems[END_REF], [START_REF] Zervogiannis | Total pressure losses minimization in turbomachinery cascades using the exact Hessian[END_REF]), while the corresponding continuous approach has been presented in [START_REF] Papadimitriou | Computation of the Hessian matrix in aerodynamic inverse design using continuous adjoint formulations[END_REF] for inviscid and in ([338], [339]) for A.5. DISCRETE ADJOINT METHOD (DAM) viscous flows. An extension to third-order sensitivity analysis for robust design is presented in ( [START_REF] Papoutsis-Kiachagias | Robust design in aerodynamics using third-order sensitivity analysis based on discrete adjoint. Application to quasi-1d flows[END_REF], [START_REF] Papadimitriou | Third-order sensitivity analysis for robust aerodynamic design using continuous adjoint[END_REF]).

Finally, in the last years the adjoint method has also been used to speed-up some studies for the assessment of performance in problems governed by uncertainty ( [START_REF] Giebmanns | Compressor Leading Edge Sensitivities and Analysis With an Adjoint Flow Solver[END_REF], [START_REF] Zamboni | Gradient-Based Adjoint and Design of Experiment CFD Methodologies to Improve the Manufacturability of High Pressure Turbine Blades[END_REF], [START_REF] Xiong | Statistical Evaluation of the Performance Impact of Manufacturing Variations for Steam Turbines[END_REF], [START_REF] Yang | Performance Impact of Manufacturing Variations for Multistage Steam Turbines[END_REF], [START_REF] Luo | Statistical evaluation of performance impact of manufacturing variability by an adjoint method[END_REF], [START_REF] Papadimitriou | Bayesian uncertainty quantification of turbulence models based on high-order adjoint[END_REF], [START_REF] Marta | Assessing Turbomachinery Performance Sensitivity to Boundary Conditions Using Control Theory[END_REF]) and for robust design and optimization of components in turbomachinery ( [START_REF] Cinnella | Hybrid Adjoint-based Robust Optimization Approach for Fluid-Dynamics Problems[END_REF], [START_REF] Papadimitriou | Aerodynamic shape optimization for minimum robust drag and lift reliability constraint[END_REF]).

A.5 Discrete Adjoint Method (DAM)

For the discrete adjoint approach, both J and R must be considered in discrete form. Therefore, as the differential form of the governing equation R Equation (A.2) is equal to zero, it is possible to plug it in the gradient equation Equation (A.3), obtaining Equation (A.10).

dJ dα = [ ∂J T ∂ω ] dω dα + [ ∂J T ∂α ] -Ψ T ([ ∂R ∂ω ] dω dα + [ ∂R ∂α ]) (A.10)
where Ψ is an arbitrary line vector with dimension N ω containing the discrete adjoint variables [15]. The first term of Equation (A.11) can be eliminated by choosing the values for the components of vector Ψ to satisfy the adjoint Equation (A.12).

∂J T ∂ω -Ψ T [ ∂R ∂ω ] = 0 ⇒ [ ∂R ∂ω ] T Ψ = ∂J ∂ω (A.12)
Thus, gradient dJ dα can be easily calculated as in Equation (A. 13).

dJ dα = ∂J T ∂α -Ψ T [ ∂R ∂α ] (A.13)
For the works here presented, the direct and the discrete adjoint codes presented in [START_REF] Cinnella | Hybrid Adjoint-based Robust Optimization Approach for Fluid-Dynamics Problems[END_REF] were used; the latter was implemented by means of the algorithmic differentiation software Tapenade [START_REF] Hascoët | The Tapenade Automatic Differentiation tool: Principles, Model, and Specification[END_REF].

A.6 Continuous Adjoint Method (CAM)

For the continuous adjoint, the adjoint equations with respect to J are derived from the continuous form of the governing equations R; afterwards, they are discretized and solved. where θ T is an arbitrary line vector with dimension N ω containing the continuous adjoint variables.

As the objective is to obtain the gradient dJ dα , Equation (A.15) can be differentiated, and with the assumption that θ T is also differentiable, it is possible to integrate the integrals containing θ T by parts, obtaining the formulation as in Equation (A.16). where K, S, and f were already defined in Equation (5.3). At this point, since θ T is an arbitrary differentiable function, it can be chosen in order to remove the dependence of dJ from the variation of the state vector dω, obtaining the differential adjoint problem defined in Equation (A.17).

] L 0 + [ θ T K dS dα ] L 0 (A.20)
In the present work, the system of linear equations defined in Equation (A.18) is discretized by a cell-centered finite volume formulation, using Rusanov's first-order scheme for space integration, while the boundary conditions in Equation (A. [START_REF] Guercio | 15 -biomass-fired organic rankine cycle combined heat and power systems[END_REF]) with a second-order backward finite differences scheme; this set of equation was solved with Gaussian elimination with partial pivoting.

A.6.1 Discrete and Continuous Adjoint Validation

In order to validate both the discrete and continuous adjoint, a comparison with gradients calculated with finite differences was performed for the baseline geometry defined in Table 5.2. For finite differences, a second-order central scheme was employed with a discretization step of 10 -5 .

In Table A.1, we compare the gradients from the finite differences method vs. the one calculated with the discrete adjoint, while Table A.2 presents the comparison with the continuous adjoint. Abstract: Robust design optimization (RDO) is an important tool for the design of industrial products under uncertainty. It combines optimization algorithms and uncertainty quantification (UQ) techniques. Quantification of uncertainties is generally too expensive for complex numerical models of engineering systems. With the aim of developing efficient RDO strategies designed for industrial applications, the coupling of parsimonious UQ techniques with a multi-objective genetic algorithm based on surrogate models (SMOGA) was studied. In this regard, a promising RDO technique was used, based on the coupling of two nested surrogate models: the first is used for UQ, while the response surface of the second is used to accelerate optimization; an infill criterion is used to update the surrogate model during optimizer convergence. Several UQ methods using information on the gradients of the solution with respect to the uncertain variables were implemented and compared in terms of precision and computational cost. We then selected a so-called "low fidelity" UQ method, i.e. inexpensive but not very accurate, and a "high fidelity" method in order to build a multi-fidelity surrogate model for robust optimization. This model allows to have an accuracy close to the high fidelity model for a much lower computation cost. The methods under investigation were applied to the RDO of organic Rankine cycles (ORC) and to the shape optimization of an ORC turbine blade grid, with very promising results.

Keywords: Robust design optimization, uncertainty quantification, organic Rankine cycles, turbomachinery Résumé: L'optimisation robuste (RDO) est un outil important pour la conception de produits industriels sous incertitude. Elle combine des algorithmes d'optimisation et des techniques de quantification de l'incertitude (UQ). La quantification d'incertitudes est généralement trop coûteuse pour des modèles numériques complexes de systèmes en ingénierie. Dans le but de développer des stratégies de RDO efficaces conçues pour des applications industrielles, le couplage de techniques UQ parcimonieuses avec un algorithme génétique multi-objectif basé sur des modèles substituts (SMOGA) a été étudié. A cet égard, une technique RDO prometteuse a été utilisée, basée sur le couplage de deux modèles substituts imbriqués: le premier est utilisé pour l'UQ, tandis que la surface de réponse du second est utilisée pour accélérer l'optimisation; un critère d'enrichissement est utilisé pour actualiser le modèle substitut pendant la convergence de l'optimiseur. Plusieurs méthodes d'UQ utilisant des informations sur les gradients de la solution par rapport aux variables incertaines ont été mises en oeuvre et comparées en termes de précision et coût de calcul. Nous avons ensuite sélectionné une méthode UQ dite «basse fidélité», c'est-à-dire peu coûteuse mais pas très précise, et une méthode «haute fidélité» afin de construire un modèle substitut multi-fidélité pour l'optimisation robuste. Ce modèle permet d'avoir une précision proche du modèle haute fidélité pour un coût de calcul bien moindre. Les méthodes étudiées ont été appliquées à la RDO de cycles thermodynamiques de Rankine à fluide organique (ORC) et à l'optimisation de forme d'une grille d'aubes de turbines ORC avec des résultats très prometteurs.

Mots clés: Optimisation robuste, quantification de l'incertitude, cycles organiques de Rankine, turbomachines

  Rankine à fluide organique (ORC) et à l'optimisation de forme d'une grille d'aubes de turbines ORC avec des résultats très prometteurs. Mots-clés : Optimisation robuste, quantification de l'incertitude, cycles organiques de Rankine, turbomachines RESUME CONTENTS A.2 Finite Differences . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . A.3 Discrete Direct Method . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . A.4 Adjoint Methods . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . A.5 Discrete Adjoint Method (DAM) . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . A.6 Continuous Adjoint Method (CAM) . . . . . . . . . . . . . . . . . . . . . . . . . . . . A.6.1 Discrete and Continuous Adjoint Validation . . . . . . . . . . . . . . . . . . . . B Global sensitivity analysis with Sobol indexes LIST OF TABLES 5.9 Cost for the computation of the optimization DOE S DOE for the LF and the HF UQ method . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 5.10 Statistical moments of the QoI PDF for the baseline geometry and the RDO optimum 5.11 Performance parameters of the optimized cascade . . . . . . . . . . . . . . . . . . . . . 6.1 Caractéristiques des pdfs utilisés pour modéliser les paramètres incertains. . . . . . . . 6.2 Paramètres de conception de la tuyère quasi 1D utilisée pour l'évaluation des techniques UQ. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 6.3 Précision et efficacité de calcul des solveurs UQ (Adj. = Adjoint, DA = Discrete adjoint, CA = Continous adjoint). . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 6.4 Solutions optimales de la RDO en fonction des différentes stratégies (DA = Discrete adjoint, CA = Continous adjoint). . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 6.5 Conditions de fonctionnement incertaines . . . . . . . . . . . . . . . . . . . . . . . . . 6.6 Données sur les grilles de calcul adoptées . . . . . . . . . . . . . . . . . . . . . . . . . . 6.7 Définition de la condition limite pour le point de simulation déterministe . . . . . . . . 6.8 Paramètres de performance de la grille . . . . . . . . . . . . . . . . . . . . . . . . . . . 6.9 Coût pour le calcul de l'optimisation DOE S DOE pour la méthode LF et la méthode HF UQ . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 6.10 Moments statistiques de la PDF QoI pour la géométrie de base et l'optimum RDO . . 6.11 Paramètres de performance de la grille optimisée . . . . . . . . . . . . . . . . . . . . . A.1 Gradients from finite differences vs. gradients from the discrete adjoint for the baseline geometry. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . A.2 Gradients from finite differences vs. gradients from the continuous adjoint for the baseline geometry. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . Chapter 1

Figure 1 . 1 :

 11 Figure 1.1: Installed ORC capacity per year (adapted from [1])

Figure 2 . 1 :

 21 Figure 2.1: Phases of modeling and simulation[START_REF] Schlesinger | Terminology for Model Credibility[END_REF]

Figure 2 . 2 :

 22 Figure 2.2: Verification process [5]

Figure 2 . 3 :Figure 2 . 4 :

 2324 Figure 2.3: Validation process [5]

Figure 2 . 6 :

 26 Figure 2.6: Classification of uncertainty and major sources of uncertainty in computerized solutions [7]

Figure 2 . 7 :

 27 Figure 2.7: General sketch of the UQ framework[START_REF] Marelli | Ingredients for an innovative uncertainty quantification platform in matlab[END_REF]

Figure 2 . 8 :

 28 Figure 2.8: Monte Carlo process

2. 4 .Figure 2 . 9 :

 429 Figure 2.9: LHS, size N = 5 and two random variables

2. 4 .Figure 2 . 10 :

 4210 Figure 2.10: Lattice based sampling technique for a two variable problem

2. 4 .Figure 2 . 11 :

 4211 Figure 2.11: General use of surrogates for problem analysis or optimization

Figure 2 .

 2 Figure 2.12: BK: from the prior to prediction passing through the posterior

Fig. 2 .

 2 Fig.2.13b depicts an example of prior sample paths generated by this covariance function.

2. 5 .

 5 UNCERTAINTY QUANTIFICATION METHODOLOGY (a) Dependence from the θ hyperparameter. (b) Example of prior pathlines (θ = 1).

Figure 2 .

 2 Figure 2.13: SE Covariance

Figure 2 . 14 :

 214 Figure 2.14: Dependence of the Matérn class of covariance functions from the ν hyperparameter

  (a) Dependence from the l hyperparameter (ν = 5/2). (b) Example of prior pathlines (l = 1, ν = 5/2).

Figure 2 .

 2 Figure 2.15: Matérn Covariance

Figure 2 . 16 :

 216 Figure 2.16: Kriging approximation of the test function with 3 samples

Figure 2 . 17 :

 217 Figure 2.17: Kriging approximation of the test function with 4, 5 and 6 samples

Figure 2 . 18 :

 218 Figure 2.18: GEK approximation of the test function with 4 samples

Fig. 2 .

 2 Fig.2.19, Fig.2.20, Fig.2.21, Fig.2.22 and Fig.2.23: they always first depict the exact analytical solution and then they compares BK and GEK predictions, showing both a scatter-plots and as contour-plots.

Figure 2 . 19 :Figure 2 . 20 :Figure 2 . 21 :Figure 2 . 23 :

 219220221223 Figure 2.19: Comparison of Kriging and GEK on 2D 2nd order test function

  [START_REF] Inc | Genetron® 245fa (r-245fa) -european refrigerants[END_REF] as a working fluid. Thermodynamic properties of the working fluid are computed by using an advanced equation of state based on the Helmholtz free energy formulation, available through the thermodynamic library Coolprop[132]. A schematic representation is provided in Fig.2.24.

Figure 2 .

 2 Figure 2.24: ORC schematic representation

3. 2 .

 2 Figure 3.1: A typical Pareto front for a 2 objective optimization problem

3. 2 .

 2 DESIGN AND OPTIMIZATION is provided in Fig. 3.2.

Figure 3 . 2 :

 32 Figure 3.2: Schematic representation of the design process of a component

3. 4 .Figure 3 . 3 :

 433 Figure 3.3: High-fidelity vs low-fidelity models [9]

3. 5 .Figure 3 . 4 :

 534 Figure 3.4: Nested sampling design for MF models[START_REF] Fernández-Godino | Issues in deciding whether to use multifidelity surrogates[END_REF]

3. 5 ,

 5 LF sampling points (in orange) and the HF ones (in blue) are first generated independently, and then the LF nearest point to each HF point is moved on top of their corresponding nearest neighbor.

Figure 3 . 5 :

 35 Figure 3.5: Nearest neighbor sampling for MF models[START_REF] Fernández-Godino | Issues in deciding whether to use multifidelity surrogates[END_REF]

3. 5 .

 5 ROBUST DESIGN OPTIMIZATION STRATEGIESsamples approximately equal to 10 times the cardinality of the parameter space to build a reasonably accurate approximation[START_REF] Loeppky | Choosing the sample size of a computer experiment. a practical guide[END_REF], the nested BK RDO strategy needs O(100 × n unc × n des ) function evaluations (with n unc the number of uncertain parameters and n des the number of design variables) in the first generation of the MOGA to generate the initial kriging surrogates for the statistical moments. Additional O(10 × n unc ) evaluations are required for each update of the external kriging surrogate. Such a number of function calls is still too expensive for complex industrial problems, even for uncertain or design spaces of low to moderate dimensionality (up to about eight uncertain or design parameters). If gradients of the cost functions with respect to the uncertain variables are

3. 5 .Figure 3 . 6 :

 536 Figure 3.6: Flowchart of the BK-based SMOGA RDO process.

3. 5 .Figure 3 . 7 :

 537 Figure 3.7: Single-objective expected improvement

  HF in S DOE,HFRun HF solver to compute the statistics µ(s(i) HF ), σ 2 (s (i)HF ) Use Equation (3.7) to construct MF surrogates: S M F,µ (S DOE,HF ), S M F,σ 2 (S DOE,HF ) 2. MF-SMOGA loop While n < n gen If MOEI infill prescribed: Solve MOEI problem on MF surrogate to choose a new sample s opt If EI>tol: Run LF solver to compute µ sopt ,σ 2 sopt Re-train S M F,µ ,S M F,σ 2 by adding s opt Else: Run HF solver to compute µ sopt ,σ 2 sopt Re-train S M F,µ ,S M F,σ 2 by adding s opt 3.6. SUMMARY OF THE CHAPTER

4. 1 . 1 Figure 4 . 1 :

 1141 Figure 4.1: Scheme of the considered ORC for a WHR application

Figure 4 . 2 :

 42 Figure 4.2: Temperature-Entropy chart of the ORC

Figure 4 .

 4 Figure 4.3: T-Q% Chart for the heat exchange in thermal loop

P

  Figure 4.4: η turb AS A FUNCTION OF PRESSURE RATIO P ev /P cond

4. 1 .Figure 4 . 5 :

 145 Figure 4.5: Validation of BK with a DOE of 120 samples

4. 1 .Figure 4 . 6 :Figure 4 . 7 :Figure 4 . 8 :

 1464748 Figure 4.6: Pareto front from the "brute-force" approach for RDO

Figure 4 . 11 :

 411 Figure 4.11: Typical uncertainty and expenditure profiles for a geothermal project[START_REF]Geothermal power: Technology brief[END_REF]

Figure 4 . 13 :

 413 Figure 4.13: Flowchart of the ORC model in design mode (nominal conditions).

Figure 4 . 14 :

 414 Figure 4.14: Input data for the design of the main components of the ORC at nominal conditions

Figure 4 . 15 :Figure 4 . 16 :

 415416 Figure 4.15: Algorithm of the ORC simulation model for off design

1 .

 1 CF1 is the mean value of the PDF of the QoI (E[W ORC,net ]), 2. CF2 is the variance of the PDF of the QoI (var[W ORC,net ]),3. CF3 is the global area of HEXsA ORC .The latter has been chosen as an indicator of the cost of the ORC: in fact, as a first approximation the cost of HEXs is proportional to their area, while the cost of components like the turbine and the pump is less sensitive to the design parameters as they benefit from a kind of economy of scale, due to the fact that the fix part of their cost is usually quite important and the variable part is usually proportional to their rated power. Consequently, the optimization problem has 3 objectives, which are the maximization of CF1, the minimization of CF2 and the minimization of CF3.The design variables of the optimization process are the source characteristics m ̇geo and T geo and the condensation temperature T cd . Such parameters are uncertain, as discussed in the preceding sections. An additional design parameter, the superheating ∆T SH , is also considered, and it is treated as deterministic. As an outcome of the UQ analysis performed (see results in Section 4.2.4), any other 4.2. ROBUST DESIGN OPTIMIZATION OF AN ORC FOR GEOTHERMAL APPLICATIONS source of uncertainty, and namely c geo , is neglected and the corresponding input parameters are set to their average values and assumed as deterministic.

Figure 4 . 17 :

 417 Figure 4.17: RDO simplified scheme

  2.2.2 is first investigated by means of the UQ algorithm described above. The ORC model is first run in the design mode, with input parameters corresponding to the average values of Tab. 4.4. Afterwards, the ORC power output under random variation of the uncertain parameters is determined from a BK surrogate of the ORC model in the simulation mode. More specifically, 50 LHS are drawn from the distribution of the uncertain parameters, and the resulting 50 sets of inputs are supplied to the ORC model to obtain the output W ORC,net . A Bayesian Kriging surrogate of the output is then used to approximate the ORC model response for any value of the uncertain parameters. A set of 10 6 Monte Carlo samples is drawn from the input distributions and an approximation of W ORC,net is evaluated from the surrogate for each of them. These values are used to build a normalized histogram, reported in Fig. 4.18. It is possible to observe that the histogram follows closely a Gaussian distribution (also reported in the figure), with a standard deviation equal to approximately 14% of the mean value.

Figure 4 . 18 :

 418 Figure 4.18: Histogram of the QoI W ORC,net under uncertainty

Figure 4 . 19 :Figure 4 . 20 :

 419420 Figure 4.19: Total Sobol Index from Global sensitivity analysis

1

 1 Assessment of gradient-based RDO strategies for an idealized turbine nozzle . . . . . 128 5.1.1 Problem definition . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 128 5.1.2 Uncertainty quantification and RDO methods . . . . . . . . . . . . . . . . . . . 130 5.1.3 Results . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 131 5.2 RDO of an ORC turbine blade . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 139 5.2.1 Problem definition . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 142 5.2.2 UQ techniques and RDO strategy . . . . . . . . . . . . . . . . . . . . . . . . . 149 5.2.3 Validation of the CFD model . . . . . . . . . . . . . . . . . . . . . . . . . . . . 150 5.2.4 Analysis of the baseline geometry . . . . . . . . . . . . . . . . . . . . . . . . . . 152 5.2.5 Robust design optimization of the cascade . . . . . . . . . . . . . . . . . . . . . 158 5.3 Summary of the chapter . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 166

. 1 )

 1 with a, b, c and d coefficients defining the geometry. The nozzle length is set to L = 10. A typical nozzle geometry is depicted in Fig. 5.1.

Figure 5 . 1 :

 51 Figure 5.1: Baseline quasi-1D nozzle

Figure 5 . 2 :

 52 Figure 5.2: Assessment of UQ methods: approximate pdf of the QoI J for MC, BK and GEK.

Figure 5 . 3 :

 53 Figure 5.3: Global sensitivity analysis with Sobol indexes

Figure 5 . 4 :Figure 5 . 5 :

 5455 Figure 5.4: Objective space

  Figure 5.6: PDF of the QoI J for the optimal solutions in Tab.5.4

Figure 5 . 7 :

 57 Figure 5.7: Baseline geometry and optimized geometries (with 3σ uncertainty intervals) for various RDO strategies.

  turbomachinery flow field in three 2D surfaces, which are usually referred to as S1, S2 and S3, as depicted in Fig. 5.8.

Figure 5 . 8 :

 58 Figure 5.8: Decomposition of the 3D flow field in three 2D surfaces, S1, S2 and S3 [12]

5. 2 .

 2 Figure 5.9: h-s chart with uncertainty at the inlet (blue hatched area) and at the outlet (red hatched area) of the cascade

Figure 5 . 10 :

 510 Figure 5.10: Blade profile modelling with FFD

Figure 5 . 11 :Figure 5 . 13 :

 511513 Figure 5.11: Coarse mesh, with the identification of boundary surfaces

5. 2 .Figure 5 . 15 :Figure 5 . 16 :Figure 5 . 18 :Figure 5 . 19 :

 2515516518519 Figure 5.15: Computational grid used for the validation of the CFD model

Figure 5 . 20 :

 520 Figure 5.20: Pressure distribution in the baseline geometry

Figure 5 . 21 :Figure 5 . 22 :

 521522 Figure 5.21: Convergence to solution of QoI J vs. mesh refinement

Figure 5 .

 5 Figure 5.23: Y plus in the simulation of the baseline geometry on the coarse mesh

  Both UQ methods have been used in the MF-based RDO strategy described in Chapter 3: the MoM is used as low-fidelity UQ model and BK as the high-fidelity one. As reported in Tab. 5.9, the LF optimization DOE S DOE,LF is composed by n LF,init = 80 samples requiring 560 CFD simulations, while the HF optimization DOE S DOE,HF contains n HF,init = 20 samples, with an equivalent computational cost of 640 CFD simulations; the computational time has been respectively 28 and 35 hours on a workstation equipped with a 16 cores Intel(R) Xeon(R) CPU E5-2637 v3 at 3.50 GHz processor. Table5.9: Cost for the computation of the optimization DOE S DOE for the LF and the HF UQ method UQ Method UQ cost for a sample Size of S DOE Total cost of S DOE Time [h] the LF and the HF optimization DOEs have been built, the MOGA initializes randomly a first population composed by 40 individuals in the design space which is let to evolve over 20 generations. During this process, an adaptive infill every 5 generations is calculated with5.2. RDO OF AN ORC TURBINE BLADEthe strategy already explained in Section 3.5.2. On the whole, 3 MOEI are performed, at generation 6, generation 11 and generation 16: the first one, adopting the UQ MoM technique, has taken 0.8 hours, while the others have required 3.2 hours together. Therefore, 1271 CFD simulations are solved for the whole RDO process, corresponding to a computational time of 67 hours, which includes the construction of both LF and HF optimization DOEs, three MOEI and the evolution of the initial population through 20 generation. This is a great improvement comparing with the TNBK approach, since the latter for the same RDO would have required 2560 CFD simulations to build the initial optimization DOE and further 96 for the three MOEI, with an estimated total optimization time of 146.5 hours approximately.The objective space containing all solutions calculated in the optimizer evolution on the response surface of the surrogate coupled with the MOGA is plotted in Fig.5.26: from the analysis of the Pareto front, depicted as a black dashed line, it is possible to identify a negative correlation between the RDO objects, since the MOGA finds that solutions minimizing the mean value of the PDF of the QoI E[J], maximize its variance var[J]. As a consequence, a unique optimum solution does not exist, as it is not possible to minimize simultaneously both E[J] and var[J] and different compromise solutions among the cost functions are possible. Among all the optimal solutions lying on the Pareto front, the one offering the best trade-off among all optimization targets have been selected; this point is represented as a red symbol in Fig 5.26.The evolution of the SMOGA seeking the optimal solutions in the objective space is further shown in Figs. 5.27 -5.30: the solutions belonging to the first generation and to the generations after each MOEI infill are highlighted in red. In the background of each plot, the reader can also find all points calculated during the SMOGA optimization, colored in grey. It appears that the SMOGA quickly reaches the region near the optimum after only two MOEI infills, as a considerable part of the generation calculated after the third MOEI lies already on the Pareto front.

Figure 5 .Figure 5 . 27 :Figure 5 . 28 :Figure 5 . 29 :Figure 5 . 30 :

 5527528529530 Figure 5.26: Objective space

5 . 7 )Figure 5 . 31 :Figure 5 . 32 :Figure 5 . 33 :

 57531532533 Figure 5.31: Comparison of the baseline and the optimal solution found

Figure 5 . 34 :

 534 Figure 5.34: Mach number distribution in the RDO optimal geometry

  6.1 et de la Fig. 6.2 montrent une explosion du nombre d'ORC installées depuis le début des années 2000 combinée simultanément à une augmentation rapide de la taille moyenne de ces installations. De la figure 6.1, il ressort que le principal domaine d'application des ORC est la géothermie, même si cette technologie se caractérise par un faible nombre d'installations : les centrales géothermiques

Figure 6 . 1 :Figure 6 . 2 :

 6162 Figure 6.1: Capacité ORC installée par an (adapté de [1])

Figure 6 . 3 :

 63 Figure 6.3: Exemple d'application employant un cycle non régénératif à un niveau de pression [3]

Figure 6 . 4 :

 64 Figure 6.4: Exemple d'une installation de 3,2 MW pour la récupération de la chaleur résiduelle d'une aciérie[START_REF] Council | Portugal: Feature on pico alto geothermal orc power plant[END_REF]

5 . 8 .

 58 La disposition de la centrale, en fonction de la disponibilité de la place pour tous les composants. R ÉSUM É ÉTENDU 6. Le fluide de travail. 7. Le type d'environnement (intérieur ou extérieur) et la présence d'acidité dans le cas des centrales géothermiques. Les restrictions de sécurité, en fonction du pays et des applications. 9. La réglementation sur les fluides et leurs interactions avec la pollution (c'est-à-dire le protocole de Montréal ou la réglementation européenne sur les gaz fluorés).

  tant que le R ÉSUM É ÉTENDU manque de connaissance peut se produire dans les modèles physiques ou dans les paramètres d'entrée, avec un risque potentiel d'affecter la fiabilité de la simulation. L'incertitude peut être classée en deux catégories : l'incertitude aléatoire et l'incertitude épistémique. La première est liée à la variabilité physique du système ou de son environnement et ne peut être réduite : la seule façon de la traiter est de la caractériser en réalisant des expériences supplémentaires, afin d'obtenir davantage de données modélisant les variables qui peuvent être utilisées en conséquence dans une approche probabiliste. D'autre part, l'incertitude épistémique découle des hypothèses et des simplifications apportées à la formulation physique : ainsi, elle peut être considérée comme réductible en réalisant davantage d'expériences et en utilisant les informations pour améliorer les modèles physiques.

Figure 6 . 5 :

 65 Figure 6.5: Classification de l'incertitude et des principales sources d'incertitude dans les solutions informatisées [7]

Figure 6 . 6 :

 66 Figure 6.6: Esquisse générale du cadre des UQ[START_REF] Marelli | Ingredients for an innovative uncertainty quantification platform in matlab[END_REF]

Figure 6 . 7 :

 67 Figure 6.7: Représentation schématique du processus de conception d'un composant

  Figure 6.13: Évolution de la convergence pour BK, GEK, MoM et MF

Figure 6 . 14 : 4 Figure 6 . 15 :

 6144615 Figure 6.14: PDF de la QoI J pour les solutions optimales en Tab.5.4

Figure 6 . 18 :Figure 6 . 19 :Figure 6 . 20 :Figure 6 . 21 :

 618619620621 Figure 6.18: Maillage grossier, avec identification des surfaces limites

Figure 6 . 22 :Figure 6 . 23 :

 622623 Figure 6.22: Grille de calcul utilisée pour la validation du modèle CFD

Figure 6 . 25 :

 625 Figure 6.25: Champ d'entropie dans la géométrie de base

Figure 6 . 26 :

 626 Figure 6.26: Distribution des nombres de Mach dans la géométrie de base

Figure 6 . 28 :Figure 6 . 29 :

 628629 Figure 6.28: Convergence vers la solution du QoI J vs. raffinement du maillage

Figure 6 .

 6 Figure 6.30: Y plus dans la simulation de la géométrie de base sur le maillage grossier

LFigure 6 .Figure 6 . 34 :

 6634 figure 6.33.L'évolution de la SMOGA à la recherche des solutions optimales dans l'espace des objectifs est illustrée plus en détail dans les figures 6.34 -6.37 : les solutions appartenant à la première génération et aux générations suivantes sont mises en évidence en rouge. En arrière-plan de chaque tracé, le lecteur peut également trouver tous les points calculés lors de l'optimisation SMOGA, colorés en gris.Il semble que la SMOGA atteigne rapidement la région proche de l'optimum après seulement deux remplissages MOEI, car une partie considérable de la génération calculée après la troisième MOEI se trouve déjà sur le front de Pareto.La solution optimale trouvée a donc été comparée à la géométrie de base : La figure6.38 montre les deux profils. Il apparaît que le RDO a conduit à un'aube qui est plus chargée à l'avant que la ligne de base, avec une charge plus élevée dans la partie initiale du profil et plus faible à l'extrémité. L'aube est sévèrement transformée par l'optimiseur, tout en conservant une forme régulière et lisse. Cependant,

Figure 6 . 35 :Figure 6 . 36 :

 635636 Figure 6.35: Génération après le premier remplissage MOEI dans la convergence SMOGA

Figure 6 . 37 :

 637 Figure 6.37: Génération après le troisième remplissage du MOEI dans la convergence SMOGA

Figure 6 . 38 :Figure 6 . 39 :Figure 6 . 40 :

 638639640 Figure 6.38: Comparaison de la ligne de base et de la solution optimale trouvée

Figure 6 . 41 :

 641 Figure 6.41: Distribution des nombres de Mach dans la géométrie optimale des RDO

Figure 6 . 42 :

 642 Figure 6.42: Répartition de la pression sur l'aube. Comparaison entre la ligne de base et l'optimum RDO

Equation (A. 10 )

 10 can be arranged to obtain Equation (A.11).

A. 6 .

 6 CONTINUOUS ADJOINT METHOD (CAM)Considering the quasi-1D nozzle geometry described in Section 5.1.1 and the definition of J provided in Equation (5.2), one can write the differential form of J as in Equation (A.14). form of the governing equation R Equation (A.2) is equal to zero, it is possible to multiply it by a vector θ T , then integrate it over the domain, and finally, subtract it from the variation of the cost function dJ, obtaining the relation in Equation (A.15). J = J -

  Summary of the chapter . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 74

	Chapter 2
	Uncertainty quantification

Contents 2.1 Verification and Validation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 32 2.2 Error and Uncertainty . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 36 2.3 UQ Framework . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 37 2.4 Overview about uncertainty propagation . . . . . . . . . . . . . . . . . . . . . . . . . . 40 2.4.1 Monte Carlo method and improved sampling techniques . . . . . . . . . . . . . 40 2.4.2 Polynomial Chaos Expansion (PCE) . . . . . . . . . . . . . . . . . . . . . . . . 43 2.4.3 Surrogate Models . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 44 2.5 Uncertainty quantification methodology . . . . . . . . . . . . . . . . . . . . . . . . . . 48 2.5.1 Methods of Moments (MoM) . . . . . . . . . . . . . . . . . . . . . . . . . . . . 49 2.5.2 Bayesian Kriging (BK) . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 49 2.5.3 Gradient Enhanced Kriging (GEK) . . . . . . . . . . . . . . . . . . . . . . . . . 57 2.6 Preliminary validations . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 58 2.6.1 Simple analytic test functions . . . . . . . . . . . . . . . . . . . . . . . . . . . . 59 2.6.2 UQ of a (simple) ORC . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 68 2.7

  Despite the multiple simplifications already presented in the employed ORC model, the dominant influence of turbine isentropic efficiency η t on the overall ORC efficiency is an important conclusion of this study, that justifies all the efforts that the scientific community has been making in last years in order to perform the robust design of turbines.

	2.6. PRELIMINARY VALIDATIONS 2.6. PRELIMINARY VALIDATIONS
	In fact, assuming that the time necessary to train the surrogate models, to get the prediction
	and to derive the gradient of the QoI (the latter just in case of GEK) is negligible if compared with
	the computational time of the black-box simulation model (which here is the non expensive ORC
	model), the most influential parameter on the overall time required for BK and GEK is represented
	a) Kriging 16 samples by the number of samples used to build the DOE. This assumption can be considered true when the b) GEK 4 samples Chapter 3
	computational time of the black-box model is of the order of seconds (or longer). Finally, a global
	sensitivity analysis is carried out, by mean of an ANOVA decomposition aiming to identify the greatest contribution to the variance of the QoI. In order to get this information, 1 st order Sobol Indexes and Robust design optimization
	c) Kriging 20 samples Total Sobol Indexes are calculated by performing a MC on the response surface of BK and GEK; the d) GEK 5 samples
	GEK results depicted in Tab. 2.1 and Fig. 2.27 are obtained.
	Contents
	Table 2.1: 1 st order and Total Sobol Indexes for the 5 uncertain parameters 3.1 Introductory notions in optimization . . . . . . . . . . . . . . . . . . . . . . . . . . . . 76
	e) Kriging 25 samples 3.2 Design and optimization . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 78 f) GEK 7 samples Sobol 1 st Order Sobol Total Index
	Figure 2.25: Benchmarch for Kriging and GEK on ORC test case T cond 74 % 74 % P ev 4 % 4 % 3.3 Robust Design Optimization . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 81
	η t 3.4 Multi-fidelity surrogate models . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 82 21 % 22 %
	η p	0	0 %
	∆T sh 3.5 Robust Design Optimization Strategies . . . . . . . . . . . . . . . . . . . . . . . . . . 85 0 0 %
	3.5.1 BK-Based Robust Design Optimization . . . . . . . . . . . . . . . . . . . . . . 87
	3.5.2 Multi-Fidelity RDO . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 91
	3.6 Summary of the chapter . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 94
	.26.
	Fig. 2.25 shows how GEK outperforms Kriging method: in fact, it needs just a 7 samples DOE
	to get an accurate prediction of the test function where Kriging requires at least 20 samples. Benefits Figure 2.26: Coefficient of variation vs. N samples (BK=Kriging CK=GEK)
	appear to be obvious by considering the computational time necessary for both surrogates (Eq. 2.28 Figure 2.27: Total Sobol Indexes for the 5 uncertain parameters
	and 2.29).
	GEK time : t = N samples DOE t CP U ORC model +
	Kriging time : t = N samples DOE t CP U ORC model + t gradient deriv +
				(2.29)
		t Surrogate training + t Surrogate training +	(2.28)
		t Surrogate prediction t Surrogate prediction

  Algorithm 2 MF RDO loop with the MOEI adaptation. 1. Initialization: LHS LF DOE with n LF,init samples S DOE,LF =

	{ s (1) LF , s (2) LF , ..., s (n LF,init ) LF

} LHS HF DOE with n HF,init << n LF,init samples S DOE,HF = { s (1)

  Preliminary assessment of the TNBK RDO strategy for a simplified ORC . . . . . . . 96 Summary of the chapter . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 126

	Chapter 4
	RDO of ORC systems
	Contents
	4.1

4.1.1 ORC specification . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 96 4.1.2 Methodology . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 101 4.1.3 Results . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 102 4.2 Robust design optimization of an ORC for geothermal applications . . . . . . . . . . . 107 4.2.1 ORC model . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 109 4.2.2 ORC Uncertainty quantification . . . . . . . . . . . . . . . . . . . . . . . . . . 114 4.2.3 ORC optimization strategy . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 118 4.2.4 UQ results and global sensitivity analysis . . . . . . . . . . . . . . . . . . . . . 120 4.2.5 RDO results . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 122 4.3

  4.1. PRELIMINARY ASSESSMENT OF THE TNBK RDO STRATEGY FOR A SIMPLIFIED ORCexplained by Gaussian pdf while Uniform pdf is used for c p gas . A summary of all the uncertain variables just presented is provided in Tab. 4.1 together with the whole amount of information necessary to model them as PDF.

		Table 4.1: List of uncertain variables
	# Uncertain Variable	PDF
	1	m gas

Table 4 .

 4

		2: Number of evaluations of the ORC model
	Method	DOE UQ DOE Optimizzation Nr. Evaluations
	Brute Force	10000	30 × 60	18000000
	TNBK	120	15	1800

  2.1 while Sections 4.2.2 and 4.2.3 describe the employment of BK for the UQ and TNBK for the RDO, respectively. Numerical results are presented in Section 4.2.4 and Section 4.2.5.

Table 4 .

 4 3: Some major results from the ORC design module

	#	Results from the ORC design module
	1	ORC mass flow rate
	2	Pressure at evaporator (P eva )
	3	∆T SH
	4	Turbine Inlet Pressure
	5	Turbine Outlet Pressure
	6	Pressure at condenser (P cd )
	7 Thermodynamic points of the ORC @ nominal point
	8	Pump mechanical power
	9	Cooling water mass flow rate
		Thermal power input
		Turbine mechanical power
		Turbine electric power
		ORC net power (W ORC,net )
		ORC efficiency
		Exchange area of PRE
		Exchange area of EVA
		Exchange area of SH
		Exchange area of CD
		Diameter of the pipe between EVA and turbine

-6 m 2 , a new guess value of T eva and T cd is considered and a new run is performed. To iterate through T eva and T cd the multi-objective

4.2. ROBUST DESIGN OPTIMIZATION OF AN ORC FOR GEOTHERMAL APPLICATIONS

Table 4 .

 4 4: Parameters of ORC uncertain variables, modelled by Gaussian PDFs

	Parameter Distribution	Mean	Variance
	m ̇geo	Gaussian	162.0 kgs -1	72.25 kg 2 s -2
	T geo	Gaussian	388.15 K	6.25 K 2
	c geo	Gaussian	4.2 kJkg -1 K -1 0.0025 kJ 2 kg -2 K -2
	T cw in	Gaussian	297.15 K	4 K 2

Table 4 .

 4

		Sobol 1st order Sobol Total
	m ̇geo	17.6%	18.0%
	T geo	64.4%	64.7%
	c geo	0.40%	0.50%
	T cd	17.1%	17.4%

5: Global sensitivity analysis with Sobol indexes

Table 5 .

 5 2: Design parameters for the quasi 1D nozzle used for the assessment of the UQ techniques.

	Design parameter	PDF
	a	∼ N (1.75, 0.0175)
	b	∼ N (0.699, 0.00699)
	c	∼ N (0.8, 0.008)
	d	∼ N (4.00, 0.004)

Table 5 .

 5 3: Accuracy and computational efficiency of the UQ solvers (Adj. = Adjoint, DA = Discrete adjoint, CA = Continous adjoint).

	MC	1.0E+05	0	3.508 0.6724	-	0.0%	3.0E+06
	BK	60	0	3.507 0.6724	-0.03%	0.0%	1.8E+03
	GEK (DA)	15	15	3.506 0.6724	-0.05%	0.0%	9.0E+02
	GEK (CA)	15	15	3.467 0.6972	-1.16%	3.69%	6.8E+02
	MoM (DA)	1	1	3.451 0.6464	-1.63%	-3.87%	6.0E+01
	MoM (CA)	1	1	3.451 0.6757	-1.63%	0.49%	4.5E+01

UQ method CFD solves Adj. solves E[J] var[J] err % E[J] err % var[J] time [s]

Table 5 .

 5

		5: Uncertain operating conditions
	Variable	PDF Type	PDF mean value	CoV
	P T,in	Gaussian	28.9 bar	4.5%
	T T,in	Gaussian	277.3 degC	0.6%
	P out	Gaussian	10.05 bar	2.0%

Table 5 . 6 :

 56 Data about the adopted computational grids

	Mesh	Number of elements	Average element quality
	Coarse	4772	0.971
	Medium	16453	0.982
	Fine	56885	0.987

Table 5 . 8 :

 58 Performance parameters of the cascade

	Performance parameter	∆s	∆P t	γ

P 6.582 Jkg -1 K -1 173029 P a 9.38%

Table 5 .

 5 10: Statistical moments of the QoI PDF for the baseline geometry and the RDO optimum

	Geometry E[J]	√	var[J]
	baseline	7.3275	2.8680
	optimum 5.5498	2.3339
	Finally, a CFD simulation of the optimized geometry is performed at similar boundary conditions
	than the baseline case (see Tab.

Table 5 .

 5

11: Performance parameters of the optimized cascade Performance parameter s P t γ P 4.618 Jkg -1 K -1 118300 P a 6.28%

  . Pour réduire le nombre d'appels de fonction coûteux, il est crucial de sélectionner des méthodes et des optimiseurs d'UQ parcimonieux, le coût global de la RDO étant généralement le produit du coût des deux approches[START_REF] Congedo | Shape optimization of an airfoil in a BZT flow with multiple-source uncertainties[END_REF]. Les exemples passés de RDO sur les simulations CFD coûteuses incluent diverses formes de solveurs UQ basés sur le chaos polynomial non intrusif expansion[START_REF] Cinnella | Robust optimization of dense gas flows under uncertain operating conditions[END_REF][START_REF] Congedo | Shape optimization of an airfoil in a BZT flow with multiple-source uncertainties[END_REF] ou des modèles substitut tels que la collocation stochastique simplex[START_REF] Congedo | A simplex-based numerical framework for simple and efficient robust design optimization[END_REF] ou kriging[117]. Toutes nécessitent un nombre de solutions CFD qui augmente rapidement avec le nombre de paramètres incertains, et leur couplage direct avec les optimiseurs MOGA n'est pas abordable sur le plan des calculs pour les applications industrielles, surtout si des ordinateurs massivement parallèles ne sont pas disponibles.Une option intéressante pour réduire le coût des solutions UQ consiste à utiliser des informations de gradient. On a déjà vu qu'une méthode simple d'approximation des moments statistiques de la fonction de coût par des expansions en série de Taylor est le MoM de premier ordre. Une alternative à conception déterministe x ∈ R n des et sur un vecteur de paramètres incertains ξ ∈ R nunc . Notez que certains des paramètres de conception peuvent également être incertains. Le problème des RDO est donc formulé en utilisant l'espérance et la variance de J comme mesures de robustesse, ce qui conduit

	à la solution du problème d'optimisation déterministe à deux objectifs dans l'équation (6.1).
	min x E[J ], min x var[J ]	(6.1)

  est utilisé : il fournit un front de Pareto approximatif des solutions optimales correspondant à différents compromis entre performance moyenne et robustesse pour les différentes QoI en présence. Par souci de simplicité, dans ce qui suit, seul le cas d'un seul QoI, m = 1, est pris en compte, mais l'approche peut être étendue à des QoI multiples. Les statistiques requises sur les conflits d'intérêts sont calculées au moyen de certaines méthodes non intrusives de QI, déjà présentées au chapitre 2. Ainsi, le premier ingrédient du processus RDO est une approche UQ efficace, qui fournit des approximations précises de E[J ] et var[J ] sur la base d'un ensemble d'échantillons déterministes N de la solution. × n pop × n gen . Le coût de calcul peut être considérablement réduit en effectuant des analyses déterministes N × n pop en parallèle à chaque génération de la NSGA[START_REF] Hercus | Robust shape optimization of uncertain dense gas flows through a plane turbine cascade[END_REF][START_REF] Congedo | Efficient robust optimization techniques for uncertain dense gas flows[END_REF], mais : (i) le nombre requis de noyaux de calcul peut dépasser les ressources de calcul disponibles, et (ii) même avec une mise à l'échelle parfaite en parallèle à chaque génération, le délai d'exécution de la RDO est au moins égal au coût moyen d'une simulation multiplié par n gen . Pour réduire le coût de calcul, un second modèle substitut (externe) est introduit pour prédire la réponse des fonctions de coût aux paramètres de conception.Dans ce qui suit, deux méthodologies RDO sont présentées : toutes deux sont construites sur une approche RDO basée sur SMOGA introduite dans[START_REF] Cinnella | Robust optimization using nested Kriging surrogates: application to supersonic ORC nozzle guide vanes[END_REF], reposant sur le couplage de deux métamodèles de krigeage bayésien (BK) emboîtés : le premier est utilisé pour calculer les statistiques requises des fonctions objectives dans l'espace des paramètres incertains, tandis que le second est utilisé pour modéliser la réponse de ces statistiques aux variables de conception. Une telle approche est parfois aussi appelée "combined kriging"[START_REF] Keane | Robust design optimization using surrogate models[END_REF]. Un critère de remplissage est utilisé pour mettre à jour le second métamodèle de krigeage lors de la convergence vers l'optimum. Cette technique a été appliquée avec succès à la conception des aubes de turbine pour les cycles de Rankine organiques [117] et au RDO du cycle thermodynamique[START_REF] Serafino | Assessment of an Innovative Technique for the Robust Optimization of Organic Rankine Cycles[END_REF]. En supposant que chaque métamodèle de krigeage nécessite un nombre d'échantillons approximativement égal à 10 fois la cardinalité de l'espace des paramètres pour construire une approximation raisonnablement précise[START_REF] Loeppky | Choosing the sample size of a computer experiment. a practical guide[END_REF], la stratégie BK RDO imbriquée nécessite O(100 × n unc × n des ) des évaluations de fonctions (avec n unc le nombre de paramètres incertains et n des le nombre de variables de conception) dans la première génération du MOGA pour générer les krigeages initiaux de métamodèle pour les moments statistiques. Des évaluations supplémentaires de O(10 × n unc ) sont nécessaires pour chaque mise à jour du métamodèle de krigeage externe. Un tel nombre d'appels de fonction est encore trop coûteux pour des problèmes industriels complexes, même pour des espaces incertains ou de conception de dimension faible à modérée (jusqu'à environ huit paramètres incertains ou de conception). Si les gradients des fonctions de coût par rapport aux variables incertaines sont disponibles, les métamodèles GEK peuvent être utilisés pour réduire le , un métamodèle BK externe a été construit pour décrire les variations de la moyenne et de la variance des QoI dans l'espace de conception, en utilisant des échantillons n init choisis selon un plan d'expériences préliminaire (DOE). Dans ces références, un UQ BK de premier niveau utilisant des échantillons N a été utilisé pour évaluer les fonctions de coût. Par conséquent, comme dans cette configuration la méthode est basée sur deux métamodèles du krigeage bayésien emboîtés, la méthode a été nommée "TNBK" dans[START_REF] Serafino | Assessment of an Innovative Technique for the Robust Optimization of Organic Rankine Cycles[END_REF].Pour réduire davantage le coût de calcul, une approche similaire peut être utilisée, selon laquelle un métamodèle GEK ou le MoM peut remplacer le BK de premier niveau pour l'étape UQ ; un organigramme du processus RDO résultant est fourni dans la Figure6.8. Flowchart of the BK-based SMOGA RDO process. La construction du métamodèle de deuxième niveau nécessite des exécutions déterministes N ×n init du solveur de simulation directe dans le cas de BK et des solveurs direct et adjoint pour GEK (voir l'annexe A pour un aperçu des techniques de calcul des gradients), qui peuvent être parallélisées de manière simple. Le coût supplémentaire modéré des solveurs adjoints pour la GEK est en général contrebalancé par un coût moindre de N nécessaire pour obtenir une précision donnée. Si la MoM est utilisée comme solveur UQ, le coût du métamodèle initial est réduit à n init pour les exécutions déterministes de la simulation directe et des solveurs adjoints. Les tests numériques montrent qu'un métamodèle BK suffisamment précis peut être obtenu en fixant n init = 10n des . Une fois que les estimations de la moyenne et de la variance du QoI ont été obtenues, un BK de deuxième niveau est construit et couplé avec le MOGA. Afin de contrôler la précision des fonctions de coût approximatives, une stratégie de remplissage adaptative est adoptée pour enrichir le métamodèle BK externe au cours de l'évolution. À cette fin, à la fin de l'itération MOGA, un nouvel échantillon est sélectionné, évalué et ajouté à la DOE initiale : ensuite, le modèle BK est retravaillé.

	R ÉSUM É ÉTENDU C'est pourquoi une nouvelle stratégie multifidélité pour la RDO est ici aussi proposée : elle peut réduire considérablement le nombre d'appels de fonction requis en tirant parti d'une MoM de premier ordre peu coûteuse (mais de faible précision) avec une BK ou GEK de plus haute fidélité. Par conséquent, deux modèles de fidélité différents sont fusionnés en utilisant une méthodologie similaire à celle de [183] pour générer un métamodèle pour l'optimisation des MOGA. Les détails des deux stratégies sont donnés ci-après. Dans [117, 40]R ÉSUM É ÉTENDU
	Le couplage direct du MOGA avec le solveur de l'UQ est trop coûteux pour un code de simula-tion complexe et coûteux (c'est-à-dire CFD), en raison du nombre élevé d'évaluations de fonctions. Par exemple, l'exécution du MOGA avec une population de n pop individus sur n gen générations et l'utilisation d'échantillons N pour UQ conduisent à un nombre global d'évaluations de la qualité de le contexte des problèmes de RDO, car il nécessite également le gradient des moments statistiques de la pdf du QoI par rapport aux variables de conception. Obtenir une telle information en utilisant des méthodes adjointes efficaces n'est pas une tâche triviale ; d'autre part, les approximations par différence finie sont facilement applicables, mais au prix d'un coût de calcul considérable pour les espaces de conception à haute dimension. l'information d'environ N nombre d'échantillons pour l'étape UQ ; cependant, un MOGA basé sur GEK n'est pas simple dans Figure 6.8:

  MOEI << n gen . Une mise à jour du MOEI implique l'exécution de l'algorithme UQ pour l'échantillon supplémentaire à intégrer au BK de deuxième niveau. La simulation déterministe directe N et les exécutions adjointes (si nécessaire) requises pour les solveurs UQ BK et GEK peuvent être effectuées R ÉSUM É ÉTENDU en parallèle. Enfin, le coût de la RDO basée sur la stratégie de remplissage adaptative en termes d'évaluation coût-fonction est donné par N ×(n init +n M OEI ) calculs CFD directs pour BK, N ×(n init + n M OEI ) calculs CFD directs et adjoints pour GEK, et (n init + n M OEI ) calculs directs et adjoints pour MoM. Le délai d'exécution, dans le cas d'une mise en oeuvre parallèle parfaitement échelonnée, est approximativement égal à n M OEI +1 runs. Sur la base d'expériences numériques [122], une adaptation du MOEI toutes les trois à cinq générations de MOGA est généralement suffisante pour obtenir une approximation précise de l'optimum. La boucle RDO complète basée sur BK est décrite dans le pseudo-code fourni ci-après dans Algorithme 1.Bien que le SMOGA basé sur BK permette une réduction considérable des appels de fonction pendant la boucle RDO, le coût de calcul global reste important pour la conception en haute dimension et les espaces incertains. Si le coût de calcul est plus élevé pour les solveurs BK ou GEK que pour la MoM, le premier peut fournir une estimation précise des statistiques de QoI, ce qui n'est pas le cas pour la MoM de premier ordre. Dans le but de parvenir à un compromis entre la précision des LF , σ 2 LF , µ δ , et σ 2 δ les moyennes et les variances du modèle LF et de la fonction de divergence δ, respectivement. Afin d'améliorer la précision du métamodèle MF lors de la convergence SMOGA, un remplissage adaptatif basé sur le critère MOEI est utilisé. Dans ce cas, cependant, le modèle LF ou HF peut être utilisé pour le remplissage. Dans les calculs suivants, la stratégie de[START_REF] Meliani | Multi-fidelity efficient global optimization: Methodology and application to airfoil shape design[END_REF] est adoptée : pour le remplissage, la priorité est donnée au modèle LF le moins cher, et le modèle HF n'est utilisé que lorsque l'amélioration obtenue avec le niveau de fidélité le plus faible est inférieure à une tolérance donnée, tol = 10 -4 dans les présents calculs. Dans les deux cas, le rééchantillonnage au même endroit est évité.

	R ÉSUM É ÉTENDU
	sont données par :
	{ µ HF = ρµ LF + µ δ σ 2 HF = ρ 2 σ 2 δ LF + σ 2	(6.3)
	avec µ
	Pour un aperçu de toutes ces techniques, le lecteur est prié de s'adresser à [184, 185]. Parmi elles, le critère d'amélioration attendue (EI) [186] a été sélectionné parce qu'il offre un bon compromis entre l'exploitation et l'exploration. Comme la RDO est intrinsèquement un processus d'optimisation multi-objectifs, la fonction EI est une surface d'un espace paramétrique hyperdimensionnel et une formulation plus complexe est développée, appelée approche de l'amélioration attendue multi-objectifs (MOEI) (voir [188]). La précision du métamodèle est rapidement améliorée avec peu de remplissages MOEI, de sorte que avancé, basé sur un modèle substitut multifidélité de l'espace de conception remplaçant le précédent métamodèle BK externe (MF-SMOGA). Dans le métamodèle MF, le solveur MoM UQ est utilisé comme modèle basse fidélité (LF), et le solveur GEK (ou BK, si les gradients des fonctions de coût par rapport aux variables incertaines ne sont pas disponibles) UQ est le solveur haute fidélité (HF). Dans la présente implémentation, un DoE initial est exécuté à la première itération du SMOGA, où les points d'échantillonnage HF et LF sont générés indépendamment. La corrélation auto-régressive dans l'équation 3.7 est utilisée pour corriger la différence entre les modèles LF et HF. Plus précisément, l'implémentation proposée dans [183] est ici adoptée, et le modèle LF est utilisé comme fonction de base pour l'expression du terme de régression d'un modèle de krigeage universel ; par conséquent, le terme m(ξ) dans Eq. 2.6 devient : i n métamodèles de Krigeage UQ et l'efficacité de calcul de la MoM, cette section présente un SMOGA m(ξ) = ∑ f i (ξ)β i + f LF β ρ (6.2)
	où β ρ est une estimation du coefficient ρ de l'équation (3.7) au moyen d'une régression GP. supposant que
	les modèles LF et HF sont indépendants, la moyenne µ HF et la variance σ 2 HF du modèle haute fidélité

  Le système d'équations est complété par l'équation d'état pour les gaz thermiquement et caloriquement parfaits, P = (γ -1)ρ(e T -v 2 /2). Les équations gouvernantes sont discrétisées par une formulation en volume fini centrée sur la cellule, en utilisant le schéma d'intégration spatiale du premier ordre en amont de Rusanov et le pas de temps explicite Runge-Kutta en quatre étapes[START_REF] Hirsch | Numerical Computation of Internal and External Flows[END_REF].Des conditions limites caractéristiques basées sur les invariants de Riemann sont imposées à l'entrée et à la sortie de la tuyère. Des conditions d'écoulement sonore sont prescrites à l'entrée, de sorte que tous les invariants de Riemann entrent dans le domaine. La plage de variation de la pression totale est telle qu'un choc est toujours créé dans le divergent. En conséquence, les conditions d'écoulement à la sortie sont toujours subsoniques. Dans ce cas, nous imposons la pression statique de sortie, qui est traitée comme déterministe, et fixée à P out = 0, 6 bar. Sur la base d'une étude préliminaire de maillage, une grille de calcul de 300 cellules uniformément espacées est utilisée dans tous les calculs suivants.

Table 6 .

 6 1: Caractéristiques des pdfs utilisés pour modéliser les paramètres incertains. Pour le problème actuel, deux méthodes d'UQ basées sur les gradients ont été utilisées : la GEK et la MoM. Comme le cas test n'est pas particulièrement coûteux à simuler, un solveur UQ basé sur le métamodèle BK est également envisagé : cette technique, qui n'utilise pas d'informations sur les gradients, a été considérée ici comme une référence.BK est également utilisé comme modèle substitut externe pour le SMOGA, qui utilise l'algorithme NSGA-II comme optimiseur. Le SMOGA nécessite des évaluations des fonctions de coût des RDO par la méthode UQ pour chaque individu et pour chaque génération. Pour réduire le coût de calcul de la RDO, une possibilité est de construire un métamodèle BK ou GEK unique sur un espace de paramètres étendu aux variables de conception. En raison de la malédiction de la dimensionnalité, la construction d'un métamodèle fiable sur un si grand espace de paramètres peut entraîner un coût inacceptable de l'approximation krigeante. Pour éviter ce problème et pour faciliter l'utilisation de l'UQ basée sur les gradients, y compris la MdM, on choisit de construire à la place un modèle substitut séparé faisant correspondre une variable dans l'espace de conception à l'espace de la fonction de coût. les méthodes BK, GEK et MoM UQ par rapport à l'échantillonnage MC de référence. Ensuite, les méthodes sont appliquées à un échantillon de géométries de tuyères et utilisées pour construire des métamodèles mono ou multi-fidélité utilisés dans la boucle SMOGA RDO.Une des géométries de la tuyère est sélectionnée dans l'espace de conception en assignant les coefficients géométriques fixés pdfs normal avec un écart-type égal à 1% de la moyenne. Les paramètres pdf pour les variables géométriques sont fournis dans l'onglet. 6.2. Les paramètres de fonctionnement et thermodynamiques se voient attribuer la même fdp que dans le Tab. 6.1. Les pdfs sont échantillonnées afin de construire des métamodèles BK et GEK, tandis que la MoM est appliquée en calculant le QoI et son gradient à la valeur attendue des paramètres d'entrée.

	Quantity PDF Distribution	Distribution parameters
	a	Gaussian	µ ∈ [1.5 -2.0], CoV = 0.01
	b	Gaussian	µ ∈ [0.6 -0.8], CoV = 0.01
	c	Gaussian	µ ∈ [0.7 -0.9], CoV = 0.01
	d	Gaussian	µ ∈ [3.9 -4.1], CoV = 0.01
	P T,in	Uniform	[0.90-1.10]
	γ	Uniform	[1.39-1.41]

Table 6 .

 6 2: Paramètres de conception de la tuyère quasi 1D utilisée pour l'évaluation des techniques UQ.Un résumé des résultats de l'UQ est donné dans l'onglet. 6.3, où sont indiqués la moyenne et la variance approximatives du QoI J selon les différentes méthodes UQ. Un calcul de référence basé sur l'intégration des MC sur des échantillons de 10 5 est effectué pour fournir une solution de référence. gradient, le coût de calcul est réduit de plus d'un tiers par rapport à BK. La précision de la GEK est confirmée par l'examen de la figure 6.10, qui montre le pdfs empirique complet de J calculé avec PC, BK et GEK. Un très bon accord entre MC, BK et GEK est observé.La précision est moins satisfaisante pour l'approche adjointe continue, en raison du calcul moins précis des gradients. Ceci reflète les erreurs numériques introduites par la discrétisation des équations adjointes et le traitement des conditions aux limites, qui ne sont pas complètement cohérentes avec les erreurs de discrétisation introduites par le solveur CFD direct. D'autre part, le solveur continu d'équations adjointes est développé indépendamment du solveur direct et, en ce sens, il est non intrusif.Le coût de calcul des échantillons GEK utilisant l'adjoint direct ou continu est approximativement le même. La méthode MoM est évidemment moins précise que les autres solveurs UQ. Néanmoins, les erreurs sur la moyenne et la variance calculées restent très raisonnables, malgré la présence d'un choc dans la divergence (où l'expansion de la série de Taylor n'est pas définie) et la plage d'incertitude relativement large sur la pression du réservoir (environ 20 %). Le choc est en pratique régularisé par

	Design parameter	PDF
	a	∼ N (1.75, 0.0175)
	b	∼ N (0.699, 0.00699)
	c	∼ N (0.8, 0.008)
	d	∼ N (4.00, 0.004)
	Pour le présent problème de test bon marché, le temps CPU requis pour l'échantillonnage MC est
	d'environ 3 ×10 6 secondes sur un ordinateur personnel ayant un CPU Intel(R) Xeon(R) E5-1620 v3
	à 3,50 GHz. Les résultats sont également présentés pour l'échantillonnage MC sur le métamodèle
	BK. Ce dernier utilise 60 évaluations de fonctions à des points sélectionnés selon un échantillonnage
	hypercube latin (LHS) de l'espace incertain à six dimensions. La taille de l'échantillon correspond	à
	la règle empirique N = 10n unc . Cela est déjà suffisant pour obtenir des erreurs extrêmement faibles
	par rapport à la moyenne et à la variance MC, tout en réduisant le coût de calcul global de l'UQ de
	trois ordres de grandeur. L'UQ GEK basé sur le solveur adjoint discret offre une précision similaire	à
	celle de BK en utilisant seulement 15 échantillons. Malgré les solvants adjoints supplémentaires pour
	les calculs de

diffusion numérique dans les solveurs direct et adjoint, tandis que les plages d'incertitude de tous les autres paramètres sont raisonnablement faibles. Alors que les calculs de MoM discrets et continusadjoints donnent strictement les mêmes résultats pour E[J] (qui n'utilise pas d'informations sur le gradient), ils prévoient des résultats différents pour l'écart type. Les erreurs légèrement plus faibles obtenues pour la méthode continue-adjointe sont l'effet d'une compensation des erreurs pour le cas présent. Comme la méthode MoM ne nécessite qu'un calcul CFD direct et un calcul CFD adjoint, son coût de calcul est essentiellement de 1/15 de l'échantillonnage GEK. Dans l'ensemble, la MoM de premier ordre fournit une estimation raisonnablement précise des statistiques d'ordre inférieur de la QoI et un très bon compromis entre coût et précision pour le problème actuel. Sa précision devrait toutefois diminuer pour les plages d'incertitude plus importantes. Pour cette raison, le modèle MoM est classé comme un modèle de fidélité inférieur à BK ou GEK.

Table 6 .

 6 3: Précision et efficacité de calcul des solveurs UQ (Adj. = Adjoint, DA = Discrete adjoint, CA = Continous adjoint). Cela est dû à la nature hautement non linéaire du problème de la CFD compressible. En raison de cette forte interdépendance de la qualité de l'information sur les paramètres d'entrée, il n'est pas possible de négliger aucun d'entre eux. En conséquence, les calculs de la RDO de la section suivante traitent l'ensemble des six paramètres d'entrée comme des variables aléatoires.

	UQ method CFD solves Adj. solves E[J] var[J] err % E[J] err % var[J] time [s]
	MC	1.0E+05	0	3.508 0.6724	-	0.0%	3.0E+06
	BK	60	0	3.507 0.6724	-0.03%	0.0%	1.8E+03
	GEK (DA)	15	15	3.506 0.6724	-0.05%	0.0%	9.0E+02
	GEK (CA)	15	15	3.467 0.6972	-1.16%	3.69%	6.8E+02
	MoM (DA)	1	1	3.451 0.6464	-1.63%	-3.87%	6.0E+01
	MoM (CA)	1	1	3.451 0.6757	-1.63%	0.49%	4.5E+01

Figure 6.10: Évaluation des méthodes UQ : pdf approximatif du QoI J pour MC, BK et GEK. Les résultats de l'UQ pour MC, BK et GEK peuvent être utilisés pour effectuer une analyse de sensibilité globale et pour identifier les paramètres aléatoires contribuant le plus à la variance du QoI au moyen d'une décomposition ANOVA. Plus précisément, des échantillons MC exacts ou des métamodèles sont utilisés pour calculer les indices de Sobol dans l'espace des paramètres défini par toutes les entrées incertaines. À cette fin, la bibliothèque Python SALib [229] a été adoptée. Seuls les résultats GEK discrets (présentés dans la figure 6.11) sont pris en compte ici, les méthodes MC, BK et GEK continu conduisant à des conclusions similaires. La figure présente les indices de Sobol du premier ordre pour chaque paramètre incertain, et la somme des indices d'ordre supérieur, correspondant à l'interaction entre les paramètres lorsque ceux-ci sont modifiés simultanément. La pression totale à l'entrée P T,in semble être le paramètre le plus influent lorsqu'il est pris seul. Ceci est cohérent avec la plage beaucoup plus large de sa pdf. Néanmoins, les termes d'interaction sont très significatifs. les résultats sont comparés avec les RDO de simple fidélité. Pour les métamodèles monofidélité, une DOE LHS initiale de 10 × n des = 40 UQ a été effectuée avant la première génération. Ensuite, la Figure 6.11: Analyse de sensibilité globale avec les indices de Sobol première population a été initialisée de manière aléatoire dans l'espace de conception et l'AG a été laissée évoluer sur 80 générations. Les remplissages MOEI ont été effectués toutes les 5 générations. Pour le métamodèle MF, le DOE initiale contenait 40 MoM continues LF et 10 UQ GEK discrètes HF. Les remplissages MOEI basse fidélité sont ajoutés toutes les 5 générations, sauf pour les derniers remplissages, pour lesquels le critère de tolérance suggère un remplissage HF. L'évolution de la recherche SMOGA dans l'espace des objectifs E[J] vs. var[J] est représentée sur la figure 6.12 pour divers optimiseurs RDO. Pour le présent problème de conception inverse, le front de Pareto s'effondre sur un point unique correspondant à la minimisation simultanée des deux objectifs E[J] et var[J].

le premier remplissage de la MOEI, à la génération après le deuxième remplissage de la MOEI et à la dernière génération. En outre, en gris, en arrière-plan de chaque sous-placette, le lecteur peut Figure 6.12: Objective space trouver tous les points calculés avec chaque technique lors de l'optimisation SMOGA. Il semble que le Krigeage, la GEK, la MoM et la MF évoluent de manière similaire, car ils atteignent rapidement la zone proche de l'optimum après seulement deux remplissages MOEI. Les solutions optimales calculées avec chaque technique sont présentées dans Tab 6.4. Malgré la précision différente des techniques UQ utilisées, toutes les boucles RDO convergent vers des solutions très similaires identifiées par les valeurs attendues des quatre paramètres de conception [a, b, c, d], qui ne diffèrent que sur le 6ème chiffre décimal ou moins pour BK, GEK et MF, et sur le faire preuve de prudence en préférant systématiquement cette méthode pour les problèmes de RDO, en raison de sa précision moindre. Enfin, la MF SMOGA RDO a un temps de calcul 8 fois plus petit

Table 6 .

 6 4: Solutions optimales de la RDO en fonction des différentes stratégies (DA = Discrete adjoint, CA = Continous adjoint). Comme Enfin, la géométrie optimisée est représentée pour toutes les méthodes employées dans la figure 6.15 et elle est comparée à celle de base. Le tracé est obtenu en propageant la fdp des coefficients géométriques par Eq. (6.4), avec une moyenne égale aux valeurs optimales de RDO et CoV = 0, 01

	UQ method	E[a]	E[b]	E[c]	E[d]	E[J]	var[J]	Optimization time [h]
	BK	1.500000 0.800000 0.900000 3.900000 0.28509 0.08545	∼ 20
	GEK (DA) 1.500002 0.799998 0.900000 3.900000 0.28771 0.08368	∼ 10
	GEK (CA) 1.500000 0.799999 0.900000 3.900000 0.28743 0.08416	∼ 8
	MoM (DA) 1.500042 0.799961 0.898917 3.940293 0.28514 0.08465	∼ 0.7
	MoM (CA) 1.500012 0.799993 0.899999 3.900004 0.28815 0.08502	∼ 0.5
	MF model 1.500000 0.799205 0.899942 3.900006 0.28501 0.08471	∼ 2.5

vérification supplémentaire des résultats des RDO, le plan optimal est recalculé en utilisant l'échantillonnage MC, ce qui donne E[J] = 0, 32325 et var[J] = 0, 09487, ce qui est plutôt en accord avec les estimations de la SMOGA, qui semblent être légèrement trop optimistes dans la prédiction des fonctions objectives. En effet, la fdp complète de la QoI calculée en propageant les paramètres incertains à travers la solution CFD pour la géométrie optimale de la tuyère en utilisant les méthodes BK et GEK (illustrée dans la Fig.

6

.12 à côté de la distribution MC) présente des écarts modérés par rapport à celle de la MC. Globalement, les distributions optimales sont beaucoup plus proches de 0 en moyenne et elles présentent une variance plus faible que la géométrie de base étudiée dans la section précédente, ce qui montre que la RDO améliore effectivement les deux critères. L'objectif de cette RDO est l'optimisation de la forme du profil de l'aube pour la première grille

  grille est supposée fonctionner dans des conditions de fonctionnement incertaines: la pression totale d'entrée P T,in , la température totale d'entrée T T,in et la pression statique de sortie P out sont des sources d'incertitude aléatoire et elles sont modélisées sous forme de PDF gaussiennes. Tab. 6.5 indique la valeur moyenne de ces trois PDF et leur coefficient de variation CoV , défini comme le rapport entre la valeur moyenne et l'écart-type. Ces incertitudes sont également représentées dans le graphique h-s de la figure6.16 : le manque de connaissances sur P T,in et T T,in détermine une incertitude à l'entrée de la grille qui est identifiée avec la zone hachurée en bleu. Dans le même graphique, l'incertitude sur P out est représentée par la zone hachurée rouge. Dans la marge inférieure gauche du graphique, il est possible de trouver la courbe de saturation, ce qui signifie qu'à l'entrée de la grille, le fluide de travail est fortement surchauffé. ÉTENDU il est possible de trouver la surface de l'aube, qui est appelée E sur la Fig.6.18.Dans chaque maille, pour résoudre le sous-couche visqueuse proche des parois de l'aube, un remplissage de la couche limite a été adopté ; on peut en voir une vue rapprochée sur la figure6.19, qui est tirée de la grille moyenne. De plus, comme les effets de bord de fuite sont importants dans le compte des pertes totales, cette région a été particulièrement affinée dans chaque grille de calcul adoptée. La figure6.20 montre un détail du raffinement dans la zone du bord de fuite pour le maillage fin. La qualité des mailles a été évaluée au moyen du paramètre de qualité gamma GMSH, défini comme le rapport entre le rayon inscrit de l'élément et le rayon circonscrit de l'élément [259] : toutes les mailles présentent une distribution similaire de cet estimateur de qualité. À titre d'exemple, la figure6.21

	Table 6.5: Conditions de fonctionnement incertaines
	Variable	PDF Type	PDF mean value	CoV
	P T,in	Gaussian	28.9 bar	4.5%
	T T,in	Gaussian	277.3 degC	0.6%
	P out	Gaussian	10.05 bar	2.0%
	pertes en bout de paroi sont presque totalement négligés dans la simulation de l'écoulement qui circule
	dans la grille. Malgré cette limitation, cette simplification permet d'éviter une augmentation du
	nombre de paramètres de conception et peut être utilisée de manière fiable pour comprendre de
	nombreux phénomènes impliqués dans l'écoulement des turbomachines [257].
	Pour étudier le domaine, trois grilles de calcul non structurées ont été générées avec le code open
	source GMSH : elles ont été construites avec la même stratégie, qui est illustrée dans la Fig. 6.18, et la
	seule différence entre elles est le raffinement. Ainsi, les mailles grossières, moyennes et fines présentent

Figure 6.16: Diagramme h-s avec incertitude à l'entrée (zone hachurée bleue) et à la sortie (zone hachurée rouge) de l'aube Figure 6.17: Modélisation du profil des pales avec FFD déformant le profil de base pour générer une pale optimale avec des performances améliorées. Comme l'optimisation est effectuée sur un plan 2D qui coïncide avec la surface de l'aube, on s'attend à une approximation dans l'évaluation des performances du système : les écoulements secondaires et les une surface d'entrée (identifiée comme A sur la Fig. 6.18), une surface de sortie (indiquée comme B) et deux surfaces périodiques (appelées C et D) ; au milieu de la région délimitée par ces quatre surfaces, R ÉSUM É représente la distribution de la qualité des éléments pour la maille moyenne ; de plus, le tableau 6.6 indique la taille et la valeur moyenne de la distribution de la qualité des éléments pour chaque grille. La valeur maximale de Y + est toujours inférieure à 1.2.

  la région de la paroi proche de la couche limite et la région du courant libre. Les équations linéarisées ont été résolues au moyen de la méthode FGMRES et d'une factorisation incomplète de la partie inférieure supérieure avec un modèle épars basé sur la connectivité comme pré-conditionneur linéaire. Les propriétés des fluides sont estimées à l'aide de l'équation d'état de Peng-Robinson[START_REF] Peng | A new two-constant equation of state[END_REF], tandis que des valeurs constantes sont supposées pour les propriétés de transport.Dans le présent travail, l'accent est mis sur l'entropie, car il s'agit d'une quantité commode dans l'évaluation des performances en grille. En fait, d'après la thermodynamique, toute irréversibilité crée inévitablement une quantité d'entropie bien définie, ce qui réduit l'efficacité du processus. Comme le débit dans les turbomachines générales est étroitement adiabatique, l'entropie ne peut pas diminuer, mais elle peut seulement augmenter si une certaine dissipation est présente. De plus, contrairement à d'autres quantités comme l'enthalpie totale, la pression totale ou l'énergie cinétique, l'entropie ne change pas si le cadre de référence est stationnaire ou en rotation. Pour ces raisons, dans[START_REF] Denton | The 1993 IGTI Scholar Lecture: Loss Mechanisms in Turbomachines[END_REF] l'auteur conclut que "la seule mesure rationnelle de la perte dans une machine adiabatique est la création d'entropie" et il suggère de considérer l'entropie comme "la fumée" qui est créée dans le flux chaque fois qu'il se produit quelque chose de nuisible à l'efficacité, comme par exemple dans les couches limites des pales ou dans les ondes de choc, et dont la concentration à la sortie de la machine est proportionnelle

	J = s ¯out -s ¯in	(6.13)
	Comme la RDO exploite un optimiseur MOGA, la déformation maximale autorisée pour les points de
	contrôle FFD doit être spécifiée : dans l'optimisation actuelle, elle a été fixée à 20% dans le système
	de référence du treillis.
	Avant d'effectuer le RDO, une validation en deux étapes du modèle CFD a été effectuée. Tout
	d'abord, les données expérimentales des cas tests fiables de turbomachines [13] pour l'air sec dans
	des conditions de gaz idéales ont été comparées aux résultats d'une simulation avec le code SU2 : le
	Ces dernières années, plusieurs travaux ont proposé des analyses intéressantes de l'UQ dans les cas test sélectionné est le scénario d'écoulement subsonique MUR-43, où le nombre de Mach de sortie
	turbomachines pour les applications de la ORC [264, 265, 266]. Dans ces travaux, l'expansion du chaos isentropique est de 0,84, qui a été modélisé avec la même stratégie de modélisation présentée dans la
	polynomial ou le krigeage bayésien ont été utilisés pour évaluer l'influence des incertitudes aléatoires section 5.2.1. Le maillage adopté est représenté sur la figure 6.22.

Jameson-Schmidt-Turkel

[START_REF] Jameson | Origins and further development of the jameson-schmidt-turkel scheme[END_REF]

a été utilisé pour la discrétisation des flux inviscides convectifs, tandis qu'une méthode implicite d'Euler a été employée pour la discrétisation temporelle. Les contraintes turbulentes sont modélisées avec le modèle SST k-omega RANS

[START_REF] Menter | Two-equation eddy-viscosity turbulence models for engineering applications[END_REF]

, couramment utilisé dans les applications de turbomachines car il est capable de résoudre à la fois et épistémiques sur des quantités telles que le nombre de Mach de sortie, le débit massique, l'efficacité isentropique ou la perte de pression totale. à la perte rencontrée à l'intérieur de celle-ci. Par conséquent, ici, la différence entre l'entropie moyenne à la sortie du domaine et la valeur moyenne de l'entropie à l'entrée est considérée comme le QoI J pour l'UQ effectuée dans le RDO. La définition formelle de J est fournie dans l'équation 6.13. L'objectif est d'identifier les meilleures valeurs de déplacement pour les points de contrôle FFD en paramétrant le profil de l'aube qui minimise à la fois la valeur moyenne et la variance du QoI J.

Table 6 . 7 :

 67 6.22, tandis que les propriétés thermodynamiques réelles du gaz sont modélisées avec l'équation d'état de Peng-Robinson.Les résultats de cette deuxième étape de la validation sont représentés sur la figure6.24 : les valeurs du coefficient de pression C P = (P t,in -P )/(P t,in -P out ) dans le sens du courant calculées avec SU2 sont représentées par des carrés rouges et montrent un bon accord avec les résultats de[START_REF] Harinck | The influence of molecular complexity on expanding flows of ideal and dense gases[END_REF], Définition de la condition limite pour le point de simulation déterministe P T,in T T,in P out Simulation point 28.9 bar 277.3 degC 10.05 bar

	indiqués par des points bleus.
	Une première simulation CFD déterministe de la géométrie de base a été réalisée au point de
	fonctionnement défini par les conditions énumérées dans le Tab 6.7, correspondant aux moyennes
	des distributions indiquées dans le Tab. 6.5: le taux d'expansion de 2,752 et le nombre de Mach

isentropique attendu à la sortie est de 1,446. Les principaux résultats de cette grille fortement chargée sont présentés dans la figure 6.25, la figure 6.26 et la figure

6

.27.

R ÉSUM É ÉTENDU

Table 6 .

 6 9: Coût pour le calcul de l'optimisation DOE S DOE pour la méthode LF et la méthode HF UQ UQ Method UQ cost for a sample Size of S DOE Total cost of S DOE Time [h] Une fois que les DOE d'optimisation LF et HF ont été construites, le MOGA initialise aléatoirement une première population composée de 40 individus dans l'espace de conception qui est laissé à évoluer sur 20 générations. Au cours de ce processus, un remplissage adaptatif toutes les 5 générations est calculé avec la stratégie déjà expliquée dans la section 3.5.2. Au total, 3 MOEI sont réalisées, à la génération 6, à la génération 11 et à la génération 16 : la première, adoptant la technique MoM de l'UQ, a pris 0.8 heures, tandis que les autres ont nécessité 3.2 heures ensemble. Par conséquent, 1271 simulations CFD sont résolues pour l'ensemble du processus RDO, ce qui correspond à un temps de calcul de 67 heures, qui comprend la construction des DOE d'optimisation LF et HF, trois MOEI et l'évolution de la population initiale jusqu'à la génération 20. Il s'agit d'une grande amélioration par rapport à l'approche TNBK, puisque cette dernière aurait nécessité, pour la même RDO, 2 560 simulations CFD pour construire la DOE d'optimisation initiale et 96 autres pour les trois MOEI, avec un temps d'optimisation total estimé à 146.5 heures environ.

	LF	7 CFD	80 samples	560 CFD	28
	HF	32 CFD	20 samples	640 CFD	35

Table 6 .

 6 10: Moments statistiques de la PDF QoI pour la géométrie de base et l'optimum RDO

	Geometry E[J]	√	var[J]
	baseline	7.3275	2.8680
	optimum 5.5498	2.3339

Enfin, une simulation CFD de la géométrie optimisée est réalisée dans des conditions aux limites

  As all thermodynamic conditions are fixed at the nozzle inlet, so the boundary conditions ; consequently, solving θ T becomesthe differential problem defined by the linear system in Equation (A.18) with the boundary conditions provided in Equation (A.[START_REF] Guercio | 15 -biomass-fired organic rankine cycle combined heat and power systems[END_REF]).Once θ T has been determined, the variation of the cost function dJ with respect to the design parameters α can be easily calculated by considering just the variation of dS, as stated in Equation (A.20).

	A.6. CONTINUOUS ADJOINT METHOD (CAM)
	dJ dα	= -	∫ L 0	dθ T dx	f	dS dα	dx -	∫ L 0	dθ T dx	K	dS dα	dx +	[ θ T f	dS dα
	∫ L 0	(P -P des )	δP δα	dx -	∫ L 0	dθ T dx	S	df dα	dx +	∫ L 0	θ T ds dx	dK dα	dx +	[ θ T S	df dα	] L 0	= 0	(A.17)
																			[	df dα	]	are
																			0
	equal to zero∫ L 0	(P -P des )	δP δα	dx -	∫ L 0	dθ T dx	S	df dα	dx +	∫ L 0	θ T ds dx	dK dα	dx = 0	(A.18)
														[ θ T S	df dα	] x=L	= 0	(A.19)

Table A .

 A 1: Gradients from finite differences vs. gradients from the discrete adjoint for the baseline geometry.

	Quantity Finite Differences Discrete Adjoint Error %
	∂J ∂a ∂J ∂b ∂J ∂c ∂J ∂d	341296 -854462 -2473058 497319	341302 -854564 -2473429 497371	0.002% 0.012% 0.015% 0.010%

Table A .

 A 2: Gradients from finite differences vs. gradients from the continuous adjoint for the baseline geometry.

		Aldo SERAFINO
	Multi-fidelity Robust Design Optimization Methods for
		Organic Rankine Cycles
	Quantity Finite Differences Continuous Adjoint Error %
	∂J ∂a ∂J ∂b ∂J ∂c ∂J ∂d	341296 -854462 -2473058 497319	341313 -854609 -2473571 497401	0.005% 0.017% 0.021% 0.016%

(a) BK, 1 st generation (b) BK, generation after

Toutefois, comme nous l'avons déjà vu au chapitre 2, l'ensemble du processus de conception décrit dans la figure 6.7 peut être affecté par plusieurs sources d'incertitude, sur lesquelles il faut agir :• les variables environnementales du problème,• les paramètres de conception,• le modèle numérique approximatif utilisé pour simuler le système réel.En conséquence, le processus de conception déterministe qui vient d'être présenté n'est plus applicable. Pour surmonter cette difficulté, une technique classique consiste à considérer les valeurs moyennes pour les conditions limites et à appliquer ensuite lors de la conception certaines marges de sécurité ; une autre approche plus moderne consiste à passer du point de vue déterministe au point de vue stochastique et à effectuer la quantification de l'incertitude.En outre, la méthodologie standard de l'optimisation déterministe présente d'autres limites dont il faut tenir compte, car les systèmes optimisés au sens classique peuvent généralement être très sensibles à de petits changements[START_REF] Bufi | Robust optimization of ORC turbine expanders[END_REF]. Par conséquent, en raison de l'effet des erreurs de modélisation, il peut arriver que la solution optimale trouvée diffère de la vraie. Enfin, même dans le cas d'un modèle parfait dépourvu d'erreurs, la réalisation de la conception optimale prédite pourrait être affectée par toutes les sources d'incertitudes, tant aléatoires (c'est-à-dire les conditions environnementales et les tolérances géométriques) qu'épistémiques (par exemple, une connaissance imparfaite de la distribution

la MoM, mieux adaptée aux plages d'incertitude élevées et à la pdf générique, consiste à exploiter les informations de gradient pour construire un métamodèle de haute qualité à partir d'un nombre réduit d'échantillons. Une telle approche est utilisée par exemple dans les métamodèles GEK.La parallélisation massive est d'une grande aide pour accélérer le processus RDO[START_REF] Congedo | Shape optimization of an airfoil in a BZT flow with multiple-source uncertainties[END_REF][START_REF] Hercus | Robust shape optimization of uncertain dense gas flows through a plane turbine cascade[END_REF][START_REF] Cinnella | Robust optimization using nested Kriging surrogates: application to supersonic ORC nozzle guide vanes[END_REF], mais elle n'est pas applicable rapidement pour une utilisation industrielle de routine. Un moyen de réduire considérablement le nombre requis d'appels de fonctions consiste à remplacer les coûteux solveurs CFD ou UQ par des modèles substitut, tels que les fonctions de base radiales[START_REF] Cinnella | Optimal airfoil shapes for viscous transonic flows of bethe-zel'dovich-thompson fluids[END_REF], les réseaux neuronaux artificiels [158], et le krigeage[117,[START_REF] Cinnella | Robust optimization using nested Kriging surrogates: application to supersonic ORC nozzle guide vanes[END_REF], en approximant les variations des fonctions de coût à travers l'espace de conception. Une telle approche est appelée algorithme génétique multi-objectif basé sur des métamodèles (SMOGA).Enfin, des réductions supplémentaires du temps de calcul peuvent être obtenues en combinant des modèles avec différents niveaux de fidélité pendant l'optimisation.En suivant la méthode RDO de Taguchi, l'objectif est de rechercher une méthodologie qui permette d'optimiser un ensemble de QoI, J = J (x, ξ), J ∈ R m en fonction d'un vecteur de paramètres de

Pour le présent travail, la méthode de déformation libre (FFD)[255] a été utilisée comme stratégie de paramétrage de la géométrie : cette méthode permet de conserver un faible nombre de paramètres de conception tout en assurant la réussite de la modélisation des géométries complexes. L'inconvénient est qu'elle ne permet pas de contrôler certaines grandeurs techniques, comme par exemple l'épaisseur des pales, qui sont généralement traitées comme des contraintes d'optimisation. Malgré ce point, il est largement adopté dans l'optimisation de forme CFD[START_REF] Vitale | Multistage turbomachinery design using the discrete adjoint method within the open-source software su2[END_REF], en raison de sa flexibilité et de sa facilité d'utilisation. L'application de la FFD au problème actuel est présentée dans la Fig.6.17 : un réseau de 12 points de contrôle, à savoir les noeuds de A à L, est utilisé pour déformer la géométrie de base, le profil VKI LS89, au moyen d'un produit tenseur bivarié (puisque le problème est en 2D) du polynôme de Bernstein[255]. Pour respecter la contrainte géométrique autour de la corde axiale et pour éviter des déformations incontrôlables du bord de fuite (son épaisseur est considérée comme fixe, tandis que l'angle du cône est maintenu libre de se déplacer), les quatre coins (points A, F, G et L, qui sont colorés en noir) sont contraints et seuls les déplacements verticaux des 8 autres points de contrôle (colorés en rouge) sont autorisés.Cette paramétrisation conduit donc à un problème d'optimisation de la forme avec huit points de conception, qui sont les déplacements des huit points de contrôle FFD libres de se déplacer en

Acknowledgements

to Prof. Paola Cinnella for her guidance, patience and encouragement throughout the development of this project. of ENERTIME and the French organisation for research and technology ANRT.

ABSTRACT

List of Tables

List of Figures

SUMMARY OF THE CHAPTER

Summary of the chapter

In the present Chapter, a quick overview of methods for uncertainty quantification has been provided. Among all techniques introduced, the deterministic method of moments (MoM) and the stochastic Bayesian Kriging (BK) have been chosen for the analysis carried out in the following. Furthermore, a gradient-enhanced Kriging (GEK) has been proposed and assessed vs. BK on some analytical test functions and in the UQ of a simple ORC. SIMPLIFIED ORC Figure 4.9: Comparison of optimal solutions from both the strategies as a function of P ev The results are reported in Tab.4.5 and Fig. 4. [START_REF] Guercio | 15 -biomass-fired organic rankine cycle combined heat and power systems[END_REF]. The table shows that the first-order and total indexes are very close to each other, indicating that coupled effects due to the simultaneous variation of multiple inputs are small and that the system output is, to a good approximation, a linear function of the inputs. This is consistent with the Gaussian shape of the output PDF in Fig. 4.18 (a linear transformation of a Gaussian is a Gaussian). Fig. 4.19 highlights the relative influence of the input parameters on the variance of the QoI W ORC,net . The results show that the source temperature T geo is by far the most important influential parameter. The brine mass flow rate m ̇geo and the condenser temperature T cd also play a significant role, while the effect of c geo can be neglected.

RDO results

The NSGA II algorithm allows to find an approximation of the Pareto front of sub-optimal solution, corresponding to different compromise solutions among the various objective. For the present RDO problem, the objective space has three dimensions and a graphical representation of the Pareto-optimal solutions is difficult. the first two cost functions, as the optimizer finds some design maximizing the mean value of the PDF of W ORC,net , maximizing at the same time its variance. Such a trend is also observed in Fig. 4. [START_REF] Campana | Orc waste heat recovery in european energy intensive industries: Energy and ghg savings[END_REF], where it appears that the maximization of the expected net power output corresponds to higher values of A ORC . As a consequence, a unique optimum solution does not exist and different compromise solutions among the cost functions are possible. Among all the designs lying on the Pareto front, the one offering the best trade-off among all optimization targets is selected; this design is represented as a red symbol in Fig. 4.20 and Fig. 4.21 and it is identified by the parameters reported in Tab.4.6. The main characteristics of the final RDO design are listed in Tab.4.7.

The selected sub-optimal designs is input to the UQ solver, to verify the accuracy of the Kriging estimates of the cost functions by comparison with the full model. The distance between this solution and the line of the optimal solutions calculated with the TNBK approach, i.e. the surrogate approximation error, is found to be always lower than 5%. Despite this deviation, the selected design point APPLICATIONS • uncertainties in inlet total pressure P T,in described as a uniformly-distributed random variable with imposed lower and upper bounds;

• uncertainties in the gas properties, here represented by the specific heat ratio γ, which is also assumed as uniformly-distributed.

The characteristics of the random parameters are listed in Table 5.1. In the inverse design process, the uncertain geometric parameters are also (uncertain) design variables: for this reason, their mean is not fixed, but varies within ranges corresponding to the bounds of the design space. This means that, even for designs corresponding to the upper/lower bounds, a realization of the nozzle geometry may lie outside the prescribed limits, due to geometric tolerances.

Uncertainty quantification and RDO methods

For the present problem, two gradient-based UQ methods have been employed: the GEK and the MoM. As the test case is not particularly expensive to be simulated, a UQ solver based on the BK surrogate is also considered: this technique, which does not use gradient information, has been here considered as a reference.

BK is also used as the external surrogate model for the SMOGA, which employs the NSGA-II algorithm as the optimizer. The MOGA requires evaluations of the RDO cost functions through the UQ method for each individual and for each generation. To reduce the computational cost of the RDO, one possibility is to construct a unique BK or GEK surrogate over a parameter space extended to the design variables. Due to the curse of dimensionality, constructing a reliable surrogate over such a large parameter space may result in an unacceptable cost of the kriging approximation. To avoid

• τ is the viscous stress tensor that can be written in vector notation as:

• c p is the heat capacity at constant pressure, depending on ρ and e by means of an equation of state c p (ρ, e).

• µ v1 = µ dyn + µ turb .

• µ dyn is the fluid dynamic viscosity.

• µ turb is the turbulent viscosity, computed with a RANS turbulence model.

• P r d and P r t are the dynamic and turbulent Prandtl number, which is here considered equal to 0.9.

The governing equations are solved by using the SU2 code [START_REF] Economon | SU2: An Open-Source Suite for Multiphysics Simulation and Design[END_REF]: for the analyses presented in the following, the Jameson-Schmidt-Turkel central scheme [START_REF] Jameson | Origins and further development of the jameson-schmidt-turkel scheme[END_REF] has been used for the discretization of convective inviscid fluxes, while an Euler implicit method has been employed for time discretization. Turbulent stresses are modelled with the SST k-omega RANS model [START_REF] Menter | Two-equation eddy-viscosity turbulence models for engineering applications[END_REF], commonly used in turbomachinery application since it is able to solve both the near-wall region of the boundary layer and the free-stream region. The linearized equations have been solved by means of the FGMRES method and an incomplete lower upper factorization with connectivity-based sparse pattern as a linear pre-conditioner. Fluid properties are estimated with Peng-Robinson [START_REF] Peng | A new two-constant equation of state[END_REF] equation of state, while constant values are assumed for transport properties.

UQ techniques and RDO strategy

In the last years, several works have proposed interesting UQ analysis in turbomachinery for ORC applications [START_REF]Sensitivity of Supersonic ORC Turbine Injector Designs to Fluctuating Operating Conditions[END_REF][START_REF] Razaaly | Impact of geometric, operational, and model uncertainties on the non-ideal flow through a supersonic orc turbine cascade[END_REF][START_REF] Romei | Amplification of Operational Uncertainty Induced by Nonideal Flows in Supersonic Turbine Cascades[END_REF]. In these works, polynomial chaos expansion or Bayesian Kriging have been used to evaluate the influence of aleatory and epistemic uncertainties on quantities like the exit Mach number, mass flow rate, isentropic efficiency or total pressure loss.

RDO OF AN ORC TURBINE BLADE

Figure 5.17: Validation of the CFD model with test case [START_REF] Harinck | The influence of molecular complexity on expanding flows of ideal and dense gases[END_REF] for propane blade axial chord calculated with SU2 are depicted as a dashed red line and they fit quite well the experimental values, reported as blue dots.

After this first step of the validation, the CFD modelling strategy has been further tested on a case study which is closer to the problem object of the RDO that is here carried out: the inviscid flow of gaseous propane through the LS89 cascade. For this case, the RANS SU2 results can be compared with inviscid flow simulations available from [START_REF] Harinck | The influence of molecular complexity on expanding flows of ideal and dense gases[END_REF]. For the SU2 calculation, computational grid is the mesh of Fig. 5.15, while the real gas thermodynamic properties are modelled with the Peng-Robinson equation of state.

The results of this second step of the validation are plotted in Fig. 5.17: the values of the pressure coefficient C P = (P t,in -P )/(P t,in -P out ) along the streamwise calculated with SU2 are depicted as red squares and they show a good agreement with the results from [START_REF] Harinck | The influence of molecular complexity on expanding flows of ideal and dense gases[END_REF], reported as blue dots. La figure 6.9 illustre une géométrie de tuyère typique. L'objectif de l'optimisation est de déterminer les paramètres de conception a, b, c, d en Eq. 6.4

Analysis of the baseline geometry

fournissant le meilleur ajustement à la distribution de pression cible sous des incertitudes multiples, au sens de l' Eq.(6.1).

On suppose que le flux est régi par les équations d'Euler pour les flux de quasi-1D (Eq. 6.6) :

Gradient calculation from codes

As already explained in Ch.2 , GEK and MoM are interesting and powerful techniques for UQ, but they require as input also the gradient of the QoIs with respect to the uncertain parameters. Because of the fact that in the present work most of the QoIs are outcomes of simulation solvers, it can be worth providing an overview of the main methodologies that nowadays are available to calculate derivatives of quantities from numerical codes.

This topic has a long history, dating back to the end of the '80s of the last century [START_REF] Sobieszczanski-Sobieski | The case for aerodynamic sensitivity analysis[END_REF] [268]

[269]: in fact, as soon as CFD started to be considered a reliable tool for engineering design process, the whole CFD community invested a lot of efforts in efficient and accurate evaluation of sensitivities because of the huge interest in using them for optimization purposes.

Despite all the remarkable progress over almost 30 years in scientific calculation, the implementation of techniques allowing the calculation of gradients in CFD codes still remain an open question and the reasons of this have probably to be found in two major causes [START_REF] Peter | Numerical sensitivity analysis for aerodynamic optimization: A survey of approaches[END_REF].

The first is likely to be the questionable suitability of gradient-based optimization methods for many engineering problems, as by nature they are all local methods; thus, in a design space with many local optima, instead of searching the best global optimum, they will tend to find the nearest local optimum with respect to the starting point. This argument is directly connected to another discussion that is still ongoing in the scientific community for decades about the performance of gradient-based methods and non-deterministic ones in optimization for engineering purposes. The second critical issue is the lack of cheap and reliable methods for an accurate evaluation of gradients.

As a consequence of this, just in the last couple of years adjoint methods for the calculation of A.2. FINITE DIFFERENCES of the sensitivities of J after that the governing equations have been linearised; thus, the dominant cost in this case is the solution of N α linear equations. On the other hand, the AM for each element of vector J a dual linear problem is constructed: by solving the original primal problem and the dual linear one (both with dimension N ω ), it is possible to obtain a cheap and accurate evaluation of sensitivities of J with respect to any design variable. The cost of evaluating the full derivative matrix dJ dα is therefore dominated by the cost of N j dual linear problems. As problems concerning aerodynamics are typically characterized by a low number of cost functions J and a large number of design variables α, the adjoint approach is almost always preferable.

Another classification of sensitivity evaluation techniques is possible: there are continuous methods, where the governing equations are linearized and, if it is the case, adjointed before being discretized, and discrete methods, where the linearization and adjointing is performed on the discretization of the non-linear system [START_REF] Peter | Numerical sensitivity analysis for aerodynamic optimization: A survey of approaches[END_REF].

In the following sections of the present chapter an overview of some techniques for the evaluation of gradients dJ dα is provided; these are:

• finite differences

• discrete adjoint method (DAM)

• continuous adjoint method (CAM).

The continuous direct method, even if it has been used for evaluating gradients for the application of GEK in the test case about the UQ of the ORC (Par. 2.6.2), is not considered useful for QoIs dealing with CFD and it is therefore here neglected.

A.2 Finite Differences

Finite differences is probably the easiest method to be applied because it is not intrusive as it does not require any change in the solver code. The basic idea relies on the approximation of derivatives with the difference quotients (see for instance [START_REF] Quarteroni | Numerical Mathematics[END_REF] or any book about calculus and numerical methods).

Several schemes having different orders are available [START_REF] Jovanović | Analysis of Finite Difference Schemes: For Linear Partial Differential Equations with Generalized Solutions[END_REF] [276]: usually the higher is the order of the scheme, the larger is the number of solver code solutions that it is required.

A.2. FINITE DIFFERENCES

For instance, considering a central finite difference second order scheme, to calculate the derivatives dJ dα it is necessary to perturb the mesh geometry twice:

1. the first time is to calculate the flow solution in the perturbed states J(ω(α + δα), α + δα)

satisfying the governing equations Eq.A. [START_REF] Schlesinger | Terminology for Model Credibility[END_REF].

2. the second time is to calculate the flow solution in the perturbed states J(ω(α -δα), α -δα)

satisfying the governing equations Eq.A.5.

Once that the two perturbed states have been calculated, the derivatives of the objective function J can be obtained by means of Eq. A.6

Without considering the evaluation of the objective function J(ω(α), α), the entire matrix of the gradient dJ dα may be evaluated at a cost of 2×N α , making the method not feasible for high-dimensional design spaces.

Finite Differences method allows also to get an estimation of the approximation error of the derivatives; for instance, the error ε between the real derivatives of J(ω(α), α, ) and the approximation with the second order scheme presented in Eq.A.6 is provided by Eq.A.7.

Eq.A.7 shows a serious problem of the Finite Difference Method: in order to get a good approximation of the derivatives, it is desiderable to choose a small value of for the step size ∥δα∥, to minimize truncation errors; unfortunately, as arithmetics in computers has a finite precision, too small values of ∥δα∥ should not be used in order to avoid that rounding errors become significant, leading to subtractive cancellation errors. Thus, step size ∥δα∥ is critical to the accuracy of the result and the only way A.3. DISCRETE DIRECT METHOD to guarantee accuracy is to perform a convergence study on ∥δα∥ for each parameter, at considerable cost [START_REF] Peter | Numerical sensitivity analysis for aerodynamic optimization: A survey of approaches[END_REF]. As most issues come from the operation of adding and subtracting in Eq.A.6 numbers whose values are too close, an improvement can be obtained by employing the Complex Step Derivation method [277] [278]. Applications of this method can be found in [START_REF] Squire | Using Complex Variables to Estimate Derivatives of Real Functions[END_REF], [START_REF] Anderson | Sensitivity Analysis for Navier-Stokes Equations on Unstructured Meshes Using Complex Variables[END_REF], [START_REF] Martins | The connection between the complex-step derivative approximation and algorithmic differentiation[END_REF],[282], [START_REF] Martins | A Coupled-Adjoint Sensitivity Analysis Method for High-Fidelity Aero-Structural Design[END_REF].

However this method cannot be considered not intrusive as usually for studies related to CFD the solver must be modified to deal with complex variables.

A.3 Discrete Direct Method

The Discrete Direct Methods differs from the Finite Differences one: starting from problem stated in Eq.A.3, it aims to calculate the gradient dJ dα directly. In order to do that, partial derivatives ∂J ∂α and ∂J ∂ω can be usually obtained by hand on by mean of Algorithmic Differentiation (AD) [START_REF] Griewank | Evaluating Derivatives: Principles and Techniques of Algorithmic Differentiation[END_REF], [START_REF] Forth | Recent Advances in Algorithmic Differentiation[END_REF], [START_REF] Hascoët | The Tapenade Automatic Differentiation tool: Principles, Model, and Specification[END_REF].

In Eq.A.3 the most critical element to be evaluated is the term dω dα ; the Discrete Direct Method propose a solution to this problem that uses the governing equation (Eq.A.2) as an ancillary equation:

under the assumption that the discrete form of R is at least once continuously differentiable with respect to the flow field and the geometry in a neighbourhood of (ω(α), α), it is possible to differentiate it with respect to α (Eq.A.8).

dR(ω(α), α) = [ ∂R ∂ω

]

where notations I and II are still used with the same purposes already explained in section A.1.

By solving the linear system in Eq.A.8, it is possible to calculate dω dα as in Eq.A.9.

Usually it is rare that R is differentiable everywhere [START_REF] Peter | Numerical sensitivity analysis for aerodynamic optimization: A survey of approaches[END_REF], but in practice this rarely causes difficulties [START_REF] Matsuzawa | Optimum shape design using adjoint equations for compressible flows with shock waves[END_REF]. Concerning computational efforts, this technique can be considered as an improvement, if compared with the Finite Differences Method; however, the calculation of Eq.A.9 still requires considerable resources, as a series of N α decoupled linear systems with N ω dimension has to be solved.

For this reason, the Discrete Direct method for CFD applications was mainly pursued in the past A.4. ADJOINT METHODS ( [START_REF] Pini | Turbomachinery design optimization using adjoint method and accurate equations of state[END_REF], [START_REF] Bristow | Subsonic panel method for the efficient analysis of multiple geometry perturbations[END_REF], [START_REF] Bristow | Subsonic panel method for designing wing surfaces from pressure distribution[END_REF], [START_REF] Elbanna | Determination of aerodynamic sensitivity coefficients in the transonic and supersonic regimes[END_REF], [START_REF] Baysal | Aerodynamic design optimization using sensitivity analysis and computational fluid dynamics[END_REF], [START_REF] Taylor | Aerodynamic Shape Sensitivity Analysis and Design Optimization of Complex Configurations Using Unstructured Grids[END_REF], [START_REF] Newman | Overview of Sensitivity Analysis and Shape Optimization for Complex Aerodynamic Configurations[END_REF]) before the development of the adjoint methods.

A.4 Adjoint Methods

As aerodynamic problems usually involve a great number of design variables, the adjoint approach has already proved that it is more efficient than other strategies to compute the gradient of the objective function. In fact, in this case the total cost for the gradient computation does not depend on the number of design variables and it is approximately equal to that of solving the flow equations and then they are discretized; at the end, the outcome of both the processes is a linear system of equations that can be easily numerically solved to compute the field of the adjoint variables. Even if both approaches have their advantages and disadvantages in terms of requirements at the development stage, they are able to reach the same level of accuracy in the calculation of sensitivity derivatives [START_REF] Mohammadi | Shape Optimization in Fluid Mechanics[END_REF], [START_REF] Nadarajah | A Comparison of the Continuous and Discrete Adjoint Approach to Automatic Aerodynamic Optimization[END_REF].

Appendix B

Global sensitivity analysis with Sobol indexes

In this work glabal sensitivity analysis has been performed by means of the calculation of Sobol Indexes. More specifically, for each input parameter the first-order and total Sobol indexes are extracted and used to evaluate the influence of the uncertain input variables on the QoI.

One should consider the function J = J(x), where x = x 1 , x 2 , ..., x N is a vector containing N random variables described by known PDF. The main idea of this method is to decompose the output variance V into the contributions associated with each input factor.

For a single variable x i , the first-order Sobol index S i is calculated as follows:

This index measures the effect on the variance of the output function J when perturbing each input parameter individually.

The total Sobol index S T i is calculated as below:

where x ˜i indicates the set of all random variables except x i .

The total Sobol index S T i represents the overall contribution of a given input to the total variance of the output, when such a parameter is perturbed individually or in conjunction to other parameters. By construction, the Sobol indexes are numbers comprised between 0 and 1. Higher values indicate a stronger sensitivity of the QoI to a given parameter.

In the present work Sobol indexes have been calculated by means of a Monte Carlo sampling.

The reader in referred to [START_REF] Sobol | Global sensitivity indices for nonlinear mathematical models and their monte carlo estimates[END_REF][START_REF] Saltelli | Global Sensitivity Analysis: The Primer[END_REF] for more details on Sobol global sensitivity analysis.