key: cord-0809617-twlt1t80
authors: Barros, Gabriel F.; Grave, Malú; Viguerie, Alex; Reali, Alessandro; Coutinho, Alvaro L. G. A.
title: Dynamic mode decomposition in adaptive mesh refinement and coarsening simulations
date: 2021-08-02
journal: Eng Comput
DOI: 10.1007/s00366-021-01485-6
sha: 8551c16fb47cd212593d3d61fd5831917b7a1b41
doc_id: 809617
cord_uid: twlt1t80

Dynamic mode decomposition (DMD) is a powerful data-driven method used to extract spatio-temporal coherent structures that dictate a given dynamical system. The method consists of stacking collected temporal snapshots into a matrix and mapping the nonlinear dynamics using a linear operator. The classical procedure considers that snapshots possess the same dimensionality for all the observable data. However, this often does not occur in numerical simulations with adaptive mesh refinement/coarsening schemes (AMR/C). This paper proposes a strategy to enable DMD to extract features from observations with different mesh topologies and dimensions, such as those found in AMR/C simulations. For this purpose, the adaptive snapshots are projected onto the same reference function space, enabling the use of snapshot-based methods such as DMD. The present strategy is applied to challenging AMR/C simulations: a continuous diffusion–reaction epidemiological model for COVID-19, a density-driven gravity current simulation, and a bubble rising problem. We also evaluate the DMD efficiency to reconstruct the dynamics and some relevant quantities of interest. In particular, for the SEIRD model and the bubble rising problem, we evaluate DMD’s ability to extrapolate in time (short-time future estimates).

Data-driven methods are currently revolutionizing the modeling, prediction, and control of complex systems. Increasingly, researchers are considering data-driven approaches for a diverse range of complex systems, such as turbulent flows, climate sciences, epidemiology, finance, robotics, and many other different applications [12] . Even with the availability of better hardware and advances in techniques and algorithms, numerical simulations of these systems are still resource-demanding: strong nonlinearities, multiple scales, and large dimensionalities are typical examples of complexities found in modern applications. With the assembly of modern mathematical methods, unprecedented data availability, and increasing computational resources, previously complex, challenging problems can now be tackled within the new research field entitled scientific machine learning (SciML).

SciML is a core component of artificial intelligence and computational technology that can be trained, with scientific data, to augment or automate human skills [6] . This emerging research area aims at the opportunities and challenges in the context of complex applications across science and engineering, and other interdisciplinary fields. A wide range of SciML methods can be categorized regarding the type of information available, and their intended use [13] . For instance, many state-of-the-art methods and applications within the SciML scope can be found in [13, 36, 45, 46, 52, 60] . In this study, we focus on Dynamic Mode Decomposition (DMD), an unsupervised SciML method that can extract the most dynamically relevant low-rank structures from large-dimensional data observed in dynamical systems. DMD can be applied to both numerical [56] and experimental data [57] .

The standard DMD procedure for numerical simulations consists of stacking the m snapshots (discrete solutions in space for a given time step) in columns to create a matrix and map the dynamics using a linear operator. The DMD procedure assumes that the m collected snapshots have spatial dimensionality n, where n ≫ m , such that the snapshots matrix has dimension n × m . This occurs when snapshots are obtained from experimental data with sensors in fixed positions or numerical simulations considering a fixed mesh. However, in many situations, this is not always achievable. For numerical simulations using Adaptive Mesh Refinement/ Coarsening (AMR/C), for instance, spatial adaptivity leads to solutions computed in meshes that constantly change in time. Adaptive meshes lead to a different number of nodes, nodal coordinates and numbering, and mesh topologies. Although there are strategies regarding missing spatial data or subsampling for the DMD [14, 41] , specific treatments related to AMR/C approaches in the finite element context are not explored in the literature. In the present study, we develop a strategy to project all the snapshots of a given simulation with different dimensionality onto a reference target mesh with minor accuracy loss, enabling the use of any SVD-based data-driven technique such as DMD.

This paper is structured as follows: Sect. 2 describes the relation between the discretization of PDEs in space and time and a dynamical system. This section introduces Dynamic Mode Decomposition, our method of choice for short-time future estimates and extrapolation. Section 3 describes our strategy to deal with simulations that consider AMR/C in their evolution (e.g., in the case of dimensionality of the output vector, as well as of mesh topology and/or node numbering, varying in time). In Sect. 4 we describe the numerical applications in this study: the use of DMD on a continuous SEIRD model for COVID-19 and two fluid dynamics problems, a density-driven gravity current, and a bubble rising problem. We show efficiency and accuracy results for the signal reconstruction. Moreover, future time step predictions using DMD are evaluated for the SEIRD model and the bubble rising problem. In Sect. 5, we draw our final remarks and conclusions.

Solving partial differential equations (PDEs) using fast, accurate, reliable, and robust methods is crucial for many industrial and scientific applications. Several methods (such as finite elements, finite differences, and many others) are responsible for approximating the infinite-dimensional PDEs into finite-dimensional spaces. The discretization of these equations allows the process to be automated. In the present study, we focus on using the finite element method for spatial discretization of the PDEs. That is, consider a generic transient parametric PDE such as where N is a nonlinear operator, is a vector of parameters (e.g., diffusion, density, viscosity, etc.), = (t; ) is a given function and the solution is a function of spatial coordinates , temporal coordinates t, and parameters such that = ( , t; ) . The equation is equipped with boundary and initial conditions where and are the Dirichlet and Neumann boundary conditions, respectively, T is the final time, and 0 is the initial condition for . The domain Ω ⊂ ℝ nsd is bounded by Lipschitz continuous boundaries Γ N ∪ Γ D = Γ ⊂ ℝ nsd−1 , and is the unit outward normal to Γ N . The union of boundaries and domain is represented as Ω . The standard finite element method consists of discretizing Ω into a mesh composed of nodes and elements. Each element has its domain Ω e ⊂ ℝ nsd and boundary Γ e ⊂ ℝ nsd−1 . The weak form of the system can be obtained by integrating Eq. (1) in its strong form against a weighting function ∈ H 1 (Ω) , where H 1 (Ω) is the Sobolev space of the square-integrable functions with an integrable first weak derivative, and applying the divergence theorem. Being P k (Ω e ) the space of polynomials of degree equal or less than k over Ω e , the function spaces are defined as Therefore, the semi-discrete finite element formulation for Eq. (1) is: find h (t) ∈ S h t such that ∀ h ∈ W h :

= on Γ D × (0, T],

where the L 2 inner product over the domain Ω is indicated by (⋅, ⋅) . The weak form given by Eq. (5) naturally accommodates several finite element formulations, from Galerkin to Variational Multiscale methods [1, 8, 19, 39, 53] . After the temporal discretization of Eq. (5), the equation can be translated into a discrete-time dynamic system. In this system, the state vector h at the time instant k + 1 can be written such that:

where F represents the discrete-time flow map of the system and incorporates information regarding the parameters , mesh size, solver tolerances, etc. In the present study, we consider that the measurements of the system are the state vectors themselves, that is, h k . Analyzing the evolution in time of a discretized PDE as a dynamical system is a key concept for introducing Dynamic Mode Decomposition (DMD).

DMD is an equation-free, data-driven method that provides accurate assessments of the dominant structures in a given complex system [43] . DMD provides a decomposition of data into spatio-temporal modes that correlates the data across spatial features and also associates them to unique temporal Fourier modes. The main idea of the method is to efficiently compute the regression of linear/nonlinear terms to a least-square linear dynamics approximation from experimental or numerical observable data. Despite its first appearance in the fluid dynamics context [55, 56] , DMD has been used in many other applications such as epidemiology [51] , biomechanics [16] , urban mobility [3] , climate [44] and aeroelasticity [28] , especially in structure extraction from data and control-oriented methods.

We can now apply DMD on the dynamical system described in Eq. (6) . Consider a dataset h containing the observations in time of the dynamical system h k for k = 0, 1, … , m , where m + 1 is the total number of observations. The dataset can be split into two datasets h is the Moore-Penrose pseudoinverse of h 1 . However, we avoid the computation of the full matrix since is a n × n matrix. Also, the computation of the full Moore-Penrose

pseudoinverse is not advisable due to its ill-conditioning. Instead, we can compute the SVD of h 1 as where ∈ ℝ n×m and ∈ ℝ m×m are the left and right singular vectors and ∈ ℝ m×m is a diagonal matrix with real, non-negative, and decreasing entries named singular values.

The singular values 0 ≥ 1 ≥ 2 ≥ ⋯ ≥ m−1 are hierarchical and can be interpreted in terms of how much the singular vectors influence the original matrix h 1 . For the DMD procedure, considering the Eckart-Young Theorem [23] , the optimal low-rank update approximation matrix h 1 , when subjected to a truncation rank r, can be written as where r ∈ ℝ n×r is a matrix containing the first r columns of , r ∈ ℝ m×r contains the first r columns of , and r ∈ ℝ r×r is the diagonal matrix containing the first r singular values. The pseudoinverse can be approximated as and, instead of computing ∈ ℝ n×n , we can obtain ̃ , a r × r projection of as Note that ̃ is unitarily similar to . Further mathematical details regarding the optimization problem (the best-fitting matrix ̃ ) and the influence of the Eckart-Young Theorem on constraints of the problem can be found in [40] . Now we can compute the eigendecomposition of ̃ :

where is a diagonal matrix containing the discrete eigenvalues j and the matrix contains the eigenvectors j of ̃ . The DMD basis can be written as and the signal reconstruction as being the vector containing the projected initial conditions such that = † h 0 , and eig is a diagonal matrix whose entries are the continuous eigenvalues i = ln( i )∕Δt o , where Δt o is the time step size between the outputs. In the present study, Δt o = jΔt , where j = 1, 2, … , m , being Δt the time step size used in the temporal integration of the PDEs. For instance, if one chooses to output the solutions once every two time steps, the time step size between the two observations will be two times larger than the time step size

used to compute time integration, that is, Δt o = 2Δt . It is important to mention that the snapshots sampling frequency affects directly the DMD's ability to capture the dynamics. For lower dominant frequencies, a larger Δt o is more adequate, while smaller values of Δt o are required for capturing rapid dynamics [56] . The form of (15) can be regarded as a generalization of the Sturm-Liouville expansion for a differential problem:

where i and i are the ith Sturm-Liouville eigenfunctions and eigenvalues for a given differential operator.

DMD can be seen as a dimensionality reduction method due to its inherent ability to extract the r most relevant dynamical modes, where r is often much smaller than the snapshot matrix rank m. However, a strategy to determine the number of relevant modes is not straightforward, and is an active topic in DMD research [59] . Even though the choice of r for DMD may require some trial and error, some techniques can be used to find a good starting point. A hard threshold technique [43] consists of choosing r such that where is a tolerance threshold, set, e.g., to 10 −6 . This method implies that more than 100(1 − )% of the variance in the data is retained by the approximation. For the case where DMD is used on experimental (or numerical but noisy) data, more sophisticated solutions are presented in the literature [22, 29] .

Another important consideration when using DMD is the SVD algorithm. The SVD can represent a significant part of the computational effort, meaning that improvements in the SVD performance lead to significant CPU time gains. For many SVD-based methods (such as DMD), there is no need to compute the SVD for the whole matrix, since the method aims at extracting the first r dominant structures in the matrix. Many algorithms are designed in this direction to make this computation more efficient. One important contribution in this direction is seen in [58] , where the method of snapshots was proposed, paving the way to more algorithms and techniques. In this paper, we employ the randomized SVD (rSVD) algorithm [25, 38] , a non-deterministic algorithm able to compute the near-optimal low-rank approximation of a given large dataset with good efficiency.

Regarding numerical methods to approximate PDEs, the use of finer meshes in finite element discretizations usually leads to more accurate solutions. On the other hand, reducing the number of equations n eq in the nonlinear system is crucial for efficiency, especially considering that the optimal computational complexity of a single physics transient finite element simulation is O(n 4 3 eq ) [15] . The duality between the two statements describes a well-known trade-off between accuracy and efficiency in the finite element context. Despite a considerable research effort in the past decades, strategies to generate tailored meshes to maximize the accuracy while minimizing the computational effort are still an open research topic. Milestones addressing this subject are finite element a posteriori error estimators/indicators [2] , techniques such as adaptation, interpolation [17, 47] and projection [18, 27, 50] , as well as libraries and frameworks containing automated versions of AMR/C techniques [4, 42] .

The general structure of the AMR/C scheme is given in Algorithm 1 and illustrated in Fig. 1 . Three criteria are fundamental in an AMR/C algorithm: remeshing, flagging, and stopping. The remeshing criterion defines whether the computed solution at a given time step requires remeshing driven by global a posteriori error estimators (or indicators) and/or by calling the refinement/coarsening procedure at every j time steps. Next, all mesh elements are visited and flagged for refinement or coarsening. The element flagging criterion is often represented by local a posteriori error estimators or indicators, i.e. flux-jumps of the solution gradient. In possession of the flagged elements, the remeshing algorithm is invoked. Finally, the previous mesh solution must be projected or interpolated into the target mesh. This is done by the projection/interpolation algorithm. The whole process is repeated until the stopping criterion is achieved. This criterion could be error equidistribution (until a certain threshold), a given element size, the maximum number of Fig. 1 Illustration of a mesh refinement procedure. A local a posteriori error estimator or indicator flags an element for refinement (in green) using the solution computed in the mesh on the left. The mesh is refined (or coarsened) according to the flagged elements, and the process can be restarted until a given criterion is met (error level, element size, maximum number of elements, etc.). Note that the initial and the final mesh differ in the number of degrees of freedom and topology elements in the mesh, or the number of refinement levels. Note that the meshes generated in Algorithm 1, guided by the error estimation procedure, have different dimensions (number of degrees of freedom) and different topologies, which poses difficulties for DMD (or, in fact, for any snapshot-based method).

The structure of the exact DMD algorithm relies on the fact that the measurements of the state vectors, i.e., h 0 … h m have the same dimensionality. This structure could be exemplified in numerical experiments as fixed discretizations in space, i.e., fixed meshes or static sensors in experimental data. However, finite element simulations equipped with AMR/C strategies provide solutions in different function spaces, depending on the mesh used to compute the solution on a given time step. Spatial adaptivity on transient finite element simulations leads to meshes with a different number of nodes, numbering, and nodal coordinates. It can also lead to different mesh topologies and structures. For that reasons, snapshots obtained by AMR/C simulations cannot be stacked in columns to construct the snapshot matrix. Even if one considers an AMR/C strategy that restricts the adaptive meshes to preserve the dimensionality of the snapshots, the difference between the nodal coordinates of the various meshes will lead to misleading dynamics captured by DMD.

In the present study, we circumvent this issue by considering the L 2 -projection [50] of the numerical simulation results for different meshes into a reference target mesh. Recent work on reducing these projection costs for related reduced-order modeling techniques (though not DMD) may be found in [26, 32, 37] . The projection or interpolation of numerical solutions between finite element meshes is a well-known computational mechanics subject. Many issues regarding boundary conditions, data visualization, or coupling arise from this kind of problem. The choice of a proper method to successfully project functions in different finite element spaces is not a trivial task since conservation may not be satisfied [18, 27, 50] . This issue is not addressed in the present study since the mesh projection occurs as a post-processing phase after computing the AMR/C solutions. Therefore, it cannot lead to cumulative errors. Furthermore, we are also careful to choose reference target meshes with characteristic length equivalent with the existent in the AMR/C meshes to avoid any major accuracy losses. For the sake of generality, we consider the L 2 -projection as our method of choice, and it can be defined as follows [61] . Assuming that the solution on the donor mesh h ( ) = ∑ n j=1 j j ( ) to be projected onto the target mesh must satisfy the orthogonality condition, where V is a finite-dimensional subspace of L 2 (Ω) defined by the target mesh and the interpolant proj is the optimal interpolant in the L 2 -norm for V . The orthogonal projection can be defined in terms of the following linear system where is the mass matrix and is the projection matrix. The mass matrix is usual in finite element computations. The matrix, however, can present technical difficulties (see "Appendix").

Therefore, our strategy consists of applying the L 2 -projection onto a tailored reference target mesh capable of representing the many scales in time and space of all the snapshots. This routine is inserted in the code and invoked after the adaptive procedure in Algorithm 1 at every output

time step Δt o . This tailored reference mesh is described in this work as target reference mesh and should not be confused with the target mesh generated during the AMR/C procedure. The new solution (computed on the adaptive donor mesh) is projected onto the reference target mesh and exported as a simulation output. In this work, we export as output files for visualization purposes, although there is no restriction on stacking the snapshots on the snapshot matrix during the simulation runtime or dumping only proj on disk. We considered Gmsh [30] , an open-source robust mesh generator, as our software of choice for defining and creating the target meshes for this study. The output for each time step is the snapshot with constant dimensionality n such that all nodes in space are correctly mapped and capturing the dynamics existent in the system. The procedure is summarized in Algorithm 2. This strategy is relatively simple since the L 2 -projection consists of solving a linear system where the generated matrix is a mass matrix and requires no extra outputs for storing the projected solutions since the projection can be applied right after the AMR/C code. The mass matrix is generated in the finite element context by a selfadjoint operator, enabling more efficient solvers. In terms of versatility, the L 2 -projection method is flexible because a solution obtained for a given mesh can be naturally projected onto reference target meshes with different topologies and dimensionalities. Also, since the mesh projection is a vital part of AMR/C algorithms, finite element algorithms frequently present efficient implementations of interpolation or projection techniques. Figure 2 shows an example where a solution obtained by an adaptive mesh simulation is projected onto two meshes with different topologies.

For this example, we consider the square domain Ω = [−1, 1] × [−1, 1] and the function u defined in Ω such that This function is approximated on a structured finite element mesh discretized into 10 × 10 cells where each cell is divided into two triangular elements. An AMR/C procedure is invoked to refine three times the transition between = 1 and = 0 , creating a new mesh containing 1672 elements and 857 nodes. Figure 2a shows the new mesh generated after the AMR/C procedure and the function approximated by the resulting function space. Two meshes-one structured and another unstructured-are considered for projection. The structured mesh contains 80 × 80 cells, resulting in 12,800 triangular elements and 6561 nodes. The element sizes of the structured mesh are similar to the smallest elements in the adaptive mesh. The unstructured mesh presents a smaller characteristic length than the structured mesh and contains 17,088 elements and 8705 nodes. 

Remark In this study, the projection is carried out inside the finite element simulation since the projection computational cost is practically negligible in comparison with the overall time required for solving nonlinear systems of a complex numerical simulation. However, if one does not have access to the finite element simulation codes used, the reconstruction of the solution of each time step can still be done offline. Output files of various formats contain information regarding the mesh used (such as nodal coordinates and connectivities) for visualization purposes. By properly reading these files, the solutions can be reconstructed under a finite element framework (such as FEniCS [4] or libMesh [42] ) or on a code developed by the user. However, this totally non-intrusive approach can significantly increase the computational cost due to several I/O operations that are often extremely low compared to computational intensive operations. Since the mesh constantly varies in an AMR/C finite element simulation, each solution obtained in the simulation has to be imported, reconstructed under its original mesh, and projected. After the projection, the user can choose to dump the solution in the disk or stack the snapshots on the snapshot matrix. For the first case, another I/O operation would be invoked. This non-intrusive strategy is especially suitable (and restricted) to be used with DMD, which is also a non-intrusive algorithm. For intrusive snapshot-based methods such as Proper Orthogonal Decomposition (POD), one needs to project finite element matrices onto the computed basis and, therefore, requires access to the code [62] .

In this section, for the sake of generality, we apply our method in several applications, with different systems of equations, numerical formulations, mesh topologies, spatial dimensions, refinement criteria, and finite element libraries.

We compare the results between the AMR/C solution, the fixed mesh solution, and the DMD results for all cases. First, we test the DMD short-time future prediction capabilities on a continuous SEIRD model for COVID-19, a nonlinear system of diffusion-reaction equations. The equations are considered using a Galerkin finite element discretization and are solved using libMesh [42] , a high-performance C++ finite element library. The error estimators, refinement/ coarsening strategies built-in on libMesh can be seen on [34, 63] . Also, the L 2 -projection algorithm is embedded in libMesh. We explore the results in one and two spatial dimensions, where the 1D case is a hypothetical example, and the 2D case describes the COVID-19 evolution in the Lombardy region in Italy [34, 63, 64] . The simulation obtains the results for 44 and 60 days for the 1D and 2D cases, respectively. Since we want to predict 14 days in the future, we only feed DMD with snapshot matrices containing the first 30 and 46 days in the 1D and 2D cases, respectively, such that the predicted results can be compared with the results obtained in the simulations. To apply DMD to the SEIRD data obtained on an AMR/C simulation, we project the adaptive solution onto a mesh with characteristic length compatible with the obtained simulation results. The discretization of the reference mesh used for the projection is as fine as the final refinement level on the adaptive meshes. Next, we consider two fluid dynamics applications where the governing equations for both cases are advection-dominated.

To circumvent the LBB condition and spurious oscillations regarding dominant advection, these equations consider the residual-based variational multiscale (RBVMS) formulation [1, 8, 19, 35, 39 , 53] on a finite element discretization. We consider the use of DMD on a 2D density-driven gravity flow and a 3D bubble rising simulation. The density-driven gravity current is modeled by the coupling of the incompressible Navier-Stokes equation and the advection-diffusion equation. We consider a lock-exchange problem, where a tank is filled with two fluids of different densities, separated by a lock. The simulation starts when the lock is removed, and the difference in the densities of the fluids generates the driving forces responsible for the motion of the fluids. Unlike the other examples in this work, the implementation of this numerical test is made on the FEniCS v.2019.1 framework [4] , a high-performance Python/C++ finite element library. The refinement/coarsening algorithm and L 2 -projection algorithm for this example are part of the framework. We consider an interface-tracking error indicator for the AMR/C simulations, and the mesh is refined following the bisection method [54] . DMD is considered to reconstruct the solution, and the results are compared to the fixed mesh results and the results obtained by AMR/C simulations. Finally, we extend our analysis to a 3D bubble rising case [33] , a two-phase incompressible flow problem where the interface is captured by the convected level-set method. This model is implemented on libMesh, taking advantage of the same refinement/coarsening strategies as well as the projection algorithm used in the SEIRD numerical tests. In this example, we test the projection of the adaptive solutions onto three different meshes containing the different scales existent on the AMR/C simulation and evaluate the results. We consider DMD to predict the bubble geometry and dynamics for a short time in the future. For all the numerical tests proposed, we evaluate the results in terms of efficiency and accuracy. For efficiency purposes, we compute the ratio between the computational time required to run the finite element simulations and the time required to run DMD separately. We refer to this quotient as speedup. The finite element code is responsible for computing the snapshots and projecting the results onto the reference target meshes proposed, while the DMD code imports the output files, extract the snapshots, computes the approximation, and outputs the results. Also, we provide a table describing how the projection routine affects the overall computational time required for the simulations for all examples. In terms of accuracy, we evaluate the results in

DMD is a matrix comprising the approximations obtained by DMD and || ⋅ || F denotes the Frobenius norm. A more detailed analysis is also done in terms of relative error in time . For that, we plot the curves of the relative errors (in terms of L 2 -norm) of each snapshot.

where || ⋅ || 2 denotes the L 2 -norm. Also, to avoid unphysical results, some quantities of interest are evaluated and compared using both simulation and approximation results. For the SEIRD model, we plot the results regarding the total population. This quantity should be constant in time according to the hypothesis of the model. For the lockexchange simulation, we compute the mass during the simulation and the front position. Since the simulation is considered on a closed tank, the mass must be kept constant during the simulation. For the 3D bubble rising problem, we plot the quantities of interest related to the geometry (volume and sphericity) and dynamics (center of mass and rise velocity) of the bubble. In this case, we test meshes with different minimum characteristic lengths. The results are shown and discussed below.

The outbreak of COVID-19 in 2020 has led to a surge in interest in the mathematical modeling of infectious diseases. This new virus is responsible for infecting millions of people worldwide and impacting the economy in an unprecedented way. Therefore, numerical simulation of the virus' dynamics may help provide short-term prediction models for forecasting the number of future cases. In this perspective, it is possible to develop strategic planning in the public health system to avoid deaths and manage patients. Disease transmission may be modeled as compartmental models, in which the population under study is divided into compartments and has assumptions about the nature and time rate of transfer from one compartment to another [11] . Here, we work with a spatio-temporal SEIRD model, presented in [34, 63, 64] , and given by,

, and d( , t) denote the densities of the susceptible, exposed, infected, recovered, and deceased populations, respectively. The sum of all the compartments with the exception of d( , t) is represented by n pop which is the total living population. A e characterizes the Allee effect (persons), that takes into account the tendency of outbreaks to cluster around large populations, i and e denote the transmission rates between symptomatic and susceptible individuals and asymptomatic and susceptible

individuals, respectively (units days −1 ), denotes the incubation period (units days −1 ), e corresponds to the asymptomatic recovery rate (units days −1 ), i the symptomatic recovery rate (units days −1 ), represents the mortality rate (units days −1 ), and s , e , i , r are the diffusion parameters of the different population groups as denoted by the sub-scripted letters (units km 2 persons −1 days −1 ). Note that all these parameters can be considered time and space-dependent.

We also compute the compartment c, the cumulative field of the i compartment.

For the numerical solution of (21)-(25), we discretize in space using a Galerkin finite element variational formulation. The resulting systems of equations are stiff, leading us to employ implicit methods for time integration. We apply the Backward Differentiation Formula (BDF2), which offers second-order accuracy while remaining unconditionally stable. We implement the whole model in libMesh [42] . We additionally make use of AMR/C, allowing us to resolve multiple scales. One may find more details about the methods in [34, 63] .

First, we use a simple 1D continuous SEIRD model for COVID-19 with adaptive mesh refinement to validate the L 2 -projection and DMD. This example was introduced in [63] and reproduced in [34] . Basically, we consider a 1D region Ω = [0, 1] with initial conditions that represents a large population centered around x = 0.35 with no exposed persons and a small population centered around x = 0.75 with some exposed individuals, as shown in Fig. 3 . Thus, we set s = s 0 and e = e 0 as follows: We use an adapted mesh with initially 125 elements, and after the refinement, the smallest element has a size 0.002. At the beginning of the simulation, we refine uniformly the whole domain into two levels and, after that, we apply the adaptive mesh refinement every 4 time steps. The idea is that the AMR/C strategy will keep this spatial resolution on more dynamically relevant regions while coarsening other regions in the domain. As a target reference mesh, we consider the uniformly refined mesh, such that all the domain contains the minimal spatial resolution obtained by the AMR/C simulation.

The results of this simulation are validated against the results from [34, 63] . Figure 4 show the solution of the 1D SEIRD example at t = 30 days for the fixed mesh simulation, the adaptive mesh simulation, and the projection of the adaptive solution onto the target reference mesh used in the fixed mesh simulation. We observe a good agreement between the solutions. This agreement is a positive indicator that the DMD can be used on the projected meshes with little to no error compared to the DMD on a fixed mesh simulation. In terms of efficiency, Table 1 shows the computational effort required for the reference mesh projection embedded on the adaptive finite element code in comparison with the time used for the simulation code itself (AMR/C FEM).

After some initial experiments, we set r = 15 for all compartments. In this example, the finite element simulation computes the results for 44 days, that is, 176 time steps. However, we only consider the first 30 days in the snapshot matrix for reconstruction. DMD is set to approximate the results for a further 14 days so the results can be compared with the 44 days results obtained in the finite element simulation. In other words, we want to evaluate the DMD ability to predict the COVID-19 scenario two weeks in the future given the data of the last 30 days. For numerical reasons, we considered an initial 3 days shift in the snapshots since some compartments are initialized with zeros, affecting how DMD captures the dynamics. The 44th day solution for the adaptive simulation and the 44th day prediction considering DMD are seen in Fig. 5 . We observe good agreement between predictions and the numerical solutions of the simulations for most compartments. The speedup and overall relative error between the DMD approximation and the snapshots are seen in Table 2 .

We extend our analysis by solving the continuous SEIRD model and applying DMD to a 2D real world domain that is the Lombardy region in Italy. The spread of the COVID-19 has been studied in this region using the continuous SEIRD model with accurate results [63, 64] . Here, we reproduce this simulation with the solver developed in [34] which invokes adaptive mesh refinement every 4 time steps. We use the same parameters as shown in [63, 64] . It is important to point that, in this simulation, the transmission rates and diffusion parameters vary with time in order to reproduce the effects of restrictions during the simulated period.

For this simulation, an unstructured mesh is considered due to the complex geometry imposed by the domain. The mesh is generated using Gmsh and is uniformly refined as the simulation starts. After refining the whole mesh in one level, the mesh presents a minimum spatial resolution of approximately 1 km. This procedure allows the solver to coarsen the regions where no significant dynamics are observed while preserving the scales of the regions of interest. Figure 6 shows the variation of the number of nodes in time for the AMR/C strategy. The coarsening approach improves the simulation performance significantly since the average number of nodes (and respectively, the number of equations n eq ) used in the adaptive simulation is approximately half with respect to the case of a fixed mesh with the same spatial resolution considered for the entire domain. In this example, the reference target mesh is the uniformly refined unstructured mesh considered in the early stages of the simulation, presenting 13,158 nodes and 25,340 elements. The simulation considers a time step size of Δt = 0.25 days for the numerical integration and Δt o = Δt = 0.25 days for the observations. Initial conditions for the Lombardy domain are the same presented on [63, 64] and are seen in Fig. 7 , while compartments r and d are initialized to zero.

We then proceed to run both adaptive and fixed mesh simulations and, to apply DMD in the adaptive mesh results, we consider the proposed projection scheme. Figure 8 shows the s compartment solution at t = 46 days for both simulations and the projected adaptive solution onto the reference mesh, revealing that the results are in good agreement. In terms of efficiency, Table 3 shows results for the computational time required for the projection compared to the finite element code. That said, the adaptive snapshots can now be assembled into a snapshot matrix for the DMD reconstruction and prediction.

The DMD analysis is made in the same way as presented in the 1D case: the simulation outputs the projected snapshots for the first 60 days (240 snapshots). The snapshot matrix assembles the information regarding 46 days of simulations, while DMD approximates the results for 60 days. The idea is to predict two weeks in the future, given the data observed in the past 46 days. For this example, we consider the SVD truncation at r = 20 for all compartments. Again, an initial 3 days shift in the data is considered to avoid issues with the compartments initialized to zero. Results for the 60th day comparing the computed numerical solutions and the DMD predictions are seen in Figs. 9 and 10. We present the projected solutions, the DMD prediction, and the relative error in space between the two results from left to right. We can note that most compartments show results in agreement with the simulations, while the exposed compartment reveals more pronounced differences than the other compartments. Figure 11 shows the relative error in time between the DMD results and the projected snapshots. The first thing we notice is that the curves are different for each compartment. This discrepancy occurs due to the different parameters for each equation in the SEIRD model, which largely affects the dynamics of the system. The dynamics for each compartment are different since each compartment presents different coupling, diffusion, and reaction parameters. Also, regarding this issue, since the parameters are time and space-dependent, sudden changes in their values can affect the dynamics of the system as well as DMD's dynamics mapping ability. Some sudden changes in the s and e compartments related to stricter public policies considered to reduce the transmission rates (parameters i and e ) are incorporated into the model. Since the variation in the parameters is not introduced smoothly, DMD's ability to map sudden changes in the dynamics of the system is reflected by the existence of some spikes on the curves of the relative errors in time. Comparing the reconstruction and prediction stages, we observe that the errors tend to grow as soon as the prediction stage starts (dashed line). The exposed compartment, which yielded most of the oscillations due to parameter changing on the reconstruction stage, presented the same behavior on the prediction phase around day 49. We also note that the exposed compartment yields a larger relative error for the 60th day in comparison with the other compartments. Table 4 shows the overall relative error and the speedups for the six compartments approximations. Comparing these results with the results presented in Figs. 9, 10 and 11, we can conclude that the predictions are reasonably accurate in comparison with the numerical solutions, specially when considering the time required for calculation.

Another important analysis to be done is the conservation property of the continuous SEIRD model. As mentioned before, the standard L 2 -projection does not guarantee conservation among the projections. Figure 12 shows the total population during the simulation, normalized by the total population modeled in the initial conditions. The total population is computed as the sum of the integral of the compartments (excluding c) divided by the sum of the integral of the elements of the mesh. Since the SEIRD model does not consider any population growth, the value must be theoretically constant for all the simulations. From the figure, we observe that the population is kept constant during all the adaptive simulation, and this was preserved by the projected solutions and the DMD reconstruction stage. That is, we can note that the L 2 -projection does not yield conservation issues in this example. For the prediction phase, DMD preserves the total population for several days in the future. However, it presents a slight increase (around 0.1% ) for predictions over 10 days, which does not affect the results significantly. This increase can be explained by the relative errors behavior, observed in Fig. 11 , as DMD computes future estimates.

This section evaluates the DMD use on two cases involving AMR/C in computational fluid dynamics: the reconstruction of a 2D density-driven gravity current simulation and the temporal prediction on a 3D rising bubble. Different from the previous cases, the test cases presented in this section are advection-dominant. To approximate the governing equations, we use a finite element RBVMS formulation [1, 19, 35, 39, 53] . In the first case, we reconstruct the solution and evaluate important quantities of interest regarding densitydriven gravity flows. On the bubble rising simulation, we show how our strategy works on a 3D mesh, and we evaluate the DMD ability to predict the quantities of interest of the rising bubble in time. 

In this section, we consider a long numerical simulation that consists of a lock-exchange between two fluids, the heavy fluid, A, and the lighter fluid, B, based on the numerical example in [49] . The difference between their densities is such that the Boussinesq hypothesis is considered valid. Moreover, particles in the heavy fluid have negligible inertia and are much smaller than the smallest length scales of . 9 Comparison between computed and predicted solutions at t = 60 days for the susceptible, exposed, and infected compartments the buoyancy-induced fluid motion. Thus, the dimensionless governing equations are

where is the fluid velocity, is the concentration field, p is the pressure, g = (0, −1) is the vector pointing in the direction of gravity, Sc = 1.0 is the Schmidt number and Gr = 5 × 10 6 is the Grashof number, two dimensionless numbers that relate viscous effects with diffusion and buoyancy effects, respectively. A Grashof number of this magnitude indicates a turbulent flow. The field = A ∕ B is the concentration and is responsible for mapping the evolution of fluid interactions. The time step size considered for this simulation consists on Δt = 0.01 s for a total simulation time Figure 13 illustrates the domain and the initial conditions.

To solve the governing equations, we implement the RBVMS formulation [35] for Eq. (28) using the FEniCS 2019.1 [4] framework to generate the snapshots for this example. The adaptive mesh refinement procedure returns the solution for , p, and . For this example, we only consider the snapshots of for our calculations, that is, an overall data reduction of 75% . Details of the formulation of the problem can be found in [7, 35] . We consider a fixed mesh simulation with 701 × 101 nodes and 700 × 100 cells, where each cell is divided into two linear triangles. We consider an interface-tracking adaptive mesh error indicator that flags and refines the mesh where the two fluids interact for the adaptive mesh simulation. The error indicator for the mesh refinement is |∇ | being larger than a given tolerance. For this purpose, a mesh containing 175 × 25 cells is considered, the interfaces are refined considering two levels of refinement, and the mesh refinement is invoked at every time step. Figure 14 shows the results at t = 10 s for both fixed and adaptive mesh simulations and the projection of the adaptive solution onto the same mesh used in the fixed mesh simulation.

AMR/C in this problem is advisable since the dynamics are predominant on the interface between the fluids. Most of the domain is not affected in the early stages of the simulation, and the use of fine meshes outside these regions may represent unnecessary computational effort. Figure 15 shows the number of nodes in the mesh during the simulation time for the adaptive mesh compared with the fixed 701 × 101 nodes mesh.

We observe that the number of nodes in the adaptive simulation is smaller than the fixed mesh for all the simulation time. Table 5 shows the time required for the simulation to run in comparison with the time spent on projecting the solutions onto the reference target mesh. Now we proceed applying DMD to the projected solution and reconstructing the solutions. Figure 16 shows the relative error for the reconstruction using different values of the rank r. The results are also confirmed in Table 6 , where the relative error between the reconstructed and projected snapshot matrices is shown and the speedup computed for each case. We observe that, for increasing values of r, the relative error decreases for all steps and affects the overall relative error of the matrices. That is, inserting more structures in the DMD basis yields better accuracy in terms of overall relative error. As for the speedup, we notice that it decreases for Table 6 Relative error between reconstructed data and the snapshots and speedup between DMD and the numerical simulation increasing values of r. Such a decrease occurs because larger values of r directly affect the rSVD algorithm performance [7, 25] and increase the dimensions of the matrices for the computation of the DMD basis. Most importantly, we can evaluate the quantities of interest common to density-driven gravity currents. Figure 17 shows the front position and mass conservation for both simulations and the respective reconstructions. We note that, for these quantities all reconstructions are in extremely good agreement with those computed with the fixed and adaptive meshes, approaching them, as expected, for higher values r.

We now study a bubble rising 3D benchmark, whose task is to track the evolution of a three-dimensional bubble rising in a liquid column. The initial configuration is described in Fig. 18 .

For this problem, we couple the Navier-Stokes equations with an interface capturing method called convected levelset [20, 33, 65] . The convected level-set is a method used to represent the interface between two phases and, by a convection equation, to move the interface as the flow evolves. A force that has an important role in bubble problems is the surface tension , which is applied using the Continuum Surface Model (CSF) [10] .

We write the governing equations in their dimensional form as where is the density, is the dynamic viscosity, is the acceleration of gravity vector, is the level-set function, , and S a function related to the level-set signed distance function. For the temporal integration, we apply the Backward-Euler method to the Navier-Stokes equations, while for the convected level-set, we use the BDF2 method. One may find more details about the governing equations and methods in [33] .

The initial configuration consists of a spherical bubble of radius R = 0.25 m centered at [0.5, 0.5, 0.5] m in a [1 × 1 × 2] m domain. The no-slip boundary condition is Table 7 lists the parameters used for this simulation.

We use an adapted mesh, initially with 10 × 10 × 20 cells, with each cell divided into 6 linear tetrahedra. We refine the initial region where the bubble is located into two levels and, after the refinement, the smallest element has a size of 0.025 m. The adaptive mesh refinement is based on the flux jump of the level-set function error, in which h max = 2 . We apply the adaptive mesh refinement every four time steps. The interface is modeled with E = 0.05 , and the time step size is defined as Δt = 0.0125 s. We output the projected solutions at every 2 time steps such that Δt o = 0.025 s. In this example, we consider three tetrahedral meshes for our projection strategy, presented in Fig. 19 . The three meshes named coarse, intermediate and fine, present characteristic lengths similar to the three scales existent in the refinement levels of the adaptive mesh simulation. The coarse mesh represents the initial mesh on the adaptive simulation, with 12,000 elements and 2541 nodes. The intermediate mesh presents smaller elements equivalent to the generated elements after the first refinement level on the AMR/C simulation, totalizing 96,000 elements and 18,081 nodes and the fine mesh contains the smallest scales presented on the adaptive numerical solution with 768,000 elements and 136,161 nodes. The figures show the level-set solution on half of the domain and at t = 3 s for the adaptive mesh solution and the respective projections. We also present on Table 8 the time required for the projection of the solutions of on the three To verify if the bubble geometry and dynamics are being preserved during the mesh projections, we evaluate the quantities of interest such as the bubble volume, sphericity, center of mass, and rise velocity. Bubble volume and sphericity are related to the bubble geometry, while the position of the bubble center of mass and rise velocity regard the bubble dynamics. The results are compared for the AMR/C simulation output and the projections on Fig. 20 . The quantities of interest are computed using the adaptive snapshots as well as the projected snapshots and compared with each other. Regarding the geometry quantities of interest, we observe that the coarse mesh projection solution affects the bubble's geometry, leading to bad results regarding the bubble volume and sphericity compared with the AMR/C simulation. For the intermediate and the fine meshes, the bubble's geometry is not largely affected, and the results are compatible with the AMR/C results. As for the quantities of interest related to the bubble dynamics, no significant difference is observed on the projection of the three meshes. We observe that the center of mass is not affected by the projections compared with the simulation results. As for the rise velocity, we observe some minor differences in all projection cases.

The simulation proceeds until T = 3.0 s, yielding a dataset containing 240 snapshots regarding the solution of for each target mesh from the projections. We do not consider the use of velocities and pressure in the DMD analysis, reducing the required data in 80% . We consider the results for the first 2.75s to construct the basis and predict the last 0.25 s. We compare the DMD results in terms of relative error between the snapshot matrix and the obtained solutions for each case. The DMD solution for the coarse mesh is compared to the projected adaptive solution onto the coarse mesh and so forth. Results are evaluated for multiple values of r, such that r = {5, 10, 15, 30, 45, 60} . The results regarding accuracy and performance for multiples values of r and the three meshes are presented in Table 9 .

The results for r = 60 are shown in Fig. 21 . We observe that the errors are stable for the reconstruction case, that is, the DMD solution before t = 2.75 s. From that point, shown as a dashed line on the figure, the errors begin to grow exponentially for each predicted time step, while still remaining below 1% until around 2.9 s. We observe that the errors in the reconstruction case are different regarding the mesh used. That is, the errors are larger with respect to the minimum characteristic length of the projection meshes. However, we observe that the errors grow at the same rate on the prediction phase independently of spatial discretization.

We also show the results considering r = 60 for reconstruction at t = 2.75 s and prediction at t = 3.0 s in Fig. 22 . We compare the DMD results for the three projected meshes, the adaptive mesh solution, and a fixed mesh solution. The fixed mesh solution is obtained by running the simulation with a mesh of 768,000 elements and 136,161 nodes, as the fine mesh used in the projection. We observe initially that the bubble geometry is better defined on the reconstruction than on the prediction figure. This better definition is directly related to the errors observed in Fig. 21 . We observe that the coarse mesh results do not capture the bubble geometry with We now proceed comparing the results in terms of quantities of interest for the DMD results. Figure 23 shows the bubble volume, sphericity, center of mass, and rise velocity for the DMD results compared to the adaptive solution results. The same issue regarding sphericity on the coarse mesh projection is observed on the coarse mesh DMD results. However, for the intermediate and fine meshes, the values match the results observed for the projection. We observe results in conformity for the center of mass and rise velocity as well. In terms of prediction, we observe that DMD accurately predicts the volume and the center of mass evolution. As for the other quantities of interest, the increasing exponential errors in the DMD prediction structure affect the quantities of interest for long time future predictions.

In this work, we propose a strategy to enable data-driven snapshot-based methods on finite element solutions obtained by AMR/C simulations. The use of AMR/C algorithms in finite element approximations of PDEs is known to reduce memory usage and to increase the efficiency of the simulations without compromising the accuracy. The simulation adapts the mesh with the evolution of the solution, refining regions of interest and coarsening regions that are not of interest. The adaptation process leads to different meshes during the simulations, and the solutions have different dimensions and topologies (or different connectivities and nodal coordinates) whenever the AMR/C algorithm is invoked. Snapshots with different dimensions prevent the use of snapshot-based algorithms such as DMD, where the snapshot matrix is built by stacking the solutions for each observation in columns. In this study, we considered a strategy to project the adaptive solutions on reference meshes such that all snapshots present the same dimensions and nodal indices, enabling the construction of snapshot matrices. The method employed to project the solutions onto the target reference mesh is the L 2 -projection, a simple, fast and versatile approach. The L 2 -projection is a common strategy present in several finite element libraries and frameworks and consists basically of solving a linear system, where the matrix is obtained by a self-adjoint operator enabling the use of efficient solvers. Despite presenting drawbacks, especially regarding properties' conservation, we investigate the use of the L 2 -projection as a postprocessing tool without any relevant issues. By postprocessing, we mean that the L 2 -projection algorithm is invoked only to insert the solutions on a reference function space for each time step and output the files. This strategy does not yield significant additional computational effort as we observe that the L 2 -projection routine required around 1% of the computational time required for the adaptive finite element code to run in most cases and 5% in the worst case. When the source code is not available to invoke the L 2 -projection, one can construct the solutions from the output files and project them onto the target reference mesh in a complete non-intrusive workflow.

We test the algorithm on several models presenting different dynamics and underlying physics. First, we present the results for DMD on a continuous SEIRD model for COVID-19 for fictional data in 1D and real data for Lombardy, Italy, in 2D. The idea of considering short-time future estimates on COVID-19 models could improve the decision-making of public policies to avoid further contamination, and the use of AMR/C in the simulations, coupled with the presented projection scheme and DMD, can lead to fast and reliable predictions. The simulations are implemented using the libMesh library using its refinement built-in functions. We compare the solution, its projection in the reference mesh, and the adaptive solution and notice that no relevant errors are found for the projection strategy. In this example, the L 2 -projection consists of 1.47% and 1.48% of the total time required for the code to run for the 1D and 2D cases, respectively. As for the use of short-time future estimates, we consider DMD for predicting two weeks in the future given a set of snapshots. For the 1D case, we feed DMD with snapshots covering 30 days, while for the 2D case, the snapshot matrix comprises 46 days of simulation results. The DMD results are presented in terms of efficiency and accuracy compared to the projected solutions, and we observe that the predictions are mostly in good agreement except for the exposed compartment. To test the approach's agnosticism regarding the dynamics of the systems, we consider two different fluid dynamics problems: a 2D density-driven gravity current and a 3D bubble rising problem. We use a 3000-snapshot simulation to generate a basis containing the spatio-temporal coherent structures for the density-driven gravity flow. The lock-exchange simulation is solved using the FEniCS framework. The computational effort regarding the projection, in this case, corresponds to 5.06% of the overall computational time. As for the DMD analysis, we evaluate the results for different values of the rank r and observe that, as we increase it, in this case, the overall relative errors decrease significantly. However, when evaluating the quantities of interest such as mass conservation and front position, even the worst case yields good approximations. For the 3D example, we use the libMesh library, testing three different reference meshes considering the three element sizes existent in the adaptive simulation. We investigate the occurrence of projecting the adaptive solution onto meshes that do not necessarily guarantee that all the scales obtained in the solution are preserved in terms of quantities of interest related to the bubble geometry and dynamics. We observe that the use of L 2 -projection on coarser meshes leads to issues regarding the shape of the bubble, affecting quantities of interest such as bubble volume and sphericity. When considering the L 2 -projection on meshes with approximately the same spatial resolution as the finest scale in the adaptive meshes, these effects are mitigated. As for the bubble dynamics, even coarser meshes reveal similar results for the bubble center of mass and rise velocity in comparison with finer mesh projections and the adaptive solution itself. In terms of computational effort, the L 2 -projection corresponds to 0.03% , 0.14% and 1.09% of the total simulation time for the coarse, intermediate, and fine meshes, respectively. Proceeding to the DMD analysis, we present the results for a 220 time step reconstructions and a 20 time step prediction for the bubble rising problem. We test the approximations for multiple values of r being r = 60 our best approximations. The results computed by DMD are presented in terms of efficiency and accuracy compared to the projected solutions. We also reevaluate the quantities of interest, now with the DMD approximation. For the reconstruction, results are practically identical to those observed for the projection, while for the prediction, the center of mass and bubble volume yield good results while sphericity and rise velocity are not in agreement with those obtained in the simulation.

Summarizing our findings, we highlight that the present approach, despite its simplicity:

• is versatile since every finite element library that presents AMR/C is often equipped with efficient projection algorithms; • is fast in the sense that the L 2 -projection only requires solving a sparse, well-conditioned linear system = where is the mass matrix, obtained by a self-adjoint operator, enabling the use of efficient solvers; • is framework-agnostic, as noted in the results of our simulations implemented on libMesh and FEniCS, requiring only the standard procedures already in place in both libraries; • is problem-agnostic, as observed on meshes with different spatial dimensions (1D, 2D, and 3D), topologies (structured and unstructured), and equations (convectiondominated and reaction-diffusion systems); • preserves the dynamics computed on the adaptive snapshots for all the tested examples.

The results observed in this study reveal that the use of mesh projections for adaptive solutions preserves the dynamics inherent to the finite element solutions. However, this approach on large problems could be inefficient or prohibitive since the mesh projection on a sufficiently fine target reference mesh can lead to unfeasible DMD storage and computation times. As a future work, aiming to circumvent this issue, we plan to explore the use of data compression on DMD, once used for streaming purposes [24] , in the mesh projection context. Also, it is known that DMD approximates the nonlinear dynamics existent in the dataset using a linear operator, being fragile to problems presenting a strong nonlinearity [43] . To deal with this, several Koopman-based frameworks, such as the Extended DMD [66] , the Kernel DMD [67] , and LANDO [5] , use nonlinear observables to improve temporal predictions. In the future, we plan to evaluate our projection approach on nonlinear observables for further investigations.

In the current work, we have performed this by using L 2 -projection of all snapshots onto a reference (or projected) mesh. Let h denote the finite element solution on the original mesh. Then we define the L 2 -projection onto the reference mesh as for all proj in V proj . One may quickly verify that this satisfies the Galerkin orthogonality property on the reference mesh, i.e., for all proj in V proj , the relation: is satisfied. At the algebraic level, the problem (30) may be expressed as:

where is the mass matrix on the FEM space V where −T denotes the left-sided pseudoinverse of (see e.g. [31] ). Thus, rank( ) = dim( h ), DMD may properly represent all relevant dynamics on the reference (projected) mesh.

Theorem If every degree of freedom of the original mesh is also a degree of freedom of the reference mesh, then the matrix is full-rank.

Proof By assumption, for every N h i , the corresponding N proj i (assumed to be more fine and hence greater in number) is defined such that This idea may be alternatively expressed as

where the different V i are finite element spaces corresponding to the different adapted meshes. Such a condition was also used in [62] . We, therefore, have that where 0 < ij ≤ 1 for all i, j is a volumetric constant relating the measures of supp(N proj j ) and supp(N h j ) . By (38) , the matrix ̃ , defined such that which has the same rank (that is, dim(V h )) as h by construction. We then observe that, for i, j ≤ dim(V h ):

Thus, we can express as Assume now for the sake of contradiction that is not fullrank. This implies that there exists a ∈ ℝ dim(V h ) such that = . Then implying that ̃ = 0 , contradicting the full-rank of ̃ . Hence, we can conclude is full-rank. The above theorem requires an assumption that, while often the case in practice, is not strictly necessary in our 

experience to obtain acceptable results using DMD. As mentioned, a similar condition was required for the related POD method in [62] , in which a common finite element space V n was implicitly constructed from the different spaces from each adapted mesh. Issues regarding projection between meshes were also studied in [9, 21, 48] , in which the interested reader may find more theoretical results regarding projection between meshes. ◻

A review of variational multiscale methods for the simulation of turbulent incompressible flows

A posteriori error estimation in finite element analysis

Understanding mass transfer directions via data-driven models with application to mobile phone data

The Fenics project version 15

Kernel learning for robust dynamic mode decomposition: linear and nonlinear disambiguation optimization (LANDO)

Workshop Report on Basic Research Needs for Scientific Machine Learning: Core Technologies for Artificial Intelligence. United States

Dynamic mode decomposition for density-driven gravity current simulations

Computational fluid-structure interaction: methods and applications

Generalized grid transfer operators for multigrid methods applied on toeplitz matrices

A continuum method for modeling surface tension

Mathematical models in epidemiology

Data-driven science and engineering: machine learning, dynamical systems, and control

Machine learning for fluid mechanics

Compressed sensing and dynamic mode decomposition

Extreme-Scale AMR, SC '10: Proceedings of the 2010 ACM/IEEE international conference for high performance computing, networking, storage and analysis

Dynamic mode decomposition analysis of high-fidelity CFD simulations of the sinus ventilation

Computational grids: generations, adaptation and solution strategies series in computational and physical processes in mechanics

Locally constrained projections on grids

Variational multiscale methods in computational fluid dynamics

Convection of local level set function for moving surfaces and interfaces in forming flow

A study of prolongation operators between non-nested meshes. domain decomposition methods in science and engineering XIX

De-noising by soft-thresholding

The approximation of one matrix by another of lower rank

Compressed dynamic mode decomposition for background modeling

Randomized matrix decompositions using R

Reduced-order modelling of an adaptive mesh ocean model

Conservative interpolation between volume meshes by local Galerkin projection

Fasel U (2020) Data-driven nonlinear aeroelastic models of morphing wings for control: data-driven nonlinear aeroelastic models

The optimal hard threshold for singular values is 4∕ √ 3

Gmsh: a 3-d finite element mesh generator with built-in pre-and post-processing facilities

POD reduced-order modeling for evolution equations utilizing arbitrary finite element discretizations

A new convected levelset method for gas bubble dynamics

Adaptive mesh refinement and coarsening for diffusion-reaction epidemiological models

Numerical simulation of particle-laden flows by the residual-based variational multiscale method

A deep collocation method for the bending analysis of Kirchhoff plate

A hyper-reduction method using adaptivity to cut the assembly costs of reduced order models

Finding structure with randomness: probabilistic algorithms for constructing approximate matrix decompositions

Multiscale and stabilized methods

Low-rank dynamic mode decomposition: optimal solution in polynomial-time

Sparsity-promoting dynamic mode decomposition

libmesh: a c++ library for parallel adaptive mesh refinement/coarsening simulations

Dynamic mode decomposition: data-driven modeling of complex systems

Multiresolution dynamic mode decomposition

Deep learning of material transport in complex neurite networks

Reaction diffusion system prediction based on convolutional neural network

Applied computational fluid dynamics techniques: an introduction based on finite element methods

Optimal constrained interpolation in mesh-adaptive finite element modeling

High-resolution simulations of particle-driven gravity currents

Interpolation with restrictions between finite element meshes for flow problems in an ale setting

Discovering dynamic patterns from infectious disease data using dynamic mode decomposition

Physics-informed neural networks: a deep learning framework for solving forward and inverse problems involving nonlinear partial differential equations

Recent developments in variational multiscale methods for large-eddy simulation of turbulent flow

Mesh refinement processes based on the generalized bisection of simplices

Spectral analysis of nonlinear flows

Dynamic mode decomposition of numerical and experimental data

Application of the dynamic mode decomposition to experimental data

Turbulence and the dynamics of coherent structures I, II and III

Modal analysis of fluid flows: an overview

Image-based modelling for adolescent idiopathic scoliosis: mechanistic machine learning analysis and prediction

Galerkin projections between finite element spaces

Pod-Galerkin reduced-order modeling with adaptive finite element snapshots

Simulating the spread of COVID-19 via a spatially-resolved susceptible-exposedinfected-recovered-deceased (SEIRD) model with heterogeneous diffusion

Diffusion-reaction compartmental models formulated in a continuum mechanics framework: application to covid-19, mathematical analysis, and numerical study

Convected level set method for the numerical simulation of fluid buckling

A data-driven approximation of the Koopman operator: extending dynamic mode decomposition

A kernel-based method for data-driven Koopman spectral analysis

Publisher's Note Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations

Funding Open access funding provided by Università degli Studi di Pavia within the CRUI-CARE Agreement.

In this appendix, we provide a more detailed discussion of several points raised in Sect. 3 regarding the use of DMD on adapted meshes. As discussed in the main text, given the way DMD works by acting on the snapshot matrix, it is, therefore, essential for its proper application that all snapshots are of equal dimension; if the original dataset has been produced using different meshes, this may not be the case. Thus, it is crucial in such a case to use some sort of post-processing technique to provide a uniform snapshot size across the time series. For practical purposes, however, ensuring simply that all snapshots of uniform dimension is a necessary but not sufficient condition to usefully apply DMD. As illustrated in Fig. 24 , we show an example of two meshes (labelled as 1 and 2) that both an have equal number of nodes. Assuming the same order of finite element approximation is used on each mesh, the resulting vectors h 1 and h 2 will have the same dimension. However, while the number of degrees of freedom is the same on each mesh, they are topologically distinct. Hence, to properly perform DMD in this instance, we must have a third mesh, which we call a reference mesh (labeled as mesh 3), on which we project both h 1 and h 2 . We note that this mesh is refined everywhere, compared to the locally refined mesh 1 and mesh 2. In this way, we guarantee that the fine-scale frequencies captured by all the different mesh levels are properly represented by the DMD modes. Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http:// creat iveco mmons. org/ licen ses/ by/4. 0/.