1 Introduction
Timedependent nonlinear problems arise in many important disciplines such as engineering, science, and technologies. They are numerically solved if it is not possible to solve them analytically. Depending on the complexity and size of the governing equations and problem domains, the problems can be computationally expensive to solve. It may take a long time to run one forward simulation even with high performance computing. For example, a simulation of the powder bed fusion additive manufacturing procedure shown in [27] takes a week to finish with cores. Other computationally expensive simulations include the 3D shocked spherical Helium bubble simulation appeared in [4] and the inertial confinement fusion implosion dynamics simulations appeared in [1]. The computationally expensive simulations are not desirable in the context of parameter study, design optimization, uncertainty quantification, and inverse problems where several forward simulations are needed. A Reduced Order Model (ROM) can be useful in this context to accelerate the computationally expensive forward simulations with good enough approximate solutions.
We consider projectionbased ROMs for nonlinear dynamical systems. Such ROMs include the Empirical Interpolation method (EIM)
[6, 22], the Discrete Empirical Interpolation Method (DEIM) [12, 15]and the Gauss–Newton with Approximated Tensors
[10, 11], the Best Point Interpolation Method (BPIM) [31], the Missing Point Estimation (MPE)
[5], and Cubature Methods [3, 19, 20, 24]. There a hyperreduction (The term first used in the context of ROMs by Ryckelynck in [33]) is necessary to efficiently reduce the complexity of nonlinear terms for a considerable speedup compared to a corresponding full order model. The DEIM and GNAT approaches take discrete nonlinear term snapshots to build a nonlinear term basis. Then, they select a subset of each nonlinear term basis vector to either interpolate or datafit in a leastsquares sense. In this way, they reduce the computational complexity of updating nonlinear terms in an iterative solver for nonlinear problems. The EIM, BPIM, and MPE approaches take the similar hyperreduction to the ones in DEIM and GNAT except for the fact that they build nonlinear term basis functions in a continuous space. Whether they work on a discrete or continuous space, they follow the framework of reconstructing “gappy” data, first introduced in the context of reduced order models by
[18]. The cubature methods, in contrast, takes a different approach. They approximate nonlinear integral as a finite sum of positive scalar weights multiplied by the integrand evaluated at sampled elements. The cubature methods developed in [3, 19, 20] do not require building nonlinear term basis. They solve the NonNegative LeastSquares (NNLS) problem directly with the nonlinear term snapshots. Recently, Hernandez, et al., in [24] developed the Empirical Cubature Method (ECM) approach that builds nonlinear term basis to solve a smaller NNLS problem. The requirement of building nonlinear term basis in hyperreduction results in computational cost and storage in addition to solution basis construction. The cost is significant if the corresponding Full Order Models (FOMs) are largescale. The largescale problem requires large additional storage for nonlinear term snapshots and largescale compression techniques to build a basis. In particular, the cost of the hyperreduction in the recentlydeveloped space–time ROM (i.e., STGNAT) [13] is even more significant than aforementioned spatial ROMs (e.g., EIM, DEIM, GNAT, BPIM, MPE, and ECM). The best nonlinear term snapshots of the STGNAT method are obtained from the corresponding space–time ROMs without hyperreduction (i.e., STLSPG) [13], which is not practical for a largescale problem due to the cumbersome size.^{1}^{1}1The full size, implicitly handled by the STLSPG approach due to the nonlinear term and corresponding Jacobian updates, is the FOM spatial degrees of freedom multiplied by the number of time steps.
Therefore, a novel and efficient way of constructing nonlinear term basis needs to be developed.This paper shows a practical way of avoiding nonlinear term snapshots for the construction of nonlinear term basis. The idea comes from a simple fact that the nonlinear terms are related with solution snapshots through underlying time integrators. In fact, many time integrators approximate time derivative terms as a linear combination of the solution snapshots. It implies that the nonlinear term snapshots belong to the subspace spanned by the solution snapshots. Furthermore, a subspace needed for nonlinear terms in a hyperreduction is determined by the range space of the solution basis matrix possibly multiplied by a nonsingular matrix (e.g., the volume matrix). Therefore, the solution snapshots can be used to construct nonlinear term basis. This leads to our proposed method, the Solutionbased Nonlinear Subspace (SNS) method, that provides two savings for constructing nonlinear term basis because of

no additional collection of nonlinear term snapshots (i.e., storage saving).

no additional compression of snapshots (e.g., no additional singular value decomposition of nonlinear term snapshots, implying computational cost saving).
The first saving is especially big for GNAT and STGNAT becuase they involve expensive collection procedure for their best performace.
1.1 Organization of the paper
We start our discussion by describing the timecontinuous representation of the FOM in Section 2. Section 2 also describes the timediscrete representation of the FOM with onestep Euler time integrators. The subspace inclusion relation between the subspaces spanned by the solution and nonlinear term snapshots is described for the Euler time integrators. Several projectionbased ROMs (i.e., the DEIM, GNAT, and STGNAT models) are considered in Section 3 for the SNS method to be applied. The SNS method is described in Section 4 and applied to those ROMs. Section 4 also introduces the conforming subspace condition to justify the SNS method. Section 5 shows an error analysis for the SNS method regarding the oblique projection error bound. It analyzes the effect of the nonsingular matrix that is used to form a nonlinear term basis. Section 6 reports numerical results that support benefits of the SNS method and Section 7 concludes the paper with summary and future work. Although the SNS method is mainly illustrated with the Euler time integrators throughout the paper, it is applicable to other time integrators. Appendix A considers several other time integrators and the subspace inclusion relation for each time integrator. The following time integrators are included: the Adams–Bashforth, Adams–Moulton, BDF, and midpoint Runge–Kutta methods.
2 Full order model
A parameterized nonlinear dynamical system is considered, characterized by a system of nonlinear ordinary differential equations (ODEs), which can be considered as a resultant system from semidiscretization of Partial Differential Equations (PDEs) in space domains
(2.1) 
where denotes a nonsingular matrix, denotes time with the final time , and denotes the timedependent, parameterized state implicitly defined as the solution to problem (2.1) with . Further, with denotes the scaled velocity of , which we assume to be nonlinear in at least its first argument. The initial state is denoted by , and denotes the parameters with parameter domain .
A uniform time discretization is assumed throughout the paper, characterized by time step and time instances for with , , and . To avoid notational clutter, we introduce the following time discretizationrelated notations: , , , and .
Two different types of time discretization methods are considered: explicit and implicit time integrators. As an illustration purpose, we mainly consider the forward Euler time integrator for an explicit scheme and the backward Euler time integrator for an implicit scheme. Several other time integrators are shown in Appendix A.
The explicit Forward Euler (FE) method numerically solves Eq. (2.1), by timemarching with the following update:
(2.2) 
Eq. (2.2) implies the following subspace inclusion:
(2.3) 
By induction, we conclude the following subspace inclusion relation:
(2.4) 
which shows that the span of nonlinear term snapshots is included in the span of scaled solution snapshots. The residual function with the forward Euler time integrator is defined as
(2.5) 
The implicit Backward Euler (BE) method numerically solves Eq. (2.1), by solving the following nonlinear system of equations for at th time step:
(2.6) 
Eq. (2.6) implies the following subspace inclusion:
(2.7) 
By induction, we conclude the following subspace inclusion relation:
(2.8) 
which shows that the span of nonlinear term snapshots is included in the span of scaled solution snapshots. The residual function with the backward Euler time integrator is defined as
(2.9) 
3 Projectionbased reduced order models
Projectionbased reduced order models are considered for nonlinear dynamical systems. Especially, the ones that require building a nonlinear term basis are of our interest: the DEIM, GNAT, and STGNAT approaches.
3.1 Deim
The DEIM approach applies spatial projection using a subspace with . Using this subspace, it approximates the solution as (i.e., in a trial subspace) or equivalently
(3.1) 
where denotes a basis matrix and with denotes the generalized coordinates. Replacing with in Eq. (2.1) leads to the following system of equations with reduced number of unknowns:
(3.2) 
For constructing , Proper Orthogonal Decomposition (POD) is commonly used. POD [7] obtains
from a truncated Singular Value Decomposition (SVD) approximation to a FOM solution snapshot matrix. It is related to principal component analysis in statistical analysis
[26] and Karhunen–Loève expansion [29] in stochastic analysis. POD forms a solution snapshot matrix, , where is a solution state at th time step with parameter for and . Then, POD computes its thin SVD:(3.3) 
where and are orthogonal matrices and is a diagonal matrix with singular values on its diagonals. Then POD chooses the leading columns of to set (i.e., in MATLAB notation). The POD basis minimizes over all with orthonormal columns, where denotes the Frobenius norm of a matrix , defined as with being an th element of . Since the objective function does not change if is postmultiplied by an arbitrary orthogonal matrix, the POD procedure seeks the optimal –dimensional subspace that captures the snapshots in the leastsquares sense. For more details on POD, we refer to [25, 28].
Note that Eq. (3.2) has more equations than unknowns (i.e., an overdetermined system). It is likely that there is no solution satisfying Eq. (3.2). In order to close the system, the Galerkin projection solves the following reduced system with :
(3.4) 
Applying a time integrator to Eq. (3.4) leads to a fully discretized reduced system, denoted as the reduced OE. Note that the reduced OE has unknowns and equations. If an implicit time integrator is applied, a Newton–type method can be applied to solve for unknown generalized coordinates each time step. If an explicit time integrator is applied, time marching updates will solve the system. However, we cannot expect any speedup because the size of the nonlinear term and its Jacobian, which need to be updated for every Newton step, scales with the FOM size. Thus, to overcome this issue, the DEIM approach applies a hyperreduction technique. That is, it projects onto a subspace and approximates as
(3.5) 
where , , denotes the nonlinear term basis matrix and denotes the generalized coordinates of the nonlinear term. The generalized coordinates, , can be determined by the following interpolation:
(3.6) 
where is the sampling matrix and is the
th column of the identity matrix
. Therefore, Eq. (3.5) becomes(3.7) 
where is the DEIM oblique projection matrix. The DEIM approach does not construct the sampling matrix . Instead, it maintains the sampling indices and corresponding rows of and . This enables DEIM to achieve a speedup when it is applied to nonlinear problems.
The original DEIM paper [12] constructs the nonlinear term basis by applying another POD on the nonlinear term snapshots^{2}^{2}2the nonlinear term snapshots are with the backward Euler time integrator and with the forward Euler time integrator obtained from the FOM simulation at every time step. This implies that DEIM requires two separate SVDs and storage for two different snapshots (i.e., solution state and nonlinear term snapshots). Section 4 discusses how to avoid the collection of nonlinear term snapshots and an extra SVD without losing the quality of the hyperreduction.
The sampling indices (i.e., ) can be found either by a row pivoted LU decomposition [12] or the strong column pivoted rankrevealing QR (sRRQR) decomposition [15, 16]. Depending on the algorithm of selecting the sampling indices, the DEIM projection error bound is determined. For example, the row pivoted LU decomposition in [12] results in the following error bound:
(3.8) 
where is the condition number of and it is bounded by
(3.9) 
On the other hand, the sRRQR factorization in [16] reveals tighter bound than (3.9):
(3.10) 
where is a tuning parameter in the sRRQR factorization (i.e., in Algorithm 4 of [23]).
3.2 Gnat
In contrast to DEIM, the GNAT method takes the LeastSquares Petrov–Galerkin (LSPG) approach. The LSPG method projects a fully discretized solution space onto a trial subspace. That is, it discretizes Eq. (2.1) in time domain and replaces with for in residual functions defined in Section 2 and Appendix A. Here, we consider only implicit time integrators because the LSPG projection is equivalent to the Galerkin projection when an explicit time integrator is used as shown in Section 5.1 in [9]. The residual functions for implicit time integrators are defined in (2.9), (A.8), and (A.12) for various time integrators. For example, the residual function with the backward Euler time integrator^{3}^{3}3Although the backward Euler time integrator is used extensively in the paper as an illustration purpose, many other time integrators introduced in Appendix A can be applied to all the ROM methods dealt in the paper in a straight forward way. after the trial subspace projection becomes
(3.11) 
The basis matrix can be found by the POD as in the DEIM approach. Note that Eq. (3.11) is an overdetermined system. To close the system and solve for the unknown generalized coordinates, , the LSPG ROM takes the square norm of the residual vector function and minimize it at every time step:
(3.12) 
The Gauss–Newton method with the starting point is applied to solve the minimization problem (3.12) in GNAT. However, as in the DEIM approach, a hyperreduction is required for a speedup due to the presence of the nonlinear residual vector function. The GNAT method approximates the nonlinear residual term with gappy POD [18], whose procedure is similar to DEIM, as
(3.13) 
where , , denotes the residual basis matrix and denotes the generalized coordinates of the nonlinear residual term. In contrast to DEIM, the GNAT method solves the following leastsquares problem to obtain the generalized coordinates at nth time step:
(3.14) 
where , , is the sampling matrix and is the th column of the identity matrix . The solution to Eq. (3.14) is given as
(3.15) 
where the Moore–Penrose inverse of a matrix with full column rank is defined as . Therefore, Eq. (3.13) becomes
(3.16) 
where is the GNAT oblique projection matrix. Note that it has a similar structure to . The GNAT projection matrix has a pseudoinverse instead of the inverse becuase it allows . The GNAT method does not construct the sampling matrix . Instead, it maintains the sampling indices and corresponding rows of and . This enables GNAT to achieve a speedup when it is applied to nonlinear problems.
The sampling indices (i.e., ) can be determined by Algorithm 3 of [11] for computational fluid dynamics problems and Algorithm 5 of [10] for other problems. These two algorithms take greedy procedure to minimize the error in the gappy reconstruction of the POD basis vectors . The major difference between these sampling algorithms and the ones for the DEIM method is that these algorithms for the GNAT method allows oversampling (i.e., ), resulting in solving leastsquares problems in the greedy procedure. These selection algorithms can be viewed as the extension of Algorithm 1 in [12] (i.e., a row pivoted LU decomposition) to the oversampling case. The nonlinaer residual term projection error associated with these sampling algorithms is presented in Appendix D of [11]. That is,
(3.17) 
where is a triangle matrix from QR factorization of (i.e., ).
We present another sampling selection algorithm that was not considered in any GNAT papers (e.g., [10, 11]). It is to use the sRRQR factorization developed originally in [23] and further utilized in the WDEIM method of [16] for the case of . That is, applying Algorithm 4 of [23] to produces an index selection operator whose projection error satisfies
(3.18) 
This error bound can be obtained by setting the identity matrix as a weight matrix in Theorem 4.8 of [16].
Finally, the GNAT method minimizes the following leastsquares problem at every time step, for example, with the backward Euler time integrator:
(3.19) 
The GNAT method applies another POD to a nonlinear residual term snapshots to construct . The original GNAT paper [11] collects residual snapshots at every Newton iteration from the LSPG simulations^{4}^{4}4We denote the LSPG simulation to be the ROM simulation without hyperreduction. That is, LSPG solves Eq. (3.12) without any hyperreduction if the backward Euler time integrator is used. for several reasons:

The GNAT method takes LSPG as a reference model (i.e., denoted as Model II in [11]). Its ultimate goal is to achieve the accuracy of Model II.

The residual snapshots taken every time step (i.e., at the end of Newton iterations at every time step) of the FOM are small in magnitude.

The residual snapshots taken from every Newton step gives information about the path that the Newton iterations in Model II take. GNAT tries to mimic the Newton path that LSPG takes.
Some disadvantages of the original GNAT approach include:

The GNAT method requires more storage than DEIM to store residual snapshots from every Newton iteration (cf., The DEIM approach stores only one nonlinear term snapshot per each time step).

The GNAT method requires more expensive SVD for nonlinear residual basis construction than the DEIM approach because the number of nonlinear residual snapshots in the GNAT method is larger than the number of nonlinear term snapshots in DEIM.

The GNAT method requires the simulations of Model II which are computationally expensive. For a parametric global ROM, it is computationally expensive because it requires as many Model II simulations as there are training points in a given parameter domain.
Section 4 discusses how to avoid the collection of nonlinear term snapshots and the extra SVD without losing the quality of the hyperreduction.
3.3 Space–time GNAT
The STGNAT method takes the space–time LSPG (STLSPG) approach. Given a time integrator, one can rewrite a fully discretized system of equations to Eq. (2.1) in a space–time form. For example, if the backward Euler time integrator^{5}^{5}5 Here we only consider the backward Euler time integrator for the simplicity of illustration. However, the extension to other linear multistep implicit time integrators is straight forward. is used for time domain discretization, then the corresponding space–time formulation to Eq. (2.6) becomes
(3.20) 
where
(3.21) 
Note that denotes the parameterized space–time state implicitly defined as the solution to the problem (3.20) with and . Further, denotes the space–time nonlinear term that is nonlinear in at least its first argument with . The space–time residual function corresponding to Eq. (3.20) is defined as
(3.22) 
To reduce both the spatial and temporal dimensions of the fullorder model, the STLSPG method enforces the approximated numerical solution to reside in an affine ‘space–time trial subspace’
(3.23) 
where with and whose elements are all one. Here, denotes the Kronecker product of two matrices and the product of two matrices and is defined as
(3.24) 
Further, denotes a space–time basis vector. Thus, the STLSPG method approximates the numerical solution as
(3.25) 
where , denotes the generalized coordinate of the STLSPG solution.
A space–time residual vector function can now be defined with the generalized coordinates as an argument. Replacing in Eq. (3.22) with gives the residual vector function
(3.26) 
where and denotes a space–time basis matrix that is defined as and denotes the generalized coordinate vector that is defined as
(3.27) 
Note that vanishes. Eq. (3.26) becomes
(3.28) 
Reduced space–time residual vector functions for other time integrators can be defined similarly. We denote as a reduced spate–time residual vector function with a generic time integrator.
The space–time basis matrix can be obtained by tensor product of spatial and temporal basis vectors. The spatial basis vectors can be obtained via POD as in the DEIM and GNAT approaches. The temporal basis vectors can be obtained via the following three tensor decompositions described in [13]:

Fixed temporal subspace via THOSVD

Fixed temporal subspace via STHOSVD

Tailored temporal subspace via STHOSVD
The STHOSVD method is a more efficient version of THOSVD. Thus, we will not consider THOSVD. The tailored temporal subspace via STHOSVD has a LL1 form that has appeared, for example, in [34, 14]. Therefore, we will refer to it as the LL1 decomposition.
Because of the reduction in spatial and temporal dimension, the space–time residual vector cannot achieve zero most likely. Thus, the STLSPG method minimizes the square norm of and computes the STLSPG solution:
(3.29) 
The STLSPG ROM solves Eq. (3.29) without any hyperreduction. As in the DEIM and GNAT approaches, a hyperreduction is required for a considerable speedup due to the presence of the space–time nonlinear residual vector. Therefore, the STGNAT method approximates the space–time nonlinear residual terms with gappy POD [18], which in turn requires construction of a space–time residual basis. Similar to the GNAT method, the STGNAT method approximates the space–time nonlinear residual term as
(3.30) 
where , , denotes the space–time residual basis matrix and denotes the generalized coordinates of the nonlinear residual term. The STGNAT solves the following space–time leastsquares problem to obtain the generalized coordinates, :
(3.31) 
where , , is the sampling matrix and is the th column of the identity matrix . The solution to Eq. (3.31) is given by
(3.32) 
Therefore, Eq. (3.30) becomes
(3.33) 
where is the STGNAT oblique projection matrix. Note that has the same structure as . The STGNAT method does not construct the sampling matrix . Instead, it maintains the sampling indices and corresponding rows of and . This enables the STGNAT to achieve a speedup when it is applied to nonlinear problems.
Section 5.3 in [13] discusses three different options to determine the sampling indices (i.e., ). However, all these three options are simple variations of Algorithm 3 in [11] and Algorithm 5 in [10]. They all minimize the error in the gappy reconstruction of the POD basis vectors for nonlinear space and time residuals. Therefore, the space–time nonlinear residual projection error due to (3.33) is similar to the ones in the GNAT method. That is,
(3.34) 
where is a triangle matrix from QR factorization of (i.e., ). On the other hand, one can also apply the sRRQR factorization in Algorithm 4 with a tuning parameter of [23] to to obtain that is associated with a tighter error bound for the projection error:
(3.35) 
This error bound can be obtained by setting the identity matrix as a weight matrix in Theorem 4.8 of [16].
Finally, the STGNAT solves the following leastsquares problem at every time step, for example, with the backward Euler time integrator:
(3.36) 
The original STGNAT paper introduces three different ways of collecting space–time residual snapshots, that are in turn used for the space–time residual basis construction (see Section 5.2 in [13]). Below is a list of the approaches introduced in [13] and explains advantages and disadvantages of each:

STLSPG ROM training iterations. This approach takes the space–time residual snapshot from every Newton iteration of the STLSPG simulations at training points in . This case leads to the number of residual snapshots, , where is the number of Newton iterations taken to solve the STLSPG simulation for the training point . This approach is not realistic for a largescale problem because it requires training simulations of the computationally expensive STLSPG ROM, where denotes the cardinality of the set . Furthermore, it requires the extra SVD on the residual snapshots, which is not necessary for our proposed method.

Projection of FOM training solutions. This approach takes the following steps:

take FOM state solution at every Newton iteration.

rearrange them in the space–time form (i.e., in Eq. (3.21)).^{6}^{6}6
Note that these are FOM solutions from time marching algorithms in which each time step results in different number of Newton iterations if implicit time integrators are used. Some time steps take a smaller number of Newton iterations than other time steps. However, in order to rearrange each Newton iterate state in the space–time form, we must have the same number of Newton iterations at each time step. Therefore, for the time steps that have converged with a smaller number of Newton iterations than other time steps, we pad the solution state vectors of the Newton iterations beyond the convergence with the converged solution. This only applies to an implicit time integrator because an explicit time integrator does not require any Newton solve.

project them onto the space–time subspace, (i.e., ).

compute the corresponding space–time residual (e.g., in Eq. (3.22) in the case of the backward Euler time integrator).

use those residuals as residual snapshots.
This approach simply requires projections and evaluations of the space–time residual. However, it requires the extra SVD on the residual snapshots, which is not necessary for our proposed method.


Random samples.
This approach generates a random space–time solution state samples (e.g., via Latin hypercube sampling or random sampling from uniform distribution) and computes the corresponding space–time residual (e.g.,
in Eq. (3.22) in the case of the backward Euler time integrator). This approach simply requires random sample generations and evaluations of the space–time residual. However, random samples are hardly correlated with actual data. Therefore, it is likely to generate poor space–time residual subspace. Furthermore, it requires the extra SVD on the residual snapshots, which is not necessary for our proposed method.
4 Solutionbased Nonlinear Subspace (SNS) method
Finally, we state our proposed method that avoids collecting nonlinear term snapshots and additional POD for the DEIM and GNAT approaches or additional tensor decomposition for the STGNAT method. We propose to use solution snapshots to construct nonlinear term basis in the DEIM, GNAT, and STGNAT approaches. A justification for using the solution snapshots comes from the subspace inclusion relation between the subspace spanned by the solution snapshots and the subspace spanned by the nonlinear term snapshots as shown in Eqs. (2.4) and (2.8) with the forward and backward Euler time integrators.^{7}^{7}7Subspace inclusion relations for other time integrators are shown in Appendix A.
4.1 DeimSns
We are going back to Eq. (3.2) and replace the nonlinear term with the approximation in Eq. (3.7):
(4.1) 
Eq. (4.1) is an overdetermined system, so it is likely that it will not have a solution. However, if there is a solution, then necessary conditions for Eq. (4.1) to have a nontrivial solution (i.e., ) are and
(4.2) 
The second condition says that the intersection of and needs to be nontrivial if there is a nontrivial solution to Eq. (4.1). Typically, we build first, using the POD approach explained in Section 3.1 and is set by . Therefore, the intersection of and can be controlled by the choice of we made. The larger the intersection of those two range spaces are, the more likely it is that there is a solution to Eq. (4.1). Given , the largest subspace intersection it can be is , i.e.,
(4.3) 
We call this condition as the conforming subspace condition. The conforming subspace condition leads to two obvious choices for :

The first choice is to ensure . If , then the range space of the left and righthand sides of Eq. (4.1) are the same. This leads Eq. (4.1) to become
(4.4) Because Eq. (4.4) is an overdetermined system and unlikely to have a solution, applying the Galerkin projection to Eq. (4.4) becomes:
(4.5) Assuming is invertible, Eq. (4.5) becomes:
(4.6) For the special case of being an identity matrix, Eq. (4.6) becomes:
(4.7) 
The second choice is to ensure . This can be achieved by taking an extended solution basis, with and . Then we set , which leads to . The obvious choice for is to take a larger truncation of the left singular matrix from SVD of the solution snapshot matrix than . Note that this particular choice of results in the first columns of being the same as (i.e., where with and ). By setting , Eq. (4.1) becomes
(4.8) Because Eq. (4.8) is unlikely to have a solution, applying the Galerkin projection to Eq. (4.8) becomes:
(4.9) (4.10) Assuming is invertible, Eq. (4.9) becomes:
(4.11) For the special case of being an identity matrix, Eq. (4.11) becomes
(4.12)
The DEIMSNS approach solves either Eq. (4.6) or Eq. (4.11) depending on the choice of above. Applying a time integrator to Eq. (4.6) or Eq. (4.12) leads to a reduced OE, whose solution, , can be lifted to find the full order approximate solution via .
Additionally, the subspace inclusion relations^{8}^{8}8 Eq. (2.4) of the forward Euler time integrator, Eq. (A.3) of the Adams–Moulton time integrator, and Eq. (A.16) of the midpoint RungeKutta method show the subspace inclusion relations for explicit time integrators. show that the subspace spanned by the solution snapshots include the subspace spanned by the nonlinear term snapshots. This fact further motivates the use of solution snapshots to build a nonlinear term basis. Indeed, numerical experiments show that the solution accuracy obtained by the DEIMSNS approach is comparable to the one obtained by the traditional DEIM approach. For example, see Figs. 3, and 4 in Section 6.1.1.
The obvious advantage of the DEIMSNS approach over DEIM is that no additional SVD or eigenvalue decomposition is required, which can save the computational cost of the offline phase.
Remark 4.1
Although the numerical experiments show that the two choices of above give promising results, the error analysis in Section 5 shows that the nonlinear term projection error bound increases by the condition number of (see Theorem 5). This is mainly because the orthogonality of is lost when or . This issue is resolved by orthogonalizing (e.g., apply economic QR factorization ) before using it as nonlinear term basis. Then apply, for example, Algorithm 4 of [23] to the transpose of orthogonalized one (e.g., ) to generate a sampling matrix . This procedure eliminates the condition number of in the error bound because (3.10) is valid.
Remark 4.2
4.2 GnatSns
The GNAT method needs to build a nonlinear residual term basis, , as described in Section 3.2. The nonlinear residual term is nothing more than linear combinations of the nonlinear term and the time derivative approximation as in Eq. (3.11). Thus, the subspace spanned by the nonlinear term residual snapshots are included in the subspace spanned by the solution snapshots. This motivates the use of the solution snapshots for the construction of a nonlinear residual term basis. Therefore, the same type of the nonlinear term basis in the DEIMSNS approach can be used to construct the nonlinear residual term basis in the GNATSNS method:

The first choice is to set . Thus, for example, the GNATSNS method solves the following leastsquares problem with the backward Euler time integrator:
(4.15) which can be manipulated to
(4.16) Note that the terms in norm in Eq. (4.16) is very similar to the discretized DEIM residual before Galerkin projection (i.e., applying the backward Euler time integrator to Eq. (4.4) gives the terms in norm in Eq. (4.16)). They are only different by the fact that one is an inverse and the other one is the Moore–Penrose inverse. In fact, if , then Eq. (4.16) is equivalent to applying the backward Euler time integrator to Eq. (4.4) and minimize the norm of the corresponding residual.
For the special case of being an identity matrix, Eq. (4.16) becomes
(4.17)
Comments
There are no comments yet.