M Allabou

R Bouclier
email: bouclier@insa-toulouse.fr

P A Garambois

J Monnier
email: monnier@insa-toulouse.fr

Reduction of the Shallow Water System by an Error Aware POD-Neural Network Method: Application to Floodplain Dynamics

Keywords:

In this study, we elaborate on and evaluate a new reduced basis method for model reduction of the shallow water equations using Proper Orthogonal Decomposition (POD) and artificial Neural Networks (NNs). The method begins with the POD technique to construct reduced bases from high-resolution solutions, followed by training two deep NNs to learn associated coefficients in the reduced bases. The approach follows an offline-online strategy: the POD reduced basis, along with the training of the NNs, is performed in an offline stage, and then the surrogate model can be used in an online stage for real-time predictions. The method takes into account the POD-based projection error, enabling the attainment of higher accuracy while preserving a limited number of POD modes, even in the delicate situation of convection-dominated flow problems. This point is crucial in our approach since it enables to limit the output dimension of the NNs, thus providing the opportunity to employ smaller NNs (with less parameters), and therefore smaller dataset (the snapshots) to train the NNs. The process is non-intrusive: it does not require opening the highresolution model code. The method is evaluated on a real-world test case aimed at simulating inundation of the Aude river (Southern France). The results show that the proposed method provides satisfying accuracy for the hydraulic variables (water elevation, discharge) compared to the reference high-resolution 2D shallow water model. Overall, the method is promising, particularly for performing real-time simulations of flood plain dynamics.

Introduction

The Reduced Order Model (ROM) theory (see, e.g., [START_REF] Quarteroni | Reduced basis methods for partial differential equations: an introduction[END_REF]13,[START_REF] Benner | Model Order Reduction. Volume 2: Snapshot-Based Methods and Algorithms[END_REF]) has gained significant attention in the field of computational sciences, offering techniques to tackle the computational challenges posed by complex systems. The goal is to simplify high-dimensional models while preserving their essential features. Model reduction methods can be broadly categorized into intrusive and non-intrusive approaches, each with its advantages and limitations, see e.g. [START_REF] Quarteroni | Reduced basis methods for partial differential equations: an introduction[END_REF]13,[START_REF] Benner | Model Order Reduction. Volume 2: Snapshot-Based Methods and Algorithms[END_REF] and references therein. In this paper, we adopt an offline-online non-intrusive strategy, applied to non-linear hyperbolic problems, specifically to the 2D Sallow Water Equations (SWEs). This PDEs system describes well the behavior of flood plain dynamics. Capturing complex features such as non linear travelling waves, over long time integration periods, using e.g. Finite Volume (FV) schemes is CPU-time consuming, which is impractical for real-time simulations or optimization and uncertainty quantification. To overcome these challenges, model reduction techniques offer an alternative approach by reducing the dimensionality of the problem while preserving the crucial dynamics.

However, ROMs are particularly challenging in this situation of convection-dominated problems, see e.g. [START_REF] Quarteroni | Reduced basis methods for partial differential equations: an introduction[END_REF][START_REF] Shady E Ahmed | Nonlinear proper orthogonal decomposition for convection-dominated flows[END_REF][START_REF] Barnett | Neural-network-augmented projection-based model order reduction for mitigating the kolmogorov barrier to reducibility[END_REF] and references therein. Applying ROM techniques aims to strike a balance between computational efficiency and accuracy.

The reduction of these complex systems enables faster simulations, facilitates real-time decision-making processes, and opens avenues for optimization and uncertainty quantification tasks.

The literature addressing ROMs for non linear hyperbolic systems like SWEs is already quite large. For intrusive ROM approaches, let us mention for example [START_REF] Steinstraesser | Modified parareal method for solving the two-dimensional nonlinear shallow water equations using finite volumes[END_REF][START_REF] Ştefănescu | Comparison of POD reduced order strategies for the nonlinear 2D shallow water equations[END_REF][START_REF] Mocayd | Polynomial surrogates for open-channel flows in random steady state[END_REF][START_REF] Strazzullo | POD-Galerkin model order reduction for parametrized nonlinear time-dependent optimal flow control: an application to shallow water equations[END_REF] and for non-intrusive ones [START_REF] Dutta | A greedy non-intrusive reduced order model for shallow water equations[END_REF][START_REF] Wang | Non-intrusive reduced order modeling of unsteady flows using artificial neural networks with application to a combustion problem[END_REF]12], see Table 1.

To reduce the 2D SWEs model (considered as non-parametrized), [START_REF] Ştefănescu | Comparison of POD reduced order strategies for the nonlinear 2D shallow water equations[END_REF][START_REF] Steinstraesser | Modified parareal method for solving the two-dimensional nonlinear shallow water equations using finite volumes[END_REF] propose POD-DEIM (Proper Orthogonal Decomposition -Discrete Empirical Interpolation Method) approaches which are intrusive. [START_REF] Dutta | A greedy non-intrusive reduced order model for shallow water equations[END_REF] proposes a nonintrusive method based on POD and Radial Basis Functions (RBF) interpolation. Considering the 2D SWEs as a parametrized model, [START_REF] Mocayd | Polynomial surrogates for open-channel flows in random steady state[END_REF][START_REF] Strazzullo | POD-Galerkin model order reduction for parametrized nonlinear time-dependent optimal flow control: an application to shallow water equations[END_REF] propose intrusive approaches based on generalized polynomial chaos and POD-Galerkin, respectively. [START_REF] Wang | Non-intrusive reduced order modeling of unsteady flows using artificial neural networks with application to a combustion problem[END_REF]12] address non-intrusive POD-NN (Neural Network) based reduction methods to reduce parametrized SWEs. These POD-NN based methods rely on the usual offline-online strategy: a computationally expensive offline phase is performed beforehand to construct the reduced-order model, while the online phase efficiently computes the reduced system's response given a new parameter value, without further involving the fullorder model (non-intrusive aspect). In these POD-NN methods, the POD enables to identify the dominant modes of the system, capturing the most significant features, while the NN serves as approximating the coefficients in the POD-based reduced basis. Such a combination of POD and NN seems to be first introduced in [START_REF] Hesthaven | Non-intrusive reduced order modeling of nonlinear problems using neural networks[END_REF], and later in [START_REF] Wang | Non-intrusive reduced order modeling of unsteady flows using artificial neural networks with application to a combustion problem[END_REF].

In these pioneering studies, the method showed promising efficiency for relatively simple PDE-based models, but also for the steady incompressible Navier-Stokes equations modeling an academic driven cavity viscous flows. However, such POD-NN based methods turned out to be insufficiently accurate for complex situations such as those presented in the present study (a real-life flooding event).

In this work, we extend the aforementioned POD-NN method by learning, in addition to the reduced (projected) solution, the projection error, through an additional NN. For non-linear convection-dominated or hyperbolic systems such as the one considered here, this enhancement turns out to be essential for constructing an accurate and robust surrogate model. Indeed, it allows for the correction of the reduced solution obtained from the standard POD-NN method while limiting the output dimension of the NNs, i.e., the dimension of the reduced spaces. This limitation enables to maintain reasonably sized datasets (snapshots) for training. Without this, achieving an accurate machinelearning-type method may be infeasible. Ultimately, the proposed method preserves twofold advantages of the original POD-NN method: (i) it is generic, meaning it can be applied to any PDEs based problems solved by any discretization schemes (FV or FE for example); (ii) it is non-intrusive as it only requires to run, in a black-box fashion, the highresolution model code multiple times in the offline phase. Due to the incorporation of this error-learning step, we refer to our approach as an Error-Aware POD-NN (EA-POD-NN) reduction method.

The paper is organized as follows. After this introduction, Section 2 introduces the considered parametrized 2D SWEs model, and its numerical resolution. Then, Section 3 presents, in a general manner, the developed EA-POD-NN methodology. It outlines the steps involved in the offline and online phases of the method and describes the corresponding algorithms. Particular care is given in this section to highlight the limitations of the existing POD-NN approach, thereby motivating the integration of the projection error within the process. Section 4 presents the numerical experiments conducted to evaluate the performance and accuracy of the EA-POD-NN method. It discusses the selected SWEs problems, the setup of the simulations, and the comparison of results between the initial highresolution model and the ROM using the proposed EA-POD-NN approach. Eventually, Section 5 provides a brief conclusion, summarizing the key findings and main contributions of this study, and identifying potential areas for future research based on the proposed approach.

The -parametrized model and the reference high-resolution solutions

In this section, first, the general concept of High-Resolution (HR) solutions corresponding to a given set of input parameters is recalled. Next, the flow model considered in our study is presented, followed by a discussion of the numerical scheme and computational software employed to solve it.

Basic principle of the high-resolution solutions generation

The reduction method developed in this study can be formally applied to any -parametrized non-linear PDE-based model with the input parameter = ( 1 , … , ) ∈ ℝ , which typically represents a set of physical parameters in the bulk and/or boundary condition parameters. From a reduction point of view, the important characteristic of is its dimension : has to be quite small, say = (10) at most, otherwise an additional reduction method of should be considered in a pre-processing step.

More precisely, we set  such that ∈  with dim() = . A set of parameters are fixed by sampling in some way the parameter space . We obtain the reference parameter set = { } =1 , ∈ () . Next, the -parametrized model, denoted as  , is employed to generate the corresponding vector solutions ℎ :

this is the classically called HR solutions set, also called the snapshots set. These snapshots ℎ constitute the set

= ℎ ( 1 )| … | ℎ ( ) .

The -parametrized mathematical flow model

In the developed application, the considered model  relies on the 2D SWEs, employed in particular to simulate river and floodplain flows dynamics. In this context, the parameter could be related to boundary conditions, to the initial conditions, or to the friction coefficient, etc.

For a given computational domain Ω ⊂ ℝ 2 and a time interval [0, ], the considered 2D SWEs are those written in conservative form as:

( ) ⎧ ⎪ ⎪ ⎨ ⎪ ⎪ ⎩ ℎ + ( ) = 0 in Ω × ]0, ] + ⊗ ℎ + ℎ 2 2 = ( ; ℎ) + ( ; ℎ, ) in Ω × ]0, ] Plus Initial Conditions, plus Boundary Conditions (B.C.)( ) . ( 1
)
∀( , ) ∈ Ω×[0, ], the variable ℎ( , ) represents the water depth (in ) and ( , ) stands for the discharge (in 3 -1 ).

Then, we have = ℎ , where = ( , ) denotes the depth-averaged velocity (in -1 ). Moreover:

• ( ; ℎ) is the gravity source term: ( ; ℎ) = -ℎ∇ , with the gravity magnitude and the bed elevation;

• ( ; ℎ, ) is the friction source term:

( ; ℎ, ) = - 2 ‖ ‖ ℎ 7∕3
, with the Manning-Strickler roughness coefficient.

The B.C. can be a mix of conditions necessary for real-world applications. At outflow boundaries, different types of flow conditions can be imposed, depending on the available information, see [START_REF]Open source computational software DassFlow (Data Assimilation for Free Surface Flows)[END_REF] for details. At the inflow boundaries, discharge time series are imposed. In the case of inflow boundaries, the inflow discharge time series , ( ; ), ∈ 1.. , are imposed. For the presentation of our method in a general manner in Section 3, let us denote by ℎ = (ℎ, ) = (ℎ, ℎ ) = (ℎ, , ) the HR solution field of [START_REF] Shady E Ahmed | Nonlinear proper orthogonal decomposition for convection-dominated flows[END_REF].

In this work, our primary focus is on the most important input parameter for operational users; that is, the inflow discharge signal at upstream (inflow B.C.). Consequently, in the numerical applications, the parameter will be related to the discharge functions , ( , ). Moreover, for simplicity in notation in the following, we will incorporate the time parameter into the parameter . This is not only convenient but also scientifically sound, considering that represents a dimension of  for constructing the snapshot set.

Numerical solver

Given a parameter value , the system  is numerically solved by the first-order FV method implemented into the open-source computational software DassFlow 2D [START_REF]Open source computational software DassFlow (Data Assimilation for Free Surface Flows)[END_REF][START_REF] Monnier | Inverse algorithms for 2D shallow water equations in presence of wet dry fronts: Application to flood plain dynamics[END_REF]. The solver relies on a well-balanced Godunov-type scheme and the explicit Euler time scheme, with additional treatment to be stable and robust at wet-dry front dynamics.

The boundary conditions required for real-world flows are available. The mesh is a mix of triangles and quadrangles.

The DassFlow 2D kernel codes are written in Fortran 90, using MPI library. Moreover, the Fortran computational kernels are wrapped in Python. This allows for seamless integration of these physics-based computations with other Python libraries. These libraries can serve various purposes, including simple pre or post-processing of computed fields or data, as well as incorporating Python libraries like PyTorch for deep learning, and so forth.

The Error-Aware POD-NN reduction method

Let us now detail the proposed reduction method that allows to construct an efficient and accurate surrogate model, here applied to the 2D SWEs-based model  . The presentation in this section is crafted in a general manner as the proposed approach is generic and, consequently, can be applied in various contexts. Specifically, the presentation follows the different aspects of the method, namely (i) the projection onto a RB, (ii) the POD technique for constructing the RB, (iii) the learning of coefficients of the reduced solution using a deep NN, and (iv) the involvement of the projection error.

Reduced basis and projection

Our method being based on a RB approach, let us start by properly defining the complete space, the reduced space through its RB, and the relations between both. For illustration purpose, one can refer to Fig. 1 which will be detailed in the following.

Let us denote ℎ as the initial complete functional space that contains the discrete solution computed using the HR model. Additionally, we introduce the basis of shape functions Φ( ) = { ( )} ℎ =1 that generates ℎ , meaning that ℎ = span Φ( ) with ( ℎ ) = ℎ . Then, let us define the reduced space ⊂ ℎ such that = span Ξ( ), where Ξ( ) = { ( )} =1 constitutes the RB. Therefore, we have ( ) = , and we expect << ℎ . Ultimately, the initial complete solution ℎ ( ; ) and the reduced one ( ; ) can be written in algebraic form as follows:

ℎ ( ; ) = ( ) ℎ ( ) ; ( ; ) = ( ) ( ). (2)
Here, ( ) (respectively, ( )) and ℎ ( ) (respectively, ( )) are vectors in ℝ ℎ (resp. in ℝ ) that collect the shape functions and associated coefficients. In the remainder of this Section, since the spaces ℎ and are the same for any ∈ , we will omit the dependence of the solutions on these parameters for more simplicity in the notations.

Remark 1. Note that the definition of ℎ encompasses all common discretization methods encountered in scientific computing, such as FE, FV, Discontinuous Galerkin, and so on. For instance, Φ( ) simply represents the nodal Lagrange shape functions in the FE context or corresponds to constant functions per cell in a FV framework. In our application, we will focus on the FV method, but it is worth noting that the proposed approach is applicable to other contexts, as will be underlined later.

As a first relation between the two spaces, we can construct the matrix ∈ ℝ ℎ × , which encodes the change of variables from the RB Ξ( ) to the complete one Φ( ). More precisely, it is defined as:

= 1 | … | , ( 3
)
where denotes the vector gathering the coefficients of function ( ) in the complete basis Φ( ). Thus, we can write:

( ) = ( ), (4)
and it follows that any reduced solution ( ) ∈ can be represented in the complete basis Φ( ) as:

( ) = ( ) = ( ) = ( ) ℎ . ( 5
)
where

ℎ = ∈ ℝ ℎ (6)
refers to the vector that contains the coefficients of the reduced solution ( ) in the complete basis Φ( ).

Next, as a second relation, we can express the orthogonal projection from the complete space ℎ onto the reduced space . In the continuous setting, the projection is defined by the projection operator  as follows:

 ∶ ℎ → ℎ ( ) ↦ ∑ =1 ( ℎ ( ), ( )) ℎ ( ), (7)
where (⋅, ⋅) ℎ represents a scalar product on ℎ . In our application, since the numerical solution is discretized using piecewise constant functions (see Remark 1), we will simply consider the 2 -scalar product. With this in mind and transitioning to the algebraic setting, the projection matrix ∈ ℝ ℎ × ℎ reads:

∀ ℎ ∈ ℝ ℎ , ℎ = ℎ = ℎ , (8)
where

ℎ = ℎ ∈ ℝ (9)
refers to the vector that collects the coefficients in the RB Ξ( ) of the projection of ℎ ( ) onto the reduced space (see Fig. 1). With all this in hand, we can eventually define the projection error of ℎ ( ) onto in the algebraic form as follows:

ℎ = ℎ - ℎ = ℎ - ℎ = ℎ - ℎ , ( 10
)
where ℎ ∈ ℝ ℎ (see Fig. 1 again).

The key point now is to build a RB, or equivalently the matrix , such that is small and ℎ ≈ ℎ , meaning ℎ ≈ 0 for any ∈ . This will be performed with the classical POD method, which is reviewed in next Section.

The POD reduced basis

The POD method is far from recent (see, e.g., [START_REF] Leask | The structure of inhomogeneous turbulent flows[END_REF] for an application in turbulent fluid mechanics) and is now widely used in many fields of scientific computing (refer to the following books [START_REF] Quarteroni | Reduced basis methods for partial differential equations: an introduction[END_REF]13,[START_REF] Benner | Model Order Reduction. Volume 2: Snapshot-Based Methods and Algorithms[END_REF] to name a few). The interest of the is the change of variables matrix; ℎ is the projection vector in ( ) of ℎ onto ; and ℎ is the corresponding projection error (expressed in basis ( )).

POD strategy is that it leads to a RB that is optimal in a chosen norm with respect to a collection of complete solutions.

In this regard, it may be noted that the Principal Component Analysis (PCA) is essentially the same as the POD method but applied to problems in statistics (the aim of the latter being to convert primitive variables into uncorrelated ones).

To briefly present the construction of the POD RB, let us start by considering the set = { } =1 ∈ () , which comprises parameter samples obtained by sampling the hypercube . These parameters are used to generate a set of complete numerical solutions with the HR model. These solutions are referred to as snapshots and their corresponding coefficient vectors are stored in the so-called snapshot matrix ∈ ℝ ℎ × , defined as:

= ℎ ( 1 )| … | ℎ ( ) . ( 11
)
The principle of the POD method is then to perform a Singular Value Decomposition (SVD) of the rectangular matrix and to retain the most significant left singular vectors, following the Schmidt-Eckart-Young theorem. These retained singular vectors constitute the POD modes that form the RB.

In our application, we will be in the case of > ℎ since numerous snapshots will be necessary to accurately capture the underlying physics of the model (see Section 4). For the SVD, we therefore build the correlation matrix as follows:

= . ( 12
)
Note that, again, we use here the 2 -scalar product due to the regularity of our discrete fields (see Remark 1). The resulting matrix ∈ ℝ ℎ × ℎ is the smallest matrix that represents the correlation between the snapshots, and it is obviously dense, symmetric, and positive. The next step is to compute the eigenvalues { } ℎ =1 ( ∈ ℝ+) and eigenvectors { } ℎ =1 ( ∈ ℝ ℎ and ‖ ‖ 2 = 1) of , which correspond, respectively, to the squares of the singular values and the left singular vectors of : = , 1 ≤ ≤ ℎ .Finally, as mentioned above, (introduced in Eq. ( 3)) is simply constructed with the vectors { } =1 associated with the largest { } =1 (and, again, we expect << ℎ ).

By construction, is an orthonormal basis of , meaning that = , where is the identity operator (of dimension ). However, ≠ ℎ since ℎ = ℎ only approximates ℎ due to the fact that ⊂ ℎ (remind Eq. ( 8) and Fig. 1). An error estimation is actually available in this context. More precisely, let us define the space of semi-orthonormal 1 matrices of dimension :

 ⊥ = { ∈ ℝ ℎ × , = }.
The interest of the POD matrix is that it minimizes, over all possible -dimensional semi-orthonormal matrices ∈  ⊥ , the projection error [START_REF] Mocayd | Polynomial surrogates for open-channel flows in random steady state[END_REF] onto the whole set of snapshots; that is, the sum of the errors between all the snapshots { ℎ ( )} =1 and their projection onto the associated subspace spanned by . Mathematically, this means that among the semi-orthonormal matrices of dimension , the POD one is the optimal one and it satisfies:

∑ =1 ‖ ℎ ( ) - ℎ ( )‖ 2 2 = min ∈ ⊥ ∑ =1 ‖ ℎ ( ) - ℎ ( )‖ 2 2 = ℎ ∑ = +1 (13)
This result derives from the Schmidt-Eckart-Young theorem, see, e.g., [START_REF] Quarteroni | Reduced basis methods for partial differential equations: an introduction[END_REF].

Benefiting from the above error estimation, the dimension of the RB can be fairly defined. In practice, we can compute the following ratio:

( ) = ∑ =1 ∑ ℎ =1 , (14)
and choose a target expressed as 1 -2 . is then defined as the minimal value such that:

( )) ≥ 1 -2 , ( 15
)
meaning that the signal of the snapshots captured by the most significant POD modes is conserved (in the 2 -norm here) within (100 -)%. The chosen value for will be specified and accounted for in Section 4.

The construction of the RB is obviously performed during the offline phase. This actually constitutes the initial step of the present reduction method (refer to Algorithm 1 for an overview of the complete offline phase carried out in this work).

1 Semi-orthonormal property denotes here the orthonormal property for rectangular matrices.

Once again, it is worth noting that the strategy employed so far can be applied to any -parametrized model, solved by any discretization method. It only requires the HR solutions to form the snapshot matrix (see Eq. ( 11)).

Remark 2. Let us mention that if ℎ >

, the alternative correlation matrix = can be built to ensure minimal computational cost and memory storage. In this case, the left singular vectors of have to be computed from the eigenvectors of , see, e.g., [START_REF] Quarteroni | Reduced basis methods for partial differential equations: an introduction[END_REF] for details.

Remark 3. It may also be noticed that other scalar products and associated norms like, e.g., the 1 (Ω)-norm, can be considered depending on the initial numerical model considered. To do so, one only has to modify the projection [START_REF] Cohen | Approximation of high-dimensional parametric pdes[END_REF] and the correlation matrix (12) accordingly.

The POD-NN method

The goal now is to compute, for a new set of parameters ∈  which does not belong to , a good approximation of the complete solution ℎ ( ; ) (or equivalently of the complete solution vector ℎ ( )) in the reduced space . To achieve this, the novel method proposed in the present study relies on an enriched version of the POD-NN method first introduced in [START_REF] Hesthaven | Non-intrusive reduced order modeling of nonlinear problems using neural networks[END_REF]. In the present section, we detail the original POD-NN method before enriching it in the next section.

The POD-NN method basic idea is to start from the projection defined in Eq. ( 8) (see also Fig. 1) and to learn, using a deep NN, the coefficients in the RB Ξ( ) of the projection of ℎ ( ) onto ; that is, to learn ℎ ( ) ∈ ℝ , see Eq. ( 9). More precisely, let us introduce the following application:

1 ∶  ⊂ ℝ → ℝ ↦ ℎ ( ) = ℎ ( ). (16)
It maps each set of input parameters ∈  to the coefficients ℎ ( ) in the RB Ξ( ) of the projection of the complete solution ℎ ( ) onto .

In the offline phase, this non-linear mapping 1 is approximated by employing an artificial NN, particularly a multilayer perceptron (see Fig. 2). Artificial NNs have proved to be universal approximators, due to their rich interpolation space [START_REF] Hornik | Multilayer feedforward networks are universal approximators[END_REF][START_REF] Chen | Universal approximation to nonlinear operators by neural networks with arbitrary activation functions and its application to dynamical systems[END_REF]. In practice, training the aforementioned NN can be performed using supervised learning based on inputoutput pairs obtained from the snapshots (see Eq. ( 11)). Denoting by  1 ( 1 )( ) the considered NN with associated parameters (weights and biases) gathered in vector 1 , this leads to minimizing the following loss function for the training:

 1 ( 1 ) = 1 ∑ =1 ‖ 1 ( 1 )( ) - ℎ ( )‖ 2 2 , ( 17
)
where, once again, the the same reason, the number of input neurons of the network must be low, i.e., the number of physical parameters ( ∈ ℝ ) must be small, which is also a limitation of more conventional reduced methods based on POD.

Next, during the online phase, all that is required is to perform a forward pass of the NN  1 given a new (unseen) set of physical parameters and the optimized 1 obtained after training. This process computes the reduced coefficients ( ) and is known to be very fast, as the forward evaluation of a NN simply involves affine and element-wise activation transformations. This makes the online phase achievable in real-time. Finally, the reduced solution can be written in the original basis Φ( ) using Eq. ( 6): ℎ ( ) = ( ) ∈ ℝ ℎ (which is also very fast to perform). In particular, this allows the utilization of existing subroutines from the initial computational code for post-processing. The entire online phase performed in this work is summarized in Algorithm 2. The procedures mentioned above correspond to Steps 1 and 2 of the online phase.

Remark 5. In the proposed method, we ultimately only need to utilize the initial computational code associated with the HR model  in a black-box fashion to run it multiple times ( times) with different parameters. Subsequently, the method is entirely independent of the considered HR model  and underlying code since, in particular, the online phase is decoupled from  . Consequently, the method has the twofold advantages of being generic (it can be applied to any type of problem solved by any discretization method), and non-intrusive (we never need to delve into the code itself). Here, we leverage on these properties to apply the strategy to a non-linear problem with a non-affine dependence on the parameters and where a Galerkin approach is not used at the discretization level, which is a well-known scenario in which traditional reduced-order methods fail.

The Error-Aware POD-NN method

The POD-NN method considered so far requires the ability to strongly reduce the dimension of the problem. More precisely, we need to have =  [START_REF] Mocayd | Polynomial surrogates for open-channel flows in random steady state[END_REF], so that it is possible to employ a deep NN with a number of parameters

1
reasonable compared to the number snapshots (see Fig. 3 that illustrates the different constraints). Therefore, in some situations, especially in the case of convection-dominated flow problems, it is requested to further enhance the accuracy of the method to mitigate the Kolmogorov barrier to reducibility (see, e.g., [START_REF] Barnett | Neural-network-augmented projection-based model order reduction for mitigating the kolmogorov barrier to reducibility[END_REF] for the concept of reduction complexity). To achieve this in a straightforward manner with minimal computational overhead, we propose to integrate the ROM error estimation. Note that this generic idea has been developed, e.g., in [START_REF] Dunhui | Error estimation of the parametric non-intrusive reduced order model using machine learning[END_REF] within the framework of Gaussian processes. The idea is here to learn, in addition to the projection of ℎ ( ) onto , the projection error introduced in Eq. ( 10) (see also Fig. 1):

ℎ ( ) = ℎ ( ) - ℎ ( ), (18)
and to benefit from the latter to correct the reduced solution obtained so far. For this purpose, we follow the same procedure used to construct the reduced solution, but this time starting with the snapshot error matrix ∈ ℝ ℎ × :

= ℎ ( 1 )| … | ℎ ( ) , (19)
which collects the projection errors of the already computed snapshots { ℎ ( )} =1 (see Eq. ( 11)). Consequently, we first apply the POD strategy (same operations as in Section 3) to compute the corresponding error RB { ( )} =1

(where we again anticipate << ℎ ) or, in practice, the error RB matrix ∈ ℝ ℎ × :

= 1 | … | . ( 20
)
Note that also serves as the change-of-variable matrix from basis { ( )} =1 to the initial complete one Φ( ). At this stage, similarly as in Eq. ( 16), we introduce the following application:

2 ∶  ⊂ ℝ → ℝ ↦ ℎ ( ) = ℎ ( ), (21)
which maps each set of input parameters ∈  to the coefficients ℎ ( ) of the projection of the error ℎ ( ) onto the subspace spanned by { ( )} =1 . The remaining task of the extended offline phase is therefore to train a new artificial NN, denoted by  2 ( 2 )( ), to approximate the non-linear mapping 2 (refer to Fig. 2 once again). As for the reduced solution, the latter is performed by minimizing the following loss function:

 2 ( 2 ) = 1 ∑ =1 ‖ 2 ( 2 )( ) - ℎ ( )‖ 2 2 . ( 22
)
For a better understanding and conciseness, the additional operations to learn the projection error in the offline phase are outlined in Algorithm 1 (see Steps 6, 7 and 9).

Next, during the online phase, it is possible to estimate the projection error at a new set of physical parameters by simply performing a forward pass of the NN  2 using the optimized 2 obtained after training. The reduced projection error can also be written in the original basis Φ( ) such as: ℎ ( ) = ( ) ∈ ℝ ℎ . In this basis, we eventually complement the former reduced solution ℎ ( ) of Section 3 with the present projection error, which yields the final corrected reduced solution as follows:

ℎ ( ) + ℎ ( ) = ̃ ℎ ( ) ≈ ℎ ( ). (23)
The overall online phase is outlined in Algorithm 2; more precisely, Steps 3 to 5 constitute the additional operations to perform to correct the initial reduced solution with a learning of the projection error. This phase can be carried out very quickly since it merely consists in performing forward passes of trained NNs.

As mentioned in the introduction, because our approach incorporates NNs with POD and takes into account the projection error, we refer to it as an Error-Aware POD-NN reduction method. Complementing the POD-NN strategy with the projection error enables to attain higher accuracy while preserving a limited number of POD modes. This is not equivalent as setting the basic POD-NN method with ( + ) modes. Indeed, the two NN outputs are uncorrelated; more importantly, the Error Aware POD-NN method enables to limit the output dimension of the NNs, thus providing the opportunity to employ smaller NNs (i.e. with less parameters), and therefore smaller datasets (the snapshots) to train the NNs (see Fig. 3).

Algorithm 1:

The Error-Aware POD-NN reduction method: the offline phase function [ ,  1 ( 1 ), ,  2 ( 2 ) ] = POD-NN-OFFLINE( , ).

1. Generate the reference parameters values by sampling the hypercube : = { } =1 ∈ () .

2. Compute the associated complete numerical solutions with the HR model and form the snapshot matrix:

= ℎ ( 1 )| … | ℎ ( ) ∈ ℝ ℎ × .

Build up the correlation matrix

= ∈ ℝ ℎ × ℎ . 4. Compute its eigenelements ( , ) ℎ =1 , with ∈ ℝ + and ∈ ℝ ℎ ‖ ‖ 2 = 1 : = , 1 ≤ ≤ ℎ .
5. Take the vectors { } =1 associated with the largest { } =1 (by using criterion [START_REF] Hornik | Multilayer feedforward networks are universal approximators[END_REF] based on )

to form the RB matrix:

= 1 | ⋯ | ∈ ℝ ℎ × with << ℎ .
6. By making use of the POD projection error defined as ℎ ( ) = ℎ ( )ℎ ( ), form the snapshot error

matrix: = ℎ ( 1 )| … | ℎ ( ) ∈ ℝ ℎ × .
7. Apply the POD procedure to (see Steps 3 to 5 above) to compute the error RB matrix:

= 1 | … | ∈ ℝ ℎ × with << ℎ .
8. Train artificial NN  1 which approximates application 1 defined in [START_REF] Diederik | Adam: A method for stochastic optimization[END_REF], by minimizing the loss function

 2 ( 1 ) = 1 ∑ =1 ‖ 1 ( 1 )( ) - ℎ ( )‖ 2
2 . 9. Train artificial NN  2 which approximates application 2 defined in [START_REF] Rasheed | Digital twin: Values, challenges and enablers from a modeling perspective[END_REF], by minimizing the loss function

 2 ( 2 ) = 1 ∑ =1 ‖ 2 ( 2 )( ) - ℎ ( )‖ 2 2 .
Algorithm 2: The Error-Aware POD-NN reduction method: the online phase (real-time computations)

function [ ℎ ( )] = POD-NN-ONLINE( ; B ,  1 ( 1 ), B ,  1 ( 2 )).
1. Given a new set of input physical parameters gathered in , and taking the optimal 1 obtained from Algorithm 1 (Step 9), evaluate the output ( ) of artificial NN  1 .

Deduce the reduced solution in basis

Φ( ): ℎ ( ) = ( ) ∈ ℝ ℎ .
3. Given the same new set of input parameters , and taking the optimal 2 obtained from Algorithm 1 (Step 10), evaluate the output ( ) of artificial NN  2 .

4. Deduce the reduced projection error in basis Φ( ): ℎ ( ) = ( ) ∈ ℝ ℎ .

Deduce the final corrected reduced solution: ̃

ℎ ( ) = ℎ ( ) + ℎ ( ).

Numerical experiments and discussions

In this section, we numerically evaluate the proposed EA-POD-NN method on a real-world flooding event of relatively high magnitude occurred on the Aude River in the southeastern region of France in 2018. The flood hydraulic model was built in [START_REF] Garambois | Coupled 2D hydrologic-hydraulic catchment scale flood modeling with data assimilation capabilities: the DassHydro platform[END_REF] based on HR terrain elevation data provided by IGN (French geographical national institute) and hydrometeorological data from SCHAPI (French national flood forecasting center) and Météo France. The DassFlow 2D software [START_REF]Open source computational software DassFlow (Data Assimilation for Free Surface Flows)[END_REF][START_REF] Monnier | Inverse algorithms for 2D shallow water equations in presence of wet dry fronts: Application to flood plain dynamics[END_REF] is used for numerically solving the 2D SWEs (see Section 2.2) over an unstructured mesh covering the spatial domain Ω that contains a floodplain located at the confluence of the Aude and Fresquel rivers. The latter is inflowed by two real hydrographs , ( , ) at = 2 upstream points on Ω the border of Ω (see Fig. 4).

Again, the primary focus is on the parameter related to the inflow discharge functions which is the most important factor controlling flooding dynamics, particularly crucial in operational context where fast and accurate hydraulic simulations of flooding are needed after hydro-meteorological forecasts of discharge. In this context, two test cases are carried out, representative of real complexity of concomitant or not flood inflows:

1. Case with simultaneous inflows where the classic POD-NN method only is sufficient to reduce the problem.

2. Case with non-simultaneous inflows, hence which introduces more complex non-linear waves interactions, where the POD-NN method struggles to accurately reproduce the solution. This is where the proposed projection error correction approach becomes essential to enhance the accuracy of the ROM.

Therefore, the investigated test cases demonstrate the limitations of the classic POD-NN method for reducing a SW model and showcase the significant improvement in the ROM's quality achieved through the proposed EA-POD-NN approach.

Let us introduce at this stage some notations for the physical and numerical quantities that will be used throughout this section, especially in the majority of the figures:

• = √ 2 + 2
denotes the velocity norm in ( -1 ).

• ℎ in ( ) (resp. in ( -1 )) denotes the reference 2D hydraulic model water depth (resp. the reference 2D hydraulic model velocity norm).

• ℎ □ in ( ) (resp. □ in ( -1 )) denotes the approximated water depth (resp. the approximated velocity norm),

where □ stands for POD-NN or EA-POD-NN depending on the considered method for the computation.

• ℎ ( , ) = ℎℎ □ (resp. ( , ) = -□ ) denotes the difference between ℎ (resp.

) and ℎ □ (resp. □ ).

• Finally, ℎ ( ) (resp. ( )) denotes the spatial average over Ω of the error on water depth (resp. of the error on velocity norm).

we have , ( peak,k , max,k ; ) for = 1, 2. The concept of utilizing max and peak as representative parameters is rooted in their ability to characterize the fundamental aspects of flood initiation. max quantifies the peak discharge intensity, reflecting the maximum flow rate reached during the flood. On the other hand, peak encapsulates the time it takes for the river to transition from its initial state to the point of maximum discharge, offering insights into the flood's temporal evolution. This simple parameterization enables to represent the inflow hydrographs while preserving their essential features, i.e. peak flow magnitude and timing, therefore ensuring consistent flood inundation modeling as shown by the limited errors in space and time on submersion depth and velocity norm between the reference 2D SW model run with real or triangular hydrographs, see Fig. 4(Bottom). In particular, the errors are relatively lower at peak time and maximum flooding rather than during more dynamic phases of submersion and recessions i.e. corresponding to hydrograph rising limb and also during recession.

As mentioned previously, two test cases are considered with different input parameters, corresponding to increasingly complex inflow signals and resulting SW model responses over the floodplain:

1. Simultaneous inflow hydrographs: = ( peak , max,1 , max,2 ; ), so = 4.

2. Non-simultaneous inflow hydrographs: = ( peak,1 , max,1 , peak,2 , max,2 ; ), so = 5.

For both cases, the different steps of Algorithm 1 are performed to achieve the offline phase of the reduction method (see Section 3). Case-specific parameters values are detailed for each experiment in the following paragraphs.

Snapshots matrix construction

The starting point of the learning phase (offline phase) is the set of parameters and the corresponding snapshots for the inflows. This selection allows us to concentrate on the period where the flood dynamics exhibit significant changes, aligning with the goals of SW model reduction. Also, this approach ensures that the NNs focus on capturing meaningful flow variations that occur starting from the rising phase.

Next, for each pair of inflow hydrographs ; =1,2 ( , ) resulting from the parameter set { } =1 , theparameterized hydraulic model  is run to generate the HR solution ℎ ( ). In the discrete setting, the dimension of the output vectors ℎ ( ) = [ , , ] is ℎ = 3 × = 21801 where = 7267 counts the number of cells of the mesh covering Ω. Note eventually that each model output , or consists in temporal snapshots of model outputs spatial fields. This enables to obtain the ( ℎ × ) snapshot matrix as follows:

= [ , , ] = ⎡ ⎢ ⎢ ⎢ ⎢ ⎣ ( 1 ; )| ⋯ | ( ; )| ( 1 ; )| ⋯ | ( ; ) ( 1 ; )| ⋯ | ( ; )| ⎤ ⎥ ⎥ ⎥ ⎥ ⎦ .
Once is built, one can perform the other steps of the offline phase for the POD-NN and EA-POD-NN algorithms (see Algorithm 1).

RB matrix construction

The second step involves constructing the reduced matrix ∈ × ℎ by applying the POD method to the snapshots matrix . Indeed, since we have ℎ ≤ , we compute the POD modes of the correlation matrix = , see Eq. ( 12). To truncate the first modes, we ask to conserve 99.99% of the energy of the system. In other words, given a tolerance 2 = 10 -2 and by using Eq. 15, we obtain the reduced dimension = 5. Following the same approach for the error matrix , by taking 2 = 10 -1 , the reduced dimension for errors is = 17.

The plot in log-scale of the eigenvalues of the correlation matrix is given in Fig. 5. Note that the behaviour of the eigenvalues of the correlation matrix of the errors ( = ⋅ ) is similar for both tested cases with simultaneous and non-simultaneous inflows.

Remark 7. The Kolmogorov -width [START_REF] Cohen | Approximation of high-dimensional parametric pdes[END_REF][START_REF] Quarteroni | Reduced basis methods for partial differential equations: an introduction[END_REF] enables to measure the degree of reductibility of model. A link between the Kolmogorov -width and the POD model construction is presented in [START_REF] Quarteroni | Reduced basis methods for partial differential equations: an introduction[END_REF]. This link relies on the decay of the eigenvalues. In Fig. 5, we observe a slow decay of the eigenvalues, which is typical in the case of non-linear convectiondominated problems [START_REF] Barnett | Neural-network-augmented projection-based model order reduction for mitigating the kolmogorov barrier to reducibility[END_REF]. That is, the error in the right hand side of Eq.( 13) is large. This turns the SWEs reduction problem to be challenging.

Remark 8.

For the EA-/POD-NN approach, certain numerical constraints must be adhered to. These constraints are summarized in Tab 2.

Case with simultaneous inflows

Let us now focus on the results obtained for the first test case involving simultaneous inflow hydrographs. The conditions that must be adhered to in the the offline phase for both numerical test cases. These conditions concern the following parameters: representing the number of snapshots, param denoting the number of parameters of the NN, and indicating the number of neurons in the last layer of the NN.

The parameter space

As mentioned above, the same peak for both inflows is considered here, so that the parameter reads = ( peak , max,1 , max,2 , ) and is of dimension The snapshot matrix ∈ ℝ ℎ × is constructed with a number of snapshots fixed as = × = 189 000 (with = 1 000 simulations using the DassFlow software and = 189 writing times, i.e. a = 18min writing time step). Note that quasi similar results could be obtained with = 36min, as will be chosen for the second case with non-simultaneous hydrographs because of memory limitations. Note finally that such writing time steps remain significantly lower than the flood propagation time that is around 3ℎ over the studied zone. The sampled parameter sets , blue points for offline phase, red and brown points being random parameter sets that will be used for online prediction; (left) For the first hydrograph ,1 on the Aude River and (right) for the second hydrograph ,2 on the Fresquel River (corresponding to red and blue hydrographs, respectively, in Fig. 4). (Bottom) Input inflow hydrographs for (left) red and (right) brown points that will be used for online prediction (vertical dashed lines correspond to the time of maximum flooding at which spatial error analysis will be performed).

The training phase (during the offline phase)

For the training of the NNs, the data of size = 189000 is classically partitioned into 80% for the actual training, 10% for the test and 10% for the validation. To ensure that the number of NNs parameters remains small relative to the data dimension, we determine the number of parameters of the NN  1 as 3% (resp. of the NN  2 as 5%) of the training data dimension , yielding 1 = 4405 for  1 with = 5 for the output layer, and 2 = 9977 with = 17 for the output layer of  2 in the case of simultaneous inflows (resp. = 16 in the case of non simultaneous inflows).

The number of hidden layers for  1 is 10 with 20 neurons per layer and for  2 it is 10 with 30 neurons per layer. The chosen activation function is 'RELU'. We train the NNs over 10000 epochs which leads to = 10 -7 in terms of loss function. The loss function on the validation data set is observed to well decrease too. An Adam optimizer [START_REF] Diederik | Adam: A method for stochastic optimization[END_REF] is used for the gradient descent to minimize the loss function. The latter computes the loss between the HR projected output

ℎ ( ) = ℎ ( ) (resp. ℎ ( ) = ℎ ( )), 1 ≤ ≤ , and the NN output  1 ( 1 )( ) (resp.  2 ( 2 )( ))
using the Mean Square Error (MSE), in line with Eqs. ( 17) and [START_REF]Open source computational software DassFlow (Data Assimilation for Free Surface Flows)[END_REF].

The online prediction phase

After successfully completing the training phase and optimizing the NNs parameters, we move to the prediction stage, where the model demonstrates its ability to rapidly forecast hydraulic variables (ℎ, , )( , ) for new parameter values, denoted as new . These new parameter values are in the parameter space range ( new ∈ ) but distinct from all the points in (let us recall that these new parameter values correspond to the red and brown points in Fig. 6 (Top), which are different from the blue points that form ). In contrast, the chosen parameter defines the times at which the snapshots are performed in the offline phase and the prediction time in online phase. These two points new are carefully selected: the first (brown point in Fig. 6(Top)) is chosen to be in close proximity to the training phase data points, while the second (red point in Fig. 6(Top)) is taken at a considerable distance from the training data points.

This deliberate selection allows us to rigorously evaluate the predictive capabilities and precision of our approximation model. The performances obtained in reducing the 2D SW model are depicted in their spatio-temporal dimensions with several graphes (see Figs. 7 to 10), especially at maximal flooding which is a crucial instant for operational forecasting, with related spatial statistics of ℎ ( , ) and ( , ), but also in terms of temporal variability of ℎ ( ) and ( ).

More precisely, spatial performances of the online solution □ ( ) at maximum flooding are first shown in Figs. 7 and8 for ,1 = (11.11, 800, 380; 24) (red point in Fig. 6) and ,2 = (12, 1000, 340; 25.2) (brown point in Next, the temporal variation of errors spatial averages | ℎ |( ) (or | |( )) in Fig. 10 shows that submersion depth (or flow velocity magnitude) stays below 0.1 (or 0.1 -1 ) over the whole spatio-temporal domain for both reduction approaches. For instance, the median error value is 0.006 (or 0.005 -1 ) with EA-POD-NN and 0.0252 (or 0.016 -1 ) with POD-NN for ,1 and 0.007 (or 0.007 -1 ) with EA-POD-NN and 0.025 (or 0.019 -1) with POD-NN for ,2 . Also, the maximum value of the temporal variation of errors spatial averaged

| ℎ |( ) (or | |( ))
is 0.04 (or 0.03 -1 ) with EA-POD-NN and 0.12 (or 0.08 -1 ) with POD-NN for ,1 and 0.04 (or 0.03 -1 ) with EA-POD-NN and 0.10 (or 0.08 -1 ) with POD-NN for ,2 .

Overall, both the classic POD-NN and proposed EA-POD-NN demonstrate good performance in reproducing flow depth and velocity in this moderately complex dynamic flooding scenario. The accuracy achieved in both spatial and temporal aspects of flow variables highlights the applicability of both methods for reducing the 2D SW model in this case with simultaneous inflows. Nevertheless, a slight improvement in accuracy is already noticeable with the proposed EA-POD-NN, suggesting that this approach may also prove high-performance in more complex flooding scenarios, as demonstrated below.

Case with non-simultaneous inflows

The parameter space and training phase

This case is for non-simultaneous inflow hydrographs and the parameter reads = ( peak,1 , max,1 , peak,2 , max,2 ; ) numerical memory capacity constraints but this remains significantly lower than the flood propagation time that is around 3ℎ over the studied zone. Note that with = 18min as in the previous case (results not presented): (i) the POD-NN method leads to an inaccurate reduction in terms of predicted physical variables, (ii) the EA-POD-NN gives model reduction performances that are comparable to those obtained with = 36min.

and
Note that the training phase is performed in this case with the same tolerances and reduced dimension = 5 as in first case (cf. section 4.1.3). The only difference is that the reduced dimension for errors reduced matrix ∈ ℝ ℎ × is = 16.

The prediction phase (online phase)

Evaluation against reference SW model with triangular inflows. For the prediction stage, again, we perform a prediction for a new parameter value, denoted as new . This new parameter value is in the parameter space range, i.e.

new ∈ , but is distinct from all the calibration parameter sets that compose . This new value new corresponds to the red point in Fig. 11 and is selected at a considerable distance from the training data points. Fig. 12 shows a The prediction performances with this online phase are analyzed, in terms of ℎ ( , ) and ( , ), for both POD-NN and EA-POD-NN methods, with similar graphes as before: spatial fields of errors at maximum flooding are shown in Fig. 12 and the corresponding values are depicted by histograms of spatial values of at maximum flooding time in function of physical variable values in Fig. 13. Again, reduction error is significant with the proposed EA-POD-NN method over this more complex case with non-simultaneous hydrographs. This is confirmed by the temporal variation and the Cumulative Distribution Function (CDF) of absolute error spatial average shown in Fig. 13 which clearly shows that the error with EA-POD-NN is 10 times lower than with POD-NN. Namely, the median value of error on ℎ (or ) is 0.02 (or 0.012 -1) with POD-NN and 1.6 × 10 -3 (or 9 × 10 -4 -1) with EA-POD-NN. Also, the temporal times those one obtained by the EA-POD-NN (0.20 (or 0.11 -1)). Again, the EA-POD-NN is able to reduced the SW model in this complex case with non-simultaneous inflows where the parameter = (11.94, 922, 9.16, 198, 25) is taken out of the learning set ( ∉ ). While, the POD-NN is no more applicable to reduce the SW model.

Conclusion

From the obtained results, and in line with the objectives of this research, the following conclusions can be made:

• The basic POD-NN method as proposed in [START_REF] Hesthaven | Non-intrusive reduced order modeling of nonlinear problems using neural networks[END_REF][START_REF] Wang | Non-intrusive reduced order modeling of unsteady flows using artificial neural networks with application to a combustion problem[END_REF], see also Section 3, enables to reduce the 2D SW model on a moderatly complex case in the sense of the tested dynamics with simultaneous inflow hydrographs, parameterized in low dimension.

• In the case of more complex non-linear waves, here produced by the two non-simultaneous inflows hydrographs, the basic POD-NN method fails to be sufficiently accurate. • The proposed EA-POD-NN method, that additionally learns the projection error onto the RB, enables to reduce much more accurately the 2D SW model (comparisons based on the same parameter sampling). In the present application, this enables to well simulate the complex flood case with the two non simultaneous hydrographs.

The spatio-temporal hydraulic variables are quite accurately predicted.

• The proposed EA-POD-NN method, is able to predict with a good accuracy ROM solution of SW model in extrapolation case.

The EA-POD-NN method paves the way for immediate application to highly non linear PDEs systems, like the 2D SWEs in complex configurations. The obtained surrogate model could then be used as a digital twin [START_REF] Chinesta | Virtual, digital and hybrid twins: a new paradigm in data-based engineering and engineered data[END_REF][START_REF] Rasheed | Digital twin: Values, challenges and enablers from a modeling perspective[END_REF],

allowing for rapid flood forecasts given the registered scenario. Immediate research work will pertain in improving the capability of the method to be applicable in higher dimensional input parameters. In a flood modeling context, this may be more complex inflow hydrographs or spatially distributed friction fields. This should enable to generalize its applicability to the variety of real hydrological signals and rivers-floodplains connectivities, as well as to optimization or uncertainty quantification tasks. Moreover, the required snapshots number may be reduced by enforcing the model residual to vanish at some points, as in the Physics Informed Neural Networks methods for instance [START_REF] Raissi | Physics-informed neural networks: A deep learning framework for solving forward and inverse problems involving nonlinear partial differential equations[END_REF][START_REF] Chen | Physics-informed machine learning for reduced-order modeling of nonlinear problems[END_REF]. acknowledge SCHAPI-DGPR, Météo-France and IGN for providing data used in this work.

Figure 1 :

 1 Figure 1: Schematic representation of the different spaces: ℎ is the initial complete space (associated basis ( ));⊂ ℎ is the reduced space (associated basis ( )); ℎ and are the dimensions of ℎ and , respectively; is the change of variables matrix; ℎ is the projection vector in ( ) of ℎ onto ; and ℎ is the corresponding projection error (expressed in basis ( )).

Figure 2 :

 2 Figure 2: NN approximation of the map 1 , see Eq. (16) (or 2 , see Eq.(21)). Note that , and are on the order of 10. In practice the hidden layers are taken with the same number of neurons and the latter is larger than the numbers of input and output neurons (see Section 4 for the exact chosen values).

  Figure 2: NN approximation of the map 1 , see Eq. (16) (or 2 , see Eq.(21)). Note that , and are on the order of 10. In practice the hidden layers are taken with the same number of neurons and the latter is larger than the numbers of input and output neurons (see Section 4 for the exact chosen values).

Remark 4 .

 4 It is important to notice at this stage that the reduction step, which transforms ℎ ∈ ℝ ℎ into ℎ ∈ ℝ with << ℎ is crucial in the process. This step allows to obtain only a few ouput neurons for the NN (in our application = [START_REF] Mocayd | Polynomial surrogates for open-channel flows in random steady state[END_REF], see Section 4), making it possible to consider a deep NN with a number of parameters 1 (i.e., 1 ∈ ℝ 1 ) reasonable compared to the number of data (i.e., the snapshots) during training. Of course, for

Figure 3 :

 3 Figure3: Due to a high reduction complexity for convection dominated problems, see, e.g.,[START_REF] Barnett | Neural-network-augmented projection-based model order reduction for mitigating the kolmogorov barrier to reducibility[END_REF] for the concept of the Kolmogorov barrier to reducibility, it is requested to complement the basic POD-NN method with the projection error learning. This allows to limit the sizes of the reduced bases (here = (10) and = (10)), thus to use smaller NNs, and therefore smaller datasets (here snapshots).

Remark 6 .

 6 matrix . Here, the array = { } =1 of parameters sets is obtained by uniformly sampling times the discharge related parameters in considering upper and lower bounds around the reference values of the 2018 real flood hydrographs (the detailed values will be specified later). The last parameter defines the physical times at which the snapshots of the reference 2D SW model are computed. It is fixed a priori in function of the response dynamics of the physical model to reduce. Overall, we have = × . It is important to note that the snapshot time parameter defines the temporal repartition of the used physical model outputs which is therefore conditioning the capture of information in the reduction process of the physical model. Here hydraulic simulations were performed over a span of 60 hours. However, in the context of SW model reduction, we select the time interval starting from the point corresponding to the onset of rising water levels

Figure 5 :

 5 Figure 5: (Left) The decay of eigenvalues of the correlation matrix . (Right) The zoom on the first 10 eigenvalues. (Top) For the case with simultaneous hydrographs ( = 4). (Bottom) For the case with non-simultaneous hydrographs ( = 5).

= 4 .

 4 Following previous discussions, the parameters ranges are defined as follows: peak ∈ [11.2, 12.69], max,1 ∈ [600, 1200], and max,2 ∈ [100, 400]. Then, the snapshots time input parameter is chosen as ∈ [3.6, 60]. Recapitulating, one has  = [11.2, 12.69] × [600, 1200] × [100, 400] × [3.6, 60]. A visual representation of the parameter set sample = { } =1 with a focus around peak , max,1 , and max,2 is given in Fig. 6.

Figure 6 :

 6 Figure 6: Simultaneous inflows case. (Top) The sampled parameter sets , blue points for offline phase, red and brown points being random parameter sets that will be used for online prediction; (left) For the first hydrograph ,1 on the Aude River and (right) for the second hydrograph,2 on the Fresquel River (corresponding to red and blue hydrographs, respectively, in Fig.4). (Bottom) Input inflow hydrographs for (left) red and (right) brown points that will be used for online prediction (vertical dashed lines correspond to the time of maximum flooding at which spatial error analysis will be performed).

Fig. 6 )

 6 Fig. 6). Then, the spatial values of misfit to the reference 2D SW simulation ( ) are depicted at peak time in histograms in Fig 9 and in terms of Cumulative Distribution Function (CDF) of their spatial variability in Fig. 10. In particular, the median value of | ℎ ( , )| (or | ( , )|) is 10 -3 (or 0.7 × 10 -3 -1 ) with EA-POD-NN and 3 × 10 -3 (or 10 -3 -1 ) with POD-NN for ,1 . Also, the median value of | ℎ ( , )| (or | ( , )|) is 2 × 10 -3 (or 10 -3 -1 ) with EA-POD-NN and 5 × 10 -3 (or 4 × 10 -3 -1 ) with POD-NN for ,2 .

Figure 7 :

 7 Figure 7: Simultaneous inflows case. Comparison of the model outputs at maximum flooding time for = (11.11, 800, 380, 24.0) (red point in Fig. 6 (Top)), (Top) For the water depth ℎ. (Bottom) For the velocity norm . (Left) The reference solution: 2D SW model with triangular inflow, i.e. ℎ ( ) . (Middle) The difference ℎ or for the calssic POD-NN ℎ, -( ). (Right) The difference ℎ or for the proposed EA-POD-NN ℎ, - -( ).

Figure 8 :

 8 Figure 8: Simultaneous inflows case. Comparison of the model outputs at maximum flooding time for = (12, 1000, 340, 25.2) (brown point in Fig. 6 (Top)), (Top) For the water depth ℎ. (Bottom) For the velocity norm . (Left) The reference solution: 2D SW model with triangular inflow, i.e. ℎ ( ) . (Middle) The difference ℎ or for the classic POD-NN ℎ, -( ). (Right) The difference ℎ or for the proposed EA-POD-NN ℎ, - -( ).

Figure 9 :

 9 Figure 9: Simultaneous inflows case. The histogram of the mean and the StD of spatial values of | ℎ | and | | at maximum flooding time. (Left) For the classic POD-NN method. (Right) For the proposed EA-POD-NN method. (Top) For = (11.11, 800, 380, 24.0) (red point in Fig. 6. (Bottom) For = (12, 1000, 340, 25.2) (brown point in Fig.6).

Figure 10 :

 10 Figure 10: Simultaneous inflows case. Comparison between the classic POD-NN and the proposed EA-POD-NN methods performances. (Top) The CDF of the absolute difference | ℎ | or | |. (Bottom) The mean of the absolute difference | ℎ | or | | with respect to the real time. (Left) For = (11.11, 800, 380, 24.0) (red point in Fig. 6. (Right) For = (12, 1000, 340, 25.2) (brown point in Fig. 6).

Figure 11 :

 11 Figure 11: Non-simultaneous inflows case. (Top) The sampled parameter sets , blue points for offline phase, red point being random parameter sets for online prediction; (left) For the first hydrograph ,1 on the Aude River and (right) for the second hydrograph ,2 on the Fresquel River (red and blue hydrographs in Fig. 4). (Bottom left) The sampled parameter sets with focus on time of peak peak,1 for the first hydrograph on the Aude River and peak,2 for the second hydrograph on the Fresquel River. (Bottom right) Input inflow hydrographs for red point used for online prediction -vertical dashed lines correspond to the time of maximum flooding at which spatial error analysis are performed in each case.

  variation of errors spatial averages | ℎ |( ) (or | |( )) are shown in Fig. 14. In the latter, the maximum values of the temporal variation of errors spatial averaged is 0.002 (or 0.002 -1) with EA-POD-NN and 0.44 (or 0.27 -1) with POD-NN. The maximum value of the temporal variation of errors spatial averaged | ℎ |( ) (or | |( )) is 0.02 (or 0.02 -1) with EA-POD-NN and 1.11 (or 0.79 -1) with POD-NN.These relatively good performances show the capability of EA-POD-NN method to accurately reduce the 2D SW model in this case with non-simultaneous inflows. Results also show the limitations of the classical POD-NN to reduce of the SW model for this more complex case with non-simultaneous hydrographs.Evaluation against reference SW model with real inflows. Now let us go further and test the predictability of our model first when taking the real hydrographs in Fig.4(right) as input inflows for the reference 2D HR hydraulic model; second when taking triangular hydrographs corresponding to a parameter . The latter is given by = (11.94, 922, 9.[START_REF] Diederik | Adam: A method for stochastic optimization[END_REF] 198,[START_REF] Wang | Non-intrusive reduced order modeling of unsteady flows using artificial neural networks with application to a combustion problem[END_REF] and it is out of the range of the space parameter points given in Fig.11. The predicted hydraulic variables at peak flow are comparable to those of the reference HR hydraulic model as shown in spatial maps of Fig.15with real like hydrographs (resp. in 18 with triangular hydrographs corresponding to the paramerter =(11.94, 922, 9.16, 198, 25)). The associated error histogram is shown in Fig.16with real like hydrographs (resp. in 19 with triangular hydrographs). The temporal variation of spatially averaged error analysis is presented in Fig.17with real like hydrographs (resp. in 20 with triangular hydrographs). These results show the relatively good accuracy of the reduction compared to the reference model. For instance, in the case of real input hydrographs, the maximum value of the temporal variation of errors spatial averaged with POD-NN on ℎ (or ) is 0.99 (or 0.73 -1) and equal to 3 times those one obtained by the EA-POD-NN (0.34 (or 0.22 -1)). Moreover, with triangular hydrographs corresponding to the parameter =(11.94, 922, 9.16, 198, 25), the maximum value of the temporal variation of errors spatial averaged with POD-NN on ℎ (or ) is 1.03 (or 0.74 -1) and equal to 5

Figure 12 :

 12 Figure 12: Non-Simultaneous inflows case. Comparison of the model outputs at maximum flooding time for = (8.66, 1100, 8.38, 300, 22.2) (red point in Fig. 11), (Top) For the water depth ℎ. (Bottom) For the velocity norm . (Left) The reference solution: 2D SW model with triangular inflow, i.e. ℎ ( ) . (Middle) The difference ℎ or for POD-NN ℎ, -( ). (Right) The difference ℎ or for EA-POD-NN ℎ, - -( ).

Figure 13 :

 13 Figure 13: Non-simultaneous inflows case. (Top) The histogram of the mean and the StD of spatial values of the absolute difference | ℎ | and | | at maximum flooding time. (Left) with POD-NN method. (Right) with EA-POD-NN method for = (8.66, 1100, 8.38, 300, 22.2) (red point in Fig. 11. (Bottom) The CDF of the absolute absolute difference | ℎ | and | |. (Left) with POD-NN method. (Right) with EA-POD-NN method for = (8.66, 1100, 8.38, 300, 22.2) (red point in Fig. 11.

Figure 14 :

 14 Figure 14: Non-simultaneous inflows case. The mean on the absolute difference | ℎ | or | | with respect to the real time with POD-NN and EA-POD-NN methods for = (8.66, 1100, 8.38, 300, 22.2) (red point in Fig. 6. The horizontal dashed lines correspond to to the maximum value of the absolute difference | ℎ | or | |.

Figure 15 :

 15 Figure 15: Non-simultaneous inflows case (real hydrographs). Comparison of the model outputs at maximum flooding time. (Top) For the water depth ℎ. (Bottom) For the velocity norm . (Left) The reference solution: 2D SW model with real inflow hydrographs (Fig. 4 right), i.e. ℎ ( ) . (Middle) The difference ℎ or for POD-NN ℎ, -( ). (Right) The difference ℎ or for EA-POD-NN ℎ, - -( ).

Figure 16 :

 16 Figure 16: Non-simultaneous inflows case (real hydrographs). (Top) The histogram of the mean and the StD of spatial values of the absolute difference | ℎ | and | | at maximum flooding time with real inflow hydrographs (Fig. 4 right). (Left) with POD-NN method. (Right) with EA-POD-NN method. (Bottom) The CDF of the absolute absolute difference | ℎ | and | |. (Left) with POD-NN method. (Right) with EA-POD-NN method

Figure 17 :

 17 Figure 17: Non-simultaneous inflows case (real hydrographs). The mean on the absolute difference | ℎ | or | | with respect to the real time with POD-NN and EA-POD-NN methods with real inflow hydrographs (Fig. 4 right). The horizontal dashed lines correspond to to the maximum value of the absolute difference | ℎ | or | |.

Figure 18 :

 18 Figure 18: Non-simultaneous inflows case (triangle hydrographs). Comparison of the model outputs at maximum flooding time. (Top) For the water depth ℎ. (Bottom) For the velocity norm . (Left) The reference solution: 2D SW model with triangle inflow hydrographs for = (11.94, 922, 9.22, 198, 25), i.e. ℎ ( ) . (Middle) The difference ℎ or for POD-NN ℎ, -( ). (Right) The difference ℎ or for EA-POD-NN ℎ, - -( ).

Figure 19 :

 19 Figure 19: Non-simultaneous case (triangle hydrographs). (Top) The histogram of the mean and the StD of spatial values of the absolute difference | ℎ | and | | at maximum flooding time with triangle inflow hydrographs for = (11.94, 922, 9.22, 198, 25). (Left) with POD-NN method. (Right) with EA-POD-NN method. (Bottom) The CDF of the absolute absolute difference | ℎ | and | |. (Left) with POD-NN method. (Right) with EA-POD-NN method

Figure 20 :

 20 Figure 20: Non-simultaneous inflows case (triangle hydrographs). The mean on the absolute difference | ℎ | or | | with respect to the real time with POD-NN and EA-POD-NN methods with triangle inflow hydrographs for = (11.94, 922, 9.22, 198, 25). The horizontal dashed lines correspond to to the maximum value of the absolute difference | ℎ | or | |.

[ 12 ]

 12 Hojat Ghorbanidehno, Jonghyun Lee, Matthew Farthing, Tyler Hesser, Eric F Darve, and Peter K Kitanidis. Deep learning technique for fast inference of large-scale riverine bathymetry. Advances in Water Resources, 147:103715, 2021. [13] Jan S Hesthaven, Gianluigi Rozza, Benjamin Stamm, et al. Certified reduced basis methods for parametrized partial differential equations, volume 590. Springer, 2016.

Table 1 A

 1

	Paper	Math.	model,	ROM method	Input	vari-	Test case		Non-intrusive?
		num. scheme		able/dimension
	Stefanescu et al.,	2D	SWEs	Tensorial POD-	time / 1		Analytical solu-	Intrusive
	IJNMF 2014 [23]	(ℎ, , ), FD	DEIM			tion
	Steinstraesser et	2D	SWEs	POD-DEIM	time / 1		Analytical solu-	Intrusive
	al., JCM 2021 [4]	(ℎ, , ), FV				tion
	Dutta et al., JCP	2D SWEs (ℎ, ),	RBF-POD	time / 1		Kissimmee river	✓
	2021 [9]	FE					(USA)
	El Moçayd et al.,	1D SWEs (ℎ, ),	Generalized poly-	( , ) / 2 or 4	Garonne	river	Intrusive
	EMA 2017 [10]	FV		nomial chaos			(FR)
	Strazzullo et al.,	2D SWEs, FE	Galerkin POD	Affinely		Analytical solu-	Intrusive
	JNMA 2022 [24]				parametrized	tion
					(optimal control
					context) / 3
	Wang et al., JCP	1D Euler equa-	POD-NN	( ; 1 , 2 )	(IC	Shock tube		✓
	2019 [25]	tions, FV			param.) / 3
	Ghorbanidehno	2D	SWEs	POD-NN	Velocity data (In-	Real rivers por-	✓
	et al., AWR 2020	(ℎ, , ), FE		verse pb context)	tions (USA)
	[12]				/ 11

few references of ROM applied to SWEs (or Euler equations). (FV) Finite Volumes, (FE) Finite Elements, (FD) Finite Differences.

Acknowledgments

The application of this computational sciences study has been inspired by scientific discussions led in the MUFFINS ANR project, "MUlti-scale Flood Forecasting with INnovating Solutions", ANR-21-CE04-0021-01. The authors

General design of the numerical experiments

To start with, we outline here the design of the numerical experiments considered to reduce the 2D hydraulic simulation of the 2018 flood event that occurred in the southeastern region of France. This event was characterized by intense rainfall on upstream catchments, river overflow, and extensive flooding over the studied domain Ω. The

Representation of the inflows signals

A necessary condition to apply the reduction method without additional pre-processing is to have an input parameter of quite small dimension: = (10) (see sections 2 and 3). Hence, an effective representation of inflows is performed using a minimal number of key parameters: the maximum peak flow discharge max ( 3 -1 ), and the corresponding time of peak peak (ℎ), defined as peak = ( max -0 ) the duration between 0 the onset of rising water levels and max the corresponding attainment of max , see Fig. 4 (Top right) for the first inflow (red curve). Hence


D S Betsis

Anna A Fedorova

Oleg I Korablev

Jean-Loup Bertaux

Franck Montmessin

Franck Lefèvre

Atmospheric Aerosols Properties via Solar Infrared Occultation Observations by Spicam IR

à la diffusion de documents scientifiques de niveau recherche, publiés ou non, émanant des établissements d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

Data retrieval:

We measured the transmission of solar radiance and then calculated the slant optical depth for each altitude layer. Retrieval of the extinction coefficient profiles was done using "onion peeling" method for 10 reference wavelengths outside the absorbtion bands of atmospheric gases. Radiative properties of the aerosol may be described in Mie theory for the visible-to-infrared opacity ratio and depend on the particle size distribution. For each altitude we assume lognormal size distribution specified by first two moments of this function: effective radius (r eff ) and effective variance (ν eff ) [START_REF] Hansen | Light scattering in planetary atmosphere[END_REF]. With the retrieved size distribution and knowing the vertical profile of extinction coefficient we calculated the aerosol number density. Also we made a decision about water ice cloud presence for each altitude based on extinction and opacity profiles and model predictions.

Results:

The seasonal, latitudinal and altitudinal map of extinction coefficient was plotted. We observe the increase of extinction during the dusty season in the second half of the year. Also the higher values are at low and equatorial latitudes. Some profiles exhibit detached layers, where extinction locally increases, and such layers may be considered as clouds of water ice. Usually there are altitudes between 40 and 60 km. The spectral range of SPICAM IR alone allows to detect particles with sizes in the range 0.2-1.2 μm and to put constrains on particles with r eff > 2 μm. The retrieved effective radius ranges within 0.2-1.3 μm for mineral dust and 0.6-1.5 μm for water ice. Smaller particles were The results of SPICAM solar occultation observations may be useful for planning operations of ExoMars NIR spectrometer which will work at the same range, as well as for other solar occultation spectrometers on Exomars.

ATMOSPHERIC 2 Introduction:

 2 AEROSOLS PROPERTIES VIA SOLAR INFRARED OCCULTATION OBSERVATIONS BY SPICAM IR D.S. Betsis 1 , A.A. Fedorova 1 , O.I. Korablev 1 , J.L. Bertaux 2,1 , F. Montmessin 2 ,F.Lefevre Solar occultation measurements by the SPICAM on Mars-Express in the IR spectral range (wavelengths of 1-1.7 μm) allow to study the vertical distribution of aerosol particles in the Martian atmosphere [1]. The considered dataset includes over 1200 occultations (January 2005 -May 2017, MY 27-34), the profiles for opacity, extinction coefficient, particle size and number density were obtained. It allows characterization of the seasonal evolution and inter-annual comparisons of aerosol particles. The 28 Mars Year (MY) global dust storm is also covered by the dataset.

  observed at high latitudes (70-90°) in both hemispheres. Particles about 1 μm and larger were located commonly at low and middle latitudes. In the period of the MY28 global dust storm the clear increase of r eff value was detected almost for all latitudes and for altitudes 40-80 km. it may be interpreted as lifting of large dust particles (~1 µm) up to 80 km during the Ls≈265-300°, MY 28.

Figure 1 :

 1 Figure 1: Distribution of aerosol extinction with season, altitude and latitude, the average map for 7 Martian years, at 1158.2 nm.

Space Research Institute (IKI), 84/32 Profsoyuznaya, 117997 Moscow, Russia,

CNRS LATMOS, 11 bd d'Alembert, 78280 Guyancourt, France

Acknowledgements: Mars Express is a space mission of ESA. Science operations of SPICAM are funded by CNES and Roscosmos. We thank our collaborators at LATMOS/France for the support of science operations. DSB, AAF, and OIK acknowledge RFBR grant 16-52-16011. JLB acknowledges grant of Ministry for Education and Science of Russian Federation #14.W03.31.0017.