An optimized finite element extrapolating method for 2D viscoelastic wave equation

In this study, we first present a classical finite element (FE) method for a two-dimensional (2D) viscoelastic wave equation and analyze the existence, stability, and convergence of the FE solutions. Then we establish an optimized FE extrapolating (OFEE) method based on a proper orthogonal decomposition (POD) method for the 2D viscoelastic wave equation and analyze the existence, stability, and convergence of the OFEE solutions and furnish the implement procedure of the OFEE method. Finally, we furnish a numerical example to verify that the numerical computing results correspond with the theoretical ones. This signifies that the OFEE method is feasible and efficient for solving the 2D viscoelastic wave equation.


Introduction
Let ⊂ R  be a bounded convex polygonal domain with a smooth boundary ∂ . We consider the following initial-boundary value problem: u(x, y, t) = ϕ(x, y, t), (x, y, t) ∈ ∂ × (, T], u(x, y, ) = ϕ  (x, y), u t (x, y, ) = ϕ  (x, y), (x, y) ∈ , where u tt = ∂  u/∂t  , u t = ∂u/∂t, and ε and γ are two positive constants, f (x, y, t), ϕ(x, y, t), and ϕ  (x, y) and ϕ  (x, y) are, respectively, the source term, the boundary value function, and the initial value functions, sufficiently smooth to ensure the validity of the following analysis, and T is the time duration. As a matter of convenience, we assume that ϕ(x, y, t) =  and ε = γ =  in the remaining part of the article.
Problem  is referred to as a system of viscoelastic wave equations. It has some special and significant physical backgrounds. For instance, it can be used to describe the wave propagation phenomena of actual vibration through a viscoelastic medium (see, e.g., [, ]). Although the existence and uniqueness of its analytic solution have been proved (see, e.g., [-]), because the viscoelastic wave equation in the real-world engineering applications usually has complex known data or computed domains, the analytical solution cannot be generally solved, so one has to find its solutions numerically. For more than  years, it has been attentively studied and many numerical methods for the viscoelastic wave equation have been developed (see, e.g., [-]). Among all numerical methods, the finite element (FE) method is considered to be one of the calculating numerical methods with the best theory for the two-dimensional (D) viscoelastic wave equation (see [, ]). Nevertheless, the classical FE methods for the D viscoelastic wave equation are some macroscale systems of equations including lots of unknowns, i.e., degrees of freedom, so entail very large computational load in real-world engineering applications. As a consequence, an important issue is how to greatly lessen the number of unknowns of the classical FE methods to reduce the computational load, ease the truncated error amassing, and save CPU time in the numerical computation, while preserving the desired FE solution accuracy.
It has been proved by lots of numerical studies (see, e.g., [-]) that the proper orthogonal decomposition (POD) method is a very useful tool to reduce the number of unknowns for numerical models and ease the truncated error amassing in numerical calculations. But most existing reduced-order models, as mentioned, were established via the POD basis formed from the classical numerical solutions at all time nodes, before repetitively computing the reduced-order numerical solutions at the same time nodes, which were some valueless repetitive calculations. Since , some reduced-order FE extrapolating methods based on the POD method for partial differential equations have been established successively by Luo's team (see, e.g., [-]) in order to avert the valueless repeated computations.
However, as far as we know, there has not been any report that the POD method is used to reduce the number of unknowns in the classical FE method for the D viscoelastic wave equation. Therefore, in this article, we devote ourselves to building an optimized FE extrapolating (OFEE) method that includes very few unknowns but maintains desired accuracy via the POD method, analyzing the existence, stability, and convergence of the OFEE solutions and verifying the efficiency and feasibility of the OFEE method by some numerical simulations.
The main distinctions between the OFEE method and the other existing reduced-order FE extrapolating methods built on the POD method (see, e.g., [-]) consist in the fact that the viscoelastic wave equation not only contains three second-order derivative terms of time and of spatial variables but also includes two mixed derivative terms of time (first-order) and spatial variables (second-order) so that either the modeling of the OFEE method or the demonstration of the existence, stability, and convergence of the OFEE solutions faces more difficulties and requires more techniques than the existing other aforementioned reduced-order FE extrapolating methods. However, the OFEE method has some specific applications. Though an optimized splitting positive definite mixed FE extrapolation (OSPDMFEE) model based on the POD technique for the D viscoelastic wave equation is developed in [], it has three unknown functions and the OSPDMFEE model has more degrees of freedom than the current OFEE format, so that its theoretical analysis and numerical simulations have more difficulties than the current OFEE method. It is worth mentioning that we can discuss the existence, stability, and convergence of the reduced-order FE solutions by means of the classical FE theory. Especially, the OFEE method only employs the classical FE solutions at the initial very few time nodes to for-mulate the POD basis and build the OFEE format so that it does not have repetitive calculations, such as done in references [-]. Consequently, it is a development and an improvement of the existing aforementioned ones (see, e.g., [-]).
The remaining content of the article is organized as follows. In Section , we first present the classical FE method for the D viscoelastic wave equation and analyze the existence, stability, and convergence of the classical FE solutions. In Section , we develop the OFEE method via the POD method for the D viscoelastic wave equation, analyze the stability and convergence of the OFEE solutions, and furnish the implement procedure of the OFEE method. Next, in Section , we use some numerical simulations to verify the efficiency and feasibility of the OFEE method. Finally, in Section , we summarize our main conclusions.

The classical FE method for the 2D viscoelastic wave equation 2.1 Generalized solution for the 2D viscoelastic wave equation
The following arisen Sobolev spaces as well as their norms are well known (see []).
For convenience, we write U = H   ( ). Thus, by using Green's formula for the D viscoelastic wave equation, we obtain the following variational formulation: where (·, ·) denotes the inner product of L  ( ).
For U, we have the following Poincaré inequality: where β is a positive real.
For Problem , we have the following result.
where β is the constant in the Poincaré inequality.
Proof Because Problem  is a system of linear equations as regards the unknown function u, in order to prove the existence and uniqueness of solutions for Problem , it is necessary to prove that Problem  has only the zero solution when f (x, Thus, by the Hölder inequality, the Poincaré inequality, and the Cauchy-Schwarz inequality, we acquire Then Problem  has a unique solution such that inequality () holds.

Semi-discrete format as regards time for the 2D viscoelastic wave equation
Let N be a positive integer, t = T/N the time step size, and t i = i t. If we use (u n+u n )/ ( t) to approximate u t and (u n+ -u n + u n- )/ t  to approximate u tt for the D viscoelastic wave equation, we obtain the following semi-discrete formulation of time: where f n = f (t n ).
For Problem , we have the following.

Theorem  Under the assumptions of Theorem
showing that the series of solutions to Problem  is stable and continuously dependent on the source function f and the initial values ϕ  and ϕ  . When u is sufficiently smooth in t, we have the following error estimations: Proof Because Problem  is a system of linear equations as regards the unknown function u n , in order to prove the existence and uniqueness of solutions for Problem , it is necessary to prove that Problem  has only the set of zero solutions when f (x, y, t) = ϕ  (x, y) = ϕ  (x, y) = .
By taking v = u n+u n- in () and using the Hölder, Poincaré, and Cauchy-Schwarz inequalities, we have By summing () from  to n and using (), we obtain Thus, when f (x, y, t) = ϕ  (x, y) = ϕ  (x, y) = , from (), we obtain ∇u n  = , implying u n = . Hence Problem  has a unique solution series.
From (), we obtain which is just the inequality (). Let e n = u(t n )u n . By applying the Taylor expansion formula to () and then subtracting () taking t = t n , we obtain where β is the same constant as in the Poincaré inequality. Because e  = , e  =  (when t is sufficiently small), by summing () from  to n, we obtain This finishes the proof of Theorem .

Classical fully discrete FE method for the 2D viscoelastic wave equation
Let h be a regular triangulation of¯ . The FE subspace U h is taken as where P k (K) is the subspace formed by kth degree polynomials on K and k ≥  is an integer. Thus, the fully discrete FE formulation for the D viscoelastic wave equation () is as follows: where f n = f (t n ) and R h is the Ritz projection as follows: For Problem , we have the following.
Consequently, the solution sequence u n h to Problem  is stable and continuously dependent on the source function f and the initial values ϕ  and ϕ  . With h = O( t), we have the following error estimations: where C is a positive constant only dependent on u, but independent of the time step t and spatial mesh parameters h.
Proof (i) The existence and uniqueness of the solution sequence for Problem . Let Then Problem  can be rewritten as follows: It is obvious that, for given u n h and u n- h as well as f n (n = , , . . . , N -), F(v h ) is a bounded linear functional of v h and a(u, v) is a bilinear functional of u and v. Because u  ≤ u  and ∇u  ≤ u  , by using the Hölder inequality, we have By summing () from  to n and using (), again the Poincaré inequality, and the properties of the Ritz projection R h , we obtain From (), we immediately obtain ().
(iii) Convergence of the solution sequence for Problem . Letẽ n = u nu n h , E n = R h u nu n h , and ρ n = u n -R h u n . By subtracting Problem  from Problem , taking v = v h ∈ U h , we obtain the following system of the error equations:

By () and the properties of the Ritz projection R h , when h = O( t), we have
By summing () from  to n, we obtain From (), by the properties of the Ritz projection and Theorem , we immediately obtain ().
Remark  The full FE formulation Problem  is directly built from the semi-discrete formulation Problem  with respect to time such that one can bypass the semi-discrete formulation with respect to spatial variables and its theoretical analysis becomes simpler. Thus, as long as f (x, y, t), ϕ  (x, y), ϕ  (x, y), ε, γ , time step k, the spatial mesh size h, and the FE subspace U h are assigned, we attain the solution sequence {u n h } N n= ⊂ U h by solving Problem . We take the subsequence {u n h } L n= from the initial L solutions of {u n h } N n= as snapshots (in general, L N and √ L < , for example, L = , N = ).

The OFEE format for the 2D viscoelastic wave equation 3.1 Formulations of the POD basis and establishment the OFEE format
Let W n (x, y) = u n h (x, y) ( ≤ n ≤ L), at least one of which is supposed to be a non-zero function, and l = dim{W  , W  , . . . , W L }. Write A = (A ik ) L×L and A ik = (∇W i (x, y), ∇W k (x, y))/L. Since the matrix A is a non-negative Hermitian matrix with rank l, it has a complete set of orthonormal eigenvectors with corresponding eigenvalues λ  ≥ λ  ≥ · · · ≥ λ L > . Thus, the POD basis {ψ  , ψ  , . . . , ψ L } is given by (see []) holding the following property (see also []).

Proposition  The following estimation holds:
Then, by functional analysis (see []), there exists an extension R h : where u ∈ U. Due to (), the operator R h is bounded. We have Further, the following holds.

Lemma  For every d ( ≤ d ≤ l), the Ritz-operator R d in () satisfies
where u i h ∈ V (i = , , . . . , L) are the solutions to Problem . Further, if u ∈ H  ( ) is the solution to Problem , the extended Ritz-operator R h defined by () satisfies the following error estimations: Thus, by means of U d , the OFEE format for the D viscoelastic wave equation is described as follows: where u n h (n = , , . . . , L) are the first L solutions for Problem .
Remark  It is easily seen that Problem  at each time node includes N h unknowns (where N h is the number of vertices of triangles in h ), whereas Problem  at the same time node contains only d unknowns (d l ≤ L N N h ). For real-world engineering problems, the number N h of vertices of triangles in h can easily reach a few millions, while d is only the number of the major eigenvalues and is very small (for example, in Section , d = , but N h ≥  ×   ). Problem  here is the OFEE format for the D viscoelastic wave equation. In particular, Problem  employs only the initial few known L solutions of Problem  used to extrapolate other N -L solutions, and has no repetitive computations. The first L OFEE solutions are obtained by projecting the first L classical FE solutions into the POD basis, while the other remaining (N -L) OFEE solutions are obtained by extrapolating and iterating equation (). Therefore, it is completely different from the existing POD-based reduced-order formulations.

The error estimations of the OFEE solutions
In the following, we employ the classical FE method to deduce the error estimations of OFEE solutions for the D viscoelastic wave equation. We have the following main result.
Theorem  Under the same conditions as Theorem , Problem  has a unique solution sequence {u n h } N n= ⊂ U satisfying As a consequence, the sequence of solutions u n d to Problem  is stable and continuously dependent on the source function f and the initial values ϕ  and ϕ  . As h = O( t), we have the following error estimations: Proof ( Thus, () in Problem  can be rewritten as follows: Seek u n d ∈ U h (n = L + , L + , . . . , N ) satisfying  a(u, v) is a bilinear functional of u and v. Because u  ≤ u  and ∇u  ≤ u  , by using the Hölder inequality, we have By summing () from L +  to n, and by () and () in Lemma , we obtain When h = O( t), from () and by the properties of the Ritz projection and Theorem , we readily obtain the case of () when n = L + , L + , . . . , N .
Remark  We make some comments on Theorem : () It is known from Theorem  that, in order to not adversely affect accuracy, it is necessary to take L as L N , for example, we usually take L such that √ L < . Thus, it is unnecessary to extract total transient solutions at all time nodal points t n as snapshots such as done in [, ]. () The error (L l j=d+ λ j ) / in Theorem  gives some indication as to how to choose the number d of the POD basis, namely, it is only necessary to meet

The implement procedure of the OFEE format
Solving the OFEE format, i.e., Problem , requires the following seven steps: Step . For given ε and γ , boundary value function ϕ(x, y, t), initial value function ϕ  (x, y), and ϕ  (x, y), source term f (x, y, t), the time step size t, and the spatial grid measurement h satisfying h = O( t) solve the following classical FM formulation on the first L ( √ L < ) steps: where f n = f (t n ) and R h is the Ritz projection. This yields the snapshots W i = u n h (n = , , . . . , L).
Step . Formulate the snapshot matrix A = (A ij ) L×L , where A ij = (∇u i h , ∇u j h ) and (·, ·) is the L  -inner product.
Step . Find the eigenvalues λ  ≥ λ  ≥ · · · ≥ λ l >  (l = dim{u n h :  ≤ n ≤ L}) of A and the corresponding eigenvectors v j = (a  (n = , , . . . , , i.e., at time t = ., ., . . . , .) for Problem  (the classical FE formulation) to constitute a set of snapshots. By computing, with d =  and k =  - , we achieve the error estimation (  j= λ j ) / ≤  ×  - in Theorem , which shows that we only need to take six POD bases. Thus, the OFEE format (Problem ) at each time level has only  degrees of freedom, while the classical FE formulation (Problem ) contains more than  ×   degrees of freedom. Therefore, the OFEE format (Problem ) cannot only alleviate the computational load and save time-consuming calculations in the computational process, but also reduce the accumulation of truncation errors in the computational process. When we solve the OFEE format (Problem ) with six optimal POD bases, according to the seven steps of implementation of the OFEE format in Section ., we find that the OFEE format at t =  is still convergent, without the need to renew the POD basis. The OFEE solution obtained with the OFEE format (Problem ) is depicted graphically in Figures  and . The images in Figures  and  look very much alike, and so do those in Figures  and . Nevertheless, the OFEE solutions are probably better than the classical FE solutions due to the little accumulation of truncated errors of the OFEE format (Problem ) in the computational process. Figure  shows the absolute error between  solutions u n d of the OFEE format (Problem ) with  different numbers of POD bases and the solutions u n h of the classical FE formulation (Problem ) at t = . It shows that, when the numbers of the POD basis are larger than five, the error does not exceed  ×  - . Therefore, the error results in the numerical example above are concordant with those obtained with the theoretical approach.  This has shown that the OFEE format is feasible and efficient for solving the viscoelastic wave equation.

Conclusions
In this article, we use the POD technique to build the OFEE format for the D viscoelastic wave equation. We first extract snapshots from the initial few L (L N ) classical FE solutions for the D viscoelastic wave equation. Next, we constitute the POD basis of snapshots by means of the POD method. Then the FE subspaces of the classical FE format are replaced with the subspaces spanning the most main POD bases to build the OFEE formulation for the D time-dependent conduction-convection problem. Finally, we deduce the existence, uniqueness, stability, and convergence of the OFEE solutions of the D viscoelastic wave equation and furnish the implement procedure for the OFEE format.
Comparing the numerical simulation errors with the theoretical errors we have verified that the theoretical errors are concordant with the computing errors, thus validating both the feasibility and efficiency of the OFEE format.