{'fraction_non_alphanumeric': 0.029945334042861982, 'fraction_numerical': 0.013158531275700257, 'mean_word_length': 4.608111774281063, 'pattern_counts': {'":': 0, '<': 0, '<?xml version=': 0, '>': 0, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 1, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
An artificial neural network framework for reduced order modeling of transient flows

This paper proposes a supervised machine learning framework for model order reduction of unsteady fluid flows to provide accurate predictions of nonstationary state variables when the control parameter values vary. Our approach utilizes a training process from full-order scale direct numerical simulation data projected on proper orthogonal decomposition (POD) modes to achieve an artificial neural network (ANN) model with reduced memory requirements. This data-driven ANN framework allows for an implicit representation of the time evolution of the temporal POD modes without performing a Galerkin projection. Our POD framework can thus be considered an equation-free approach to the calculation of the temporal modes in nonlinear transient systems and can be applied to a wide range of physical and engineering applications. Furthermore, we perform a systematic analysis and chart a route to expanding the applicability of ANN architectures in model order reduction of nonlinear systems for solving one-dimensional Burgers and two-dimensional Boussinesq equations.


Introduction
The full resolution of the dynamics of engineering flows requires a computational expense that is as yet infeasible for any practical purpose. Due to these overwhelming computational demands, model order reduction type approaches for solving these problems are gaining popularity. They are extensively used for many optimization and control applications (Akhtar et al., 2015;Buffoni et al., 2006;Dyke et al., 1996;Fang et al., 2009;Fortuna et al., 2012;Freund, 1999;Noack et al., 2011;Lucia et al., 2004;Roychowdhury, 1999;Silveira et al., 1996). Among the multitude of reduced order modeling approaches, proper orthogonal decomposition (POD) has emerged as a popular technique for the study of dynamical systems (Amsallem and Farhat, 2012;Aubry et al., 1988;Banyay et al., 2014). To clarify, POD is also documented under different terminology such as the Karhunen-Loève expansion (Loève, 1955), principal component analysis (Hotelling, 1933) or empirical orthogonal function (Lorenz, 1956). POD extracts the most energetic modes from a collection of high fidelity numerical simulations of the governing equations of the dynamic system being studied. These bases are then used to reduce the degrees of freedom of the governing equation to scales that are computationally tractable. In particular, the fluid mechanics community has traditionally used POD as a method of extraction of the large scale coherent structures (represented by these modes) for the purpose of statistical pattern recognition (Lumley, 1967). This is because the global POD modes may be spanned through a considerably truncated number of bases to resolve attractors and transients well. However, the recent focus of reduced order model development has tilted more towards feedback flow control since a low dimensional model with a reduced degree of freedom allows the use of state of the art control algorithms.
One of the promising approaches of the reduced order modeling (ROM) analysis is based on POD combined with a Galerkin projection (POD-GP) which have been extensively used to provide fast and accurate simulations of large nonlinear systems (Borggaard et al., 2007;Kunisch and Volkwein, 2002;Weller et al., 2010). The Galerkin projection is devised so as to generate a reduced model where the truncated modes are evolved through time as a set of ordinary differential equations (ODEs). This process is denoted a projection as it is obtained by projecting the truncated POD modes onto the governing equations. The spatial fluctuations of the reduced order model are reconstructed entirely from these truncated modes and the orthogonal POD bases. The reduced model is considerably reduced in its degrees of freedom but represents the most significant characteristics of the governing equations from which it is determined. The small memory requirements and low computational expense of this framework makes it particularly suitable for post-processing analysis, optimization and control type problems. This approach is significantly effective for quasi-stationary, timeperiodic and decaying problems but problematic for highly non-stationary, convective and nonlinear problems. These drawbacks are generally due to the lack of adaptation of the POD modes as well as the loss of information due to the truncation being limited to only the most energetic of modes (Cordier et al., 2013;El Majd and Cordier, 2016). Amongst the various efforts used to mitigate these drawbacks include modeling the effect of truncated modes on the retained ones (an approach known as closure modeling) (Wang et al., 2011(Wang et al., , 2012 and developing strategies to construct a more representative basis (San and Borggaard, 2015;Bui-Thanh et al., 2007;Carlberg and Farhat, 2011;Iollo et al., 2000;Kunisch and Volkwein, 2010;Buffoni et al., 2006). Alternatively, machine learning offers promise to generate accurate parametric reduced order models. Recently, a supervised machine learning approach is used to develop feasible regions in the parameter space where the admissible target accuracy is achieved with a predefined reduced order basis (Moosavi et al., 2015).
The limits of the projection-based model reduction has been discussed by Peherstorfer and Willcox (2016) including the lack of access for the full-order model operators or the complexity of the codes which renders obtaining the full-model operators. The present study aims to bypass the Galerkin projection with a completely different approach through the use of an non-intrusive (i.e., without the need for access to full order model operators) reduced order model using a single hidden layer feedforward aryificial neural network (ANN) architecture. ANNs can be classified under the broader category of machine learning methods (i.e. systems that learn from data) and are a mathematical representation of the biological neural networks found in the human central nervous system. They have been used for a wide variety of applications such as function approximation, classification, data processing and dynamic systems control (Widrow et al., 1994). Briefly, an ANN can be used to setup a nonlinear relationship between a desired set of inputs and targets provided a large amount of benchmark data for the underlying relationship is available. This fitting to available data (also known as training) ideally generates a set of linear combinations of parameters and transfer functions which replicate the mean behavior of the underlying phenomena. This allows for the representation of subtle relationships which cannot be expressed explicitly in a functional form. The interested reader is directed to Demuth et al. (2014) which is an excellent introductory text on the development and application of a variety of ANN architectures and their underlying principles.
With regard to the ROM methodology, we are inspired by the works of Narayanan et al. (1999) and Khibnik et al. (2000) where the POD-GP approach was replaced by an ANN based evolution of the temporal modes of the ROM for the purpose of active separation control in a planar diffuser. The dataset for ANN training was obtained both from direct numerical simulations of the Navier-Stokes equations as well as experiments. In Sahan et al. (1997), an ANN was used to devise a feedback control circuit for a transitional flow in a grooved channel. The dataset for ANN training was obtained through a POD-GP simulation itself. Both works demonstrated the unique advantages of the ANN approach, for example, the reduction of a system of governing PDEs to a no-equation dynamic system through data driven learning. We note, in passing, that ANN also finds great utilization in the field of feedback flow control where they are used to generate a direct mapping of flow measurements to actuator control systems (Gillies, 1995(Gillies, , 1998(Gillies, , 2001Faller and Schreck, 1997;Hocevar et al., 2004;Efe et al., 2004Efe et al., , 2005Lee et al., 1997). General machine learning based methods (of which ANN is a subset) are also growing in popularity for fluid flow based applications and represent a computationally viable alternative to the full Navier-Stokes equations (Bright et al., 2013;Gautier et al., 2015;Müller et al., 1999). An excellent source of information about some viable machine learning classes can be found in Ling and Templeton (2015). ANN has also been used within the full Navier-Stokes framework as well for the purpose of quantifying errors in functional closure models (Singh and Duraisamy, 2016) and for data driven closure modeling (Milano and Koumoutsakos, 2002). Before proceeding with our documentation, we describe the major questions we wish to answer in this work: • Can we obtain a reliable ROM in a non-intrusive way without the need for access to full-order modeling operators for nonlinear systems?
• How does an ANN based evolution of the temporal modes of a POD-ROM (POD-ANN) compare to the POD-GP approach?
• Can a trained POD-ANN (using a particular dataset) be used to interpolate or extrapolate to flow problems with slight differences in temporal evolution and physics (i.e. for different parameter values such as Reynolds and Rayleigh numbers)?
• Can the ANN-POD approach be used for complicated dynamical systems such as nonlinear Burgers or coupled Boussinesq systems with complex dynamics?
To answer these questions, we investigate two test cases given by the nonlinear viscous Burgers equation (i.e. a nonlinear wave propagation or NWP problem) and a differentially heated cavity problem (DHC) with buoyancy driven fluid flow to compare both the POD-ANN and POD-GP approaches.

Proper orthogonal decomposition
In this section we explore the POD approach to reduced order modeling to capture unsteady, convective and non-periodic dynamics of the underlying governing partial differential equations. We elaborate the procedure for calculating the orthogonal bases and temporal modes for the NWP problem and note that an extension to higher dimensionality is straightforward. For a detailed discussion on this technique we refer the reader to San and Iliescu (2014).
A POD can be constructed from the scalar field variable f at different times (also known as snapshots). These snapshots are obtained by solving the governing equations we are attempting to model using a regular DNS approach. The numerical methods that are used to obtain this field data are explained in Section 6. In the following, we will utilize the index i to indicate a particular snapshot in time. For the POD approach we utilize a total of N snapshots for the field variable i.e. f i (x) for i = 1, 2, ..., N . The flow field data thus obtained is first decomposed into the mean and fluctuating part as follows wheref implies a temporal averaging of a particular point value in the field andf contains the fluctuating quantity. We must remark here thatf is a function of space alone whereaŝ f is a function of both space and time. A correlation matrix may be constructed using the fluctuating components of the snapshots to give where Ω is the entire spatial domain and i and j refer to the i th and j th snapshots. The correlation matrix C is a non-negative symmetric square matrix of side N . If we define the inner product of any two fields f 1 and f 2 as we may express the correlation matrix as C ij = (f i ,f j ). In this study, we use the wellknown Simpson's 3/8 integration rule for a numerical computation of the inner products.
The optimal POD basis functions are obtained on performing an eigendecomposition for the C matrix. This has been shown in detail in the POD literature (see, e.g., Sirovich (1987); Holmes et al. (1998); Ravindran (2000)). The eigenvalue problem can be written in the following form: where Λ = diag[λ 1 , λ 2 , ..., λ N ] is a diagonal matrix containing the eigenvalues of this decomposition and W =[ w 1 , w 2 , ...,w N ], λ j is an orthogonal basis consisting of the eigenvectors of this decomposition. The eigenvalues are stored in descending order, λ 1 ≥ λ 2 ≥ ... ≥ λ N .
The POD basis functions can be written as where w j i is the ith component of eigenvector w j . Therefore we emphasize that the POD modes are ranked according to the magnitude of their eigenvalue. The eigenvectors must also be normalized in order to satisfy the condition of orthonormality between bases: It can be shown that, for Eq. (6) to be true, the eigenvector w j must satisfy the following equation: In practice, most of the subroutines for solving the eigensystem given in Eq. (4) return the eigenvector matrix W having all the eigenvectors normalized to unity. In that case, the orthogonal POD bases are given by where φ j (x) is the jth POD basis function. We remark here that the POD should be carried out independently for each field variable for multidimensional governing equations.

Burgers equation
The viscous Burgers equation is used in our study as a test case to explain the proposed model reduction framework. In its original form, the Burgers equation can be written as where Re is a non dimensional Reynolds type number. We must mention here that this equation is generally considered a framework for preliminary evaluation of numerical methods with respect to fluid flow applications as it possess the hallmarks of general multidimensional advection-diffusion problem. We consider a challenging convective system in our model order reduction assessments since it characterizes the structure of localized flow structures such as shock waves. The above PDE can be solved exactly to obtain an analytical formulation for time evolution for the field variable u(x, t) given by (Maleewong and Sirisup, 2011) u(x, t) = x t+1 where t 0 = exp( Re 8 ). This exact expression is used to generate snapshot data for our POD-ROM analysis and representative space-time solutions are shown in Fig. 3.

Boussinesq equations
To further test the proposed reduced order modeling framework, we consider the dimensionless form of the two-dimensional incompressible Boussinesq equations on a domain where u = [u, v] T is the velocity with horizontal and vertical components represented by u and v, respectively, p is the pressure, and θ is the temperature. The Boussinesq equations are generally characterized by three dimensionless numbers. Re is the Reynolds number, the ratio of viscous effects to the inertial effects; P r is the Prandtl number, the ratio of the kinematic viscosity to the heat conductivity; and Ri is the Richardson number, the ratio of the buoyancy force to the inertial forces. Depending on the physics of the system, for example in natural convection heat transfer, one can introduce other relevant dimensionless numbers, such as the Rayleigh number Ra = RiRe 2 P r, the Grashof number Gr = RiRe 2 , and the Peclet number P e = ReP r.
Instead of the primitive variable formulation given by (11)-(13), we introduce the vorticity ω = (∇ × u) · k and the stream function ψ. The governing equations for two-dimensional incompressible Boussinesq equations can be written as the two coupled scalar transport along with the kinematic relationship between vorticity and stream function given by where ∆ refers to the Laplacian operator. The desired flow velocity components can be recovered using

Galerkin projection methodology
The POD basis functions account for the essential dynamics of the underlying governing equations. After empirical POD basis functions are obtained, a set of nonlinear ordinary differential equations (ODEs) can be derived by using a Galerkin projection. To build Galerkin projection based ROMs, we truncate the system by considering the first R largest POD basis functions with R N . These POD modes correspond to the R largest eigenvalues, λ 1 , λ 2 , ..., λ R . Using these first R largest POD basis functions the field variables are then approximated. construction of a POD-ROM is evident when it is observed that the 5 most energetic modes retain 97.3% of the energy of the system and the first 10 modes contribute to 99.7% of the total energy of the system. Some of the bases generated by POD method for this problem are shown in Fig. 5. We note that these bases are all orthogonal to each other. For the first test case we study, the Galerkin projection can be carried out in the following manner.

Galerkin projection for Burgers equation
where a k are the time dependent coefficients, and φ k are the space dependent modes. To derive the POD-ROM, we first rewrite the Burgers equation (i.e., Eq. (9)) in the following where L[f ] = ∂ 2 f ∂x 2 is the linear operator and N [f ; g] = −f ∂g ∂x is the nonlinear operator. We emphasize that the procedure can be easily extended to the general nonlinear partial differential equations such as those of the Boussinesq equations in the next subsection. By applying this projection to our nonlinear system (i.e., multiplying Eq. (19) with the basis functions and integrating over the domain), we obtain the Galerkin POD-ROM, denoted POD-GP: Substituting Eq. (1), Eq. (8), and Eq. (8) into Eq. (18), and simplifying the resulting equation by using the condition of orthogonality given in Eq. (6), the POD-GP implementation can be written as follows: where The POD-GP given by Eq. (21) consists of R coupled ODEs and can be solved by a standard numerical method (such as the third-order Runge-Kutta scheme that was used in this study). The number of degrees of freedom of the system is now significantly lower.
The vectors, matrices and tensors in Eqs. (22)-(26) are also precomputed quantities, which results in a dynamical system that can be solved very efficiently.
To complete the dynamical system given by Eq. (21), the initial condition is given by using the following projection: where u(x, t = 0) is the physical initial condition of the problem given in Eq. (10).

Galerkin projection for Boussinesq equations
The POD procedure is carried out independently for the temperature and the vorticity Then we can represent the scalar field variables as where we now need to determine dynamical equations for the time dependent coefficients {α k } and {β k }. To obtain these dynamical systems for the coefficients, we rearrange (14)-(15) using linear and nonlinear operators in the following where L and D are linear operators, and N is the nonlinear operator given by where f and g are arbitrary field variables.
Next, we apply Galerkin projection to the system by multiplying (30)-(31) with the basis functions and integrate over the domain. The Galerkin projection onto the kth basis function can be written as Substituting (28)-(29) into (33)-(34), and using the orthogonality condition given by (6), we obtain the following coupled POD reduced-order system for k = 1, 2, . . . , R: where The proper orthogonal decomposition reduced-order model (POD-ROM) given by (35)-(36) consists of 2R coupled ordinary differential equations and can be solved easily by a standard numerical method (e.g., we use a third-order Runge-Kutta scheme in this study). Once again, we emphasize that the degrees of freedom of the system are considerably reduced and this provides for a fast solution of the system. A point of clarification here is that using the kinematic relationship between stream function and vorticity given by Eq. (16), the jth basis function for the stream function, φ ψ j (x, y), can be obtained from the jth vorticity basis function by solving a Poisson equation

Artificial neural networks
The basic structure of the simple feed-forward artificial neural network consists of L layers with each layer possessing a predefined number of unit cells called neurons. Each of these layers has an associated transfer function and each unit cell has an associated bias.
Any input to the neuron has a bias added to it followed by activation through the transfer function. To describe this process using equations, we have for a single neuron in the l th layer receiving a set of inputs s l i from the (l − 1) th layer given by Demuth et al. (2014) where W l ij stands for a matrix of weights linking the l − 1 and l layers with X l−1 j being the output of the (l − 1) th layer. The output of the l th layer is now given by where b l i is the biasing parameter for the i th neuron in the l th layer. Every node (or unit cell) has an associated transfer function which acts on its input and bias to produce an output which is 'fed forward' the network. The nodes which take the raw input value of our training data set (i.e. the nodes of the first layer in the network) perform no computation (i.e. they do not have any biasing or activation through a transfer function). The next layers are a series of unit cells which have an associated bias and activation function which perform computation on their inputs. These are called the hidden layers with the individual unit cells known as neurons. Note that it is common in literature to consider these the only layers in the network. The final layer in the network is that of the outputs. The output layers generally have a linear activation function with a bias which implies a simple summation of inputs incident to a unit cell with its associated bias. In this investigation, we have used one hidden layer of neurons between the set of inputs and targets with a tan-sigmoid activation function. The tan-sigmoid function can be expressed as The transfer function F calculates the neuron's output given its net input. In theory, any differentiable function can qualify as an activation function (Zhang et al., 1998), however, only a small number of functions which are bounded, monotonically increasing and differentiable are used for this purpose. The choice of using ANN is motivated by its excellent performance as a forecasting tool (Dawson and Wilby, 1998;Kim and Valdés, 2003) and its general suitability in the machine learning and function estimation domain (e.g. see Haykin et al. (2009) and references therein).

ANN training
The training of a desired ANN is carried out by minimizing the error between the target and the inputs to determine a set of best fit parameters. These best fit parameters are the biases and linear weights that have captured the underlying relationship between the targets and inputs and may now be used to predict target data for inputs a posteriori. The advantage of using the ANN approach over traditional statistical regression models is that comparatively smaller data sets for training are suitable. We shall describe the individual features of the ANN architecture relevant to each test case in the sections that follow.
The training of an ANN architecture may be undertaken by a variety of optimization algorithms. For our purpose, we have employed the Levenberg-Marquardt minimization algorithm (Levenberg, 1944) which is a universally popular approach to ANN training, which is available in MATLAB's neural network toolbox. We remark that other training methods may also be used such as Bayesian regularization (Gençay and Qi, 2001) (which is particularly useful for noisy data) and biomimetic algorithms such as leapfrogging (Rhinehart et al., 2012). For the purpose of tracking the training of the ANN and evaluating its performance, we utilize three different subsets of our overall data with 70% being utilized for training, and 15% being utilized for validation and testing each. The training data set is used to adjust the weights of the ANN. The validation data set is used to minimize overfitting i.e. it is used to verify whether an increase in accuracy over the training data set also yields an increased accuracy over a data set that is not revealed to the network a priori. The validation data set is used for testing the final solution of the neural network.
Within the POD context, the trained ANN can be used as a 'black box' which takes the present temporal modes as the input along with a dimensionless parameter and the present time and uses the memory of its training to determine the evolution of the temporal mode.
Individual schematics for the two test cases are provided in their problem descriptions.

ANN architecture for model order reduction
In the present study, we devise a framework for developing an artificial neural network (ANN) in order to create a very fast (equation-free) integrator/solver for certain types of

PDE systems (examples are given for 1D Burgers equation for NWP and 2D Boussinesq
equations for DHC problems). The ANN is created by getting the POD basis, but choosing the weights for the inner neurons so that they input basis coefficients and output basis coefficient, to evolve from time n to time n + 1. Simulation data is used to train the ANN, i.e. make it choose the weights to best fit the data. This idea can be considered as a competitor for Galerkin-POD, and the advantages are that this is an equation-free method, and also it can be more accurate since it can be trained with more data.
The NWP test case utilizes an ANN architecture which is given in Fig. 1. In this architecture, several inputs are provided out of which the majority are temporal modes at the current timestep in the POD-ROM evolution and the remaining two are Reynolds number and time. The training data is also arranged in the same manner. Our outputs are given by just the temporal modes required for a field reconstruction in the next step. The figure referred to above shows the number of inputs (12 or 7 with 10 or 5 modes, Reynolds number and time) and the number of outputs (ten or four predicted modes at the next timestep). The middle layer l contains the neurons with the tan-sigmoid activation function and these are varied in our investigation. The next step is to test this system a posteriori.
We state for the sake of clarity that our training data utilized snapshots of simulations using Reynolds numbers between 100 and 1000 at increments of 100 each. Through this trained ANN we aim to interpolate a posteriori analyses for Reynolds numbers that are within the bounds of the training data set extremes and extrapolate for Reynolds numbers that are outside these bounds. Fig. 2   . .

Re
t n a n 1 a n 2 a n 3 a n

Numerical methods
The objective of the present work is to test and evaluate the characteristics of proposed ANN and GP frameworks and compare them with a solution obtained by a direct numerical simulation (DNS). Here, we provide a brief description of the numerical methods employed in this study.

Compact difference schemes
First of all, we implement compact difference schemes to approximate the differential operators (N , L, and D). Compact difference schemes have shown their ability to reach the objectives of high-accuracy and low computational cost for many problems in fluid dynamics (Lele, 1992).
In a compact difference scheme, the first order derivatives can be computed accordingly (Lele, 1992) which gives rise to an α-family of tridiagonal schemes with a = 2 3 (α + 2), and b = 1 3 (4α − 1). The subscript i represents the spatial grid index in the considered direction, and h is the  uniform grid spacing in that direction. Here, α = 0 leads to the explicit non-compact fourthorder scheme for first derivative computations. A classical compact fourth-order scheme, which is also known as Padé scheme, is obtained by setting α = 1/4, resulting in b = 0 with a compact dependence. The second derivative compact scheme is given by where a = 4 3 (1 − α), and b = 1 3 (10α − 1). For α = 1/10 the classical fourth-order Padé scheme is obtained. The high-order one-sided derivative formulas are used for the Dirichlet or Neumann boundary conditions to complete a tridiagonal system of equations for both the first and second order derivative formulas (Carpenter et al., 1993). In the vorticity-steam function formulation, Briley's formula is usually used to satisfy the no-slip boundary condition at the walls. In this formula the vorticity values at the boundary are computed from the stream function values according to the following third-order accurate formula where the subscript 0 represents the discrete point index on the no-slip boundary where we set ψ 0 = 0.

Time advancement scheme
Semi-discrete ordinary differential equations are obtained after a spatial discretization of the partial differential equations. To implement the Runge-Kutta scheme for the time integration, we cast the model equations in the following form where £(u) is the discrete operator of spatial derivatives including nonlinear convective terms and diffusive terms. It should be also noticed that the POD-ROM equations are also in this forms. We assume that the numerical approximation for time level n is known, and we seek the numerical approximation for time level n + 1, after the time step ∆t. The optimal third-order accurate total variation diminishing Runge-Kutta (TVDRK3) scheme is then given as (Gottlieb and Shu, 1998) u (1) = u n + ∆t£(u n ),

Numerical integration
To perform the inner products given by (3), we compute the integral of u(x, y) over the domain Ω by using the dual integration method with the Simpson's 3/8 rule (Hoffman and Frankel, 2001) Ω u(x, y)dxdy = ∆y 3 and N x and N y are even numbers representing the total number of grid points in the x and y directions, respectively.

Fast Poisson solver
The elliptic equations given in (16) and (45) can be written in the form ∇ 2 u = f . The compact fourth-order discretization scheme with a nine point stencil can be written as (Wang and Zhang, 2009) where the coefficients are a = −10(1+γ 2 ), b = 5−γ 2 , c = 5γ 2 −1, d = (1+γ 2 )/2, e = ∆x 2 /2, and γ is defined as the mesh aspect ratio, γ = ∆x/∆y. For γ = 1 the scheme is well known and is sometimes called the Mehrstellen scheme. The no-slip boundary condition implies an impermeability condition for the stream function. Therefore, the prescribed boundary condition for the stream function ψ| Γ = 0 results in the homogeneous Dirichlet boundary condition (i.e., u| Γ = 0) on the boundary Γ and suggests the use of a fast sine transform.
The procedure involves three steps. First, a forward sine transform for the source term is given byf where k and l are wavenumbers in Fourier space. Here, the symbol hat is used to represent the corresponding Fourier coefficient of the physical grid valued data with a subscript pair i, j, where i = 0, 1, 2, ..., N x , and j = 0, 1, 2, ..., N y . As a second step, we directly solve the subproblem in Fourier space using the following relationship u k,l = ef k,l 8 + 2 cos πk Nx + 2 cos πl Ny a + 2b cos πk Nx + 2c cos πl Ny + 4d cos πk Nx cos πl and finally, the physical values for u are found by performing an inverse sine transform as

Results
In the following we provide temporal mode evolution data for our test cases and make assessments of the relative performance of the POD-ANN and POD-GP approaches. Although measuring the temporal mode error is not equivalent to measuring the error in any standard function space norm or quantity of interest, we found that it simplifies our demonstration and presents an efficient way in our assessments between the proposed approach and standard POD-Galerkin model.

Burgers equation: nonlinear wave propagation (NWP) problem
In this subsection we shall elaborate about the performance of the POD-ANN and POD-GP approaches for the VBE test case using an both 5 and 10 temporal modes for the reduced order modeling approach. Fig. 3 shows the space-time behavior of the two representative VBE solutions we shall use for ROM with the low Reynolds number case showing a higher dissipative behavior as the solution advects from its initial condition. Figs. 4 and 5 show the eigenvalues and the POD basis functions obtained for this test case. It is evident that a small number of modes are enough to capture the spatial and temporal dynamics of the system well. We perform our analysis with two configurations using 5 modes and 10 modes each.  for a greater duration as compared to the POD-GP approach. At these low number of modes being preserved, the performance of the ANN-POD is remarkable. Indeed for the first mode, the ANN-POD approach preserves the amplitude and phase behavior very well.
Higher modes tend to show a slight variation in amplitude when compared to DNS data but this variation is quite negligible when compared along the POD-GP approach which shows extensive loss of phase and amplitude at higher modes. It is clear, that at a lower preservation of modes, the POD-ANN approach is superior. Fig. 7 shows a similar study where instead of truncating modes beyond the first 5, we have retained 5 more modes instead. One can immediately notice a significant improvement in the POD-GP approach. This is expected because the POD-GP approach converges to the underlying governing equations when all the modes are considered. The POD-ANN approach still performs well in terms of capturing the amplitude of the first five modes although slight variations in the phase for the highest modes are seen as shown in the Re = 1250 case. When examining modes 6-10 however (in Fig. 8), we can see that the POD-GP approach shows a considerable deviation from the truth as compared to the POD-ANN approach which does a decent job at capturing amplitude and phase behavior with a lower error. This examination also points to a more accurate behavior shown by the POD-ANN approach.

Boussinesq equation: differentially heated cavity (DHC) problem
The DHC test case is used to make assessments about the scalability of the POD-ANN approach for higher dimensional governing laws. Fig. 9 shows instantaneous temperature plots for the different Rayleigh numbers of 3.4 × 10 5 , 4.4 × 10 5 and 5.4 × 10 5 . Similar to the Reynolds number in the NWP problem, the Rayleigh number behaves as a control parameter and increasing values leads to finer scale structures in the flow field as seen. Fig. 10 shows the evolution of the temperature and the Nusselt number N u with time in the simulation.
The lower two values of the Rayleigh number show a periodic behavior of N u whereas at Ra = 5.4×10 5 we can observe the emergence of a quasiperiodic behavior of the system. This higher Rayleigh number thus represents a tougher challenge for our reduced order modeling approach as more modes may be needed to represent the dynamics of the system. We note that the following equation has been used for calculating the Nusselt number along the vertical wall A POD is performed on these three DNS simulations for the different Rayleigh numbers. Fig. 11 shows the eigenvalues of the DHC test case for the three Rayleigh numbers for both Rayleigh numbers corresponding to turbulence. Fig. 12 shows the first four basis functions obtained through the POD analysis for both the vorticity and temperature variables and one can notice similar structures in the basis functions for both variables. We remark that only 4 modes are chosen for the POD-ROM as they are sufficient to capture a high fraction of the energy in the system as shown in Table 1. Fig. 13 shows the comparison of the first temperature and vorticity temporal modes for Ra = 3.4 × 10 5 . It can be seen that the POD-ANN approach does an excellent job in capturing the temporal dynamics of this mode with a good agreement of amplitude and phase for both field variables. The POD-GP approach is seen to have an increased amplitude and phase for the vorticity temporal mode while it exhibits a phase difference for the temperature mode. POD-ANN however captures the vorticity mode exceptionally well and shows the slightest amplitude over-estimation for the temperature mode. These errors are far more smaller than the POD-GP approach. A similar trend is seen in Fig. 14 where the POD-ANN approach captures both field variable modes well in comparison to POD-GP. It also seen from Figs. 15 and 16 that the POD-GP shows increasing errors in capturing the behavior of the higher modes in comparison to the POD-ANN method. We analyse their comparative performances in more detail using phase diagrams as shown in showing more adherence to the DNS solution for the a 3 − a 1 and b 3 − b 1 phase diagrams.
On extending the same analysis to the higher Rayleigh number of 4.4 × 10 5 . It is evident quite immediately (from Fig. 18) that the higher Rayleigh number impedes the performance of the POD-GP approach where the modes are seen to increase in their amplitude and go well out of phase in comparison to POD-ANN which captures the dynamics perfectly. The next three modes also exhibit the same behavior with POD-GP showing even more errors at the highest modes. POD-ANN however shows a negligible error for these cases (Figs. 19-21 ).
The phase diagram for this Rayleigh number shown in Fig can also be observed in the modes of the DNS solution and the POD-ANN performance is even more striking. Only, for the fourth mode given in Fig. 26, can we see discernible errors from the POD-ANN approach. However, these errors are minimal in comparison to POD-GP which bears no resemblance to the true solution. POD-ANN is still able to capture phase and amplitude trends extremely well. The phase diagrams for this Rayleigh number ( Fig. 27) show how much better POD-ANN performs than POD-GP for this approach with a very close agreement to the DNS for the former. We can conclude that the POD-ANN methodology performs in a superior manner to POD-GP.

Conclusion
In this study, we have investigated the predictive performance of two different PODenabled model order reduction approaches for the time evolution of the temporal modes.
First, we look at the classical Galerkin projection approach where a system of ODEs is developed for the evolution of each temporal mode (corresponding to a scalar field variable).
Second, we have developed an equation-free artificial neural network (ANN) framework through data driven learning which predicts the evolution of the modes without the explicit time stepping approach used in the Galerkin projection. After the mathematical development of the ROM framework, the Galerkin projection approach and a discussion about the ANN implementation, extensive comparisons are made for two benchmark test cases which are the viscous Burgers equation and a differentially heat cavity problem in a Boussinesq flow setting. We notice that the POD-ANN approach outperforms POD-GP consistently across the different testing environments and proves excellent in both interpolating and extrapolating its learning from obtained data sets for simulations with different control parameters and underlying physics. We revisit the questions we posed in the introduction to conclude our investigation: • We can conclude with confidence that the POD-ANN approach is superior to the POD-GP approach for the Burgers equation and differentially heated cavity test cases especially for highly truncated systems (i.e. systems represented by a few modes only).
For the latter, the presence of more turbulent flow fields causes POD-GP to breakdown in comparison to POD-ANN which performs excellently.
• We can also state that POD-ANN is a viable tool for extrapolation and interpolation beyond the data sets used to train its learning. While it is expected to provide good results within the bounds of the training data set, the Burgers case detailed its capacity to predict the temporal mode evolution of different Reynolds numbers in an accurate manner with much greater accuracy than POD-GP.
• The differentially heated cavity was tested to ensure the scalability of the ANN approach and the use of two neural networks for each scalar field variable proved suitable to an excellent prediction of the temporal modes.
Furthermore, we present the following observations which were made during this study: • The POD-ANN approach provide a reliable model reduction paradigm without the need for access to the full-order model operators.
• Once trained, the POD-ANN can be considered to be an 'equation-free' ROM architecture.
• While increasing the retained modes in the POD gives us much better approximations using the POD-GP approach, POD-ANN is better suited to more radically truncated systems (an ANN with higher degrees of freedom may show statistical errors on time evolution).
• POD-ANN represents a viable framework for modeling highly nonlinear systems when data is available incrementally and constant updates are required. It is also very well suited to situation where model interpretability is not a major concern.
• While our training algorithm utilized the standard Levenberg-Marquardt algorithm, it must be noted that there exists many training methodologies which may better preserve the underlying statistical relationships of the inputs and outputs. The effect of training approaches is a future direction of study for this work.
We believe that the integration of machine learning tools such as ANN into the mainstream reduced order modeling ideology would represent a significant increase in its viability for both theoretical and application oriented goals. We would like to conclude that POD-ANN would be a viable ROM procedure for convective flows where the response of complex system dynamics is desired to be represented accurately by a few modes only.  DNS amplitudes of temporal modes provided for comparison.