domain: spi info

text:
Jiayao Li
email: jiayaol@bjtu.edu.cn

Frédérick Bénaben
email: frederick.benaben@mines-albi.fr

Juanqiong Gou
email: jqgou@bjtu.edu.cn

Wenxin Mu
email: wxmu@bjtu.edu.cn

A Proposal for Risk Identification Approach in Collaborative Networks Considering Susceptibility to Danger

Keywords: Risk identification, DRC chain, Susceptibility to danger, Cascading effect, Collaborative network

The paper proposes a research framework for risk identification approach in collaborative networks dedicated to develop a formalizing, structured reference for risk identification and risk mitigation and explore an effective mechanism that can motivate diverse partners to manage risks collaboratively. The approach is based on a formalized vision of Danger/Risk/Consequence chain that is defined as the primary schema of the proposed methodology. The DRC chain indicates five risk-related concepts and their interrelationships, which is able to well describe risk-related collaborative contexts. Cascading effect in the concept chain are presented for further interpreting. Furthermore, a supply chain scenario of three use cases is given to illustrate the proposed framework.

Introduction

A collaborative network is an alliance constituted by a variety of entities (e.g. organizations and people) that are largely autonomous, geographically distributed, and heterogeneous in terms of their operating environment, culture, social capital and goals, but that collaborate to better achieve common or compatible goals, and whose interactions are supported by a computer network [START_REF] Camarinha-Matos | On reference models for collaborative networked organizations[END_REF]. Collaborative networks such as virtual organizations, dynamic supply chains, professional virtual communities, collaborative virtual laboratories, etc. are complex systems associated with uncertainties in dynamic business environments [START_REF] Jamshidi | A New Decision Support Tool for Dynamic Risks Analysis in Collaborative Networks[END_REF]. It is noted that the collaboration increases the dependences among enterprises, which makes enterprises more susceptible to risks. It is critical for collaborative networks to take systematic approaches to identify risks as early as possible, and implement appropriate strategies to manage the risk propagation throughout the evolution of collaboration [START_REF] Wulan | A fuzzy logic based system for risk analysis and evaluation within enterprise collaborations[END_REF].

While the above discussion shows the importance of risk identification, however, it is worthy to note that there appears to be no overarching typology to delineate exactly what constitutes risk and how to understand risk [START_REF] Rao | Supply chain risks: a review and typology[END_REF]. [START_REF] Juttner | Supply chain risk management -outlining an agenda for future research[END_REF] argues that risk management consists of four key management aspects: (i) assessing the risk sources; (ii) defining the adverse consequences; (iii) identifying the risk drivers; and (iv) mitigating risks. [START_REF] Zeng | Rethinking the role of partnerships in global supply chains: a risk-based perspective[END_REF] indicates that risk is at least made up of three essential components: (i) a driver or drivers which trigger the risk to happen; (ii) an event with probability that signifies the occurrence of the risk; and (iii) a consequence or consequences resulted from the risk. [START_REF] Bénaben | Collaborative Systems in Crisis Management: A Proposal for a Conceptual Framework[END_REF] presents a three-dimensional framework dedicated to structure the domain of crisis management based on the Danger/Risk/Consequence chain (DRC chain). In summary, the risk-related concepts could be concluded as danger, stake, risk, event and consequence.

Anyway, risk management is a very complex domain with a lot of constraints. Consequently, it is very difficult to get a global overview of such a domain. This article is mainly dedicated to present a proposal for risk identification approach to support collaborative networked organizations. The proposed approach is based on the DRC chain (as shown in Fig. 1), which is not so far from FTA (Fault Tree Analysis) principles [START_REF] Vesely | Fault Tree Handbook[END_REF]. Furthermore, danger, risk and consequence may be considered as causal sources (in a waterfall structure) that must be formalized as models to help decision makers [START_REF] Bénaben | Collaborative Systems in Crisis Management: A Proposal for a Conceptual Framework[END_REF].

More specifically, the research methodology is proposed:

• A danger typology and a stake typology are developed for collaborative networked organizations based on the review of present related literatures.

•

The interconnections rules between dangers and stakes and their impacts on risks could be summarized, which aims to build a risk typology as a static reference model of risk knowledge.

•

In addition to danger, stake and event, risk and consequence can be also considered in a cascading effect structure [START_REF] Turoff | A collaborative dynamic scenario model for the interaction of critical infrastructures[END_REF]. Their interrelationships could be analyzed and summarized.

•

A risk knowledge base with a risk typology and interrelationships is used for the deduction of risk identification rules based on previous research works.

•

A metamodel connected with a risk knowledge base will be defined based on the research of collaborative situation metamodel [START_REF] Bénaben | A Metamodel for Knowledge Management in Crisis Management[END_REF][START_REF] Lauras | A meta-ontology for knowledge acquisition and exploitation of collaborative social systems[END_REF]. The research should contribute to a deeper and broader understanding risk based on the schema of DRC chain. Besides, the risk identification research based on such an understanding will contribute to better practices by suggesting collaborative responses from multiple partners in collaborative networks.

This article is structured according to the following sections: Section 2 presents a literature review; section 3 describes the DRC chain considering susceptibility to danger and interrelationships of the five risk-related concepts; in section 4, a supply chain scenario is given to illustrate the part of proposed methodology; section 5 concludes this research work and gives some perspectives for future works.

Literature Review

According to [START_REF] Edwards | Risk management in project organisations[END_REF], risk is the combination of the frequency, or probability, of occurrence and the consequence of special hazardous events. Risk Management corresponds to a set of activities that organizations use to control the many risks, which may undermine their ability to achieve objectives. Considering international standards on entrepreneurial management process risk, [START_REF] Rosas | Exploratory Study on Risk Management in Open Innovation[END_REF] organizes a reference risk taxonomy, which shows that risk management includes two principal dimensions, namely, assessment and treatment. Risk assessment could be summarized as three phases: analysis, identification and evaluation. Risk identification can be defined as a process of identifying the dangers, events and consequences.

It is a common assumption that participation in a collaborative network has the potential of bringing benefits to the involved partners. However, the interconnections between partners of collaborative network cause numerous new risks, of which the impacting magnitude and scope are larger than ever before [START_REF] Zeng | Rethinking the role of partnerships in global supply chains: a risk-based perspective[END_REF]. In recent years, there are few studies about the risk identification of collaborative networks by considering susceptibility to danger existing in the literature of collaborative networks and further research in this field is required.

Regarding the approach of risk identification in collaborative context, it focuses on literature review, semi-structured interviews and questionnaire. They can be called as the experienced-based methods. [START_REF] Fang | Network theory-based analysis of risk interactions in large engineering projects[END_REF] directly takes advantage of the risks identified by project manager to determine the events with negative impacts. [START_REF] Yang | Stakeholder-associated risks and their interactions in complex green building projects: a social network model[END_REF] identifies the stakeholder-associated risks through the previous risk identification literatures, and classified them into seven categories. [START_REF] Zhang | Risk identification and evaluation of customer collaboration in product development[END_REF] identifies the risks caused by customer collaboration in product development through relevant literatures by domain experts and questionnaire in the enterprises. [START_REF] Yang | Modelling stakeholder-associated risk networks in green building projects[END_REF] analyzes the research paradigms regarding risk and stakeholder analysis in green buildings through literature review. [START_REF] Zhao | A Literature Review on Risk Sources and Resilience Factors in Agri-Food Supply Chains[END_REF] undertakes a systematic literature review on risks sources and resilience factors in agri-food supply chains. [START_REF] Wang | Human safety risks and their interactions in china's subways: stakeholder perspectives[END_REF] successively uses literature review, semi-structured interviews and questionnaires to identify a list of human safety risk factors and also the cause-effect relationships among those risks.

Regarding the application of risk identification results, most results of identified risks are used to further risk evaluation and propose risk response or risk mitigation strategies. In order to investigate those risk interactions, the focus of risk evaluation methods has gradually been shifting from individual risks to networks of risk [START_REF] Zeng | Rethinking the role of partnerships in global supply chains: a risk-based perspective[END_REF], [START_REF] Wang | Human safety risks and their interactions in china's subways: stakeholder perspectives[END_REF]. It considers nodes in the network and their relationships, focusing on the structure and patterning of these relationships and seeking to identify both their causes and effects [START_REF] Scott | Social network analysis: A handbook[END_REF][START_REF] Zhou | Risk evaluation of dynamic alliance based on fuzzy analytic network process and fuzzy TOPSIS[END_REF].

We conclude the following:

• A systematic approach to the identification and categorization of risks in collaborative context is lacking.

•

The current risk identification methods mainly focus on review, expert interview and questionnaire.

•

More future works are attentive to the identification of risk interconnections.

•

The present research of the application of risk identification result also lacks of a sharing-based risk response mechanism considering capabilities and resources of partners to contribute to their collaboration.

Understanding Risk

Risk can be seen as combination of the probability of an event and its consequence. However, danger and stake are also closely related to risk with the exception of the concepts of event and consequence. DRC chain is a concepts schema that is able to describe risk-related contexts.

General Illustration of DRC Chain

In this schema, the five risk-related concepts could be defined as follows:

• Danger can be defined as any specific dangerous characteristic of the environment, which is a signal word used to indicate an imminently hazardous situation [START_REF]Graphical symbols -Safety colours and safety signs[END_REF].

•

Stake or assets can be seen as item, thing or entity that has potential or actual value to an organization [23] and potential susceptibility to dangers.

•

Risk is a potential manifestation of the danger onto some concerned stakes [START_REF] Bénaben | Collaborative Systems in Crisis Management: A Proposal for a Conceptual Framework[END_REF].

•

Event is defined as a change or outcome that triggers risks. If one risk might occur it would be due to some events [START_REF] Lahmar | Towards an Integrated Model of Supply Chain Risks: An Alignment Between Supply Chain Characteristics and Risk Dimensions[END_REF].

•

Consequence generally means a set of negative impacts of the risk occurrence.

The general illustration of DRC chain could be described as follows: Each of those negative facts is due to one (or several) event(s) that trigger(s) one (or several) risk(s); This (or these) risk(s) occur(s) because the considered area/system is concerned by one (or several) danger(s) that affect(s) one (or several) stake(s) [START_REF] Bénaben | Collaborative Systems in Crisis Management: A Proposal for a Conceptual Framework[END_REF].

Furthermore, the DRC chain could indicate the susceptibility to danger, which means the state of being very likely to be influenced or affected by danger. The following example in enterprise collaborative context illustrates susceptibility to danger (see Fig. 2).

Company C is the only one able to produce Product P for the core company CC, which is a danger for the stake CC because its Product PP must be produced by using P. Consequently, the risk is that PP may not be produced. It would be triggered if C decides not to produce P anymore (an event occurs), then cause a consequence as which CC cannot sell PP to its customers. For this general illustration of DRC chain, the demand side CC is susceptible to the danger of which there is no alternative provider regarding its required product. Different stakes have different degrees of susceptibility to danger while some stakes are not susceptible to danger.

Danger:

Company C is the only one able to produce Product P.

Stake:

Core company CC (using P to produce Product PP)

Event:

C decides not to produce P anymore.

Risk: PP may not be produced.

Consequence:

CC can't sell PP to its customers.

Fig. 2.

An illustration of "susceptibility to danger" in DRC chain

Interconnections in DRC Chain

In the schema of the DRC chain, risk is created by danger and stake while consequence is created by risk and event. It can be seen that risk and consequence are the "generated" elements. Therefore, they are the ones which directly impact the considered system. Accordingly, it is worth to focus on what effects that they might create. One of perspective that can indicate it is cascading effects that could be described as multiple connections initiated by risk and consequence in the DRC chain.

A cascading effect is an inevitable and sometimes unforeseen chain of events due to an act affecting a system [START_REF]A Dictionary of Ecology[END_REF]. If there is a possibility that the cascading effect will have a negative impact on the system, it is probable to analyze the effects with a risk or consequence analysis. Fig. 3 shows six connections initiated from consequence (see connections (1)(2)(3)) and risk (see connections (4)(5)( 6)), in order to present the cascading effect as follows: consequence and risk are regarded as the causal sources that target at danger, stake and event.

Fig. 3. Interconnections in DRC chain

To illustrate the cascading effect further, three use cases are presented in Fig. 4, which are the real ones from a cosmetic French company. The interconnections between the first two use cases (see connections [START_REF] Camarinha-Matos | On reference models for collaborative networked organizations[END_REF](2)(3) of Fig. 4) could indicate consequence as the causal source. The first use case is already described in previous section, the consequence of which is that CC cannot sell PP. It is worth mentioning that this consequence could bring about the next use case.

The second use case is described as follows: CC has to contract with another company C' that can produce P, however, which is a danger for the new provider C'. Concerning this stake, the danger could be manifested as a risk that its Product P may decline in quality. Under this circumstance, if there is a big demand for P from CC, the risk would be triggered and cause a consequence as which C' produces lower quality products and its image might be degraded.

In summary, it can be seen that one consequence in a scenario could create a new danger (CC contracts with C'), also a new stake (company C') and even an event (big demand from CC), which lead to another risk-related scenario.

The third use case is described as follows: A huge workload of employees in Company C' is a danger for the stake C', because it could create a risk that the employees may go on strike for salary increase. The risk would be triggered if C' requires them to work overtime (an event), then a consequence of employees' strike might be caused.

Regarding risk as the causal source, the creation of danger, stake and event is generally due to the actions of risk prevention and mitigation. The interconnections between these three use cases (see connections (4)(5)(6) of Fig. 4) could indicate it as below:

• Risk-danger: C' is confronted with a risk that its products may decline in quality in the second use case. In order to prevent the risk, C' might increase the workload of its employees to ensure the products' quality, which is a new danger shown in the third use case.

•

Risk-stake: CC needs to find the other company that can also produce its required Product P in order to prevent the risk in the first use case (PP may not be produced). As a result, the risk creates a new stake Company C' that is shown in the second use case.

•

Risk-event: For the sake of qualified quality of Product P, C' might requires its employees to work overtime. Consequently, the risk in the second use case (P by C' may decline in quality) creates an event in the third use case.

means that the interrelationship does not exist between the concepts. To conclude, we hold that in considered scenario consequence or risk could create new danger, stake and event, and update or delete the current state or form of danger and stake. Nevertheless, event could not be updated or deleted because we cannot change what has happened before. 4)Yes ( 4)Yes ( 5)Yes ( 5)Yes ( 5)Yes ( 6)Yes ( 6)No ( 6)No

Supply Chain Scenario Illustration

Supply chain is considered as a specific form of collaborative network. It is a stable long-term network of enterprises each having clear roles in the value chain, covering all steps from initial product design and the procurement of raw materials, through production, shipping, distribution, and warehousing until a finished product is delivered to a customer [START_REF] Camarinha-Matos | Collaborative Networked Organization -Concepts and Practice in Manufacturing Enterprises[END_REF].

In order to further illustrate the proposed DRC chain, three use cases of the supply chain scenario (presented in Fig. 5) would be given to be used to perform the progress of current work. Two partners are involved in this simple scenario: the core enterprise is the demand side, which submits orders to buy its required materials from the suppliers. Furthermore, the three use cases are described as follows:

• The first use case: Labor strike of suppliers is a danger for the core enterprise (a stake), because it would create the risk of the shortage of its required products. If the production disruption (an event) happens, the risk could be triggered. Then the consequence of overdue delivery cannot be avoided.

•

The second use case: Shortage of products required by the core enterprise is a danger for the supplier who provides them (a stake). The risk is that the core enterprise might give a negative feedback in the evaluation of this supplier, and the risk occurs if overdue delivery from the supplier (an event) happens for several times. The core enterprise might consider to change the supplier to ensure its normal operation of business (a consequence).

•

The third use case: Overdue delivery of supplier for long time is a danger for the core enterprise (a stake). It would create a risk that the core enterprise has to change the supplier. If many negative comments from the evaluation of the supplier (an event) are given to the decision-maker, the risk could be triggered. Consequently, the core enterprise needs to reselect the other suppliers to replace the tasks of the original supplier as soon as possible (a consequence).

By comparing with the three use cases (presented in Fig. 5) and referring to the cascading effect in DRC chain (presented in Fig. 3), it can be seen that a target or outcome in one situation could be the causal source to connect another situation. Accordingly, several interconnections could be concluded:

Conclusion

The presented framework of risk identification approach might be considered as a formalizing reference dedicated to identify and mitigate risk in collaborative networks. The proposed approach is compliant with the schema of Danger/Risk/Consequence chain that helps to formalize the risk-related knowledge that includes five concepts (danger, stake, risk, event and consequence) and their interrelationships. Cascading effect could be indicated in DRC chain, which could contribute to a deeper understanding of risk-related collaborative contexts. Besides, the devised danger typology and stake typology can be further used to develop a risk knowledge base for risk identification of collaborative networks. Furthermore, the risk knowledge base combined with the current methodology of metamodeling could contribute to further explore the ways in which an effective mechanism that can motivate diverse partners in collaborative networks to manage risks collaboratively. The future works would use System Dynamics [START_REF] Forrester | System dynamics, systems thinking, and soft OR[END_REF] to develop the proposed approach, which focus on: (i) developing a danger typology and a stake typology by synthesizing the growing diverse literatures; (ii) the deduction of interconnection rules of dangers and stakes in order to build a risk knowledge base for collaborative networks; (iii) developing risk identification rules based on the risk knowledge base; (iv) metamodeling with risk knowledge base that can be dedicated to support collaboration of partners, and deduce the collaborative processes of risk mitigation.

Fig. 1 .

 1 Fig. 1. A framework for risk identification approach

•Fig. 5 .

 5 Fig. 5. Use cases illustration of supply chain scenario

Table 1 .

 1 Interrelationships in DRC chain

		Danger	Stake	Event
		Create Update Delete Create Update Delete Create Update Delete
	Consequence (1)Yes (1)Yes (1)Yes (2)Yes (2)Yes (2)Yes (3)Yes (3)No (3)No
	Risk	(4)Yes (

Corresponding author.

Acknowledgments. The presented research works have been supported by "the Fundamental Research Funds for the Central Universities of China". The authors would like to thank the project partners for their advice and comments.

Danger:

Company C is the only one able to produce Product P.

Stake:

Core company CC (using P to produce Product PP)

Event: C decides not to produce P anymore.

Risk: PP may not be produced.

Consequence:

CC cannot sell PP to its customers.

Danger: CC contracts with Company C' to produce P.

Stake: Company C'

Event:

Big demand for P from CC Risk: P by C' may decline in quality.

Consequence: C' produces P with lower quality.

(1) create ( 2

Risk:

The employees may go on strike for salary increase.

Consequence:

Employees strike (6) create ( 5) create (4) create

Fig. 4. An illustration of cascading effect in DRC chain

In addition to the interrelationship "create" between consequence/risk and danger/stake/event, there are some other interrelationships presented in Table 1. We define the interrelationships in DRC chain as follows:

• "Create": Consequence/risk makes new danger/stake/event happen. • "Update": Consequence/risk makes danger/stake/event from one state or form into another. • "Delete": Consequence/risk removes or makes danger/stake/event invisible. "Yes" in Table 1 means that we have found the use cases to support this type of interconnection (it cannot be all shown in this article). Accordingly, "No" in Table 1


Wouter N Edeling

Martin Schmelzer

Richard P Dwight

Paola Cinnella

Bayesian Predictions of Reynolds-Averaged Navier-Stokes Uncertainties Using Maximum a Posteriori Estimates

published or not. The documents may come   L'archive ouverte pluridisciplinaire

I. Introduction

Direct numerical computation of turbulent flow fields is computationally intractable for most applications. Several levels of approximation are possible, according to flow scales that are resolved/modelled. In all cases, scale separation introduces unclosed terms that need to be modelled.

The choice of the appropriate modelling level remains essentially a matter of expert judgement.

In particular, it will always depend on cost versus accuracy considerations. On the other hand, even once a given level has been selected, the simplifying assumptions introduce unclosed terms for which several possible sub-models may be designed. These may differ both in their mathematical structure and their closure parameters. The common practice in turbulence modelling is often to leave the choice of a specific model structure to expert judgement, while treating model constants as adjustable parameters that are calibrated in such a way as to reproduce simple, well-documented flows. Both of the preceding aspects, however, represent sources of uncertainty in the prediction of a new flow.

While more general modelling approaches like Large Eddy Simulation (LES) have made considerable progress in the last decade, turbulence models based on Reynolds Averaged Navier Stokes (RANS) will remain the work-horse tool for engineering design and optimization in the decade to come. An extremely large variety of RANS models have been proposed in the past, ranging from simple algebraic models to sophisticated Reynolds Stress models (see [START_REF] Wilcox | Turbulence Modeling for CFD[END_REF] for a review). However, in order to arrive at any given closure model, numerous assumptions in the mathematical derivation and subsequent calibration must be made. No universally accepted and valid model has been identified in CFD literature, meaning that the appropriateness of a given model structure is highly uncertain when predicting a new configuration. Additionally, RANS models are developed and calibrated on simple flows using noisy experimental data, and applied to radically different flows with the same coefficients. Furthermore, since these calibrations are deterministic (and as experimental uncertainty is generally ignored), one obtains point estimates for the closure coefficients. This ultimately leads to predictions which are subject to unknown degrees of uncertainty, especially if the flow configuration of interest is far removed from the flow scenarios under which the model was calibrated.

The above-mentioned uncertainties in RANS models can be addressed by means of statistical tools. All such tools are meant to quantify the uncertainty in RANS closures, and should not be considered as a replacement for deterministic simulations. In the context of uncertainty quantification (UQ), uncertainties associated with the closure coefficients are referred to as parametric uncertainties. These can be effectively addressed by replacing standard deterministic calibration with statistical calibration, see e.g. [START_REF] Cheung | Bayesian uncertainty analysis with applications to turbulence modeling[END_REF][START_REF] Guillas | Bayesian calibration of the constants of the k-ε turbulence model for a CFD model of street canyon flow[END_REF][START_REF] Kato | Approach for uncertainty of turbulence modeling based on data assimilation technique[END_REF]. Specifically, Bayesian calibration allows one to infer posterior probability distributions of the closure parameters given some available set of experimental data. However, it is important to stress that a posterior distribution might only display a predictive capability for the flow case for which it was calibrated. In [START_REF] Edeling | Bayesian estimates of parameter variability in the k-ε turbulence model[END_REF], we performed multiple Bayesian calibrations for the coefficients of the well-known kε model [START_REF] Launder | Application of the energy-dissipation model of turbulence to the calculation of flow near a spinning disc[END_REF], using experimental velocity data from boundary-layers subject to a wide range of different pressure gradients. Even though the flow topology was the same for all calibration cases (i.e. boundary-layers), the resulting posterior distributions were significantly different. This variability of posterior distributions is a measure for the lack of predictive capability of a given model due to the simplifying assumptions in its mathematical structure. The model needs adjusted coefficients for each new case to counterbalance deficiencies inherent to neglecting part of the flow physics.

The second class of uncertainty, namely model-form uncertainty, arises through the assumptions that are made in the mathematical form of the closure model. Many different models exist, and their performance is known to be flow dependent, see e.g. [START_REF] Zhang | Evaluation of various turbulence models in predicting airflow and turbulence in enclosed environments by CFD: Part 2 Comparison with experimental data from literature[END_REF][START_REF] Harrison | Comparison of RANS turbulence models for prediction of film cooling performance[END_REF]. Previous attempts to tackle model-form uncertainty include calibrations of a modified eddy viscosity field [START_REF] Dow | Quantification of structural uncertainties in the k-ω turbulence model[END_REF], the introduction of statistical error terms in the calibration procedure [START_REF] Cheung | Bayesian uncertainty analysis with applications to turbulence modeling[END_REF][START_REF] Kennedy | Bayesian calibration of computer models[END_REF] or in the model transport equations [START_REF] Singh | Using field inversion to quantify functional errors in turbulence closures[END_REF], and the use of a stochastic counterpart of the Reynolds stress tensor [START_REF] Xiao | Quantifying and reducing model-form uncertainties in Reynolds-averaged Navier-Stokes simulations: A data-driven, physics-informed Bayesian approach[END_REF]. Even UQ techniques for the RANS equations with no modelling, where the transport terms are represented with direct numerical simulation data, have been developed [START_REF] Poroseva | RANS simulations of a channel flow with a new velocity/pressure-gradient model[END_REF][START_REF] Poroseva | On the accuracy of RANS simulations with DNS data[END_REF]. All these methods rely on detailed data from direct numerical simulations to infer about stochastic fields. Although promising, the above-mentioned methodologies suffer from the following drawbacks: i) DNS data are available only for relatively simple, low-Reynolds configurations; ii) inferring the stochastic fields involves costly optimization algorithms, e.g. to maximize the likelihood of observing the data; iii) extrapolating the stochastic correction terms to the prediction of unseen configurations is still an open problem, especially for realistic flow configurations.

Other means for representing model-form uncertainty include the use of a multi-model framework. A relevant example can be found in Poroseva et al. [START_REF] Poroseva | Improving the predictive capability of turbulence models using evidence theory[END_REF], where Dempster-Shafer evidence theory is used to obtain a measure of the total uncertainty in the mean velocity predictions around a RAE 2822 airfoil. Dempster' rule is used to fuse the predictions of two closure models, and using only sparse experimental data, the authors of [START_REF] Poroseva | Improving the predictive capability of turbulence models using evidence theory[END_REF] predict the velocity with quantified uncertainty at locations in the same flow where no data is available. Multi-model approaches have also found application in meteorology [START_REF] Poroseva | Multimodel approach based on evidence theory for forecasting hurricane/typhoon tracks: Further improvements[END_REF][START_REF] Diomede | Discharge prediction based on multi-model precipitation forecasts[END_REF], hydrology [START_REF] Duan | Multi-model ensemble hydrologic prediction using Bayesian model averaging[END_REF] and climate science [START_REF] Tebaldi | The use of the multi-model ensemble in probabilistic climate projections[END_REF], to name just a few.

We will focus on Bayesian Model-Scenario Averaging [START_REF] Draper | Assessment and propagation of model uncertainty[END_REF] (BMSA), which provides a coherent framework to address both parametric and model-form uncertainties, providing a robust predictive method (see for instance [START_REF] Rojas | Application of a multimodel approach to account for conceptual model and scenario uncertainties in groundwater modelling[END_REF] for an application to groundwater modelling uncertainties). Like Bayesian Model Averaging (BMA), BMSA combines the predictions from multiple models, thereby providing a measure for (closure) model uncertainty [START_REF] Hoeting | Bayesian model averaging: a tutorial[END_REF]. In addition, BMSA can combine the posterior distributions from different calibration scenarios, which allows one to compute error estimates due to scenario variability. Unless this variability is negligible, one should not assume that obtained parameter estimates are extrapolative to other flow scenarios. In other words: combining multiple non-overlapping posteriors adds uncertainty for the assumption that a posterior distribution can be applied outside the flow scenario under which it was calibrated.

In this article, we rely solely on RANS models for prediction. Hence we cannot assume that the true model which generated the data is contained in our selected model set, although in principle higher fidelity models could be included in the Bayesian average, at increased computational cost. It should therefore be noted that subject-matter expertise, applied to the selection of competing models (and scenarios), is still important for maintaining predictive accuracy. As noted by Draper [START_REF] Draper | Assessment and propagation of model uncertainty[END_REF], although multi-model approaches still introduce a bias in the prediction, due to the selection of a finite set of model structures, they play a useful role in reducing the bias with respect to predictions based on a single model structure.

In previous work [START_REF] Edeling | Predictive RANS simulations via bayesian model-scenario averaging[END_REF], BMSA was used to construct a stochastic model for predicting boundary layer flows subject to an arbitrary pressure gradient, based on multiple RANS models calibrated using selected experimental boundary-layer data. Specifically, we applied BMSA using 5 closure models and 14 boundary-layer calibration scenarios (corresponding to various external pressure gradients) to predict unobserved boundary-layer flows subject to an arbitrary pressure gradient.

Each individual calibration was performed using the Bayesian method developed by the authors of [START_REF] Cheung | Bayesian uncertainty analysis with applications to turbulence modeling[END_REF][START_REF] Oliver | Bayesian uncertainty quantification applied to RANS turbulence models[END_REF]. Other notable means of stochastic parameter estimation can be found in [START_REF] Margheri | Epistemic uncertainties in RANS model free coefficients[END_REF][START_REF] Papadimitriou | Bayesian uncertainty quantification of turbulence models based on high-order adjoint[END_REF]. Both the calibration and the prediction phases relied on a fast boundary-layer code. When addressing more complex, computationally intensive, flow configurations, a clear downside of BMSA is the associated computational burden. The first cost is associated to the calibration phase, in which the Markov-Chain Monte Carlo method [START_REF] Hastings | Monte Carlo sampling methods using Markov chains and their applications[END_REF] is used to draw samples from the posterior distribution.

Furthermore, the predictive phase of the BMSA procedure dictates the need to propagate I × K (I being the number of competing models and K the number of calibration scenarios) posterior distributions through a RANS code applied to the flow scenario of interest, which is prohibitively expensive for most practical flow configurations. A common option is to replace the expensive RANS code with a much cheaper surrogate model, such as for instance a polynomial approximation created with stochastic collocation or polynomial chaos methods [START_REF] Eldred | Comparison of non-intrusive polynomial chaos and stochastic collocation methods for uncertainty quantification[END_REF]. Still, creating the surrogate model will require samples from the full RANS code, and surrogate modeling techniques are subject to the so-called curse of dimensionality. The number of required code samples therefore rises exponentially with the number of unknown parameters.

Our goal for this paper is to assess the predictive capability of the BMSA method for computationally expensive flow cases, while using the most drastic cost reduction technique possible to speed up the prediction step. We assume that calibration cost is less critical, since it is done only once in a preliminary phase. Specifically, in this work we use the posterior distributions of [START_REF] Edeling | Predictive RANS simulations via bayesian model-scenario averaging[END_REF], and propagate them through new flow configurations described by a full RANS solver instead of a boundary layer code. Note that posterior distributions resulting from additional calibration scenarios can be added as soon as they become available.

To drastically reduce the computational cost associated with the propagation step, we approximate each posterior as a Dirac delta distribution centered at the Maximum A Posteriori (MAP) estimate. Now, instead of having to propagate I × K posteriors, we are left with I × K deterministic evaluations. In order to evaluate the predictive capability of the proposed modified BMSA approach, we selected a number of flow problems which include physics that the boundary-layer trained BMSA model was not subjected to during the calibration phase. This paper is organized as follows. First, we will give a very brief overview of the competing closure models considered in this paper. In Section III, we describe our reduced BMSA method using MAP estimates. Preliminary validations of the MAP-base BMSA approach are carried out in Section IV using the boundary layer flow cases of [START_REF] Edeling | Predictive RANS simulations via bayesian model-scenario averaging[END_REF]. Section V shows the results obtained for the different prediction flow cases. Finally, we give our conclusion in Section VI.

II. RANS turbulence models

In the present work, we restrict our attention to three transport equation models largely employed in engineering applications. These are selected from the family of so-called linear eddyviscosity models, which adopt a linear representation of the Reynolds stress tensor, based on the well-known Boussinesq assumption. We use the kε, kω and Spalart-Allmaras (SA) models.

As these are very well-known models, we do not give full details for the sake of conciseness.

Rather, we refer to [START_REF] Wilcox | Turbulence Modeling for CFD[END_REF] for a complete description of their mathematical structure. Also, all values for the closure coefficients (both nominal and perturbed), can be found in Appendix A.

III. Bayesian predictive Methodology

We describe a statistical methodology based on Bayesian Model Scenario Averaging (BMSA), for making predictions of a chosen Quantity of Interest (QoI). BMSA provides a coherent framework for combining predictions from multiple competing conceptual models and calibration data to attain a more realistic and reliable description of the predictive uncertainty.

A. BMSA formulation

Consider a particular flow of interest, including boundary-conditions, material parameters, and all other physical properties needed to define it unambiguously. This we term a scenario, denoted S. Assume that there is reference data z available for S, which may be measurements of any flow quantity or data from higher fidelity simulations (e.g. DNS). Let m CFD (S; M, θ) be a governing model, here represented by a CFD code, which takes as arguments the scenario S, a turbulence model M and its closure-coefficients θ. Note that the size of θ and the meaning of its components will depend on M . Let the model return the full state of the fluid in all variables. Then we can define an operator H z (•) that maps the state to the observed quantities z (this can be, e.g., a procedure for extracting the aerodynamic coefficients from the full solution for the flow around a body or a velocity profile at a given station). The model output is then given by H z • m CFD (S; M, θ) and can be related to the observed data by means of a statistical model. One of such models, used in [START_REF] Edeling | Bayesian estimates of parameter variability in the k-ε turbulence model[END_REF][START_REF] Edeling | Predictive RANS simulations via bayesian model-scenario averaging[END_REF] is represented by:

z = η H z • m CFD (S; M, θ) + , (1)
where , η are random-variables (RVs) representing measurement noise and model error respectively.

The latter accounts for the fact that, even with the best possible parameters, the model does not predict the true value for z, due to errors intrinsic to modelling assumptions. The multiplicative model-error term of (1) was first introduced by [START_REF] Cheung | Calculation of posterior probabilities for Bayesian model class assessment and averaging from posterior samples based on dynamic system data[END_REF]. Additive models are also available, see e.g. [START_REF] Papadimitriou | Bayesian uncertainty quantification of turbulence models based on high-order adjoint[END_REF].

We assume throughout that the CFD calculations have run on sufficiently fine grids, such that errors due to grid resolution are significantly less than the uncertainty due to modelling assumptions. We also assume that the calculations have well converged to the steady state solution. The experimental noise term is modelled using a zero-mean Gaussian distribution.

The observations available for a given S can be used to infer the model parameters θ using a Bayesian calibration approach, i.e. by applying Bayes' rule:

p (θ|z, M, S) = p (z|θ, M, S) p(θ | M, S) p(z | M, S) , (2)
where p (θ|z, M, S) is the joint posterior probability distribution of the closure coefficients under the given turbulence model and the calibration data z observed for scenario S. is the evidence of the data, which normalizes the posterior distribution. Since p(z|M, S) does not depend upon θ, it does not have to be computed while calibrating θ.

Given the complexity of the likelihood function, the posterior p(θ|z, M, S) cannot be obtained analytically and is instead evaluated numerically. We employ the Markov-Chain Monte-Carlo method [START_REF] Hastings | Monte Carlo sampling methods using Markov chains and their applications[END_REF] to draw samples from (2). To reach convergence of the Markov-chain, we observed that O 10 4 code samples were required [START_REF] Edeling | Bayesian estimates of parameter variability in the k-ε turbulence model[END_REF]. Ordinarily this would constitute an excessive strain on available computational resources in a CFD context. In previous work [START_REF] Edeling | Predictive RANS simulations via bayesian model-scenario averaging[END_REF] our experimental data consisted of boundary-layer quantities, and m CFD reduced to a fast boundary-layer code. As such, no real computational bottleneck exists during the calibration phase.

The calibration outcome is a discrete approximation of the joint posterior probability distribution of the parameters, which is then available for the prediction of new (unobserved) scenarios. Let ∆ be a quantity of interest (QoI) to be predicted for a new scenario S, for which we do not have data. The QoI may be a scalar, vector, or functional quantity derived from the flow-state. Since m CFD is a model for ∆, we write

∆ H ∆ • m CFD ( S; M, θ), (3)
where, analogous to H z , H ∆ is an observation operator extracting ∆ from the full flow-state, M is the closure model adopted for the prediction, S is the prediction scenario. Equation ( 3) is a prediction based on a single choice for M and the associated parameters, i.e. a deterministic prediction. A stochastic prediction accounting for parameter uncertainty is obtained by propagating the posterior distributions for θ through (3). This however does not account for the fact that different distributions may be obtained by calibrating the model against different scenarios, nor for the fact that multiple alternative closure models M may be used.

To estimate the effect of these uncertainties, we build the full BMSA predictive distribution of ∆ for a set of alternative conceptual models M = (M 1 , M 2 , ..., M i , ..., M I ), under different calibration scenarios S = (S 1 , S 2 , ..., S k , ..., S K ) for which we have data Z = (z 1 , z 2 , ..., z k , ..., z K ) [START_REF] Draper | Assessment and propagation of model uncertainty[END_REF][START_REF] Hoeting | Bayesian model averaging: a tutorial[END_REF]:

p(∆ | S, Z, M, S) = K k=1 I i=1 p(∆ | S, z k , M i , S k ) p(M i | z k , S k ) p(S k ) (4)
Equation ( 4) is an average of the I ×K posterior predictive distributions (ppds

) p(∆ | S, z k , M i , S k ).
These are weighted by the posterior model probabilities p(M i | z k , S k ) and by scenario probabilities p(S k ). The ppds are obtained by marginalizing ∆ over the posterior distributions of the parameters:

p(∆ | S, z k , M i , S k ) = p(∆ | S, M i , S k , θ) p(θ | z k , M i , S k ) dθ (5)
In simpler terms, this amounts to propagating the

I × K posterior distributions p(θ | z k , M i , S k )
through the new scenario S by a suitable uncertainty quantification method.

The posterior model probabilities, on the other hand, reflect how well model M i fits the data over the support of the prior. They can be determined from a separate application of Bayes' rule, i.e.

p (M i | z k , S k ) = p(z k | M i , S k )p(M i | S k ) I j=1 p(z k | M j , S k )p(M j | S k ) , (6)
where p(M i | S k ) are the prior model probabilities (generally modelled as equi-probable). Note that

p(z k | M i , S k )
is the evidence for model M i and scenario S k , previously introduced in (2). Thus, in case of a model ensemble, the denominator of ( 2) is not ignored and is computed via

p(z k | M i , S k ) = p(z k | M i , S k , θ)p(θ | M i , S k )dθ. (7)
The first two moments of p(∆ | S, Z, M, S) can be derived directly from (4) and read [START_REF] Draper | Assessment and propagation of model uncertainty[END_REF][START_REF] Edeling | Predictive RANS simulations via bayesian model-scenario averaging[END_REF]:

E[∆ | S, Z, M, S] = I i=1 K k=1 E[∆ | S, z k , M i , S k ] p(M i | z k , S k ) p(S k ) (8) V[∆ | S, Z, M, S] = I i=1 K k=1 V[∆ | S, z k , M i , S k ] p(M i | z k , S k ) p(S k ) + I i=1 K k=1 (E[∆ | S, z k , M i , S k ] -E[∆ | S, z k , M, S k ]) 2 p(M i | z k , S k ) p(S k ) + K k=1 (E[∆ | S, z k , M, S k ] -E[∆ | S, Z, M, S]) 2 p(S k ), (9)
where

E[∆ | S, z k , M, S k ] := I i=1 E[∆ | S, z k , M i , S k ]p(M i | z k , S k )
is the BMA prediction within scenario S k . Equation [START_REF] Dow | Quantification of structural uncertainties in the k-ω turbulence model[END_REF] shows that the posterior predictive distribution variance for ∆ consists of three terms: the first term is called the within model, within scenario variance, and accounts for predictive uncertainties due to the residual parameter variability after calibration of a model M i under a scenario S k . The second one is called the between model, within scenario variance. It is large if the models, using the posterior parameter distributions originating from the same calibration scenario, provide different predictions for the new configuration. It can therefore be considered as a measure for model error. The last term is the between scenario variance, which accounts for the fact that different calibration scenarios lead to different estimates of the parameters. Thus, when applying the BMSA model to a prediction scenario S for which it was not calibrated, this uncertainty should be taken into account.

B. Implementation of BMSA for CFD problems

In this paper, we focus on the prediction step, and refer the reader to [START_REF] Edeling | Bayesian estimates of parameter variability in the k-ε turbulence model[END_REF][START_REF] Edeling | Predictive RANS simulations via bayesian model-scenario averaging[END_REF] for details about the calibration step. All flow predictions presented in the following are based on the posterior parameter distributions previously computed for a class of flat plate boundary layer flows subject to various external pressure gradients [START_REF] Edeling | Predictive RANS simulations via bayesian model-scenario averaging[END_REF]. These constitute a database of posterior probability distributions, which is available for propagation through new flow configurations. The interested reader may download the samples of the posterior θ distributions used in the following applications from [START_REF] Edeling | MCMC traces of posterior distributions of closure coefficients for 6 turbulence models[END_REF]. Since these parameter estimates where determined for wall-bounded flow configurations, we expect these may be used to predict other (more complex) wall bounded flows. The validity of this working assumption will be assessed in subsequent sections.

The calibration step also provides information about the posterior model probabilities [START_REF] Launder | Application of the energy-dissipation model of turbulence to the calculation of flow near a spinning disc[END_REF].

Indeed, for each model and calibration scenario, it is possible to compute the evidence by Monte-Carlo integration of [START_REF] Zhang | Evaluation of various turbulence models in predicting airflow and turbulence in enclosed environments by CFD: Part 2 Comparison with experimental data from literature[END_REF]. For instance, for the cheap flat plate boundary layer scenarios, the evidence is computed by drawing a large number of samples,

θ n , n = 1, • • • , N , from the prior p(θ | M i , S k ),
and by approximating [START_REF] Zhang | Evaluation of various turbulence models in predicting airflow and turbulence in enclosed environments by CFD: Part 2 Comparison with experimental data from literature[END_REF] as

p(z k | M i , S k ) ≈ 1 N N n=1 p(z k | M i , S k , θ n ). (10)
The posterior model weights can now be calculated through the simple algebraic relation [START_REF] Launder | Application of the energy-dissipation model of turbulence to the calculation of flow near a spinning disc[END_REF], the results of which we tabulated in Appendix A. It should be noted that computing the evidence using samples from the prior is not always feasible. Convergence of the MC estimate [START_REF] Kennedy | Bayesian calibration of computer models[END_REF] will become slow if the region where the likelihood has high probability is significantly different from the region of high prior probability content [START_REF] Beck | Bayesian system identification based on probability logic[END_REF]. In this case more advanced sampling techniques will be required, see e.g. [START_REF] Cheung | Calculation of posterior probabilities for Bayesian model class assessment and averaging from posterior samples based on dynamic system data[END_REF][START_REF] Chib | Marginal likelihood from the Metropolis-Hastings output[END_REF].

The prediction step involves:

1. Propagating posterior parameter distributions for each competing model through the prediction scenario S;

2. Assigning a suitable probability mass function to the scenarios S k ;

3. Applying the BMSA relations (4), and/or (8),( 9) to obtain the predictive posterior distribution and/or its moments for the QoI ∆.

Steps 2. and 3. only involve a negligible or small computational load. The computational bottleneck is represented by step 1., which requires running the CFD model on scenario S for all models in M, and for a large number of samples in each posterior distribution p(θ | z k , M i , S k ), leading to I × K stochastic estimates of ∆. Since each CFD simulation is likely to be expensive itself, propagation by means of Monte Carlo sampling is unacceptably costly. The computational cost can be greatly alleviated by using surrogate models (see e.g. [START_REF] Edeling | Simplex-stochastic collocation method with improved scalability[END_REF]), yet a significant number of CFD runs (O(100) or higher) is still required even when using advanced surrogate modelling techniques. For this reason, in Section III D we introduce an approximation that drastically reduces the computational cost associated to the prediction step. In the next section we complete the description of the BMSA by discussing the choice of the scenario probability mass function (step 2).

C. Weighting the scenarios

Whereas the posterior model probabilities p (M i | S k , z k ) are informed from the calibration data z k , the scenario probabilities p (S k ) are of a predictive nature, and incorporate uncertainty for the fact that the MAP estimates are applied to a predictive scenario S for which they were not calibrated. Ideally experimental data from S would be available, which we could use to inform p (S k ). However, since in many engineering applications of CFD this will not be the case, we assume we only have data coming from the K calibration scenarios S k . Instead of informing on data, we specify p(S k ) based on model agreement per S k . Specifically, following to [START_REF] Edeling | Predictive RANS simulations via bayesian model-scenario averaging[END_REF], we set:

p(S k ) := ξ -p k K j=1 ξ -p j , ξ k = I i=1 E[∆ | S, z k , M i , S k ] -E[∆ | S, z k , M, S k ] 2 ∀S k ∈ S. (11)
Here, ξ k measures model agreement by means of an L 2 -normed variation around the BMA expec-

tation E[∆ | S, z k , M, S k ].
If all turbulence models show a high level of agreement regarding the value of ∆ when using the posteriors calibrated on S k , such a scenario is given a high probability and vice versa. The rationale is that if a given S k is similar to the predictive S at hand, the models are expected to give similar predictions given that the parameter distributions were calibrated on the same data z k . Additionally, p ∈ N 0 is an integer which controls the degree to which variation is penalized. For p = 0 we obtain a uniform distribution and for higher values of p we effectively damp the results coming from those scenarios with a low level of model agreement ξ k . In section IV we carry out numerical experiments to demonstrate that ξ k indeed represents a good model for the true discrepancy between the BMSA prediction and the data and we discuss the role of the exponent p.

Finally, note that (11) requires a minimum of two models, the use of one model will always yield a uniform distribution. With just a single model, one might still use the variability over S k to inform the p(S k ), or simply rely on expert opinion.

A valid criticism of ( 11) is that the 'best' S k will not be favoured if all models in M for a given 'wrong' S k produce predictions which are incorrect in a very similar manner. Due to the fact that we do not allow ourselves any reference data for the predictive flow scenario, we cannot fully exclude the possibility of such a situation. However, we can suggest means to decrease its risk through: i)

selecting a larger set of diverse models of different fidelity, which are unlikely to agree with each other, ii) selecting a relevant set of calibration scenarios, since in previous research [START_REF] Edeling | Predictive RANS simulations via bayesian model-scenario averaging[END_REF] we found that [START_REF] Singh | Using field inversion to quantify functional errors in turbulence closures[END_REF] does favour S k similar to S if the flow topology is comparable, and iii) limit the value of p, which gives the user control over how much p(S k ) is allowed to deviate from a uniform distribution.

θ MAP i,k := argmax θ p(θ | z k , S k , M i ) (12)
so that

p(θ | z k , S k , M i ) δ θ -θ MAP i,k . ( 13
)
The consequence of this approximation is to neglect the effect of "within-model, within-scenario" variance. The ppd variance still includes the effect of multiple models and scenarios. Note that if perfectly plentiful data were available in the training scenarios (and the models were able to fit the data exactly for some values of the closure coefficients), then the p(θ | S k , M i , z k ) would tend toward Dirac δ-functions. So one way to think of this approximation is as neglecting the effect of imperfect training. The MAP estimates (see Appendix A), were estimated from kernel-density estimates constructed with McMC samples of previous research [START_REF] Edeling | MCMC traces of posterior distributions of closure coefficients for 6 turbulence models[END_REF], see also Figure 1. While this might only be an approximation to each individual MAP value, it is the spread of different MAP values that is important in our approach, see the discussion in Section IV A.

Substituting (13) into (5) leads to an approximation of the posterior predictive distribution

p(∆ | S, Z, M, S) p(∆ | S, Z, M, S) by p(∆ | S, Z, M, S) = I i=1 K k=1 δ ∆ -∆ CFD S; M i , θ MAP i,k p(M i | S k , z k ) p(S k ) (14)
where we define ∆ CFD ( S) := H ∆ • m CFD ( S) and we used the fact that propagating a δ distribution requires running m CFD only at θ = θ MAP i,k , returning a δ distribution centered at ∆ CFD ( S) for the output. The approximate ppd is then reduced to a weighted-sum of I × K δ-functions, one at each prediction of ∆ CFD ( S) for each model, and each scenario's MAP-estimate of θ. The cost of evaluating the ppd is then I × K runs of ∆ CFD ( S). Finally, the mean and variance of ( 14), are easy to evaluate using the properties of δ-functions, and contain terms from our database (p(M i | S k , z k ) and θ MAP i,k ), a user-choice (p(S k )) and the RANS output for ∆ via m CF D :

E[∆ | S, Z, M, S] = I i=1 K k=1 ∆ CFD S; M i , θ MAP i,k p(M i | S k , z k ) p(S k ) (15) V[∆ | S, Z, M, S] = K k=1 I i=1 ∆ CFD S; M i , θ MAP i,k -E[∆ | S, z k ] 2 p(M i | S k , z k ) p(S k ) + K k=1 E[∆ | S, z k ] -E[∆ | S, z] 2 p(S k ). (16)
Note that compared to (9), the ppd variance for ∆ now contains only two terms: the between model, within scenario variance, and the between scenario variance. The effect of the MAP approximation on the ppd variance is investigated in Section IV.

IV. Preliminary validation

In this section we use our numerical boundary-layer flow database [START_REF] Edeling | Predictive RANS simulations via bayesian model-scenario averaging[END_REF] along with experimental data from [START_REF] Coles | Computation of turbulent boundary layers[END_REF] for the preliminary validation of two key ingredients of the present BMSA formulation, namely, the MAP approximation and the scenario weighting. The boundary-layer flows considered in the present calculations and the corresponding nomenclature are described in Appendix A.

A. Influence of the MAP approximation

The MAP approximation ignores the effect of imperfect training within each calibration scenario (i.e. it ignores the variance of each individual posterior), As an example, consider the results of Figure 1, where all 14 posterior distributions of the coefficient C ε2 (of the k -ε model, see Appendix A) and the corresponding MAP estimates are plotted. The variance of these marginals is not considered small, due to the fact that even at the best model parameters, there is still model inadequacy which will prevent the posteriors from approaching delta distributions [START_REF] Oliver | Bayesian uncertainty quantification applied to RANS turbulence models[END_REF][START_REF] Papadimitriou | Bayesian uncertainty quantification of turbulence models based on high-order adjoint[END_REF][START_REF] Cheung | Calculation of posterior probabilities for Bayesian model class assessment and averaging from posterior samples based on dynamic system data[END_REF]. Moreover, not all parameters will be equally sensitive to the reference data. We previously observed [START_REF] Edeling | Bayesian estimates of parameter variability in the k-ε turbulence model[END_REF][START_REF] Edeling | Predictive RANS simulations via bayesian model-scenario averaging[END_REF] that parameters with low Sobol indices will yield posteriors that do not deviate much from our (uniform)

prior. As such, the MAP approximation might seem to neglect a large amount of uncertainty. However, since we consider 14 different calibration scenarios, and since the modes of the posteriors significantly disagree on the most likely value of θ, the MAP approximations still span a large section of the θ domain, see again Figure 1 for an example. Hence, the MAP approximation is of the least influence when this span is sufficient to yield an envelope on the QoI that is similar in size to the full posterior ppd variance, obtained by propagating all posterior θ distributions. It should be noted however that this will only hold for lower-order statistical moments, all information contained in the tails of the posteriors will be lost in the MAP approximation.

We can evaluate the effect of the MAP approximation on the boundary-layer flows of the calibration phase, by computing the posterior standard-deviation of the full BMSA method, and the current proposed approach of Section III D. Figure 2 shows the predictive results for an unseen boundary layer, where ∆ is a u + profile, defined as u + := U /u τ . Here, U is the mean steamwise velocity and u τ is the friction velocity. Precisely, we used posteriors from the remaining 13 boundary layer scenarios in Table 6 (Appendix A) to predict case 1400. Note that the total variability in ∆ is fairly similar in both cases. Due to the lack of a within-model, within-scenario component in the MAP approach, a larger part of the standard-deviation in ∆ is lumped into the between-scenario component. We repeated this procedure for a number of unseen boundary layers, and obtained

Var[∆

| z] Full BMSA in M i in S k σ ∆ between Mi in S k σ ∆ between S k σ ∆ (a) Full BMSA.

Var[∆

| z] MAP approach in M i in S k σ ∆ between Mi in S k σ ∆ between S k σ ∆ (b) MAP approach.
Fig. 2 The breakdown of the posterior standard-deviation of ∆ = u + , plotted versus y + for both the full BMSA method (computed with ( 9)) and the simplified MAP approach (computed with ( 16)). Both figures represent the the same predictive flow case, and were generated with uniform p(S k ).

similar results. Thus, for those cases considered where we could compute the full BMSA ppd standard-deviation, we did not observe a significant underestimation of the posterior predictive standard-deviation when the MAP approach was used.

B. Influence of the scenario weighting

Equation ( 11) raises the question: is ξ k a good measure for the true model error contained in each S k ? Let us first define the true model error as:

ξk := I i=1 H z • m CFD S; M i , θ MAP i,k -z k 2 ,
where in our case the z k are experimental u + profiles. We can then evaluate its proximity to our error-model ξ k by performing a leave-one-out validation analysis on the K calibration scenarios.

Note that a given MAP estimate θ M AP i,k

has not been calibrated on the data of the other K -1

calibration scenarios S j , j ∈ {1, • • • , K} \ {k}.
The scenarios are roughly ordered from 'easy' zeropressure gradient scenarios at k = 1, to the most difficult, highly adverse scenarios near k = K, see Appendix A for a more complete description. By using each θ M AP i,k

to compute u + in these remaining K -1 scenarios we can evaluate both ξ k and ξk due to the availability of data. The results are shown in Figure 3, where the dotted line corresponds to perfect agreement ξ k = ξk . Each subplot shows the results of a single MAP estimate, applied to all other scenarios. Note that in general the results lie close to the ξ k = ξk line. However, MAP estimates calibrated under 'easy' zero or favorable pressure gradient scenarios can yield discrepancies between ξ k and ξk . More specifically, we observe such discrepancies when using a MAP estimate calibrated under an easy scenario to predict a QoI under a difficult scenario, see e.g. the results of applying θ M AP i,1

to predict the u + profile of S 14 .

On the other hand, MAP estimates that were calibrated under the most difficult strongly-adverse scenarios perform well across all S k , see the results for θ M AP i,13

and θ M AP i,14 .

The parameter p of ( 11) is user specified. To estimate its value we use a similar leave-one-out analysis and the experimental data from the calibration scenarios. In this case we fix a scenario S k , and use the remaining θ M AP i,j

, j ∈ {1, • • • , K} \ {k} to predict the QoI for a given value of p. Let us define a relative error in the BMSA prediction as

rel k := E [∆|z] -z k 2 z k × 100%. ( 17
)
For all S k we can now examine the behavior of ( 17) as a function of p. The results are depicted in Figure 4. For most S k the relative error is either minimized or reaches a plateau for p ∈ {1, 2, 3}.

Based on these results we will set p = 2 in predictive cases where no experimental data is available. V. Assessment of BMSA for selected flow cases

In this section, the MAP-based BMSA predictive methodology is assessed against selected flow cases, namely separated incompressible flow, and 3D transonic flow past a wing. The considered applications are rather different from our calibration scenarios, although they still involve boundary layers subject to variable pressure gradients, and exhibit features beyond the limits of applicability of linear RANS models, due to the strong non-equilibrium effects. As such, they represent severe tests for the proposed approach. S 14

Fig. 4 The relative error in the BMSA prediction of u + for all K scenarios S k .

A. Flow over periodic hills at ReH = 5600

The flow over a series of periodically arranged hills deals both with flow separation on the curved surface of the hills and reattachment between the hills. The test case is known to be especially challenging for linear eddy-viscosity models, which are not able to predict the mean effect of the unsteady fluctuation of the separation and reattachment points correctly [START_REF] Jakirlic | Extended excerpt related to the test case: Flow over a periodical arrangement of 2D hills[END_REF].

The periodicity of the hill-geometry is mimicked with periodic boundary conditions at the inlet and outlet and no-slip conditions at the walls. A volume forcing is applied to each cell, which maintains a bulk velocity of U b = 1.0 between the hill's crest and the upper surface. The functional form of the lower surface is defined according to the ERCOFTAC test case description[? ]. The considered Reynolds-number Re H = 5600 based on the hill height H and bulk velocity between the hill's crest and the upper surface is maintained by a volume forcing. The CFD simulations are performed using OpenFOAM's simpleFoam-solver [START_REF]OpenFOAM User Guide[END_REF]. We use validation data of the mean-flow field from Breuer et al. for comparison [START_REF] Breuer | Flow over periodic hills -Numerical and experimental study in a wide range of Reynolds numbers[END_REF]. Three turbulence models were used, Spalart-Allmaras, The baseline simulations for the three models show the expected behavior as reported in the literature (for a larger Re-number [START_REF] Jakirlic | Extended excerpt related to the test case: Flow over a periodical arrangement of 2D hills[END_REF]): while both Spalart-Allmaras and the [START_REF] Wilcox | Turbulence Modeling for CFD[END_REF] kω over-predict the size of the recirculation zone, characterised by a zero velocity component close to the wall, the Launder-Sharma kunder-predicts this flow feature. For S 14 the simulation using [START_REF] Wilcox | Turbulence Modeling for CFD[END_REF] kω didn't converge so that this scenario was excluded from the set, leading to only 13 scenarios being used for the BMSA predictions.

The BMSA expectation E[∆| S, z] of the stream-wise velocity, as shown in Figure 5, over-predicts the recirculation zone and gives similar velocity profiles for x = 3.0 and 4.0 close to the lower surface. Especially, for x = 0.0 the expectation doesn't capture the local maximum of the velocity close to y/H = 1.0, i.e. at the hills crest. Throughout the different locations the expectation shows differences compared to the DNS for the upper part of the velocity profile, which might be due to a compensation-effect caused by the fixed mass-flow.

The 99% confidence interval captures the DNS data for x ≥ 3, but not the region close to the upper wall for x ≥ 5 and the local velocity maxima in the area of the free-shear layer in the leeward region of x = 0 to 2.0. Interestingly for x ≤ 3.0 when the range of the confidence interval shrinks locally in y/H-direction also the DNS is still inside. However, this pattern is not the same for every y/H-position, e.g. at x = 2.0 the expectation and the DNS match for 0.5 ≤ y/H ≤ 1, but the confidence interval is large. Interestingly, in the reattachment region between the two hills, where the case is similar to a flat plate, the confidence interval is able to capture the DNS data consistently.

In general, the application of BMSA to this test case demonstrates the aforementioned fact that the used linear eddy-viscosity models suffer from restrictions, which inhibit the reproduction of effects due to e.g. streamline curvature as well as Reynolds-stress anisotropy. These are essential characteristics of this challenging test case [START_REF] Jakirlic | Extended excerpt related to the test case: Flow over a periodical arrangement of 2D hills[END_REF]. Adding second-order closures to M might remedy this deficiency, at increased computational cost. However, the BMSA confidence interval captures most of the modelling error showing robustness of the method even for this challenging test case.

B. Transonic flow over the ONERA M6 swept wing

The second application of our boundary-layer trained BMSA is the transonic flow over the ONERA M6 wing. This is a swept wing based on the ONERA D airfoil with no twist. The M6 wing is known to produce complex physics at transonic Mach numbers. See for instance the pressure distribution depicted in Figure 6. Due to the configuration of local supersonic regions [START_REF] Kuzmin | On the lambda-shock formation on ONERA M6 wing[END_REF], a λ-shock is formed where the double shock merges into a single one between y/b = 0.8 and y/b = 0.9.

A well-known CFD validation case (see e.g. [START_REF] Kuzmin | On the lambda-shock formation on ONERA M6 wing[END_REF][START_REF] Lerat | High-order residual-based compact schemes for aerodynamics and aeroacoustics[END_REF][START_REF] Huang | A third-order accurate centered scheme for turbulent compressible flow calculations in aerodynamics[END_REF]) involves the ONERA M6 wing at Mach the numerical setup can be found in [START_REF]NPARC Alliance Validation Archive, ONERA M6 wing[END_REF] and [START_REF]FLUENT -3D Transonic Flow Over a Wing[END_REF]. For these conditions experimental pressure data is available from Schmitt and Charpin [START_REF] Schmitt | Pressure distributions on the ONERA-M6-Wing at transonic Mach numbers[END_REF] at several y/b stations. The error in the reported pressure data is ∆C p = ±0.02.

Before beginning the BMSA procedure, we performed a qualitative comparison of the C p distributions using the MAP estimates of the 1400 case on 4 unstructured meshes of increasing size, wing. The mean is is very similar to the results of Figure 6. More interesting is the uncertainty in the prediction, depicted in Figure 8(b). From this standard deviation is becomes clear that the uncertainty of the RANS prediction is localized in the region of the λ-shock, especially where the two shocks merge.

The corresponding p(S k ) distribution is shown in Figure 8(a). Note that the distribution favours the results of MAP estimates from a zero (1400) and a favorable (2700) calibration scenario, as well as two adverse cases (3300 and 1200). This might be explained by the fact that the transonic flow experiences a wide range of pressure gradients, at different locations of the wing. Perhaps the weights p(S k ) should ideally reflect the inhomogeneous nature of the uncertainty, rather than computing a scalar p(S k ) for each k, and assuming it is applicable throughout the entire flow field.

To investigate, we compute the BMSA prediction along different y/b stations of the wing. Fig. 8 The results for the ONERA M6 wing.

yield very similar predictions. Other computed p(S k ) distributions are comparable, up until the neighborhood of y/b = 0.95, where the results are quite different. Here, a relatively large amount of uncertainty is present, accompanied by a different p(S k ) distribution. There is a preference for the pressure computed with MAP estimates coming from a strongly-adverse calibration scenario (namely scenario 1200). These results suggest that incorporating a spatial dependence in P(S k )

does indeed reflect the local flow physics better. Specifying p(S k ) locally through [START_REF] Singh | Using field inversion to quantify functional errors in turbulence closures[END_REF] assumes that there is no spatial correlation.

VI. Conclusion

A Bayesian Model-Scenario Averaging methodology was applied to estimate turbulence modelling uncertainties in RANS computations. Specifically, the goal of this paper was to assess the predictive capability of a BMSA model, trained on experimental velocity profiles for a set of flat plate boundary layers, when predicting other, more computationally expensive flow problems. We considered two severe validation cases, namely the separated flow over periodic hills, and a shocked transonic flow over a 3D wing. These contain physics that were not present during the training phase, such as compressibility effects and separation. The full BMSA approach requires the propagation of a set of posterior distributions, which proves to be computationally intractable for industrialrelevant test cases. We therefore proposed to approximate the posterior distributions with Dirac δ functions centered at the Maximum A Posteriori (MAP) estimates of the posterior distributions. This is the most drastic cost reduction possible, with just one code evaluation per closure modelcalibration scenario pair.

The main advantage of an approach such as BMSA, is that one can make predictions with quantified uncertainty. Such information is not available with the baseline models, and provides insight into where, and to what extend, the closure models fail to be trustworthy. The confidence region generated by BMSA thus provides meaningful bounds on the QoI. In our selected cases we observed an overlap between the confidence intervals and the reference data, especially when the local flow physics resembled the attached boundary-layers of the calibration phase. Where there was no overlap the discrepancy between the propagated model-scenario pairs was small, which might be remedied by including additional closure models in the set, e.g. nonlinear eddy viscosity models.

  Furthermore, p(θ | M, S) is the joint prior probability distribution for the model coefficients, which summarizes available knowledge about θ before observing any data. The likelihood function p (z|θ, M, S), which stems from the statistical model (1), describes the probability of observing the data for scenario S, given a model and a realization of the closure coefficients. Evaluating p (z|θ, M, S) for a given realization of θ involves running m CFD and is, as such, a costly operation. Finally p(z|M, S)

D.

  Speeding up BMSA predictions: MAP approximation of the posterior parameter distributions A way to drive the estimate of the individual predictive distributions p(∆ | S, z k , M i , S k ) to an acceptable computational cost consists in approximating the marginal posterior probability distributions p(θ | z k , M i , S k ) with Dirac-δ functions at their maximum a posteriori (MAP) values:

Fig. 1

 1 Fig. 1 All marginal posterior distributions of Cε2 (obtained from kernel-density estimates using the McMC samples of [23]), and the corresponding MAP estimates.

Fig. 3

 3 Fig.3The results of the cross validation study. In each of the k subplots (k = 1, • • • , 14), the MAP estimates are fixed to those obtained from calibration scenario S k . The MAP estimates are then used to compute both the true error ξk and the modelled error ξ k (11) when applied to prediction scenarios Sj, j ∈ {1, • • • , K} \ {k}. The true errors ξk were computed using the experimental u + data from[START_REF] Coles | Computation of turbulent boundary layers[END_REF].

Launder-

  Sharma kand Wilcox (2006) kω, with y + ≤ 1.0. A mesh convergence study was conducted for each model using baseline coefficients resulting in meshes with sizes of 100 × 110 (Spalart-Allmaras and Launder-Sharma k -) and of 150 × 140 (Wilcox (2006) kω) cells in x × y direction respectively.

Fig. 5

 5 Fig. 5 BMSA for the periodic hill flow.

0. 835 ,Fig. 6 A

 8356 Fig. 6 A typical Cp distribution on the suction side of the ONERA M6 wing. The black lines correspond to the y/b = [0.2, 0.44, 0.65, 0.80, 0.90, 0.95, 0.99] stations where Cp data is available from [39].

  namely {396, 753, 901, 1180} × 10 3 elements. The results at z/b = 0.80 and z/b = 0.95 are shown in Figure 7. Note that the coarsest mesh smears the double shock at y/b = 0.80, and both mesh 1 and 2 smear the shock at y/b = 0.95. Between the two finest meshes we observed little variation, and selected therefore mesh 3 to perform the sampling with MAP estimates. Our model set for this case consists of the kε and SA model. We computed the mean C p prediction (15) and the standard deviation based on (16) for the entire suction side of the M6 The Cp distributions at z/b = 0.8. The Cp distributions at z/b = 0.95.

Fig. 7 The

 7 Fig.7The Cp distributions of the 4 meshes at two spanwise locations. All computations were done using the MAP estimates of calibration case 1400.

  The P(S k ) distribution at y/b = 0.44. | z] 99% confidence Schmitt&Charpin ± ∆C p (d) The BMSA prediction at y/b = 0.44 plus 99% confidence region. The P(S k ) distribution at y/b = 0.95. The BMSA prediction at y/b = 0.95 plus 99% confidence region.

  For each station, taking the local C p distribution as QoI, we recompute p(S k ). This allows us to build up a picture of how p(S k ) varies from the wing root to the tip. A selection of results, at y/b

		0.18
		0.16
		0.14
		0.12
		0.10
	P(S	0.08
		0.06
		0.04
		0.02
		0.00
	stations where we have experimental validation data, are displayed in Figures 8(c)-8(f). First notice that the amount of uncertainty at y/b = 0.44 is very low, a result that could already be inferred from Figure 8(b). The p(S k ) distribution at this station displays a preference for mostly those results computed using zero or favourable pressure-gradient MAP estimates. Note however, that the final C 1400 1300 2700 6300 1100 2100 2500 2400 3300 0141 1200 4400 2133 2134

p distribution is insensitive to P(S k ) in cases of very low uncertainty, seeing that all S k k ) (a) The p(S k ) distribution for the entire wing.

Table 3

 3 MAP estimates of the SA model for the different calibration scenarios. The nominal values are depicted on the first row.

	Scenario C b1	C b2	Cv1	σ	Cw2	Cw3
	nominal 0.1355 0.622 7.1	0.67 0.3	2.0
	S1	0.1354 0.5377 7.1833 0.8009 0.3664 2.0619
	S2	0.1185 0.6152 6.881 0.602 0.349 1.8144
	S3	0.1264 0.6302 7.1017 0.5711 0.2702 2.4675
	S4	0.144 0.5629 7.1384 0.7964 0.3731 2.4636
	S5	0.1299 0.7561 7.6842 0.746 0.3682 1.7774
	S6	0.1305 0.5678 6.6114 0.658 0.2931 1.8615
	S7	0.1606 0.6559 9.0895 0.8144 0.3153 2.3111
	S8	0.1048 0.519 7.7523 0.6889 0.2468 2.047
	S9	0.1353 0.5042 8.5117 0.5256 0.3033 2.054
	S10	0.1693 0.6837 9.0427 0.7917 0.2384 2.0496
	S11	0.1286 0.5441 6.8871 0.7676 0.3671 1.5532
	S12	0.1214 0.7211 7.022 0.5924 0.257 1.641
	S13	0.1062 0.6327 7.9545 0.6648 0.2601 1.5889
	S14	0.0983 0.5603 6.5696 0.5643 0.2938 2.1707

Table 4

 4 The posterior model probabilities p(Mi|z k ) when M consists of the k -ε and SA model.We utilized a uniform prior p(Mi).

	Scenario k -ε	SA
	S1	6.54763e-02 9.34524e-01
	S2	1.31256e-01 8.68744e-01
	S3	9.98800e-01 1.19985e-03
	S4	2.09486e-01 7.90514e-01
	S5	9.95836e-01 4.16372e-03
	S6	2.36294e-01 7.63706e-01
	S7	1.00000e+00 4.51629e-07
	S8	1.00000e+00 4.34707e-11
	S9	1.00000e+00 0.00000e+00
	S10	1.00000e+00 0.00000e+00
	S11	9.99867e-01 1.32751e-04
	S12	9.98324e-01 1.67614e-03
	S13	0.00000e+00 1.00000e+00
	S14	0.00000e+00 1.00000e+00

Table 5

 5 The posterior model probabilities p(Mi|z k ) when M consists of the k -ε, k -ω and SA models. We utilized a uniform prior p(Mi)

		Scenario k -ε	k -ω	SA
		S1	1.04824e-02 6.76163e-03 9.82756e-01
		S2	5.15244e-02 8.73518e-02 8.61124e-01
		S3	4.11104e-01 5.88478e-01 4.18024e-04
		S4	7.15142e-02 4.21388e-02 8.86347e-01
		S5	9.73137e-01 2.57038e-02 1.15898e-03
		S6	1.04312e-01 1.26575e-03 8.94423e-01
	.	S7	1.00000e+00 1.94074e-10 2.35604e-08
		S8	9.99999e-01 9.14473e-07 1.98968e-11
		S9	9.98710e-01 1.29035e-03 2.34332e-21
		S10	1.00000e+00 0.00000e+00 4.95032e-13
		S11	3.78159e-01 6.21680e-01 1.61533e-04
		S12	9.42438e-01 5.36331e-02 3.92870e-03
		S13	0.00000e+00 3.69311e-01 6.30689e-01
		S14	0.00000e+00 9.74174e-01 2.58255e-02

Table 6

 6 Flow descriptions of training database. Identification numbers are taken from the source[START_REF] Coles | Computation of turbulent boundary layers[END_REF]. The 'Type' column describes the pressure gradient, which ranges from favourable to strongly adverse.

	Scenario Identification Type	Description
	S1	1400	Zero	Equilibrium boundary layer at constant pressure
	S2	1300	Fav	Near-equilibrium boundary layer in moderate negative pressure gradi-
				ent
	S3	2700	Fav	Equilibrium boundary layer in mild negative pressure gradient
	S4	6300	Fav	Near-equilibrium boundary layer growing beneath potential flow on
				model spillway
	S5	1100	Mild adv Boundary layer in diverging channel
	S6	2100	Div	Boundary layer on large airfoil-like body; pressure gradient mildly neg-
				ative
	S7	2500	Mild adv Equilibrium boundary layer in mild positive pressure gradient
	S8	2400	Div	Initial equilibrium boundary layer in moderate positive pressure gradi-
				ent; pressure gradient abruptly decreases to zero, and flow relaxes to
				new equilibrium
	S9	3300	Mod adv Boundary layer, initially at constant pressure, developing into equilib-
				rium flow in moderate positive pressure gradient
	S10	0141	Str adv Boundary-layer with strong adverse pressure gradient, source [44]
	S11	1200	Str adv Boundary layer in diverging channel with eventual separation
	S12	4400	Str adv Boundary layer in strong positive pressure gradient
	S13	2133	Str adv Boundary layer on large airfoil-like body; pressure gradient strongly
				positive
	S14	2134	Str adv Boundary layer on large airfoil-like body; pressure gradient strongly
				positive, close to separation

Acknowledgment

This research was partially supported by the European Union's Seventh Framework Programme (2007-2013) under grant agreement 605036, UMRIDA project and by French ANR under grant agreement 11-MONU-0008.

APPENDIX A: TABLES OF MAP ESTIMATES, POSTERIOR MODEL PROBABILITIES, AND CALIBRATION CASES

Here we report the MAP estimates of the closure coefficients for the three turbulence models of Section 2 (tables 1-3, as well as the corresponding posterior model probabilities (Tables 4,5). We also give a description of the calibration cases.

In the case of the kε mode we use (A1)-(A2) to fix the value of σ ε and C ε1 , where we set P/ε = 2.09 for the latter, see [START_REF] Edeling | Bayesian estimates of parameter variability in the k-ε turbulence model[END_REF] for a discussion. We use other well-known algebraic relations for the kω and SA models in (A3) and (A4) respectively [START_REF] Wilcox | Turbulence Modeling for CFD[END_REF].