The calibration of a computer code is a process that reduces the uncertainty of model parameters by matching the code's predictions to experimental observations of a quantity of interest. A more faithful representation of the global uncertainty is achieved by including a model error term, a discrepancy between the physical system and the computer code. The recently proposed complete maximum a posteriori (CMP) method is able to infer both a posterior distribution of the model parameters and a model error term, improving upon traditional frameworks. On the other hand, the CMP method relies on an optimization step which increases the cost of complex calibration problems. This article proposes a surrogate-based strategy to reduce the computational cost of the CMP method. First, we build a surrogate model of the model error's hyperparameters using Gaussian processes. Second, we propose an iterative algorithm that builds a training set in regions of the parameter space that are more likely, reducing the overall cost of the algorithm and improving the accuracy of the surrogate. The proposed strategy is applied to four different examples, including a design problem in solid mechanics and a complex test case in fluid dynamics. The results show that the proposed strategy is able to accelerate the CMP method without losing accuracy, making it suitable for real-world applications. In an industrial application, we demonstrate a speed-up of almost 100 compared to the original CMP method.
Quantile regression aims to estimate the conditional quantiles of a response variable from observed data. In a Bayesian setting, Gaussian process quantile regression provides uncertainty quantification but faces significant computational challenges due to the nonconjugacy of the asymmetric Laplace likelihood and the cost of posterior inference. We develop a sparse Gaussian process framework in which the quantile function is represented through a reduced set of inducing variables and posterior inference is performed using a Laplace approximation. A decomposition of the predictive uncertainty into conditional-prior and posterior-induced variance components is then exploited to drive two complementary adaptive mechanisms: inducing-input infilling and data acquisition. These mechanisms are combined within a sequential algorithm that allocates computational effort toward the dominant source of predictive uncertainty and adaptively controls model complexity. Numerical experiments on benchmark problems demonstrate the accuracy of the Laplace approximation, the benefits of variance-based inducing-input placement, and the effectiveness of the proposed sequential enrichment strategy compared with predefined data-acquisition strategies.
This paper addresses the challenge of dimension reduction (DR) in Bayesian inference of high-resolution two-or three-dimensional fields, where a priori parametrizations require a large number of terms. The underlying idea is common to state-of-the-art methods in which the parameter space is decomposed into two subspaces, one informed by the likelihood and one constrained by the prior. DR techniques generally use gradient information from the log-likelihood to derive the corresponding subspaces. However, the gradient may be unavailable or expensive to compute accurately, for instance in the case of simulation-based inference. Inspired by approaches based on likelihood-informed subspaces, we develop a new DR method tailored for settings where gradient computation is not feasible. More specifically, we propose a gradient-free indicator for determining whether a direction is informed by the data. This indicator is derived from the posterior-to-prior covariance ratio introduced in Spantini et al. (2015). We show that, in the linear Gaussian case, this indicator combined with an approximate likelihood leads to a better posterior approximation. The method is then extended to nonlinear cases, and strategies to approximate the posterior covariance are detailed. We demonstrate the effectiveness of this DR through two high-dimensional inference problems arising from groundwater and atmospheric applications.
We propose a multifidelity formulation for generating cokriging surrogates of complex physics models. First, we show that the standard autoregressive recursive approach may be subject to substantial limitations due to possible modeler’s biases/errors. These are inherent to the process of establishing a nested hierarchy concerning the alleged fidelity of the available models. The formulation we propose mitigates this issue. At each hierarchy level, the predictor consists of a linear combination of all previous levels instead of just the underlying one. The methodology implies a slightly higher training cost for the surrogate. However, the higher training cost is acceptable, considering the effort typically required to generate data in aerospace applications. A few artificial tests, including the optimization of a two-dimensional airfoil, illustrate strengths and weaknesses of the approach.
Large cetaceans face several anthropogenic threats. Among these, collisions are a major cause of anthropogenic mortality. Assessing and limiting their impact on populations is essential, as these species play an essential ecological role. All types of vessels, including offshore racing vessels, can collide with cetaceans. When a collision occurs between an offshore racing vessel and a large cetacean, the consequences are severe for both the whale, which is often injured or even killed and the vessel, which can suffer severe damage and be forced to withdraw from the race. Our study aimed to develop an encounter model that takes the characteristics of both cetaceans and racing vessels into account to estimate the number of encounters along vessel routes. The model was applied to three different routes commonly used in offshore racing: the first between Newport, USA and Skagen, Denmark; the second between Dover, England and the Gibraltar Strait; and the third between the Gibraltar Strait and Genoa, Italy. The number of encounters was estimated to be 1.7 for Route 1, 4.1 for Route 2 and 2.6 for Route 3. The model was also used to estimate the impact of routing vessels away from any exclusion zones that may be established in areas of high cetacean abundance. This routing could significantly reduce the number of encounters and offer potential solutions to reduce collisions between cetaceans and all types of vessels. The issue of collisions is becoming increasingly important and requires the development of methods to reduce the number of collisions worldwide.
Computer models are widely used for the prediction of complex physical phenomena. Based on observations of these physical phenomena, it is possible to calibrate the model parameters. In most cases, such computer models are misspecified, and the calibration process must be improved by including a model error term. The model error hyperparameters are, however, rarely learned jointly with the model parameters to reduce the dimensionality of the problem. Sequential and nonsequential approaches have been introduced to estimate the hyperparameters. The former, such as the Kennedy and O'Hagan (KOH) framework, estimates the model error hyperparameters before calibrating the model parameters. The latter, such as the full maximum a posteriori (FMP), introduces a functional dependence between the model parameters and the model error hyperparameters. Despite being more reliable in some cases (bimodality, e.g.), the FMP method still fails to estimate correctly the posterior distribution shape. This work proposes a new methodology for treating the model error term in computer code calibration. It builds upon the KOH and FMP framework. Called the complete maximum a posteriori (CMP) method, it provides a closed-form expression for the marginalization integral over the model error hyperparameters, significantly reducing the dimensionality of the calibration problem. Such expression relies on a set of assumptions that are more general and less stringent than the ones usually employed. The CMP method is applied to four examples of increasing complexity, from elementary to real fluid dynamics problems, including or not bimodality. Compared to the true reference solution and unlike the KOH and FMP, the CMP method correctly captures the shape of the posterior distribution, including all modes and their weights. Moreover, it provides an accurate estimate of the distribution tails.
A Bayesian inference approach for inferring the source of marine pollution released from a moving source in an uncertain flow field is proposed. A Markov Chain Monte Carlo (MCMC) algorithm is developed and applied for inferring single and multiple release events from vessels moving at known velocity along a predefined path in the Mediterranean Sea. The likelihood is based on a logistic regression cost function that measures the discrepancy between the modeled spill distribution and a binary representation of the observed images. We assess the performance of the proposed methodology using a synthetic release scenario employing realistic ocean currents to drive a stochastic Lagrangian Particle Tracking (LPT) algorithm to generate a probabilistic representation of the spill distribution. The MCMC algorithm employs an adaptive scheme to robustly ensure convergence and well-mixed chains. The proposed Bayesian framework is tested by inferring the location, or injection time, and relative contributions of single and multiple moving sources, contributing to separate and common observation patches, with a focus on various scenarios that demonstrate the efficiency of our sampling algorithm. The performance of the proposed framework was further assessed by comparing the model predictions with the most probable release parameters predicted by a global optimization algorithm.
This paper proposes an effective treatment of hyperparameters in the Bayesian inference of a scalar field from indirect observations. Obtaining the joint posterior distribution of the field and its hyperparameters is challenging. The infinite dimensionality of the field requires a finite parametrization that usually involves hyperparameters to reflect the limited prior knowledge. In the present work, we consider a Karhunen-Lo{è}ve(KL) decomposition for the random field and hyperparameters to account for the lack of prior knowledge of its autocovariance function. The hyperparameters must be inferred. To efficiently sample jointly the KL coordinates of the field and the autocovariance hyperparameters, we introduce a change of measure to reformulate the joint posterior distribution into a hierarchical Bayesian form. The likelihood depends only onthe field's coordinates in a fixed KL basis, with a prior conditioned on the hyperparameters. We exploit this structure to derive an efficient Markov Chain Monte Carlo (MCMC) sampling scheme based on an adapted Metropolis-Hasting algorithm. We rely on surrogate models (Polynomial Chaos expansions) of the forward model predictions to further accelerate the MCMC sampling. A first application to a transient diffusionproblem shows that our method is consistent with other approaches based on a change of coordinates (Sraj et al., 2016). A second application to a seismic traveltime tomography highlights the importance of inferring the hyperparameters. A third application to a 2D anisotropic groundwater flow problem illustrates the method on a more complex geometry.
Oil spills at sea pose a serious threat to coastal environments. Identifying oil pollution sources could help to investigate unreported spills, and satellite imagery can be an effective tool for this purpose. We present a Bayesian approach to estimate the source parameters of a spill from contours of oil slicks detected by remotely sensed images. Five parameters of interest are estimated: the 2D coordinates of the source of release, the time and duration of the spill, and the quantity of oil released. Two synthetic experiments of a spill released from a fixed point source are investigated, where a contour is fully observed in the first case, while two contours are partially observed at two different times in the second. In both experiments, the proposed method is able to provide good estimates of the parameters along with a level of confidence reflected by the uncertainties within.
This study presents a novel approach to applying data assimilation techniques for particle-based simulations using the Ensemble Kalman Filter. While data assimilation methods have been effectively applied to Eulerian simulations, their application in Lagrangian solution discretizations has not been properly explored. We introduce two specific methodologies to address this gap. The first methodology employs an intermediary Eulerian transformation that combines a projection with a remeshing process. The second is a purely Lagrangian scheme designed for situations where remeshing is not appropriate. The second is a purely Lagrangian scheme that is applicable when remeshing is not adapted. These methods are evaluated using a one-dimensional advection-diffusion model with periodic boundaries. Performance benchmarks for the one-dimensional scenario are conducted against a grid-based assimilation filter Subsequently, assimilation schemes are applied to a non-linear two-dimensional incompressible flow problem, solved via the Vortex-In-Cell method. The results demonstrate the feasibility of applying these methods in more complex scenarios, highlighting their effectiveness in both the one-dimensional and two-dimensional contexts.
We investigate a computer model calibration technique inspired by the well-known Bayesian framework of Kennedy and O'Hagan (KOH). We tackle the full Bayesian formulation where model parameter and model discrepancy hyperparameters are estimated jointly and reduce the problem dimensionality by introducing a functional relationship that we call the full maximum a posteriori (FMP) method. This method also eliminates the need for a true value of model parameters that caused identifiability issues in the KOH formulation. When the joint posterior is approximated as a mixture of Gaussians, the FMP calibration is proven to avoid some pitfalls of the KOH calibration, namely missing some probability regions and underestimating the posterior variance. We then illustrate two numerical examples where both model error and measurement uncertainty are estimated together. Using the solution to the full Bayesian problem as a reference, we show that the FMP results are accurate and robust, and avoid the need for high-dimensional Markov chains for sampling.
A preconditioning strategy is proposed for the iterative solve of large numbers of linear systems with parameter-dependent matrix and right-hand side which arise during the computation of solution statistics of stochastic elliptic partial differential equations with random and spatially variable coefficients sampled by Monte Carlo. Building on the assumption that a truncated Karhunen-Loève expansion of a known transform of the random coefficient is available, we introduce a compact approximation of the random coefficient in the form of a Voronoi quantizer. The number of Voronoi cells, each of which is represented by a centroidal coefficient, is set to the prescribed number of preconditioners. Upon sampling the random coefficient, the linear system assembled with a given realization of the coefficient is solved using a Krylov subspace iterative solver with the preconditioner whose centroidal coefficient is the closest to the realization. We consider different ways to define and obtain the centroidal coefficients, and we investigate the properties of the induced preconditioning strategies in terms of average number of solver iterations for sequential simulations, and of load balancing for parallel simulations. Another approach, which is based on deterministic grids on the system of stochastic coordinates of the truncated representation of the random coefficient, is proposed with a stochastic dimension that increases with the number of preconditioners. This approach allows to bypass the need for preliminary computations in order to determine the optimal stochastic dimension of the truncated approximation of the random coefficient for a given number of preconditioners.
To support accidental spill rapid response efforts, oil spill simulations may generally need to account for uncertainties concerning the nature and properties of the spill, which compound those inherent in model parameterizations. A full detailed account of these sources of uncertainty would however require prohibitive resources needed to sample a large dimensional space. In this work, a variance-based sensitivity analysis is conducted to explore the possibility of restricting a priori the set of uncertain parameters, at least in the context of realistic simulations of oil spills in the Red Sea region spanning a two-week period following the oil release. The evolution of the spill is described using the simulation capabilities of Modelo Hidrodinâmico, driven by high-resolution metocean fields of the Red Sea (RS) was adopted to simulate accidental oil spills in the RS. Eight spill scenarios are considered in the analysis, which are carefully selected to account for the diversity of metocean conditions in the region. Polynomial chaos expansions are employed to propagate parametric uncertainties and efficiently estimate variance-based sensitivities. Attention is focused on integral quantities characterizing the transport, deformation, evaporation and dispersion of the spill. The analysis indicates that variability in these quantities may be suitably captured by restricting the set of uncertain inputs parameters, namely the wind coefficient, interfacial tension, API gravity, and viscosity. Thus, forecast variability and confidence intervals may be reasonably estimated in the corresponding four-dimensional input space.
A second-order accurate time-stepping scheme for solving a time-fractional Fokker–Planck equation of order $\alpha \in (0, 1)$, with a general driving force, is investigated. A stability bound for the semidiscrete solution is obtained for $\alpha \in (1/2,1)$ via a novel and concise approach. Our stability estimate is $\alpha $-robust in the sense that it remains valid in the limiting case where $\alpha $ approaches $1$ (when the model reduces to the classical Fokker–Planck equation), a limit that presents practical importance. Concerning the error analysis, we obtain an optimal second-order accurate estimate for $\alpha \in (1/2,1)$. A time-graded mesh is used to compensate for the singular behavior of the continuous solution near the origin. The time-stepping scheme scheme is associated with a standard spatial Galerkin finite element discretization to numerically support our theoretical contributions. We employ the resulting fully discrete computable numerical scheme to perform some numerical tests. These tests suggest that the imposed time-graded meshes assumption could be further relaxed, and we observe second-order accuracy even for the case $\alpha \in (0,1/2]$, that is, outside the range covered by the theory.
Static Velocity Prediction Programs (VPP) are standard tools in sailing yachts’ design and performance assessment. Predicting the maximal steady velocity of a yacht involves resolving constrained optimization problems. These problems have a prohibitive computational cost when using high-fidelity global modeling of the yacht. This difficulty has motivated the introduction of modular approaches, decomposing the global model into subsystems modeled independently and approximated by surrogate models (response surfaces). The maximum boat speed for prescribed conditions solves an optimization problem for the trimming parameters of the model constrained by compatibility conditions between the subsystems’ surrogate solution (e.g., the yacht equilibrium). The accuracy of the surrogates is then critical for the quality of the resulting VPP. This paper relies on Gaussian Process (GP) models of the subsystems and introduces an original sequential Active Learning Method (ALM) for their joint construction. Our ALM exploits the probabilistic nature of the GP models to decide the enrichment of the training sets using an infilling criterion that combines the predictive uncertainty of the surrogate models and the likelihood of equilibrium at every input point. The resulting strategy enables the concentration of the computational effort around the manifolds where equilibrium is satisfied. The results presented compare ALM with a standard (uninformed) Quasi-Monte Carlo method, which samples the input space of the subsystems uniformly. ALM surrogates have higher accuracy in the equilibrium regions for equal construction cost, with improved mean prediction and reduced prediction uncertainty. We further investigate the effect of the prediction uncertainty on the numerical VPP and in a routing problem.
This paper compares risk-averse optimization methods to address the self-scheduling and market involvement of a virtual power plant (VPP). The decision-making problem of the VPP involves uncertainty in the wind speed and electricity price forecast. We focus on two methods: risk-averse two-stage stochastic programming (SP) and two-stage adaptive robust optimization (ARO). We investigate both methods concerning formulations, uncertainty and risk, decomposition algorithms, and their computational performance. To quantify the risk in SP, we use the conditional value at risk (CVaR) because it can resemble a worst-case measure, which naturally links to ARO. We use two efficient implementations of the decomposition algorithms for SP and ARO; we assess (1) the operational results regarding first-stage decision variables, estimate of expected profit, and estimate of the CVaR of the profit and (2) their performance taking into consideration different sample sizes and risk management parameters. The results show that similar first-stage solutions are obtained depending on the risk parameterizations used in each formulation. Computationally, we identified three cases: (1) SP with a sample of 500 elements is competitive with ARO; (2) SP performance degrades comparing to the first case and ARO fails to converge in four out of five risk parameters; (3) SP fails to converge, whereas ARO converges in three out of five risk parameters. Overall, these performance cases depend on the combined effect of deterministic and uncertain data and risk parameters. Summary of Contribution: The work presented in this manuscript is at the intersection of operations research and computer science, which are intrinsically related with the scope and mission of IJOC. From the operations research perspective, two methodologies for optimization under uncertainty are studied: risk-averse stochastic programming and adaptive robust optimization. These methodologies are illustrated using an energy scheduling problem. The study includes a comparison from the point of view of uncertainty modeling, formulations, decomposition methods, and analysis of solutions. From the computer science perspective, a careful implementation of decomposition methods using parallelization techniques and a sample average approximation methodology was done . A detailed comparison of the computational performance of both methods is performed. Finally, the conclusions allow establishing links between two alternative methodologies in operations research: stochastic programming and robust optimization.
In this work, we calibrate a carbon nitridation model for a broad span of surface temperatures from existing plasma wind tunnel measurements by accounting for experimental and parametric uncertainties. A chemical non-equilibrium stagnation line model is proposed to simulate the experiments and obtain recession rates and CN densities, the measured model outputs. First, we establish the influence of the experimental boundary conditions and nitridation parameters on the simulated observations through a sensitivity analysis. Results show that such quantities are mostly affected by the efficiency of nitridation reactions at the gas-surface interface. We then perform model calibrations for each experimental condition and compare them based on the experimental data used. This allows us to check the consistency of the experimental dataset. Using only the trustworthy experimental data, we perform a calibration of Arrhenius law parameters for nitridation efficiencies considering all available experimental conditions jointly, allowing us to compute nitridation efficiencies even for surface temperatures for which there are no reliable experimental data available. The stochastic Arrhenius law agrees well with most of the data in the literature. This result constitutes the first nitridation model extracted from plasma wind tunnel experiments with accurate uncertainty estimates.
Obtaining accurate high-resolution representations of model outputs is essential to describe the system dynamics. In general, however, only spatially- and temporally-coarse observations of the system states are available. These observations can also be corrupted by noise. Downscaling is a process/scheme in which one uses coarse scale observations to reconstruct the high-resolution solution of the system states. Continuous Data Assimilation (CDA) is a recently introduced downscaling algorithm that constructs an increasingly accurate representation of the system states by continuously nudging the large scales using the coarse observations. We introduce a Discrete Data Assimilation (DDA) algorithm as a downscaling algorithm based on CDA with discrete-in-time nudging. We then investigate the performance of the CDA and DDA algorithms for downscaling noisy observations of the Rayleigh-Bénard convection system in the chaotic regime. In this computational study, a set of noisy observations was generated by perturbing a reference solution with Gaussian noise before downscaling them. The downscaled fields are then assessed using various error- and ensemble-based skill scores. The CDA solution was shown to converge towards the reference solution faster than that of DDA but at the cost of a higher asymptotic error. The numerical results also suggest a quadratic relationship between the ℓ 2 error and the noise level for both CDA and DDA. Cubic and quadratic dependences of the DDA and CDA expected errors on the spatial resolution of the observations were obtained, respectively.
The tangent linear approximation (TLA) developed in Almohammadi et al. (Combust. Flame 230, 111426) is extended to estimate the sensitivity of the ignition delay time with respect to species enthalpies and entropies. The proposed method relies on integrating the linearized system of equations governing the evolution of the state vector's partial derivatives with respect to variations in thermodynamic parameters. The sensitivity of the ignition delay time is estimated through a linearized approximation of a temperature functional. The TLA approach is applied to three gas mixtures, H-2 , n-butanol, and iso-octane, reacting in air under adiabatic, constant-volume conditions. The numerical experiments indicate that the linearized approximation of the ignition delay time's sensitivity is in excellent agreement with the finite difference estimates. This is also the case for sensitivity estimates obtained using the TLA approach. Further, significant computational speed-ups are achieved with the TLA approach, and the method scales well with the number of perturbed parameters. In the case of the H-2 mechanism, TLA is about ten times faster than finite differences, and this enhancement becomes even more substantial when more complex mechanisms are considered. (C) 2021 The Combustion Institute. Published by Elsevier Inc. All rights reserved.