Olivier Ezvan
email: olivier.ezvan@univ-eiffel.fr

Christian Soize
email: christiansoize@univ-eiffel.fr

Christophe Desceliers
email: christophe.desceliers@univ-eiffel.fr

Roger Ghanem
email: ghanem@usc.edu

Updating an uncertain and expensive computational model in structural dynamics based on one single target FRF using a probabilistic learning tool

Keywords: updating, computational dynamics, uncertainties, machine learning, probabilistic learning, PLoM

The paper presents an appropriate and efficient methodology for updating the control parameters of very large uncertain computational models, which are used for analyzing the linear vibrations in the frequency domain of highly complex structures for which there are an enormous number of intertwined local and global elastic structural modes in the broad frequency band of analysis. Moreover, the numerical cost of a single evaluation of the frequency response functions with the computational model is assumed to be very high and only one experimental frequency response function is available as a target. For decreasing the numerical cost of this challenging problem, a parameterized reduced-order model is constructed. Nevertheless, this reduction is not sufficient to be able to solve the non-convex optimization problem related to the updating. Consequently, for avoiding the call to the computational model, the probabilistic learning on manifolds is used for generating a learned set from a training set, which, coupled with conditional statistics, allows the evaluation of the cost function without calling the computational model. A numerical illustration is presented for validating the proposed methodology.

Introduction (i) Objectives of the paper.

The framework of the developments presented in this paper is that of linear vibrations of highly complex structures for which the uncertain high-dimensional model (HDM) can have 10 8 , even 10 9 degrees of freedom (dof). The number of frequency points in the frequency band to analyze the structure is of the order of several thousands, even 10 4 . The number of elastic structural modes in the frequency band of analysis is of the order of 10 6 with the simultaneous presence of local and global elastic structural modes. This dynamic system depends on a control Figure 1: Left: spent nuclear fuel canister loaded with two fuel assemblies -Right: spent nuclear fuel assembly with and without outer shell (pictures from [START_REF] Ezvan | Multiscale modal analysis of fully-loaded spent nuclear fuel canisters[END_REF]) parameter made up of several hundred parameters, mainly related to local stiffnesses of joints inside the structure. The problem consists in identifying an updated value of the control parameter on the basis of only one available experimental frequency response function (FRF) between one input dof and one output dof located on the boundary of the structure. This updating requires to solve an optimization problem on an admissible set of the control parameter and by using an adapted formulation of the cost function. The admissible set is the support of the prior probability distribution of the uncertain control parameter. This situation, which we have just described, is that of the fully-loaded spent nuclear fuel canisters studied in [START_REF] Ezvan | Dominant vibration modes for broadband frequency analysis of multiscale structures with numerous local vibration modes[END_REF][START_REF] Ezvan | Multiscale modal analysis of fully-loaded spent nuclear fuel canisters[END_REF][START_REF] Ezvan | Dominant substructural vibration modes for fully-loaded spent nuclear fuel canisters[END_REF] and depicted in Fig. 1. Those are multiscale structures made up of repeated components at various scales, which require fine mesh discretizations and which also exhibit numerous local structural elastic modes that are associated with the isolated vibrations of each of the numerous components. In [START_REF] Ezvan | Dominant substructural vibration modes for fully-loaded spent nuclear fuel canisters[END_REF], an HDM was constructed based on nearly half-a-million modes from a finite element model of about 130 million dof. Adapted algorithms were developped allowing a single run of the HDM to be performed in about 40 minutes using a cluster of 50 nodes. In this application, it is of particular interest to be able to perform the updating of physical parameters that describe the structural integrity of the internal components and structural connections, which are critical to nuclear safety.

This challenging problem is particularly difficult because the dynamical system is more than under-observed (in fact only one experimental FRF is available) and one could even think that it is badly posed. However, the FRF is known over a wide frequency band, at a large number of frequency points, which partially removes the under observability. The difficulty of the considered problem also comes from the fact that the numerical cost of one evaluation of a few FRFs with the HDM, for a single fixed value of the control parameter, and for all the frequency points, is huge, even unrealistic.

It is therefore necessary to use a methodology that makes it possible to circumvent these difficulties by considering several ingredients. The first ingredient consists in reducing the numerical cost of a single evaluation of the HDM by constructing a parameterized reduced-order model (ROM) as a function of the control parameter. Such a reduction can be performed either using the proper orthogonal decomposition (POD) of the responses of the dynamic system computed for a given small set of the control parameter values, or using the modal synthesis approach with the parameterized elastic structural modes, formulated with the classical approach or introducing a multi-frequency model reduction, with or without dynamic substructuring.

The POD approach is very well adapted to nonlinear dynamic systems formulated in the time domain and allows the construction of a basis that is independent of the control parameter, vector basis that is used for projecting the HDM in order to construct a parameterized ROM (see for instance [START_REF] Grepl | Efficient reduced-basis treatment of nonaffine and nonlinear partial differential equations[END_REF][START_REF] Nguyen | An efficient reduced-order modeling approach for non-linear parametrized partial differential equations[END_REF][START_REF] Chaturantabut | Nonlinear model reduction via discrete empirical interpolation[END_REF][START_REF] Degroote | Interpolation among reduced-order matrices to obtain parameterized models for design, optimization and probabilistic analysis[END_REF][START_REF] Carlberg | Efficient non-linear model reduction via a least-squares petrov-galerkin projection and compressive tensor approximations[END_REF][START_REF] Carlberg | A low-cost, goal-oriented compact proper orthogonal decomposition basis for model reduction of static systems[END_REF][START_REF] Amsallem | Nonlinear model order reduction based on local reduced-order bases[END_REF][START_REF] Carlberg | The gnat method for nonlinear model reduction: effective implementation and application to computational fluid dynamics and turbulent flows[END_REF][START_REF] Zahr | Progressive construction of a parametric reduced-order model for pde-constrained optimization[END_REF][START_REF] Amsallem | Design optimization using hyper-reduced-order models[END_REF]).

Since the problem under consideration is linear and is directly formulated in the frequency domain, we prefer to use the modal synthesis approach, which is based on the elastic structural modes that depend on the control parameter, and which is more accurate in the considered case. For the computational aspects related to the elastic structural modes and modal synthesis, we refer the reader, for instance, to [START_REF] Bathe | Numerical Methods in Finite Element Analysis[END_REF][START_REF] Geradin | Mechanical Vibrations[END_REF][START_REF] Ohayon | Advanced Computational Vibroacoustics -Reduced-Order Models and Uncertainty Quantification[END_REF] and for the continuum-mechanics aspects to [START_REF] Ohayon | Structural Acoustics and Vibration: Mechanical Models, Variational Formulations and Discretization[END_REF]. Due to the presence of a very large number of global and local elastic structural modes in the frequency band of analysis, a multi-frequency model reduction (see [START_REF] Ezvan | Multilevel reduced-order computational model in structural dynamics for the lowand medium-frequency ranges[END_REF][START_REF] Ezvan | Multilevel model reduction for uncertainty quantification in computational structural dynamics[END_REF][START_REF] Reyes | Multi-frequency model reduction for uncertainty quantification in computational vibroacoutics[END_REF]) or a multiscale modal analysis [START_REF] Ezvan | Multiscale modal analysis of fully-loaded spent nuclear fuel canisters[END_REF] could also be used.

For the structures considered, the techniques of dynamic substructuring are often very effective. As explained in the Review Paper [START_REF] Ohayon | Clarification about component mode synthesis methods for substructures with physical flexible interfaces[END_REF], the concept of substructures was first introduced by Argyris and Kelsey in 1959 [START_REF] Argyris | The analysis of fuselages of arbitrary cross-section and taper: A dsir sponsored reserach program on the development and application of the matrix force method and the digital computer[END_REF] and by Przemieniecki in 1963 [START_REF] Przemieniecki | Matrix structural analysis of substructures[END_REF] who introduced the static boundary functions. This work was extended by Guyan, Irons, and Hurty [START_REF] Guyan | Reduction of stiffness and mass matrices[END_REF][START_REF] Irons | Structural eigenvalue problems -elimination of unwanted variables[END_REF][START_REF] Hurty | Vibrations of structural systems by component mode synthesis[END_REF][START_REF] Hurty | Dynamic analysis of structural systems using component modes[END_REF]. Finally, Craig and Bampton in 1968 [START_REF] Craig | Coupling of substructures for dynamic analyses[END_REF] adapted the Hurty method in order to represent each substructure of the same manner consisting in using the elastic structural modes of the substructure with fixed geometrical interface and the static boundary functions on its geometrical interface. These aspects in the formulation of the continuum mechanics can be found in [START_REF] Ohayon | Structural Acoustics and Vibration: Mechanical Models, Variational Formulations and Discretization[END_REF].

For the complex structures under consideration that we have described at the beginning of this section, the simultaneous use of the substructuring technique coupled with a multi-frequency model reduction would be efficient. Nevertheless, in order to simplify the presentation of this paper, we have limited the developments to the framework of the classical modal synthesis without using substructuring technique and multi-frequency model reduction. As the reader will be able to see, the implementation of these two ingredients in the methodology proposed could be performed without significant difficulties.

Even using the parameterized ROM, the numerical cost of a single evaluation of the cost function of the optimization problem, for updating the HDM, is far too high. Indeed, any optimization algorithm for this type of nonconvex problem will have to make a large number of evaluations of the cost function or even of its gradient. This means that it is necessary to build a "surrogate model" that relates the control parameter to the frequency-sampled FRF. Such surrogate models can be chosen in the class of deterministic or statistical representations [START_REF] Bhosekar | Advances in surrogate based modeling, feasibility analysis, and optimization: A review[END_REF]. There are many methods in the class of deterministic representations based on polynomial and spline regressions or on the use of meshless representations that make it possible to calculate partial derivatives (see for instance [START_REF] Nayroles | Generalizing the finite element method: diffuse approximation and diffuse elements[END_REF][START_REF] Belytschko | Meshless methods: an overview and recent developments[END_REF][START_REF] Duarte | H-p clouds, an h-p meshless method, Numerical Methods for Partial Differential Equations[END_REF][START_REF] Breitkopf | Explicit form and efficient computation of MLS shape functions and their derivatives[END_REF][START_REF] Rassineux | Surface remeshing by local Hermite diffuse interpolation[END_REF][START_REF] Zhang | Meshless methods based on collocation with radial basis functions[END_REF]). Nevertheless, these methods are not well adapted neither to the high-dimensional problems nor to the statistical framework. For the class of statistical representations, and in particular in the framework of uncertainty quantification [START_REF] Ghanem | Handbook of Uncertainty Quantification[END_REF][START_REF] Soize | Uncertainty Quantification[END_REF], there are many possibilities. The use of Gaussian processes to build surrogate models, such as kriging methods, is very popular [START_REF] Kleijnen | Kriging metamodeling in simulation: A review[END_REF][START_REF] Dubourg | Reliability-based design optimization using kriging surrogates and subset simulation[END_REF][START_REF] Kersaudy | A new surrogate modeling technique combining kriging and polynomial chaos expansions-application to uncertainty analysis in computational dosimetry[END_REF] but these may have limitations induced by high dimensions and by the geometric complexity of the manifold that we seek to represent by a surrogate model. Certainly, another popular method, which is very efficient for many challenging problems, is the use of the polynomial chaos expansion (see for instance [START_REF] Ghanem | Stochastic Finite Elements: a Spectral Approach[END_REF][START_REF] Xiu | The Wiener-Askey polynomial chaos for stochastic differential equations[END_REF][START_REF] Soize | Physical systems with random uncertainties: chaos representations with arbitrary probability measure[END_REF][START_REF] Wan | Multi-element generalized polynomial chaos for arbitrary probability measures[END_REF][START_REF] Soize | Computational aspects for constructing realizations of polynomial chaos in high dimension[END_REF][START_REF] Blatman | Adaptive sparse polynomial chaos expansion based on least angle regression[END_REF][START_REF] Perrin | Identification of polynomial chaos representations in high dimension from a set of realizations[END_REF][START_REF] Tipireddy | Basis adaptation in homogeneous chaos spaces[END_REF][START_REF] Babaei | Robust optimization of subsurface flow using polynomial chaos and response surface surrogates[END_REF][START_REF] Soize | Polynomial chaos expansion of a multimodal random vector[END_REF][START_REF] Abraham | A robust and efficient stepwise regression method for building sparse polynomial chaos expansions[END_REF][START_REF] Luthen | Sparse polynomial chaos expansions: Literature survey and benchmark[END_REF][START_REF] Soize | Polynomial chaos representation of databases on manifolds[END_REF]). The Bayesian inference could also be considered as a powerful statistical tool for building a surrogate model defined by a statistical inverse problem (see [START_REF] Bernardo | Bayesian Theory[END_REF][START_REF] Kennedy | Bayesian calibration of computer models[END_REF][START_REF] Congdon | Bayesian Statistical Modelling[END_REF][START_REF] Carlin | Bayesian Methods for Data Analysis[END_REF][START_REF] Marin | Approximate Bayesian computational methods[END_REF][START_REF] Scott | Bayes and big data: The consensus monte carlo algorithm[END_REF][START_REF] Ghanem | Handbook of Uncertainty Quantification[END_REF][START_REF] Soize | Sampling of Bayesian posteriors with a non-Gaussian probabilistic learning on manifolds from a small dataset[END_REF] for general aspects and [START_REF] Marzouk | Stochastic spectral methods for efficient Bayesian solution of inverse problems[END_REF][START_REF] Stuart | Inverse problems: a Bayesian perspective[END_REF][START_REF] Soize | A computational inverse method for identification of non-Gaussian random fields using the Bayesian approach in very high dimension[END_REF][START_REF] Dashti | The Bayesian approach to inverse problems[END_REF][START_REF] Arnst | Itô-SDE MCMC method for Bayesian characterization of errors associated with data limitations in stochastic expansion methods for uncertainty quantification[END_REF][START_REF] Perrin | Adaptive method for indirect identification of the statistical properties of random fields in a Bayesian framework[END_REF][START_REF] Soize | Sampling of Bayesian posteriors with a non-Gaussian probabilistic learning on manifolds from a small dataset[END_REF] for more specific aspects related to statistical inverse problems). In the framework of Artificial Intelligence [START_REF] Russel | Artifical Intelligence, A Modern Approach[END_REF], machine learning tools such as the Neural Networks [START_REF] Gurney | An introduction to Neural Networks[END_REF] and Deep Learning [START_REF] Lecun | Deep learning[END_REF] are very appropriate representations for building surrogate models but generally require big data. For many problems in computational sciences and engineering based on large computational models for simulating physics as a function of a set of control parameters, big data are generally not available due to the high numerical cost for performing a large number of evaluations using the computational model (this is the framework of this paper).

(ii) Strategy used in the paper. In this work, we propose to use the probabilistic learning on manifolds (PLoM) that was proposed in [START_REF] Soize | Data-driven probability concentration and sampling on manifold[END_REF][START_REF] Soize | Probabilistic learning on manifolds[END_REF][START_REF] Soize | Probabilistic learning on manifolds (plom) with partition[END_REF] as a complementary approach to existing methods in machine learning for sampling underlying distributions on manifolds [START_REF] Talwalkar | Large-scale manifold learning[END_REF][START_REF] Gorissen | A surrogate modeling and adaptive sampling toolbox for computer based design[END_REF][START_REF] Öztireli | Spectral sampling of manifolds[END_REF][START_REF] Marzouk | Sampling via measure transport: An introduction[END_REF][START_REF] Parno | Transport map accelerated markov chain Monte Carlo[END_REF][START_REF] Perrin | Data-driven kernel representations for sampling with an unknown block dependence structure under correlation constraints[END_REF][START_REF] Tsilifis | Bayesian adaptation of chaos representations using variational inference and sampling on geodesics[END_REF][START_REF] Kevrekidis | Manifold learning for parameter reduction[END_REF]. It allows for solving unsupervised and supervised problems under uncertainty for which the training sets are small. As we have already explained, this situation is encountered in many problems of physics and engineering sciences with expensive function evaluations. The PLoM was successfully adapted to tackle these challenges for several related problems including nonconvex optimization under uncertainty [START_REF] Ghanem | Probabilistic nonconvex constrained optimization with fixed number of function evaluations[END_REF][START_REF] Ghanem | Optimal well-placement using probabilistic learning[END_REF][START_REF] Soize | Design optimization under uncertainties of a mesoscale implant in biological tissues using a probabilistic learning algorithm[END_REF][START_REF] Ghanem | Design optimization of a scramjet under uncertainty using probabilistic learning on manifolds[END_REF][START_REF] Almeida | A probabilistic learning approach applied to the optimization of wake steering in wind farms[END_REF], updating digital twin [START_REF] Ghanem | Probabilistic learning and updating of a digital twin for composite material systems[END_REF], calculation of Sobol's indices [START_REF] Arnst | Computation of sobol indices in global sensitivity analysis from small data sets by probabilistic learning on manifolds[END_REF], quantifying uncertainties with small data [START_REF] Farhat | Feasible probabilistic learning method for model-form uncertainty quantification in vibration analysis[END_REF][START_REF] Soize | Enhancing model predictability for a scramjet using probabilistic learning on manifolds[END_REF][START_REF] Soize | Probabilistic learning for modeling and quantifying model-form uncertainties in nonlinear computational mechanics[END_REF][START_REF] Guilleminot | Data-driven enhancement of fracture paths in random composites[END_REF]. More recently, extensions have been proposed for analyzing the case of high stochastic dimension by using the PLoM with partition [START_REF] Soize | Probabilistic learning on manifolds (plom) with partition[END_REF] and also for taking into account constraints in the PLoM algorithm defined by experimental statistical moments [START_REF] Soize | Physics-constrained non-Gaussian probabilistic learning on manifolds[END_REF] or by nonlinear partial differential equations [START_REF] Soize | Probabilistic learning on manifolds constrained by nonlinear partial differential equations for small datasets[END_REF], for computational fluid dynamics and nonlinear solid dynamics.

In this paper, we present a methodology for updating the uncertain computational model. As it is well known, the updating problem requires to solve an optimization problem for which several formulations are possible in the statistical context for constructing the "cost function": the least-square method, the maximum likelihood, or the Bayesian approach. The last two cannot easily be used because only one target is available. We address this challenge by using the algorithm that we have presented in [START_REF] Ghanem | Probabilistic nonconvex constrained optimization with fixed number of function evaluations[END_REF] and validated in [START_REF] Ghanem | Optimal well-placement using probabilistic learning[END_REF][START_REF] Soize | Design optimization under uncertainties of a mesoscale implant in biological tissues using a probabilistic learning algorithm[END_REF][START_REF] Ghanem | Design optimization of a scramjet under uncertainty using probabilistic learning on manifolds[END_REF], but, on the one hand by replacing the PLoM algorithm [START_REF] Soize | Data-driven probability concentration and sampling on manifold[END_REF] by the novel algorithm of PLoM [START_REF] Soize | Physics-constrained non-Gaussian probabilistic learning on manifolds[END_REF][START_REF] Soize | Probabilistic learning on manifolds (plom) with partition[END_REF] in order to apply the constraints related to the normalization and on the other hand by constructing a cost function on the basis of conditional statistics estimated with the learned set generated with PLoM. This approach allows for maintaining the number of function evaluations at the level essentially equal to that of the construction of the training set. Finally, the presented numerical application is an illustration, which is representative of all the reported difficulties of the problem posed, but which does not have the complexity of the structures mentioned at the beginning of this paragraph. This choice was made so that the reader can reproduce all the results presented from a simply and completely described dynamic system.

(iii) Novelty of the paper. This paper presents an appropriate and efficient methodology for updating high-dimensional uncertain computational models in the domain of linear vibrations, with a formulation in the frequency domain, for a broad frequency band, for which each evaluation is expensive, and for which only one experimental frequency response function is available as a target. We formulate a nonconvex optimization problem, which is solved using a probabilistic learning tool and conditional statistics, which have been specifically developed for the case of a training set constituted of a small number of points. In this framework, the methodology proposed is novel. It should be noted that we are searching for improving a computational model in the framework of an ill-posed problem, related to the nonconvexity of the optimization problem formulated for the updating, to the under observability, and to the knowledge of only one ex-perimental FRF. Nevertheless, despite this serious difficulty, we propose an approach that allows for updating the model. As explained in the description of the strategy, the updating is first performed by using the training set. Then, alternatively, the PLoM method is used as a resampling tool for trying to improve the updating with respect to the one given by the training set. Through the presented applications, it will be seen that this additional resampling may or may not improve the updating. The main explanation to this phenomena is not due to the capability of PLoM, but is due to the lack of monotonicity of the optimum solution as a function of the size of the training set. This point will be further discussed in the application.

(iv) Organization of the paper. In Section 2, we list all the hypotheses used for the proposed methodology in order to carry out the updating of the computational model. Section 3 is devoted to the definition of the parameterized computational model for the linear vibrations formulated in the frequency domain. Section 4 deals with the construction of the parameterized ROM that is obtained by projection on the parameterized elastic structural modes. The convergence of the modal synthesis is controlled for each value of the control parameter yielding parameterized ROM whose dimension depends on the control parameter. Section 5 begins with the construction of the probabilistic model of the control parameter. Then the training set is built using the stochastic ROM deduced from the parameterized ROM and the probabilistic model of the control parameter. This training set will allow the learned set to be generated, which will be used to solve the optimization problem. Section 6 deals with the formulation of the updating of the stochastic parameterized ROM. The cost function is constructed with a least-square formulation and the use of the statistical conditioning of the random parameterized frequency response function, given a value of the control parameter. For each value of the control parameter, which is proposed by the optimization algorithm, this statistical conditioning is performed with the learned set generated by PLoM and an algebraic formulation derived from the nonparametric statistics. The numerical application is presented in Section 7. All the numerical values of the parameters and data that are necessary for reproducing the results are given. The results obtained are presented, which are accompanied by a discussion. The PLoM algorithm is summarized in Appendix A.

Notations

Lower-case letters such as q or η are deterministic real variables. Boldface lower-case letters such as q or η are deterministic vectors. Upper-case letters such as X or H are real-valued random variables. Boldface upper-case letters such as X or H are vector-valued random variables.

Lower-case letters between brackets such as [x] or [η] are deterministic matrices. Boldface upper-case letters between brackets such as [X] or [H] are matrix-valued random variables.

δ αβ : Kronecker's symbol. ι: imaginary unit in set C. n u : number of dof in the computational model. n ω : number of frequency sampling points. n w : number of control parameters. C n : Hermitian vector space on C of dimension n.

[I n ]: identity matrix in M n . M n,N : set of (n × N) real matrices. M n : set of square (n × n) real matrices. M +0 n : set of positive symmetric (n × n) real matrices. M + n : set of positive-definite symmetric (n × n) matrices. R n : Euclidean vector space on R of dimension n. [x] k j : entry of matrix [x]. [x] T : transpose of matrix [x]. < x , y >: Euclidean or Hermitian inner product. ∥x∥: Euclidean or Hermitian norm < x , x > 1/2 .
2. Hypotheses used for developing the methodology (H1) We consider the linear structural vibrations of a complex heterogeneous linear mechanical system with respect to the geometry, the boundary conditions, and the materials. We are interested in characterizing the vibrations in the frequency domain Ω that is a bounded interval of R. (H2) The structure is a multi-scale system inducing the presence of a huge number of local and global elastic structural modes in the frequency band Ω.

(H3) The frequency sampling set Ω s = {ω 1 , . . . , ω n ω } ⊂ Ω is made up of a large number n ω of frequency sampling points.

(H4) Due to (H2), the HDM requires the use of a very large number n u of degrees of freedom.

(H5) The HDM depends on a control parameter w = (w 1 , . . . , w n w ) belonging to an admissible set C w ⊂ R n w in which n w is relatively small and for which a nominal value w is given in C w . The control parameter is uncertain and its nominal value must be updated by using experimental information.

(H6) For a fixed value of w in C w , we consider n o frequency response functions {ω → FRF k (ω; w), k = 1, . . . , n o } from Ω into C in which n o is relatively small (a few units). The quantity of interest is the real matrix [q(w)] ∈ M n ω ,n o depending on the control parameter w and on the frequency sampling such that [q(w)] mk = 20 log 10 (|FRF k (ω m ; w)|) with m ∈ {1, . . . , n ω } and k ∈ {1, . . . n o }.

The numerical cost of one evaluation performed with the HDM is assumed to be high.

(H7) It is assumed that only one experimental FRF (the target) is available and is sampled using the frequency points of Ω s . This single experimental FRF is associated with the last observation of the HDM whose index is n o . It means that the number n o of observations is possibly greater than the number of targets (that is only of one). This sampled experimental FRF is then denoted by {FRF targ (ω m ), m = 1, . . . , n ω }. The experimental quantity of interest (the target) is then defined as the real vector q targ = (q targ 1 , . . . ,

q targ n ω ) ∈ R n ω such that q targ m = 20 log 10 (|FRF targ (ω m )|) for m ∈ {1, . . . , n ω }.

Computational model and frequency response functions

From Hypotheses (H1), (H3), and (H4), the finite element discretization with n u degrees of freedom of the boundary value problem yields the HDM that is written, in the frequency domain and for ω ∈ Ω, as

(-ω 2 [M(w)]+ι ω [D(w)] + [K(w)]) u(ω; w) = [B] f(ω) , (1) o
(ω; w) = [A] u(ω; w) , (2)
in which w ∈ C w is the control parameter, u(ω; w) ∈ C n u is the complex vector of the response for all the degrees of freedom, [M(w)], [D(w)], and [K(w)] are the mass, damping, and stiffness matrices that are assumed to be in M + n u (set of all the positive-definite symmetric (n u × n u ) real matrices), where [B] is the controllability real matrix in M n u ,n f in which n f is the number of excited dofs, [A] is the observability real matrix in M n o ,n u , where f(ω) is the complex vector of applied loads in C n f , and where o(ω; w) is the complex vector of the observations in C n o .

Taking into account Hypothesis (H7), let us define n o frequency response functions {ω → FRF k (ω; w) : Ω → C, k = 1, . . . , n o }, which all have the same input dof, i input , and for which the

n o output dofs are {o k (ω; w), k = 1, . . . , n o }. We then have n f = 1, f(ω) = f (ω) = 1 for all ω, and [B] i,1 = δ i,i input for all i ∈ {1, . . . , n u }. Therefore, we have, FRF k (ω; w) = o k (ω; w). For ω ∈ Ω s , the quantity of interest is the real matrix [q(w)] ∈ M n ω ,n o such that [q(w)] mk = 20 log 10 (|FRF k (ω m ; w)|).

Parameterized reduced-order model

Due to Hypotheses (H3) (n ω large) and (H4) (n u very large), and as we have explained in Section 1, we have to construct a reduced-order model (ROM). However, since the HDM depends on control parameter w, a parameterized ROM has to be built. For that we use the classical approach of modal synthesis consisting in projecting Eqs. ( 1) and ( 2) on a finite family of the elastic structural modes of the associated homogeneous conservative system. As we have previously explained, a multi-frequency model reduction with or without substructuring techniques can be used for complex structural dynamics problems, but in this paper, for simplifying the presentation, we limit the development to the classical modal synthesis. For all w fixed in C w , we introduce the generalized eigenvalue problem associated with Eq. ( 1),

[K(w)] [Φ(w)] = [M(w)] [Φ(w)] [Λ(w)] , (3)
in which [Λ(w)] is the diagonal matrix of the smallest eigenvalues 0 < λ 1 (w) ≤ . . . ≤ λ N(w) (w) in which N(w) ≪ n u is defined below by using a convergence criterion. The columns of the matrix [Φ(w)] in M n u ,N(w) are the eigenvectors (the elastic structural modes) associated with the eigenvalues that are the square of the eigenfrequencies. This matrix satisfies the usual orthogonality properties

[Φ(w)] T [M(w)] [Φ(w)] = [I N(w) ] and [Φ(w)] T [K(w)] [Φ(w)] = [Λ(w)
]. Introducing the complex vector v(ω; w) in C N(w) of the generalized coordinates such that u(ω; w) = [Φ(w)] v(ω; w), the parameterized ROM can then be written as

o(ω; w, N(w)) = [A(w)] v(ω; w, N(w)) , (4)
(-ω 2 [I N(w) ]+ι ω [D(w)] + [Λ(w)]) v(ω; w, N(w)) = [B(w)] f(ω) , ( 5
) in which [A(w)] = [A][Φ(w)] ∈ M n o ,N(w) is the reduced observability matrix, where [B(w)] = [Φ(w)] T [B] ∈ M N(w)
,n u is the reduced controllability matrix, and where

[D(w)] = [Φ(w)] T [D(w)] [Φ(w)] ∈ M + N(w)
is the positive-definite reduced damping matrix. For estimating N(w), the usual L 2 criterion cannot be used because it requires the computation of o(ω m ; w) for all ω ∈ Ω s by using Eqs. ( 1) and (2) of the HDM, which is way too expensive. Two methods can be used to analyze the convergence of the ROM with respect to N(w). The first one consists in analyzing the convergence of the sequence {o(ω; w, N(w))} N(w) with respect to N(w) by using the following criterion (of the Cauchy type): for ε given independent of w, find N opt (w) such that for all N(w) ≥ N opt (w) and N ′ (w) ≥ N opt (w), we have

2 m ∥o(ω m ; w, N(w)) -o(ω m ; w, N ′ (w))∥ 2 m ∥o(ω m ; w, N(w)) + o(ω m ; w, N ′ (w))∥ 2 ≤ ε .
The second one is as follows. At the desired cutoff frequency of the modal basis, the LDL factorization of the dynamic stiffness matrix (without damping) is computed. The number of negative terms in the diagonal matrix of this factorization is equal to the number of eigenvalues below the cutoff frequency. This way, for all w, the N(w) structural elastic modes cover the same frequency band. It should thus be noted that the dimension N opt (w) of the ROM depends on w. In the application presented in Section 7, it is the latter method that will be used.

Probabilistic model of the control parameter and generation of the training set using the stochastic parameterized reduced-order model

The control parameter w = (w 1 , . . . , w n w ) with values in C w ⊂ R n w is modeled by a R n w -valued random variable W = (W 1 , . . . , W n w ) defined on a probability space (Θ, T , P), the probability distribution P W (dw) of which is given and has a support that is C w . The mean value w = E{W} = R nw w P W (dw) of W is thus known and considered as the nominal value of the control parameter. The random generator of P W (dw) allows for generating

N d independent realizations {w 1 d , . . . , w N d d } of W.
For each w j d ∈ C w with j fixed in {1, . . . , N d }, the parameterized ROM of dimension N opt (w j d ), defined in Section 4, is used for computing the observations and then deducing, using the end of Section 3, the quantity of interest [q j d ] ∈ M n ω ,n o such that, for m = 1, . . . , n ω and k = 1, . . . , n o , we have [q j d ] mk = [q(w j d )] mk = 20 log 10 (|FRF k (ω m ; w j d )|). The training set is then defined by the N d independent realizations {([q j d ], w j d ), j = 1, . . . , N d } of random variable ([Q], W) with values in M n ω ,n o × R n w . For generating the learned set using the PLoM whose algorithm is summarized in Appendix A, matrix [q j d ] is reshaped in a vector q j d ∈ R n q with n q = n ω × n o .

6. Optimization problem for model updating

Defining the cost function and formulation of the optimization problem for model updating

From Hypotheses (H6) and (H7), the quantity of interest, associated with the one experimental FRF sampled in the n ω frequency points of Ω s , is the real vector q targ = (q targ 1 , . . . , q targ n ω ) belonging to R n ω . When computed with the stochastic computational model, the corresponding quantity (see Hypothesis (H7)) is the

R n ω -valued random variable ([Q] 1n o , . . . , [Q] n ω n o ) such that [Q] mn o = [q(W)] mn o with m ∈ {1, . . . , n ω }, which depends on the R n w -valued random control parameter W.
We propose the following formulation for the model updating. The updated value w opt ∈ R n w of the nominal value w ∈ R n w of the control parameter is constructed as the solution of the optimization problem, w opt = arg min

w 0 ∈C w J(w 0 ) , (6)
in which C w ⊂ R n w is the admissible set of the control parameter and where J is the cost function.

A first formulation of the cost function could be defined for all w 0 ∈ C w by

J(w 0 ) = 1 n ω n ω m=1 (q targ m -E{[Q] mn o | W = w 0 }) 2 . (7)
In Eq. ( 7), E{[Q] mn o | W = w 0 } is the conditional mathematical expectation of the real-valued random variable [Q] mn o given W = w 0 for w 0 in C w . Since there are no uncontrolled random variables in the considered system, if W is fixed to a given value w 0 , then [Q] mn o is deterministic and consequently, the conditional probability measure of [Q] mn o given W = w 0 is a Dirac measure. Nevertheless, since the training set has a finite number N d of points, the estimation of the joint probability measure of ([Q] mn o , W) is performed using the Gaussian kernel density estimation method of the nonparametric statistics. Consequently, such a data smoothing allows for introducing a regularization of the joint probability measure whose support is the manifold defined by the graph {([q] mn o (w), w), w ∈ C w }. Therefore, the joint probability measure of random variable ([Q] mn o , W) admits a density with respect to the Lebesgue measure on R × R n w . In the framework of such a smoothing, it can then be assumed that the conditional probability measure of [Q] mn o given W = w 0 can be written as

P [Q] mno |W (dq|w 0 ) = p [Q] mno |W (q|w 0 ) dq in which q → p [Q] mno |W (q|w 0
) is the conditional probability density function of [Q] mn o with respect to the Lebesgue measure dq on R, given W = w 0 (it will be the case, see Section 6.2). Consequently,

E{[Q] mn o |W = w 0 } = R q p [Q] mno |W (q|w 0 ) dq. (8)
Note that the conditional random variable [Q] mn o |W = w 0 is a R-valued random variable. In general, the maximum of the conditional pdf, q → p [Q] mno |W (q|w 0 ), is not reached at the mean value E{[Q] mn o | W = w 0 }. Hence, there is an interest in considering a second formulation of the cost function that is defined as follows. For w 0 fixed in C w and for m fixed in {1, . . . , n ω }, let q ML m (w 0 ) be defined by

q ML m (w 0 ) = arg max q p [Q] mno |W (q|w 0 ) . (9)
Therefore, Eq. ( 7) is replaced by the following one,

J(w 0 ) = 1 n ω n ω m=1 (q targ m -q ML m (w 0 )}) 2 . ( 10
)

Estimating the conditional statistics using PLoM

As previously explained, the PLoM method is used for generating the learned set made up of

N ar ≫ N d realizations {([q ℓ ar ], w ℓ ar ), ℓ = 1, . . . , N ar } of random variable ([Q], W) with values in M n ω ,n o × R n w ,
without calling the computational model. The PLoM algorithm [START_REF] Soize | Data-driven probability concentration and sampling on manifold[END_REF][START_REF] Soize | Probabilistic learning on manifolds[END_REF][START_REF] Soize | Physics-constrained non-Gaussian probabilistic learning on manifolds[END_REF][START_REF] Soize | Probabilistic learning on manifolds (plom) with partition[END_REF] is summarized in Appendix A, in which the random matrix [Q] with values in M n ω ,n o is reshaped in the random vector Q with values in R n q with n q = n ω × n o .

For m fixed in {1, . . . , n ω } and for w 0 fixed in C w , the conditional mathematical expectation

E{[Q] mn o | W = w 0 }
is estimated using the learned set {([q ℓ ar ] mn o , w ℓ ar ), ℓ = 1, . . . , N ar }. In order to simplify the notation, in this subsection, the subscripts m and n o will be removed and consequently,

Q mn o = [Q] mn o , q mn o = [q] mn o , q mn o = [q] mn o , [q ℓ
ar ] mn o , and so on, will simply be noted as Q, q, q, q ℓ ar , and so on. Let q and σ Q be the empirical mean value and standard deviation of random variable Q estimated with {q ℓ ar , ℓ = 1, . . . , N ar }. Similarly, for i ∈ {1, . . . , n w } let w i and σ i be the empirical mean value and standard deviation of the real-valued random variable W i estimated with {w ℓ ar,i , ℓ = 1, . . . , N ar }. We introduce the normalized random variables Q and W = ( W 1 , . . . , W n w ) such that, for i = 1, . . . , n w ,

Q = q + σ Q Q , W i = w i + σ i W i .
We define w 0 = ( w 0 1 , . . . , w 0 n w ) such that

w 0 i = w i + σ i w 0 i , i = 1, . . . , n w .
Consequently, Eqs. ( 8) and ( 9) can be calculated by using the following equations,

E{Q | W = w 0 } = q + σ Q E{ Q | W = w 0 } , (11)
in which

E{ Q | W = w 0 } = 1 p W ( w 0 ) R q p Q, W ( q, w 0 ) d q , (12)
and q ML (w 0 ) = arg max q p Q| W ( q| w 0 ) .

In Eq. ( 12), p W ( w 0 ) = R p Q, W ( q, w 0 ) d q and where ( q, w) → p Q, W ( q, w) is the pdf on R × R n w of random variable ( Q, W) with respect to d q d w. The classical Gaussian KDE method [START_REF] Duong | Cross-validation bandwidth matrices for multivariate kernel density estimation[END_REF][START_REF] Duong | Feature significance for multivariate kernel density estimation[END_REF][START_REF] Filippone | Approximate inference of the bandwidth in multivariate kernel density estimation[END_REF][START_REF] Zougab | Bayesian estimation of adaptive bandwidth matrices in multivariate kernel density estimation[END_REF] is used for estimating p Q, W ( q, w) on R × R n w as

p Q, W ( q, w) = 1 N ar N ar ℓ=1 1 s SB √ 2π exp{- 1 2s 2 SB ( q -q ℓ ar ) 2 } × 1 (s SB √ 2π) n w exp{- 1 2s 2 SB ∥ w -w ℓ ar ∥ 2 } , (14)
in which s SB is the following Silverman bandwidth, s SB = (4/(N ar (3 + n w ))) 1/(5+n w ) , and where, for ℓ = 1, . . . , N ar , q ℓ ar = (q ℓ arq)/σ Q and for i = 1, . . . , n w , w ℓ ar,i = (w ℓ ar,iw i )/σ i . From Eqs. ( 12), (13),and ( 14), it can be deduced that

E{ Q | W = w 0 } = N ar ℓ=1 q ℓ ar exp{-1 2s 2 SB ∥ w 0 -w ℓ ar ∥ 2 } N ar ℓ=1 exp{-1 2s 2 SB ∥ w 0 -w ℓ ar ∥ 2 } , (15)
q ML (w 0 ) = arg max

q N ar ℓ=1 exp{- 1 2s 2 SB (( q -q ℓ ar ) 2 + ∥ w 0 -w ℓ ar ∥ 2 )} . (16)

Choice of the optimization algorithm

For each value w 0 proposed by the optimization algorithm, the evaluation of the cost function is performed using PLoM as we have explained in Sections 6.1 and 6.2. Nevertheless, an algorithm has to be chosen for solving the optimization problem defined either by Eq. ( 6) with Eqs. ( 7) and ( 8) or by Eq. ( 6) with Eqs. [START_REF] Amsallem | Nonlinear model order reduction based on local reduced-order bases[END_REF] and [START_REF] Carlberg | A low-cost, goal-oriented compact proper orthogonal decomposition basis for model reduction of static systems[END_REF]. Note that in the application of Section 7, we will consider the optimization problem defined by the second set of equations. The optimization algorithms have extensively been developed for several decades and are today really efficient (see for instance [START_REF] Zelinka | Handbook of optimization: from classical to modern approach[END_REF] for problems in a deterministic framework). Concerning the algorithms for solving optimization problems in a statistical framework, many methods have been proposed such as the gradient-based learning that is adapted to convex problems [START_REF] Lecun | Gradient-based learning applied to document recognition[END_REF][START_REF] Spall | Introduction to Stochastic Search and Optimization: Estimation[END_REF], the global search algorithms such as the stochastic algorithms, the genetic algorithm, and the evolutionary algorithms [START_REF] Coello | Evolutionary multi-objective optimization: a historical view of the field[END_REF][START_REF] Konak | Multi-objective optimization using genetic algorithms: A tutorial[END_REF]. If dimension n w is small, a grid search algorithm can also be used (see the numerical application presented in Section 7).

Optimization strategy

We assume here that the optimization problem defined by Eq. ( 6) with Eqs. ( 10) and ( 9) is considered (which will be the case in Section 7). Once the optimal control parameter w opt is found from solving this optimization problem, the HDM is rerun once, using the control parameter w opt . We then introduce an error function e dB,k associated with observation number k, which is given, using the trapezoidal rule, by:

e dB,k (w) =         1 2(n ω -1) n ω -1 m=1 f k m+1 (w)+ f k m (w)         1/2 , (17)
in which f k m (w) = [q(w)] mkq targ m 2

.

The error e dB,n o (w opt ) characterizes the quality of the model updating with respect to the target. It is recalled that the cost function in Eq. ( 6) relies on the estimation of the conditional statistics using PLoM. The PLoM method is used for generating N ar = n MC × N d realizations based on N d calls to the HDM (see Appendix A) in which n MC is an integer greater than 1. It should be noted that increasing N d may not lead to an improved updating because of the lack of monotonicity in N d of the error function defined by Eq. ( 17) (as mentioned in Section 1 (iii)). The global strategy is then to solve the optimization problem several times by considering subsets of the training set (the PLoM method being used for each subset), because a subset can potentially yield a better updating than the training set (due to the lack of monotonicity). For each optimization problem related to a given subset, we estimate the optimal value w opt of w and one call to the HDM with w = w opt is then required for estimating the error defined by Eq. ( 17). Let N c be the maximum number of calls to the HDM that we can afford. A number N d < N c of calls to the HDM are done to obtain the training set. Let N cd = N c -N d be the number of calls to the HDM that are left. This means that N cd optimization problems (with associated PLoM computation) can be carried out for improving the updating.

Numerical application

The mechanical system considered in this application has been voluntarily chosen as a simple one and is comprehensively described so that the results presented can be reproduced.

Defining the computational model

As a numerical application, a structure made up of plates in bending mode is considered, see Fig. 2. It is constituted of two structural levels, such that it exhibits both global and local Figure 2: Multiscale plate structure including the position of the input and output nodes of the considered FRF (big black dots) structural elastic modes. This structure is similar to the one analyzed in [START_REF] Soize | Stochastic reduced-order model in low-frequency dynamics in presence of numerous local elastic modes[END_REF] in which global and local modes are plotted. This way, despite its simplicity, the structure presents similar characteristics as the structures of interest as described in Section 1. The structure is rectangular with length 0.26 m (x axis) and width 0.2 m (z axis), and it includes twelve square slots with dimension 0.04 m, three of them being empty (holes) and nine of them being thin panels. These twelve slots are regularly spaced, such that the matrix frame (or frame) is of width 0.02 m. The structure is modeled using isotropic Kirchhoff plate elements. The nominal properties of the structure are now given. The plate thickness h 1 of the frame is of h 1 = 2 × 10 -3 m, its Young modulus is E 1 = 210 × 10 9 GPa, its Poisson ratio is ν 1 = 0.3, and its mass density is ρ 1 = 2600 kg/m 3 . For the nine panels, the Young modulus is E 2 = 5 × 10 9 GPa, the Poisson ratio is ν 2 = 0.3, and the mass density is ρ 2 = 5200 kg/m 3 . Regarding their thickness, it is h 2 = 1.1 × 10 -4 m for the bottom-left, bottom-right, and top-mid panels, whereas it is h ′ 2 = 9 × 10 -5 m for the six remaining panels. As boundary conditions, the four corners are fixed. The finite element model has n u = 163,032 dof. The quantity of interest is the FRF between the transverse displacement dof of the two nodes depicted in Fig. 2. Taking the bottom left corner as origin, the coordinates (m) of the two nodes are (x 1 = 0.01, z 1 = 0.08) and (x 2 = 0.23, z 2 = 0.155). We then have n o = 1. The FRF is uniformly sampled using n ω = 5000 frequency points.

Defining the probabilistic model of the control parameter W

The computational model depends on a control parameter w = (w 1 , . . . , w n w ) with dimension

n w = 4. Let k 1 = E 1 /(3(1-2ν 1 )), g 1 = E 1 /(2(1+ν 1 )), k 2 = E 2 /(3(1-2ν 2 )
), and g 2 = E 2 /(2(1+ν 2 )) denote the bulk and shear moduli of the frame and of the panels, respectively. The n w = 4 components of the control parameter w are given by w 1 = Log(k 1 ), w 2 = Log(g 1 ), w 3 = Log(k 2 ), and w 4 = Log(g 2 ) . The R n w -valued random variable W = (W 1 , . . . , W n w ) associated with w is defined as follows.

• Let K 1 denote the random bulk modulus of the frame. Its mean value is k 1 and its coefficient of variation is denoted as δ K 1 . The probability distribution of K 1 is chosen as a Gamma distribution with shape parameter α

K 1 = δ -1/2 K 1
and rate parameter β

K 1 = α K 1 /k 1 . Random variable W 1 is given by W 1 = Log(K 1 ) .
• Similarly, let G 1 denote the random shear modulus of the frame. Its mean value is g 1 and its coefficient of variation is denoted as δ G 1 . Random variable G 1 follows a Gamma distribution with shape parameter α

G 1 = δ -1/2 G 1 and rate parameter β G 1 = α G 1 /g 1 . Random variable W 2 is given by W 2 = Log(G 1 ) .
• Let K 2 denote the random bulk modulus of the panels. Its mean value is k 2 and its coefficient of variation is denoted as δ K 2 . The probability distribution of K 2 is chosen as a Gamma distribution with shape parameter α

K 2 = δ -1/2 K 2
and rate parameter

β K 2 = α K 2 /k 2 . Random variable W 3 is given by W 3 = Log(K 2 ) .
• Similarly, let G 2 denote the random shear modulus of the panels. Its mean value is g 2 and its coefficient of variation is denoted as δ G 2 . Random variable G 2 follows a Gamma distribution with shape parameter α

G 2 = δ -1/2 G 2 and rate parameter β G 2 = α G 2 /g 2 . Random variable W 4 is given by W 4 = Log(G 2 ) .

Parameterized reduced-order computational model and training set computation

The training set is obtained based on the realizations of random control parameter W. The admissible set C w of W is such that k 1 ∈ ]0, +∞], k 2 ∈ ]0, +∞], ν 1 ∈ ]0.15 , 0.45], and ν 2 ∈ ]0.15 , 0.45]. For each material, according to the paper [START_REF] Guilleminot | On the statistical dependence for the components of random elasticity tensors exhibiting material symmetry properties[END_REF], the bulk and shear moduli are statistically independent and it is shown that we have:

δ G 1 = 1 -4 + 5/δ K 1 , (18)
δ G 2 = 1 -4 + 5/δ K 2 . (19)
In addition, it is assumed that the random mechanical properties are independent from one material to the other. For the training set, we choose δ K 1 = 0.15 and δ K 2 = 0.5 . This means that the material variability is greater for the panels than for the main frame. The nominal values given in Section 7.1 are used and allow to obtain the random realizations of W through the steps detailed in Section 7.2. As explained in Section 1, for the case of an HDM, the number N c of calls to the computational model is limited by the computational resources. In this application, we then also choose a small value N c = 200 of the maximum number of calls to the HDM. The training set is then generated by N d = 175 independent realizations of W. The value of N d has been chosen to let 12 optimization problems to be solved with and without PLoM (2

× 13 + 175 = 201 ∼ 200).
For a fixed realization of W, the computational model is used for computing the corresponding realization of the FRF. As explained in Section 4, this involves computing all the structural elastic modes up to a given cutoff frequency, for which the number of elastic modes is determined by using the LDL factorization. The FRF is computed for n ω = 5000 uniformly sampled frequencies in Ω = 2π×]0, 1000] rad/s. In order to get a good convergence of the modal representation for all the values of the control parameter considered in the training set, a cutoff frequency of 1500 Hz is considered. A modal damping model is used with the same damping ratio ξ = 0.01.

Remark concerning the representativeness of the presented example

The following comments can be done concerning the representativeness of the presented example with respect to the nuclear engineering (NE) application mentioned in Section 1. Regarding the elastic modes, in the present application, there are significantly less (a few hundreds). Nevertheless, for both cases, the elastic modes include a relatively small number of global modes, which are preponderant in the structural response. This means that for both cases, the influence of local modes is rather diffuse and sometimes imperceptible for one given mode. The frequency content of the structural response is not that much more complex for the case of the NE application as for the case of the present application. Regarding the low number of parameters used (here four) compared to NE applications, for which there would be, for example, a hundred parameters, the proposed methodology would remain the same. Between the two cases, the difference would be the choice of the optimization algorithm. In either case, the number of evaluations of the cost function would be too large. We therefore need a surrogate model, which we proposed by building a statistical surrogate model using PLoM. Regarding the under-observability of the model, only one FRF is observed (only one target FRF) and for this FRF, only one "experimental" realization is considered available. This is a highly under-observed case and for the NE application, the under-observability is slightly less because several FRF are observed at several locations of the outer shell.

Model updating

In this Section, two cases are presented, corresponding to two different choices of target FRF.

Case 1. For this case, we consider a target FRF generated by the computational model with E 6), [START_REF] Amsallem | Nonlinear model order reduction based on local reduced-order bases[END_REF], and ( 9) is solved and yields w opt that is used to update the computational model and obtain the updated FRF (as explained in Section 6.4). Concerning the PLoM parameters used in the algorithm summarized in Appendix Appendix A, the value of ε PCA is fixed to 10 -3 . For the 13 values of N d belonging to S N d , which vary from 25 to 175, the reduced-order dimension ν of the PCA varies from 22 to 82, the optimal value m opt of the truncated diffusion-maps basis is given by m opt = ν + 1, and the optimal value ε opt varies from 40 to 204. As an illustration, for N d = 175, for which ν = 82, m opt = 83, and ε opt = 204, Fig. 3 displays the distribution of the eigenvalues µ α of the PCA of X and Fig. 4 displays the distribution of the eigenvalues λ α of the transition matrix. Regarding the optimization algorithm, a grid search is performed, with 35 points per axis (there are four axes, for E 1 , ν 1 , E 2 , and ν 2 ). The updated FRF is compared to the target FRF using the error defined by Eq. [START_REF] Ohayon | Structural Acoustics and Vibration: Mechanical Models, Variational Formulations and Discretization[END_REF]. In parallel, the exact same procedure is also carried out using the training set instead of the learned set. The total number of additional calls to the computational model is twice the size of set S N d , that is 26, which makes it a total of 201 calls to the computational model. It can be seen in Fig. 5 that the updating is best for N d = 125 using the learned set. The error is of 3.4 dB and the updated parameters are E The quality of the updated model is put into evidence in Fig. 6 where the nominal FRF, the target FRF, and the updated FRF using the learned set with N d = 125 and N ar = 20,000 are depicted. It can be seen that the proposed methodology is successful in updating the computational model. = 0.273. The entire procedure is kept identical and the results are shown in Figs. 7 and8. It can be seen in Fig. 7 that the updating is best for N d = 138 using the training set (the fact that the training set yields a better optimum than the learned set is discussed in Section 7.6 as a complement to the reasons already given in Section 1-(iii)). The error is of 2.2 dB and the updated parameters are E First, we recall the framework of the proposed updating methodology based on a single target FRF. For a given large computational model, only N d realizations can be computed, with N d having a small value. In this context, the objective is to improve the computational model using the single target FRF. So, the metholodogy that we have proposed is the following. For such a given value of N d , the training set and a learned set are generated. The optimization problem is solved on the one hand with the training set and on the other hand with the learned set. For this fixed value of N d , the best updated model is the one having the smallest error between the target FRF and the updated FRF (error defined by Eq. ( 17)). In this framework of the use of the methodology, the notion of convergence with respect to N d has no meaning. Figures 5 and7 effectively show that the error has large fluctuations as a function of N d . Such fluctuations are due to the small value of N d (the largest is 175), to the fact that only a single target FRF is available, that there are multiple control parameters, and that the optimization problem is nonconvex. Consequently, there is no uniqueness of the optimal control parameters. Clearly, if the training set could be constructed for a very large value N * d of N d , the error could perhaps, in this case, decrease for N d > N * d . For instance, for the application presented, the updating has also been performed using the training set with N d = 2000 without obtaining a clear convergence of the error, and the available computing resources did not permit to consider an even greater value of N d .

Conclusions

We have presented a formulation and a methodology for the updating of high-dimensional computational models with uncertainties of complex structures, which present an enormous number of global modes and local modes in the frequency band of analysis. The model updating is performed in an under-observed framework for which only one experimental frequency response function is available. The methodology we have proposed is based on the introduction of a parameterized ROM and on the use of a probabilistic learning algorithm, the PLoM that has been developed to deal with small data cases, and which allows to solve the optimization problem without using the resolution of the HDM for which a single evaluation is very expensive. The methodology presented has been validated on a representative example from the point of view of the difficulties encountered for the dynamic systems described in the introduction. Obviously the structure considered in the example presented is much simpler than those of the targeted applications, but this choice was made so that the results presented can be reproduced. Finally, as the example presented is really simple, a multi-frequency model reduction with or without substructuring techniques is not necessary for constructing the parameterized ROM. This allowed us to simplify the presentation of the methodology while remaining to the essentials. However, these methodological additions can be implemented without difficulty in the proposed method. We refer the reader to Section 7.6 for the discussion about the framework of the use of the proposed method.

Appendix A. Summary of the probabilistic learning on manifolds (PLoM) algorithm and its parameterization

The PLoM approach [START_REF] Soize | Data-driven probability concentration and sampling on manifold[END_REF][START_REF] Soize | Probabilistic learning on manifolds[END_REF][START_REF] Soize | Physics-constrained non-Gaussian probabilistic learning on manifolds[END_REF], which has specifically been developed for small data (in opposite to big data) starts from a training set D d made up of a relatively small number N d of points. It is assumed that D d is generated with an underlying stochastic manifold related to a R n -valued random variable X = (Q, W), defined on a probability space (Θ, T , P), in which Q is the quantity of interest that is a R n q -random variable, where W is the control parameter that is a R n w -random variable, and where n = n q + n w . Another R m u -valued random variable U defined on (Θ, T , P) can also be considered, which is an uncontrolled parameter and/or a noise. Random variable Q is assumed to be written as Q = f(U, W) in which the measurable mapping f is not explicitly known. The joint probability distribution P W,U (dw, du) of W and U is assumed to be given. Note that in this paper, we have not considered any uncontrolled variable. The non-Gaussian probability measure P X (x) = P Q,W (dq, dw) of X = (Q, W) is concentrated in a region of R n for which the only available information is the cloud of the points of training set D d . The PLoM method makes it possible to generate the learned set D ar for X whose n MC ≫ N d points (learned realizations) are generated by the non-Gaussian probability measure that is estimated using the training set. The concentration of the probability measure is preserved thanks to the use of a diffusion-maps basis that allows to enrich the available information from the training set. Using the learned set D ar , PLoM allows for carrying out any conditional statistics such as w → E{Q|W = w} from C w in R n q , and consequently, to directly construct metamodels in a probabilistic framework.

The training set D d is made up of the N d independent realizations x j d = (q j d , w j d ) in R n = R n q ×R n w for j ∈ {1, . . . , N d } of random variable X = (Q, W). The PLoM method allows for generating the learned set D ar made up of N ar ≫ N d learned realizations {x ℓ ar , ℓ = 1, . . . , N ar } of random vector X. As soon as the learned set has been constructed, the learned realizations for Q and W can be extracted as (q ℓ ar , w ℓ ar ) = x ℓ ar for ℓ = 1, . . . , N ar . (A.2) Construction of a reduced-order diffusion-maps basis. For preserving the concentration of the learned realizations in the region in which the points of the training set are concentrated, the PLoM relies on the diffusion-maps method [START_REF] Coifman | Diffusion maps[END_REF][START_REF] Lafon | Diffusion maps and coarse-graining: A unified framework for dimensionality reduction, graph partitioning, and data set parameterization[END_REF]. This is an algebraic basis of vector space R N d , which is constructed using the diffusion maps. (A.5) Generation of learned realizations {η ℓ ′ ar , ℓ ′ = 1, . . . , N ar } of random vector H. The MCMC generator is detailed in [START_REF] Soize | Data-driven probability concentration and sampling on manifold[END_REF]. Let {([Z(t)], [Y(t)]), t ∈ R + } be the unique asymptotic (for t → +∞) stationary diffusion stochastic process with values in M ν,m opt × M ν,m opt , of the following reducedorder ISDE (stochastic nonlinear second-order dissipative Hamiltonian dynamic system), for t > 0,

p(u j ) = 1 N d N d j ′ =1 exp{- 1 2 s 2 ν ∥ s ν s ν η j ′ -u j ∥ 2 } , ∇ u j p(u j ) = 1 s 2 ν N d N d j ′ =1 ( s ν s ν η j ′ -u j ) × exp{- 1 2 s 2 ν ∥ s ν s ν η j ′ -u j ∥ 2 } ,
in which s ν is the modified Silverman bandwidth s ν , which has been introduced in [START_REF] Soize | Polynomial chaos expansion of a multimodal random vector[END_REF],

s ν = s ν s 2 ν + N d -1 N d , s ν = 4 N d (2 + ν) 1/(ν+4)
.

(2) [W wien (t)] = [W wien (t)] [a] where {[W wien (t)], t ∈ R + } is the M ν,N d -valued normalized Wiener process.

(3) [N ] is the M ν,N d -valued normalized Gaussian random matrix that is independent of process [W wien ]. (4) The free parameter f 0 , such that 0 < f 0 < 4/ s ν , allows the dissipation term of the nonlinear second-order dynamic system (dissipative Hamiltonian system) to be controlled in order to kill the transient part induced by the initial conditions. A common value is f 0 = 4 (note that s ν < 1). [START_REF] Nguyen | An efficient reduced-order modeling approach for non-linear parametrized partial differential equations[END_REF] We then have [Z m opt ] = lim t→+∞ [Z(t)] in probability distribution. The Störmer-Verlet scheme is used for solving the reduced-order ISDE, which allows for generating the learned

2 = 3 . 7 GPa, and ν target 2 =

 2372 0.351. We consider a sampling of N d up to 175, such that N d ∈ S N d = {25, 38, 50, 63, 75, 88, 100, 113, 125, 138, 150, 163, 175}. For fixed N d , we define n MC such that N ar = n MC × N d is the closest to 20,000. For fixed (N d , n MC ) couple with N d ∈ S N d , the PLoM method is used to generate N ar = n MC × N d learned realizations, based on which the optimization problem defined by Eqs. (

Figure 3 :

 3 Figure 3: For N d = 175, ν = 82, m opt = 83, and ε opt = 204, distribution of the eigenvalues µ α of the PCA of X.

Figure 4 :

 4 Figure 4: For N d = 175, ν = 82, m opt = 83, and ε opt = 204, distribution of the eigenvalues λ α of the transition matrix in log 10 scale.

Figure 5 :Case 2 . 1 = 254 GPa, ν target 1 = 0.312, E target 2 = 5 . 8 GPa, and ν target 2

 52112582 Figure 5: Case 1: for each value of N d belonging to S N d , value of the error obtained using the learned set with N ar = N d × n MC ∼ 20 000 realizations (dashed line) and using the training set with N d realizations (solid line)

2 =Figure 6 :Figure 7 :

 267 Figure 6: Comparison between the nominal FRF (solid line), the target FRF (dash-dotted line), and the updated FRF (dashed line)

Figure 8 : 7 . 6 .

 876 Figure 8: Comparison between the nominal FRF (solid line), the target FRF (dash-dotted line), and the updated FRF (dashed line)

(A. 1 )

 1 Reduced representation. The N d independent realizations {x j d , j = 1, . . . , N d } are represented by the matrix [x d ] = [x 1 d . . . x N d d ] in M n,N d . Let [X] = [X 1 , . . . , X N d] be the random matrix with values in M n,N d , whose columns are N d independent copies of random vector X. Using the PCA of X, random matrix [X] is written as,[X] = [x] + [φ] [µ] 1/2 [H] , (A.1) in which [H] = [H 1 , . . . , H N d ] is a M ν,N d -valued random matrix, where ν ≤ n, and where [µ] is the (ν × ν) diagonal matrix of the ν positive eigenvalues of the empirical estimate of the covariance matrix of X. The (n×ν) matrix [φ] is made up of the associated eigenvectors such [φ] T [φ] = [I ν ]. The matrix [x] in M n,N d has identical columns, each one being equal to the empirical estimate x ∈ R n of the mean value of random vector X. The columns of [H] are N d independent copies of a random vector H with values in R ν . The realization [η d] = [η 1 d . . . η N d d ] ∈ M ν,N d of [H] is computed by [η d ] = [µ] -1/2 [φ] T ([x d ] -[x]). The value ν is classically calculated in order that the L 2 -error function ν → err X (ν) defined by errX (ν) = 1 -ν α=1 µ α E{∥X∥ 2 } , (A.2)be smaller than ε PCA . If ν < n, then there is a statistical reduction.

  Let [K] and[b] be the matrices such that, for all i and j in {1, . . . ,N d }, [K] i j = exp{-(4 ε DM ) -1 ∥η i d -η j d ∥ 2 } and [b] i j = δ i j b i with b i = N d j=1 [K] i j , in which ε DM > 0 isa smoothing parameter. The eigenvalues λ 1 , . . . , λ N d and the associated eigenvectors ψ 1 , . . . , ψ N d of the right-eigenvalue problem [P] ψ α = λ α ψ α are such that 1 = λ 1 > λ 2 ≥ . . . ≥ λ N d and are computed by solving the generalized eigenvalue problem [K] ψ α = [H m opt ], we have the estimate d 2 N d (m opt ) ≃ (1/n MC ) n MC ℓ=1 {∥[η ℓ ar ] -[η d ]∥ 2 }/∥[η d ]∥ 2 .

  d[Z(t)] = [Y(t)] dt , d[Y(t)] = [L([Z(t)])] dt -1 2 f 0 [Y(t)] dt + f 0 [dW wien (t)] , with [Z(0)] = [η d ] [a] and [Y(0)] = [N ] [a], in which [a] = [g m opt ] ([g m opt ] T [g m opt ]) -1 ∈ M N d ,m o pt .(1) [L([Z(t)])] = [L([Z(t)] [g m opt ] T )] [a] is a random matrix with values in M ν,m opt . For all [u] = [u 1 . . . u N d ] in M ν,N d with u j = (u j 1 , . . . , u j ν ) in R ν , the matrix [L([u])] in M ν,N d is defined, for all k = 1, . . . , ν and for all j = 1, . . . , N d , by [L([u])] k j = 1 p(u j ){∇ u j p(u j )} k , (A.6)

λ α [b] ψ α with the normalization < [b] ψ α , ψ β >= δ αβ . The eigenvector ψ 1 associated with λ 1 = 1 is a constant vector. For a given integer κ ≥ 0, the diffusion-maps basis {g 1 , . . . , g α , . . . , g N d } is a vector basis of R N d defined by g α = λ κ α ψ α . For a given integer m, the reduced-order diffusionmaps basis of order m is defined as the family {g 1 , . . . , g m } that is represented by the matrix [g m ] = [g 1 . . . g m ] ∈ M N d ,m with g α = (g α 1 , . . . , g α N d ) and [g m ] ℓα = g α ℓ . This basis depends on two parameters, ε DM and m, which have to be identified. It is proven in [START_REF] Soize | Probabilistic learning on manifolds[END_REF], that the PLoM method does not depend on κ that can therefore be chosen to 0.

We have to find the optimal value m opt ≤ N d of m and the smallest value ε opt > 0 of ε DM such that (see [START_REF] Soize | Probabilistic learning on manifolds (plom) with partition[END_REF])

with an amplitude jump equal to an order of magnitude (a factor 10 as demonstrated in [START_REF] Soize | Probabilistic learning on manifolds[END_REF]) between λ m opt (ε opt ) and λ m opt +1 (ε opt ). A further in-depth analysis makes it possible to state the following algorithm to estimate ε opt and m opt . Let ε DM → Jump(ε DM ) be the function on ]0, +∞[ defined by Jump(ε

The algorithm is the following:

-set the value of m to m opt = ν + 1; -identify the smallest possible value ε opt of ε DM in order that Jump(ε opt ) ≤ 0.1 and such that Equation (A.3) be verified. (A.4) Criterion for quantifying the concentration of the probability measure of random matrix [H m opt ]. For m ≤ N d , the concentration of the probability measure of random matrix [H m ] is defined (see [START_REF] Soize | Probabilistic learning on manifolds[END_REF]) by

Let M opt = {m opt , m opt + 1, . . . , N d } in which m opt is the optimal value of m previously defined. Theorem 7.8 of [START_REF] Soize | Probabilistic learning on manifolds[END_REF] shows that min m∈M opt d 2

, which means that the PLoM method, for m = m opt and [g m opt ] is a better method than the usual one corresponding to ar , ℓ ′ = 1, . . . , N ar } of random vector X are then calculated (see Eq. (A.1)) by x ℓ ′ ar = x + [φ] [µ] 1/2 η ℓ ′ ar . (A.6) Constraints on the second-order moments of the components of H. In general, the mean value of H estimated using the N ar learned realizations {η ℓ ′ ar , ℓ ′ = 1, . . . , N ar }, is sufficiently close to zero. Likewise, the estimate of the covariance matrix of H, which must be the identity matrix, is sufficiently close to a diagonal matrix. However, sometimes the diagonal entries of the estimated covariance matrix can be lower than 1. Normalization can be recovered by imposing constraints

in the algorithm presented in paragraph (v). For that, we use the method and the iterative algorithm presented in [START_REF] Soize | Probabilistic learning on manifolds (plom) with partition[END_REF] (that is based on Sections 5.5 and 5.6 of [START_REF] Soize | Physics-constrained non-Gaussian probabilistic learning on manifolds[END_REF]). The constraints are imposed by using the Kullback-Leibler minimum cross-entropy principle. The resulting optimization problem is formulated using a Lagrange multiplier v = (v 1 , . . . , v ν ) associated with the constraints. The optimal solution of the Lagrange multiplier is computed using an efficient iterative algorithm. At each iteration, the MCMC generator detailed in paragraph (v) is used. The constraints are rewritten as

The iteration algorithm for computing v i+1 as a function of v i is the following,

)}] (the covariance matrix), and where α i is a relaxation function (less than 1) that is introduced for controlling the convergence as a function of iteration number i. For given i 2 ≥ 2, for given β 1 and β 2 such that 0 < β 1 < β 2 ≤ 1, α i can be defined by: -for i ≤ i 2 , α i = β 1 + (β 2 -β 1 )(i -1)/(i 2 -1); -for i > i 2 , α i = β 2 .

The convergence of the iteration algorithm is controlled by the error function i → err(i) defined by err(i) = ∥b -E{h(H v i )}∥/∥b∥ . (A.8)

At each iteration i, E{h(H v i )} and [cov{h(H v i )}] are estimated by using the N ar learned realizations of H m opt (v i ) obtained by reshaping the learned realizations.