- Email: [email protected]
A Reduced-Order MFE Formulation Based on POD Method for Parabolic Equations
Acta Mathematica Scientia 2013,33B(5):1471–1484 http://actams.wipm.ac.cn
A REDUCED-ORDER MFE FORMULATION BASED ON POD METHOD FOR PARABOLIC EQUATIONS∗
Zhendong LUO (
)†
School of Mathematics and Physics, North China Electric Power University, Beijing 102206, China E-mail: [email protected]
)
Lei LI (
)
Ping SUN (
School of Mathematics and Computer Science, Guizhou Normal University, Guiyang 550001, China E-mail: [email protected]; [email protected]
Abstract In this paper, we extend the applications of proper orthogonal decomposition (POD) method, i.e., apply POD method to a mixed finite element (MFE) formulation naturally satisfied Brezz-Babu˘ska for parabolic equations, establish a reduced-order MFE formulation with lower dimensions and sufficiently high accuracy, and provide the error estimates between the reduced-order POD MFE solutions and the classical MFE solutions and the implementation of algorithm for solving reduced-order MFE formulation. Some numerical examples illustrate the fact that the results of numerical computation are consistent with theoretical conclusions. Moreover, it is shown that the new reduced-order MFE formulation based on POD method is feasible and efficient for solving MFE formulation for parabolic equations. Key words proper orthogonal decomposition method; mixed finite element formulation; parabolic equation; error estimate 2010 MR Subject Classification
1
65N30; 65M30
Introduction
Let Ω ⊂ R2 be a bounded and connected polygonal domain. Consider the following parabolic equations. Problem (I) Find u such that
∗ Received
(x, y, t) ∈ Ω × (0, T ), ut − △u = f, u(x, y, t) = 0, (x, y, t) ∈ ∂Ω × [0, T ), u(x, y, 0) = u0 (x, y), (x, y) ∈ Ω,
(1.1)
May 26, 2011; revised October 31, 2012. The research was supported by the National Science Foundation of China (11271127 and 11061009), Science Research Program of Guizhou (GJ[2011]2367), and the Co-Construction Project of Beijing Municipal Commission of Education. † Correspondence author: Zhendong LUO.
1472
ACTA MATHEMATICA SCIENTIA
Vol.33 Ser.B
where u represents unknown function, T the total time, f the given source term, and u0 (x, y) the given initial function. Let p = −∇u, then Problem (I) can be written as follows. Problem (II)
Find (u, p) such that
ut + divp = f, p + ∇u = 0, u(x, y, t) = 0, u(x, y, 0) = u0 (x, y),
(x, y, t) ∈ Ω × (0, T ), (x, y, t) ∈ Ω × (0, T ),
(x, y, t) ∈ ∂Ω × [0, T ),
(1.2)
(x, y) ∈ Ω.
The parabolic equations may be used to describe many physical phenomena of the natural environment, engineering equipment, and living organisms, such as the proliferation of gas, the infiltration of liquid, the conduction of heat, and the spread of impurities in semiconductor materials. It isn’t easy to find the exact solutions for practical engineering parabolic problems if their computational fields are usually irregular. Thus, finding their numerical solutions becomes an efficient approach. The mixed finite element (MFE) method is one of the most effective methods to find the numerical solutions for parabolic equations since it can simultaneously find unknown function (displacement) and gradient (stress), which has been widely applied in scientific computations and its theory has been established (see [1, 2]). However, the classical MFE formulations for parabolic equations include too many degrees of freedom and need to satisfy an important convergence stability condition, i.e., the Brezzi-Babuˇska (BB) inequality (see [1–3]) so that they bring many difficulties for practical applications. We adopt the same MFE formulation as that introduced in [4] which can circumvent the constraint of the BB inequality. Thus, the key problem is how to alleviate the computational load by saving time-consuming calculations in the computational process in a way that guarantees sufficiently accurate numerical solutions. It was shown that a proper orthogonal decomposition (POD) method by combining with some numerical methods for partial differential equations can provide an efficient means of generating reduced-order models and alleviating the computational load and memory requirements (see [5]). POD method was widely and successfully applied to numerous fields, including signal analysis and pattern recognition (see [6]), statistics (see [7]), geophysical fluid dynamics and meteorology (also see [7]). The POD method essentially provides an orthogonal basis for representing the given data in a certain least squares optimal sense, that is, it provides a way to find optimal lower dimensional approximations of the given data. In early times, POD method was mainly used to perform principal component analysis in computations of statistics (see [5–7]), until the method of snapshots was introduced by Sirovich (see [8]) and was then widely applied to reducing the order of the POD eigenvalue problem and searching the main behavior of a dynamic system (see [8–13]). Until ten years ago, some Galerkin POD methods for parabolic problems and a general equation in fluid dynamics were not presented (see [14, 15]). More recently, some reduced finite difference schemes and finite element formulations based on POD technique for parabolic problems were presented (see [16– 18]). Moreover, some reduced-order methods based on POD technique for four-dimensional variational data assimilation were established (see [19–21]).
No.5
Z.D. Luo et al: A REDUCED-ORDER MFE FORMULATION BASED ON POD METHOD
1473
Though some Galerkin methods, finite difference schemes, and finite element formulations for parabolic problems were reduced with POD method (see [14–18]), to the best of our knowledge, there are no published results addressing the case that the classical MFE formulation for parabolic equations is reduced with the POD method or providing error estimates between the classical MFE solutions and the reduced-order MFE solutions or the implementation of algorithm for solving the reduced-order MFE formulation. In this paper, we extend the developments in [14–18], i.e., combine the classical MFE method with POD method to establish a reduced-order MFE formulation with lower dimensions and sufficiently high accuracy for two-dimensional parabolic equations, analyze the errors between the reduced-order MFE solutions based on POD method and the classical MFE solutions, and provide the implementation of algorithm for solving reduced-order MFE formulation. Some numerical examples illustrate the fact that the results of numerical computations are consistent with theoretical conclusions. Moreover, it is shown that the reduced-order MFE formulation based on POD method is feasible and efficient for finding numerical solutions of two-dimensional parabolic equations. The main advantage in this paper consists in that it does only use the first fewer given numerical solutions on very short time span [0, T0 ] (T0 ≪ T ) as snapshots to construct POD basis and establish the reduced-order MFE formulation based on POD technique for finding the numerical solutions on total time span [0, T ]. Thus, it can sufficiently adopt the advantage of POD method, namely it utilizes the given data (on very short time span [0, T0 ] and T0 ≪ T ) to forecast (or infer) future physic phenomenons (on time span [T0 , T ]). This is the improvement and innovation for the existing reduced-order methods based on POD technique (see, e.g., [14–18]). The rest of this paper is organized as follows. Section 2 is to recall the classical fully discrete MFE formulation for two-dimensional parabolic equations introduced in [4] and to generate snapshots from the first few transient solutions computed from the equation system derived by the classical fully discrete MFE formulation. In Section 3, the optimal orthonormal POD bases are spanned by the elements of the snapshots with POD method and a reducedorder MFE formulation with lower dimensions and sufficiently high accuracy based on POD method for two-dimensional parabolic equations is established. In Section 4, the error estimates between classical MFE solutions and the reduced-order MFE solutions based on POD method and the implementation of algorithm for solving the reduced-order MFE formulation are provided. In Section 5, some numerical examples are presented for illustrating that the numerical computing errors between the reduced-order MFE solutions based on POD method and the classical MFE solutions are consistent with previously obtained theoretical results, thus validating the feasibility and efficiency of the reduced-order MFE formulation. Section 6 provides main conclusions and future tentative ideas.
2
Recall Classical MFE Formulation and Generate Snapshots
The Sobolev spaces used in this context are standard (see [22]). Let M = H01 (Ω), H = L2 (Ω)2 . For the sake of convenience, without loss of generality, we may as well assume that u0 (x, y) = 0 in the following theoretical analysis. Then a mixed variational formulation for Problem (II) is written as follows.
1474
ACTA MATHEMATICA SCIENTIA
Vol.33 Ser.B
Problem (III) Find (u, p) ∈ M × H such that, for any t ∈ (0, T ), a(p, q) + b(q, u) = 0,
b(p, v) − (ut , v) = −(f, v), u(x, y, 0) = 0,
∀q ∈ H,
(2.1)
∀v ∈ M, (x, y) ∈ Ω,
where a(p, q) = (p, q), b(q, v) = (q, ∇v), (·, ·) represents the inner product in L2 (Ω) or L2 (Ω)2 . Since a(·, ·) and b(·, ·) all are continuous bilinear functionals, a(·, ·) is coercive on H × H, and b(·, ·) satisfies continuous BB condition (as are proved in [4]), i.e., for any q ∈ H, b(·, ·) satisfies b(q, v) sup > βkvk1 , ∀v ∈ M, (2.2) q∈H kqk0 where β is a constant, if f ∈ H s (0, T ; H r (Ω)), the mixed variational formulation Problem (III) for Problem (II) exists a unique solution (u, p) ∈ M × H satisfying that: kp(i) kr+1 + ku(i+1) kr+2 6 Ckf (i) kr , ∀t ∈ (0, T ), i = 1, 2, · · · , s,
(2.3)
where C is a positive constant, f (i) represent the i-th order derivative with respect to time t. In order to obtain the MFE solutions of Problem (III), it is necessary to introduce an MFE approximation for the spatial variables of Problem (III) and an approximation of difference quotient for time derivative. ¯ (see [2, 3, 23, 24]). MFE Let {ℑh } be a uniformly regular family of triangulation of Ω spaces are taken as follows. Hh = {q h ∈ H; q h |K ∈ Pm−1 (K)2 , ∀K ∈ ℑh },
Mh = {vh ∈ M ∩ C 0 (Ω); vh |K ∈ Pm (K), ∀K ∈ ℑh }, where m > 1 and Pm (K) is the space of polynomials of degree 6 m on K. For any positive integer N , let k = T /N denote time–step and write tn = nk (0 6 n 6 N ). The time derivative for Problem (III) are approximated by one–step Euler backward difference quotient ∂¯t un = (un − un−1 )/τ . Let (unh , pnh ) be the fully discrete MFE approximation of (u, p). Then, the fully discrete MFE formulation for Problem (III) may be written as follows. Problem (IV) Find (unh , pnh ) ∈ Mh × Hh such that, for 1 6 n 6 N , n n a(ph , q h ) + b(q h , uh ) = 0, τ b(pnh , vh )
−
(unh , vh )
=
−(un−1 , vh ) h
u0h (x, y) = 0, (x, y) ∈ Ω.
n
− τ (f , vh ),
∀q h ∈ Hh ,
∀vh ∈ Mh ,
(2.4)
Put A((u, p); (v, q)) = τ a(p, q) + τ b(q, u) − τ b(p, v) + (u, v) and F (v, q) = (un−1 , v) + h τ (f n , v). Problem (IV) is equivalently rewritten as follows. Problem (V) Find (unh , pnh ) ∈ Mh × Hh such that, for 1 6 n 6 N , A((un , pn ); (v , q )) = F (v , q ), ∀(vh , q h ) ∈ Mh × Hh , h h h h h h (2.5) u0 (x, y) = 0, (x, y) ∈ Ω. h
No.5
Z.D. Luo et al: A REDUCED-ORDER MFE FORMULATION BASED ON POD METHOD
1475
It has been proved that Problem (V) has a unique set of solutions (unh , pnh ) ∈ Mh × Hh (n = 1, 2, · · · , N ) and the following theorem of error estimates (see [4]). n P Theorem 1 If f ∈ H 1 (0, T ; L2(Ω)), kf i k2m−1 is a bounded constant, and τ = O(h), i=1
then errors between the solution (u, p) to Problem II and solutions (unh , pnh ) are as follows: ku(tn ) −
unh k0
+τ
1 2
X n i=1
kp(ti ) −
pih k20
21
1
6 C(τ 1+ 2 + hm+1 ), 1 6 n 6 N.
(2.6)
If source term f (x, y, t), initial conditions u0 (x, y), triangulation parameter h, the time– step k, finite element subspaces Xh and Mh are given, we can obtain a ensemble of solutions {(unh (x, y), pnh (x, y))}N n=1 by solving Problem (IV). And then we choose the first L (in general, L ≪ N , for example, L = 20, N = 200) instantaneous solutions (uih (x, y), pih (x, y)) (1 6 i 6 L ≪ N ) from the N group of solutions {(unh (x, y), pnh (x, y))}N n=1 for Problem (IV), which are useful and of interest for us and referred to as snapshots. Remark 1 When one computes actual problems, one may obtain the ensemble of snapshots from physical system trajectories by drawing samples from experiments and interpolation (or data assimilation).
3 Generation of POD Basis and a Reduced-Order MFE Formulation Based on POD Technique for Parabolic Equations For (uih (x, y), pih (x, y)) (1 6 i 6 L ≪ N ) in Section 2, let U i (x, y) = (uih (x, y), pih (x, y)) (1 6 i 6 L) and V = span{U 1 , U 2 , · · · , U L }, (3.1) and refer to V as the space generating by the snapshots {U i }L i=1 at least one of which is assumed to be non-zero. Let {ψ j }lj=1 denote an orthonormal basis of V with l = dim V (l 6 L). Then each member of the ensemble can expressed as Ui =
l X j=1
(U i , ψ j )H×M ψ j , i = 1, 2, · · · , L,
(3.2)
where (U i , ψ j )H×M = (∇uih , ∇ψuj ) + (pih , ψ pj ), ψ j = (ψuj , ψ pj ), and ψuj and ψ pj are orthonormal basis corresponding to u and p, respectively. Definition 2 The method of POD consists in finding an orthonormal basis sequence ψ i (i = 1, 2, · · · , l) such that for every d (1 6 d 6 l) the mean square error between the elements U i (1 6 i 6 L) and corresponding d-th partial sum of (3.2) is minimized on average min
{ψ j }d j=1
such that
2 L d X
1 X
U i − (U i , ψ j )H×M ψ j
L i=1 H×M j=1
(ψ i , ψj )H×M = δij , 1 6 i 6 d, 1 6 j 6 i,
(3.3)
(3.4)
where kU i k2H×M = k∇uih k20 + kpih k20 . A solution {ψ j }dj=1 of (3.3) and (3.4) is known as a POD basis of rank d.
1476
ACTA MATHEMATICA SCIENTIA
Vol.33 Ser.B
We introduce the correlation matrix A = (Aij )L×L ∈ RL×L corresponding to the snapshots {U i }L i=1 by 1 Aij = (U i , U j )H×M . (3.5) L Since the matrix A is positive semi-definite and has rank l, the solution of (3.3) and (3.4) can be found, in addition, the following result holds (for example, see [7, 8, 14, 15]). Proposition 2 Let λ1 > λ2 > · · · > λl > 0 denote the positive eigenvalues of A and v 1 , v 2 , · · ·, v l the associated orthonormal eigenvectors. Then a POD basis of rank d 6 l is given by L 1 X √ ψi = (v i )j U j , 1 6 i 6 d, (3.6) Lλi j=1 where (v i )j ’s denote the j-th component of the eigenvector v i . Furthermore, the following error formula holds
2 L d l X X
1 X
U i −
(U , ψ ) ψ = λj . (3.7) i j H×M j L i=1 H×M j=1 j=d+1
Let M d = span {ψu1 , ψu2 , · · · , ψud } and H d = span ψ p1 , ψ p2 , · · · , ψ pd . Due to for any vd ∈ M d satisfying ∇vd ∈ H d , we define the projection P d : M h → M d denoted by, for uh ∈ Mh , (∇P d uh , q d ) = (∇uh , q d ), ∀q d ∈ H d . Then there is an extension P h : M → Mh of P d such that P h |Mh = P d : Mh → M d denoted by (see [25]), for any u ∈ M , (∇P h u, qh ) = (∇u, q h ),
∀q h ∈ Hh .
(3.8)
Set the L-projection ρd : H → H d denoted by, for any p ∈ H, (ρd p, q d ) = (p, q d ),
∀q d ∈ H d .
(3.9)
Due to (3.8) and (3.9), the operators P h and ρd are well-defined and bounded (see [15, 18]) k∇(P d u)k0 6 k∇uk1 , ku − P h uk0 6 Chk∇(u − P h u)k0 , ∀u ∈ M
(3.10)
kρd pk0 6 CkpkM , ∀p ∈ H.
(3.11)
and
Moreover, the following results hold (see [14, 15, 17, 18]). Lemma 3 For every d (1 6 d 6 l) the projection operators P d and ρd satisfy L l X 1X i [kuh − P d uih k20 + h2 k∇(uih − P d uih )k20 ] 6 Ch2 λj , L i=1
(3.12)
j=d+1
L l X 1X i kph − ρd pih k20 6 λj , L i=1 j=d+1
where (uih , pih ) ∈ V is the solution sequence of Problem (IV).
(3.13)
No.5
Z.D. Luo et al: A REDUCED-ORDER MFE FORMULATION BASED ON POD METHOD
1477
Thus, substituting H d and M d for Hh and Mh for Problem (IV), respectively, we can obtain a new reduced-order MFE formulation based on POD technique for Problem (I) as follows. Problem (VI) Find (und , pnd ) ∈ M d × H d such that n n a(pd , q d ) + b(q d , ud ) = 0, τ b(pnd , vd ) − (und , vd ) = −(un−1 , vd ) − τ (f n , vd ), d 0 ud (x, y) = 0, (x, y) ∈ Ω.
∀q d ∈ H d ,
∀vd ∈ M d ,
(3.14)
By using the same approaches as proofs of Theorem 3 in [4], we may prove that Problem (VI) exists a unique solution sequence (und , pnd ) ∈ M d × H d (n = 1, 2, · · · , N ) satisfying kund k20
+τ
n X i=1
kpid k20
n X 0 2 i 2 6 kud k0 + 2τ kf k0 exp(2nτ ),
(3.15)
i=1
which has shown that the solution sequence of Problem (VI) is stable and continuously dependent on source term f (x, y, t) and initial conditions u0 (x, y). Remark 2 If ℑh is a uniformly regular triangulation and Mh and Hh are respectively the spaces of piecewise second order polynomials and piecewise linear polynomial vectors, the number of total degrees of freedom for Problem (VI), i.e., the number of unknown quantities is 3Np +Ns ≈ 5Np (where Np and Ns are respectively the number of vertices and sides of triangles in ℑh ), while the number of total degrees of freedom for Problem (VI) is 3d (d ≪ l 6 L ≪ N ). For scientific engineering problem, the number Np of vertices of triangles in ℑh is more than ten of thousands, even more than a hundred million, while d is only the number of the first few maximal eigenvalues which is chosen the first L snapshots from the N snapshots so that it is very small (for example, in Section 5, d = 6, while Np = 2000 × 2000 = 4 × 106 ). Therefore, Problem (VI) is a reduced-order MFE formulation based on POD method for Problem (IV).
4 Error Analysis of Solutions and Implementation of Algorithm for Problem (VI) 4.1
Error Estimates of Solutions for Problem (VI) In this subsection, we recur to the classical MFE method to derive the errors between the solutions for Problem (IV) and the reduced-order MFE solutions for Problem (VI). To this end, it is necessary to introduce the following discrete Gronwall lemma (see [2, 23]). Lemma 4 (Discrete Gronwall Lemma) If {an }, {bn }, and {cn } are three positive sequences, and {cn } is monotone increase progressively, that satisfy ¯ a n + b n 6 cn + λ
n−1 X
¯ > 0, ai , λ
a0 + b 0 6 c 0 ,
i=0
then ¯ an + bn 6 cn exp(nλ), n > 0. We have the following main results for Problem (VI).
1478
ACTA MATHEMATICA SCIENTIA
Vol.33 Ser.B
Theorem 5 Under hypotheses of Theorem 1, if L2 = O(N ) and τ = O(h), then the following error estimates hold kunh
−
und k0
+k
1/2
X n i=1
kpih
−
pid )k20
1/2
1/2 l X 1/2 6C τ λj , 1 6 n 6 L.
(4.1)
j=d+1
Proof Since M d ⊂ Mh and H d ⊂ Hh , subtracting Problem (VI) from Problem (IV) respectively taking vh = vd ∈ M d and q h = q d ∈ H d yields that (pnh − pnd , q d ) + (q d , ∇(uh − und )) = 0, ∀q d ∈ H d , (uh −
und , vd )
u0h (x, y)
−
−
τ (pnh
u0d (x, y)
=
−
pnd , ∇vd )
u0h (x, y)
−
= (un−1 − un−1 , vd ), ∀vd h d d 0 P uh (x, y), (x, y) ∈ Ω.
(4.2) d
∈M ,
(4.3) (4.4)
By using (3.8), (3.9), (4.2), and (4.3), we have that τ kpnh − pnd k20 = τ (pnh − pnd , pnh − ρd pnh ) + τ (pnh − pnd , ρd pnh − pnd )
= τ (pnh − ρd pnh , pnh − ρd pnh ) − τ (ρd pnh − pnd , ∇(uh − und ))
= τ (pnh − ρd pnh , pnh − ρd pnh ) − τ (ρd pnh − pnh , ∇(uh − und ))
−τ (pnh − pnd , ∇(uh − P d unh )) − τ (pnh − pnd , ∇(P d uh − und ))
= τ (pnh − ρd pnh , pnh − ρd pnh ) − τ (ρd pnh − pnh , ∇(uh − P d unh )) −τ (pnh − ρd pnh , ∇(uh − P d unh )) − (unh − und , P d uh − und )
+(un−1 − un−1 , P d uh − und ) h d
= τ (pnh − ρd pnh , pnh − ρd pnh ) − (unh − und , P d uh − unh )
−(unh − und , uh − und ) + (un−1 − un−1 , P d uh − unh ) h d
+(un−1 − un−1 , uh − und ). h d
(4.5)
By transporting terms for (4.5) and using H¨older inequality, Cauchy inequality, and (3.10), if τ = O(h) we have that τ kpnh − pnd k20 + kunh − und k20
= τ (pnh − ρd pnh , pnh − ρd pnh ) − (unh − und , P d uh − unh )
+(un−1 − un−1 , P d uh − unh ) + (un−1 − un−1 , uh − und ) h d h d
6 τ kpnh − ρd pnh k20 + Chkunh − und k0 k∇(P d uh − unh )k0
+Chkun−1 − un−1 k0 k∇(P d uh − unh )k0 + kun−1 − un−1 k0 kuh − und k0 h d h d
6 Cτ [kpnh − ρd pnh k20 + k∇(P d uh − unh )k20 ] + Cτ kunh − und k20 1 − un−1 k20 + kuh − und k20 ]. +Cτ kun−1 − un−1 k20 + [kun−1 d h d 2 h
(4.6)
That is that τ kpnh − pnd k20 + kunh − und k20
6 2τ kpnh − pnd k20 + kunh − und k20
6 Cτ [kpnh − ρd pnh k20 + k∇(P d uh − unh )k20 ] + Cτ kunh − und k20 +Cτ kun−1 − un−1 k20 + kun−1 − un−1 k20 . h d h d
(4.7)
No.5
Z.D. Luo et al: A REDUCED-ORDER MFE FORMULATION BASED ON POD METHOD
1479
Summing from 1 to n for (4.7), noting that L2 = O(N ) (i.e., τ 1/2 = O(1/L)), and using Lemma 3 could yield that kunh − und k20 + τ
n X i=1
kpih − pid k20 6 Cτ 1/2
l X
λj + Cτ
n X i=1
j=d+1
kuih − uid k20 .
(4.8)
If τ is sufficiently small such that Cτ 6 1/2, we obtain from (4.8) that kunh − und k20 + τ
n X i=1
kpih − pid k20 6 kunh − und k20 + 2τ 6 Cτ 1/2
l X
n X i=1
λj + Cτ
kpih − pid k20 n−1 X i=0
j=d+1
kuih − uid k20 .
(4.9)
Applying Lemma 4 to (4.9) yields that kunh − und k20 + τ
n X i=1
kpih − pid k20 6 Cτ 1/2
l X
λj exp(Cτ n),
(4.10)
j=d+1
which yields (4.1) and completes the proof of Theorem 5. 2 Combining Theorems 1 and 5 yields the following result. Theorem 6 Under hypotheses of Theorem 5, the error estimates between the solution for Problem (II) and the solutions for the reduced-order Problem (VI) are X 21 n 1 ku(tn ) − und k0 + τ 2 kp(ti ) − pid k20 i=1
1
6 C(τ 1+ 2
1/2 l X m+1 1/2 +h )+C τ λj , 1 6 n 6 L. j=d+1
Remark 3 The condition N = O(L2 ) in Theorem 5 shows the relation between the number L of snapshots and the number N of all time instances. Therefore, it is unnecessary to take total transient solutions at all time instances tn as snapshots in [15, 16]. Theorems 5 and 6 have respectively provided the error estimates between the solutions of the reduced-order MFE Problem (VI) and the solutions of the classical MFE Problem (IV) and Problem (II), which have supplied a guide choosing number d of POD basis. Our method here employs the first L MFE solutions (unh , pnh ) (n = 1, 2, · · · , L) for Problem (IV) as assistant analysis. However, when one solves actual problems, he may obtain the ensemble of snapshots from physical system trajectories by drawing samples from experiments with interpolation (or data assimilation) or previous results. Therefore, the assistants (unh , pnh ) (n = 1, 2, · · · , L) could be replaced with the interpolation functions of experimental or previous results, thus rendering it unnecessary to solve Problem (IV) and requiring only to solve directly Problem (VI) such that Theorem 5 is satisfied. And then, time instances are continuously extrapolated forward and POD basis is ceaselessly renewed, the rules of future development and change of natural phenomenon would be very well simulated. 4.2 Implementation of Algorithm for Problem (VI) In the following, we give the implementation of algorithm for solving Problem (VI) which consists of seven steps.
1480
ACTA MATHEMATICA SCIENTIA
Vol.33 Ser.B
Step 1 Generate the snapshots ensemble U i (x, y) = (uih , pih ), i = 1, 2, · · · , L = O(N 1/2 ), which may be the first L solutions for Problem (IV), physical system trajectories by drawing samples from experiments with interpolation (or data assimilation), or previous results. Step 2 Generate the correlation matrix A = (Aik )L×L , where Aik = L1 (U i , U k )H×M , (U i , U k )H×M = (∇uih , ∇ukh ) +(pih , pkh ), and (·, ·) is L2 -inner product. Step 3 Solving the eigenvalue problem Av = λv, v = (a1 , a2 , · · · , aL )T
(4.11)
obtains eigenvalues λ1 > λ2 > · · · > λl > 0 (l = dim{U 1 , U 2 , · · · , U L }) and corresponding eigenvectors v k = (ak1 , ak2 , · · · , akL ) (k = 1, 2, · · · , l). Step 4 For error δ needed, decide on the triangulation parameter h, the time step increment τ , the degree m of interpolation polynomial, and the amount d of POD basis such that 1/2 l P 1 τ 1+ 2 + hm+1 + k 1/2 6 δ. λj j=d+1
Step 5
Generate POD basis ψ k (x, y) = (ψuk (x, y), ψ pk (x, y)) (k = 1, 2, · · · , d): (ψuk (x, y), ψ pk (x, y)) =
L L 1 X k i 1 X k i √ ai u h , √ ai ph . Lλk i=1 Lλk i=1
Step 6 Taking M d = span {ψu1 (x, y), ψu2 (x, y), · · · , ψud (x, y)} and H d =span {ψp1 (x, y), ψ p2 (x, y), · · · , ψ pd (x, y)} and solving Problem (VI) which only includes 3d degrees of freedom yield the solutions (und , pnd ) (n = 1, 2, · · · , L, L + 1, · · · , N ). Step 7 If kun−1 − und k0 > kund − un+1 k0 and kpn−1 − pnd k0 > kpnd − pn+1 k0 (n = L, L + d d d d n n 1, · · · , N −1), then (ud , pd ) (n = 1, 2, · · · , N ) are the solution for Problem (VI), whose errors are 1/2 l P not more than k 1+1/2 + hm+1 + k 1/2 λj . Else, i.e., if kun−1 − und k0 < kund − un+1 k0 d d j=d+1
or kpn−1 − pnd k0 < kpnd − pn+1 k0 (n = L, L + 1, · · · , N − 1), let U i = (uid , pid ) (i = n − L, n − L − d d 1, · · · , n − 1), repeat Step 1 to Step 6 until the solutions (und , pnd ) (n = 1, 2, · · · , L, L + 1, · · · , N ) are obtained.
5
Some Numerical Experiments
In this section, some numerical examples are used to validate the feasibility and efficiency of the reduced-order MFE formulation, i.e., Problem (VI) based on POD method for twodimensional parabolic equations. For the sake of convenience, without loss of generality, herein we consider a two-dimensional parabolic equations by Choosing computing field Ω = [0, 2]×[0, 2], total time T = 0.4, u0 (x, y) = 2 sin πx sin πy, and f (x, y, t) = − sin πx sin πye−2π t in Problem (I) as an example, whose ideas and approaches could directly apply to numerical computation for two-dimensional parabolic equations with real practical applied background. We first divide the field Ω into 2000×2000 small squares with side length △x = △y = 0.001, and then link the diagonal of the square to divide each square into two triangles in the same
No.5
Z.D. Luo et al: A REDUCED-ORDER MFE FORMULATION BASED ON POD METHOD
1481
√ direction which consists of triangulation ℑh . Thus h = 2 × 0.001. In order to make the condition τ = O(h) satisfied in Theorem 5, we take time step size as τ = 0.002. We first find a sequence of numerical solutions (unh , pnh ) ∈ Mh × Hh (n = 1, 2, · · · , 200, i.e., at time t = 1τ, 2τ, · · · , 200τ ) of the classical MFE formulation, i.e., Problem (IV) where Mh is the subspace of piecewise 2-th degree polynomials and Hh the subspace of piecewise 1-th degree polynomial vectors. And then, we chose the first 20 numerical solutions unh (n = 1, 2, · · · , 20) as snapshots U i = (uih , pih ) (i = 1, 2, · · · , 20). Finally, we find 20 eigenvalues which is arranged in a non-increasing order and 20 eigenvectors correspondent to the 20 eigenvalues and we construct 20 POD bases ψ j (j = 1, 2, · · · , 20) with Step 2, 3, and 5 of Subsection 4.2. Take the first 6 POD bases ψ j (j = 1, 2, · · · , 6) from 20 POD bases ψ j = (ψuk (x, y), ψ pk (x, y)) (j = 1, 2, · · · , 20) to span into subspaces M d and H d and compute a numerical solution (und , pnd ) (n = 200, i.e., at t = 200τ ) with the reduced-order MFE formulation, i.e., the reduced-order Problem (VI) according seven steps in Subsection 4.2, without renewing POD basis it is still convergent and depicted graphically on the left charts in Figures 1–3. While the classical MFE solution (unh , pnh ) ∈ Mh × Hh (where Mh is the subspace of piecewise 2-th degree polynomial and Hh the subspace piecewise 1-th degree polynomial vector) is obtained by solving Problem (IV) when n = 200, i.e., at t = 0.4 and is depicted graphically on the right charts in Figures 1–3. Every two charts in Figures 1–3 are exhibiting quasi-identical similarity, respectively, but the reduced-order MFE solutions are better than the classical MFE solutions (it is due to use more initial numerical solutions, i.e., six POD bases which is an open question).
Fig.1
Left chart is the reduced-order MFE solution pn 1d when d = 6 at t = 0.4; right chart is
the classical MFE solution pn 1h at t = 0.4
Fig.2
Left chart is the reduced-order MFE solution pn 2d when d = 6 at t = 0.4; right chart is
the classical MFE solution pn 2h at t = 0.4
1482
ACTA MATHEMATICA SCIENTIA
Fig.3
Vol.33 Ser.B
Left chart is the reduced-order MFE solution un d when d = 6 at t = 0.4; right chart is
the classical MFE solution un h at t = 0.4
Fig.4
When t = 2, the errors (log10 ) between solutions of Problem (VI) with different
number of POD bases for a group of 20 snapshots and the classical MFE formulation Problem (IV)
√ When we take 6 POD bases and τ = 0.002 and h = 2 × 0.001, by computing we h i1/2 20 P obtain that τ 1/2 λj + τ 1+1/2 + h2 6 3 × 10−4 of theoretical errors in Theorem 6. j=7
Figure 4 computationally shows the errors (log10 ) between the reduced-order MFE solutions (und , pnd ) with 20 different POD bases and the classical MFE solution (unh , pnh ) at t = 200τ (i.e., n = 200), respectively. Comparing the classical MFE formulation with the reduced-order MFE formulation containing 6 POD bases implementing the numerical simulation computations when total time t = 200τ , we find that for the classical MFE formulation with piecewise 2-th degree polynomials for unh and piecewise linear polynomials vector for pnh , which has 5 × 2000 × 2000 = 2 × 107 degrees of freedom, the required computing time is 16 minutes, while for the reduced-order MFE formulation with 6 POD bases, which has only 18 degrees of freedom, the corresponding time is only 8 seconds, i.e., the required computing time to solve the classical MFE formulation is about as 120 times as that to do the reduced-order MFE formulation with 6 POD bases, while the errors between their respective solutions do not exceed 3 × 10−4 . It is also shown that finding the approximate solutions for two-dimensional parabolic equations with the reduced-order MFE formulation is computationally very effective. And the results for numerical examples are consistent with those obtained for the theoretical case. Though our examples are in sense recomputing what we have already computed by classical MFE
No.5
Z.D. Luo et al: A REDUCED-ORDER MFE FORMULATION BASED ON POD METHOD
1483
formulation (but only use the first 20 numerical solutions (unh , pnh ) (n = 1, 2, · · · , 20) at the first 20 steps), when we compute actual problems, we may also construct the snapshots and POD basis with interpolation or data assimilation by drawing samples from experiments, then solve directly the reduced-order MFE formulation, while it is unnecessary to solve classical MFE formulation such that the computational load could be alleviated and time-consuming of calculations in the computational process saved.
6
Conclusions
In this paper, we have employed the POD method to establish a new reduced-order MFE formulation for two-dimensional parabolic equations, analyzed the errors between the solutions of its classical MFE formulation and the solutions of the reduced-order MFE formulation based on POD method, and provide the steps of the implementation of algorithm for solving the reduced-order MFE formulation, i.e., Problem (VI), which have shown that our present method has improved and innovated the existing methods (for instance, the methods in [1, 2, 4, 13– 18]). Comparing with the theoretical error estimates with those of numerical examples shows that the results for numerical examples are consistent with those obtained for the theoretical case, further validates both the feasibility and efficiency of the reduced-order MFE formulation. Though snapshots and POD basis of our numerical examples are constructed with the first fewer solutions of the classical MFE formulation, when one computes actual two-dimensional parabolic equations, this process can be omitted in actual applications and one may construct the snapshots and POD basis with interpolation or data assimilation by drawing samples from experiments, then directly solves the reduced-order MFE formulation or Problem (VI), while it is unnecessary to solve Problem (IV) such that the computational load could be alleviated and time-consuming of calculations in the computational process saved. Therefore, the method in this paper holds a good prospect of extensive applications. Future research work in this area will aim at extending the reduced-order MFE formulation, applying it to a set of more complicated PDEs such as the atmosphere quality forecast system and the ocean fluid forecast system. References [1] [2] [3] [4] [5] [6] [7] [8] [9] [10]
Thom´ ee V. Galerkin Finite Element Methods for Parabolic Problems. Berlin: Springer, 1997 Luo Z D. Mixed Finite Element Methods and Applications. Beijing: Science Press, 2006 Brezzi F, Fortin M. Mixed and Hybrid Finite Element Methods. New York: Springer-Verlag, 1991 Li L, Sun P, Luo Z D. A new mixed finite element formulation and error estimates for parabolic equations. Acta Mathematica Scientia, 2012, 32A(6): 1158–1165 Fukunaga K. Introduction to Statistical Recognition. New York: Academic Press, 1990 Jolliffe I T. Principal Component Analysis. Berlin: Springer-Verlag, 2002 Holmes P, Lumley J L, Berkooz G. Turbulence, Coherent Structures, Dynamical Systems and Symmetry. Cambridge, UK: Cambridge University Press, 1996 Sirovich L. Turbulence and the dynamics of coherent structures: Part I-III. Quarterly of Applied Mathematics, 1987, 45(3): 561–590 Lumley J L. Coherent structures in turbulence//Meyer R E, ed. Transition and Turbulence. Academic Press, 1981 Aubry N, Holmes P, Lumley J L, Stone E. The dynamics of coherent structures in the wall region of a turbulent boundary layer. J Fluid Dynamics, 1988, 192: 115–173
1484
ACTA MATHEMATICA SCIENTIA
Vol.33 Ser.B
[11] Moin P, Moser R D. Characteristic-eddy decomposition of turbulence in channel. J Fluid Mech, 1989, 200: 417–509 [12] Rajaee M, Karlsson S K F, Sirovich L. Low dimensional description of free shear flow coherent structures and their dynamical behavior. J Fluid Mech, 1994, 258: 1–29 [13] Joslin R D, Gunzburger M D, Nicolaides R A, Erlebacher G, Hussaini M Y. A self-contained automated methodology for optimal flow control validated for transition delay. AIAA J, 1997, 35(5): 816–824 [14] Kunisch K, Volkwein S. Galerkin proper orthogonal decomposition methods for parabolic problems. Numerische Mathematik, 2001, 90(1): 117–148 [15] Kunisch K, Volkwein S. Galerkin proper orthogonal decomposition methods for a general equation in fluid dynamics. SIAM J Num Anal, 2002, 40(2): 492–515 [16] Sun P, Luo Z D, Zhou Y J. Some reduced finite difference schemes based on a proper orthogonal decomposition technique for parabolic equations. Appl Num Math, 2010, 60(1/2): 15–164 [17] Luo Z D, Chen J, Sun P, et al. Finite element formulation based on proper orthogonal decomposition for parabolic equations. Science in China Series A: Mathematics, 2009, 52(3): 587–596 [18] Luo Z D, Chen J, Xie Z H, et al. A reduced second-order time accurate finite element formulation based on POD for parabolic equations (in Chinese). Science in China Series A: Mathematics, 2011, 41(5): 447–460 [19] Cao Y H, Zhu J, Navon I M, et al. A reduced-order approach to four-dimensional variational data assimilation using proper orthogonal decomposition. Int J Num Methods Fluids, 2007, 53(10): 1571–1583 [20] Tian X J, Xie Z H, Dai A G. An ensemble-based explicit four-dimensional variational assimilation method. J Geophys Res, 2008, 113(D21124), doi: 10.1029/2008JD010358 [21] Tian X J, Xie Z H, Sun Q. A POD-based ensemble four-dimensional variational assimilation method. Tellus, 2011, 63A: 805–816 [22] Adams R A. Sobolev Spaces. New York: Academic Press, 1975 [23] Girault V, Raviart P A. Finite Element Approximations of the Navier-Stokes Equations//Theorem and Algorithms. New York: Springer-Verlag, 1986 [24] Ciarlet P G. The Finite Element Method for Elliptic Problems. Amsterdam: North-Holland, 1978 [25] Rudin W. Functional and Analysis. 2nd ed. McGraw-Hill Companies, Inc, 1973
A REDUCED-ORDER MFE FORMULATION BASED ON POD METHOD FOR PARABOLIC EQUATIONS∗
Zhendong LUO (
)†
School of Mathematics and Physics, North China Electric Power University, Beijing 102206, China E-mail: [email protected]
)
Lei LI (
)
Ping SUN (
School of Mathematics and Computer Science, Guizhou Normal University, Guiyang 550001, China E-mail: [email protected]; [email protected]
Abstract In this paper, we extend the applications of proper orthogonal decomposition (POD) method, i.e., apply POD method to a mixed finite element (MFE) formulation naturally satisfied Brezz-Babu˘ska for parabolic equations, establish a reduced-order MFE formulation with lower dimensions and sufficiently high accuracy, and provide the error estimates between the reduced-order POD MFE solutions and the classical MFE solutions and the implementation of algorithm for solving reduced-order MFE formulation. Some numerical examples illustrate the fact that the results of numerical computation are consistent with theoretical conclusions. Moreover, it is shown that the new reduced-order MFE formulation based on POD method is feasible and efficient for solving MFE formulation for parabolic equations. Key words proper orthogonal decomposition method; mixed finite element formulation; parabolic equation; error estimate 2010 MR Subject Classification
1
65N30; 65M30
Introduction
Let Ω ⊂ R2 be a bounded and connected polygonal domain. Consider the following parabolic equations. Problem (I) Find u such that
∗ Received
(x, y, t) ∈ Ω × (0, T ), ut − △u = f, u(x, y, t) = 0, (x, y, t) ∈ ∂Ω × [0, T ), u(x, y, 0) = u0 (x, y), (x, y) ∈ Ω,
(1.1)
May 26, 2011; revised October 31, 2012. The research was supported by the National Science Foundation of China (11271127 and 11061009), Science Research Program of Guizhou (GJ[2011]2367), and the Co-Construction Project of Beijing Municipal Commission of Education. † Correspondence author: Zhendong LUO.
1472
ACTA MATHEMATICA SCIENTIA
Vol.33 Ser.B
where u represents unknown function, T the total time, f the given source term, and u0 (x, y) the given initial function. Let p = −∇u, then Problem (I) can be written as follows. Problem (II)
Find (u, p) such that
ut + divp = f, p + ∇u = 0, u(x, y, t) = 0, u(x, y, 0) = u0 (x, y),
(x, y, t) ∈ Ω × (0, T ), (x, y, t) ∈ Ω × (0, T ),
(x, y, t) ∈ ∂Ω × [0, T ),
(1.2)
(x, y) ∈ Ω.
The parabolic equations may be used to describe many physical phenomena of the natural environment, engineering equipment, and living organisms, such as the proliferation of gas, the infiltration of liquid, the conduction of heat, and the spread of impurities in semiconductor materials. It isn’t easy to find the exact solutions for practical engineering parabolic problems if their computational fields are usually irregular. Thus, finding their numerical solutions becomes an efficient approach. The mixed finite element (MFE) method is one of the most effective methods to find the numerical solutions for parabolic equations since it can simultaneously find unknown function (displacement) and gradient (stress), which has been widely applied in scientific computations and its theory has been established (see [1, 2]). However, the classical MFE formulations for parabolic equations include too many degrees of freedom and need to satisfy an important convergence stability condition, i.e., the Brezzi-Babuˇska (BB) inequality (see [1–3]) so that they bring many difficulties for practical applications. We adopt the same MFE formulation as that introduced in [4] which can circumvent the constraint of the BB inequality. Thus, the key problem is how to alleviate the computational load by saving time-consuming calculations in the computational process in a way that guarantees sufficiently accurate numerical solutions. It was shown that a proper orthogonal decomposition (POD) method by combining with some numerical methods for partial differential equations can provide an efficient means of generating reduced-order models and alleviating the computational load and memory requirements (see [5]). POD method was widely and successfully applied to numerous fields, including signal analysis and pattern recognition (see [6]), statistics (see [7]), geophysical fluid dynamics and meteorology (also see [7]). The POD method essentially provides an orthogonal basis for representing the given data in a certain least squares optimal sense, that is, it provides a way to find optimal lower dimensional approximations of the given data. In early times, POD method was mainly used to perform principal component analysis in computations of statistics (see [5–7]), until the method of snapshots was introduced by Sirovich (see [8]) and was then widely applied to reducing the order of the POD eigenvalue problem and searching the main behavior of a dynamic system (see [8–13]). Until ten years ago, some Galerkin POD methods for parabolic problems and a general equation in fluid dynamics were not presented (see [14, 15]). More recently, some reduced finite difference schemes and finite element formulations based on POD technique for parabolic problems were presented (see [16– 18]). Moreover, some reduced-order methods based on POD technique for four-dimensional variational data assimilation were established (see [19–21]).
No.5
Z.D. Luo et al: A REDUCED-ORDER MFE FORMULATION BASED ON POD METHOD
1473
Though some Galerkin methods, finite difference schemes, and finite element formulations for parabolic problems were reduced with POD method (see [14–18]), to the best of our knowledge, there are no published results addressing the case that the classical MFE formulation for parabolic equations is reduced with the POD method or providing error estimates between the classical MFE solutions and the reduced-order MFE solutions or the implementation of algorithm for solving the reduced-order MFE formulation. In this paper, we extend the developments in [14–18], i.e., combine the classical MFE method with POD method to establish a reduced-order MFE formulation with lower dimensions and sufficiently high accuracy for two-dimensional parabolic equations, analyze the errors between the reduced-order MFE solutions based on POD method and the classical MFE solutions, and provide the implementation of algorithm for solving reduced-order MFE formulation. Some numerical examples illustrate the fact that the results of numerical computations are consistent with theoretical conclusions. Moreover, it is shown that the reduced-order MFE formulation based on POD method is feasible and efficient for finding numerical solutions of two-dimensional parabolic equations. The main advantage in this paper consists in that it does only use the first fewer given numerical solutions on very short time span [0, T0 ] (T0 ≪ T ) as snapshots to construct POD basis and establish the reduced-order MFE formulation based on POD technique for finding the numerical solutions on total time span [0, T ]. Thus, it can sufficiently adopt the advantage of POD method, namely it utilizes the given data (on very short time span [0, T0 ] and T0 ≪ T ) to forecast (or infer) future physic phenomenons (on time span [T0 , T ]). This is the improvement and innovation for the existing reduced-order methods based on POD technique (see, e.g., [14–18]). The rest of this paper is organized as follows. Section 2 is to recall the classical fully discrete MFE formulation for two-dimensional parabolic equations introduced in [4] and to generate snapshots from the first few transient solutions computed from the equation system derived by the classical fully discrete MFE formulation. In Section 3, the optimal orthonormal POD bases are spanned by the elements of the snapshots with POD method and a reducedorder MFE formulation with lower dimensions and sufficiently high accuracy based on POD method for two-dimensional parabolic equations is established. In Section 4, the error estimates between classical MFE solutions and the reduced-order MFE solutions based on POD method and the implementation of algorithm for solving the reduced-order MFE formulation are provided. In Section 5, some numerical examples are presented for illustrating that the numerical computing errors between the reduced-order MFE solutions based on POD method and the classical MFE solutions are consistent with previously obtained theoretical results, thus validating the feasibility and efficiency of the reduced-order MFE formulation. Section 6 provides main conclusions and future tentative ideas.
2
Recall Classical MFE Formulation and Generate Snapshots
The Sobolev spaces used in this context are standard (see [22]). Let M = H01 (Ω), H = L2 (Ω)2 . For the sake of convenience, without loss of generality, we may as well assume that u0 (x, y) = 0 in the following theoretical analysis. Then a mixed variational formulation for Problem (II) is written as follows.
1474
ACTA MATHEMATICA SCIENTIA
Vol.33 Ser.B
Problem (III) Find (u, p) ∈ M × H such that, for any t ∈ (0, T ), a(p, q) + b(q, u) = 0,
b(p, v) − (ut , v) = −(f, v), u(x, y, 0) = 0,
∀q ∈ H,
(2.1)
∀v ∈ M, (x, y) ∈ Ω,
where a(p, q) = (p, q), b(q, v) = (q, ∇v), (·, ·) represents the inner product in L2 (Ω) or L2 (Ω)2 . Since a(·, ·) and b(·, ·) all are continuous bilinear functionals, a(·, ·) is coercive on H × H, and b(·, ·) satisfies continuous BB condition (as are proved in [4]), i.e., for any q ∈ H, b(·, ·) satisfies b(q, v) sup > βkvk1 , ∀v ∈ M, (2.2) q∈H kqk0 where β is a constant, if f ∈ H s (0, T ; H r (Ω)), the mixed variational formulation Problem (III) for Problem (II) exists a unique solution (u, p) ∈ M × H satisfying that: kp(i) kr+1 + ku(i+1) kr+2 6 Ckf (i) kr , ∀t ∈ (0, T ), i = 1, 2, · · · , s,
(2.3)
where C is a positive constant, f (i) represent the i-th order derivative with respect to time t. In order to obtain the MFE solutions of Problem (III), it is necessary to introduce an MFE approximation for the spatial variables of Problem (III) and an approximation of difference quotient for time derivative. ¯ (see [2, 3, 23, 24]). MFE Let {ℑh } be a uniformly regular family of triangulation of Ω spaces are taken as follows. Hh = {q h ∈ H; q h |K ∈ Pm−1 (K)2 , ∀K ∈ ℑh },
Mh = {vh ∈ M ∩ C 0 (Ω); vh |K ∈ Pm (K), ∀K ∈ ℑh }, where m > 1 and Pm (K) is the space of polynomials of degree 6 m on K. For any positive integer N , let k = T /N denote time–step and write tn = nk (0 6 n 6 N ). The time derivative for Problem (III) are approximated by one–step Euler backward difference quotient ∂¯t un = (un − un−1 )/τ . Let (unh , pnh ) be the fully discrete MFE approximation of (u, p). Then, the fully discrete MFE formulation for Problem (III) may be written as follows. Problem (IV) Find (unh , pnh ) ∈ Mh × Hh such that, for 1 6 n 6 N , n n a(ph , q h ) + b(q h , uh ) = 0, τ b(pnh , vh )
−
(unh , vh )
=
−(un−1 , vh ) h
u0h (x, y) = 0, (x, y) ∈ Ω.
n
− τ (f , vh ),
∀q h ∈ Hh ,
∀vh ∈ Mh ,
(2.4)
Put A((u, p); (v, q)) = τ a(p, q) + τ b(q, u) − τ b(p, v) + (u, v) and F (v, q) = (un−1 , v) + h τ (f n , v). Problem (IV) is equivalently rewritten as follows. Problem (V) Find (unh , pnh ) ∈ Mh × Hh such that, for 1 6 n 6 N , A((un , pn ); (v , q )) = F (v , q ), ∀(vh , q h ) ∈ Mh × Hh , h h h h h h (2.5) u0 (x, y) = 0, (x, y) ∈ Ω. h
No.5
Z.D. Luo et al: A REDUCED-ORDER MFE FORMULATION BASED ON POD METHOD
1475
It has been proved that Problem (V) has a unique set of solutions (unh , pnh ) ∈ Mh × Hh (n = 1, 2, · · · , N ) and the following theorem of error estimates (see [4]). n P Theorem 1 If f ∈ H 1 (0, T ; L2(Ω)), kf i k2m−1 is a bounded constant, and τ = O(h), i=1
then errors between the solution (u, p) to Problem II and solutions (unh , pnh ) are as follows: ku(tn ) −
unh k0
+τ
1 2
X n i=1
kp(ti ) −
pih k20
21
1
6 C(τ 1+ 2 + hm+1 ), 1 6 n 6 N.
(2.6)
If source term f (x, y, t), initial conditions u0 (x, y), triangulation parameter h, the time– step k, finite element subspaces Xh and Mh are given, we can obtain a ensemble of solutions {(unh (x, y), pnh (x, y))}N n=1 by solving Problem (IV). And then we choose the first L (in general, L ≪ N , for example, L = 20, N = 200) instantaneous solutions (uih (x, y), pih (x, y)) (1 6 i 6 L ≪ N ) from the N group of solutions {(unh (x, y), pnh (x, y))}N n=1 for Problem (IV), which are useful and of interest for us and referred to as snapshots. Remark 1 When one computes actual problems, one may obtain the ensemble of snapshots from physical system trajectories by drawing samples from experiments and interpolation (or data assimilation).
3 Generation of POD Basis and a Reduced-Order MFE Formulation Based on POD Technique for Parabolic Equations For (uih (x, y), pih (x, y)) (1 6 i 6 L ≪ N ) in Section 2, let U i (x, y) = (uih (x, y), pih (x, y)) (1 6 i 6 L) and V = span{U 1 , U 2 , · · · , U L }, (3.1) and refer to V as the space generating by the snapshots {U i }L i=1 at least one of which is assumed to be non-zero. Let {ψ j }lj=1 denote an orthonormal basis of V with l = dim V (l 6 L). Then each member of the ensemble can expressed as Ui =
l X j=1
(U i , ψ j )H×M ψ j , i = 1, 2, · · · , L,
(3.2)
where (U i , ψ j )H×M = (∇uih , ∇ψuj ) + (pih , ψ pj ), ψ j = (ψuj , ψ pj ), and ψuj and ψ pj are orthonormal basis corresponding to u and p, respectively. Definition 2 The method of POD consists in finding an orthonormal basis sequence ψ i (i = 1, 2, · · · , l) such that for every d (1 6 d 6 l) the mean square error between the elements U i (1 6 i 6 L) and corresponding d-th partial sum of (3.2) is minimized on average min
{ψ j }d j=1
such that
2 L d X
1 X
U i − (U i , ψ j )H×M ψ j
L i=1 H×M j=1
(ψ i , ψj )H×M = δij , 1 6 i 6 d, 1 6 j 6 i,
(3.3)
(3.4)
where kU i k2H×M = k∇uih k20 + kpih k20 . A solution {ψ j }dj=1 of (3.3) and (3.4) is known as a POD basis of rank d.
1476
ACTA MATHEMATICA SCIENTIA
Vol.33 Ser.B
We introduce the correlation matrix A = (Aij )L×L ∈ RL×L corresponding to the snapshots {U i }L i=1 by 1 Aij = (U i , U j )H×M . (3.5) L Since the matrix A is positive semi-definite and has rank l, the solution of (3.3) and (3.4) can be found, in addition, the following result holds (for example, see [7, 8, 14, 15]). Proposition 2 Let λ1 > λ2 > · · · > λl > 0 denote the positive eigenvalues of A and v 1 , v 2 , · · ·, v l the associated orthonormal eigenvectors. Then a POD basis of rank d 6 l is given by L 1 X √ ψi = (v i )j U j , 1 6 i 6 d, (3.6) Lλi j=1 where (v i )j ’s denote the j-th component of the eigenvector v i . Furthermore, the following error formula holds
2 L d l X X
1 X
U i −
(U , ψ ) ψ = λj . (3.7) i j H×M j L i=1 H×M j=1 j=d+1
Let M d = span {ψu1 , ψu2 , · · · , ψud } and H d = span ψ p1 , ψ p2 , · · · , ψ pd . Due to for any vd ∈ M d satisfying ∇vd ∈ H d , we define the projection P d : M h → M d denoted by, for uh ∈ Mh , (∇P d uh , q d ) = (∇uh , q d ), ∀q d ∈ H d . Then there is an extension P h : M → Mh of P d such that P h |Mh = P d : Mh → M d denoted by (see [25]), for any u ∈ M , (∇P h u, qh ) = (∇u, q h ),
∀q h ∈ Hh .
(3.8)
Set the L-projection ρd : H → H d denoted by, for any p ∈ H, (ρd p, q d ) = (p, q d ),
∀q d ∈ H d .
(3.9)
Due to (3.8) and (3.9), the operators P h and ρd are well-defined and bounded (see [15, 18]) k∇(P d u)k0 6 k∇uk1 , ku − P h uk0 6 Chk∇(u − P h u)k0 , ∀u ∈ M
(3.10)
kρd pk0 6 CkpkM , ∀p ∈ H.
(3.11)
and
Moreover, the following results hold (see [14, 15, 17, 18]). Lemma 3 For every d (1 6 d 6 l) the projection operators P d and ρd satisfy L l X 1X i [kuh − P d uih k20 + h2 k∇(uih − P d uih )k20 ] 6 Ch2 λj , L i=1
(3.12)
j=d+1
L l X 1X i kph − ρd pih k20 6 λj , L i=1 j=d+1
where (uih , pih ) ∈ V is the solution sequence of Problem (IV).
(3.13)
No.5
Z.D. Luo et al: A REDUCED-ORDER MFE FORMULATION BASED ON POD METHOD
1477
Thus, substituting H d and M d for Hh and Mh for Problem (IV), respectively, we can obtain a new reduced-order MFE formulation based on POD technique for Problem (I) as follows. Problem (VI) Find (und , pnd ) ∈ M d × H d such that n n a(pd , q d ) + b(q d , ud ) = 0, τ b(pnd , vd ) − (und , vd ) = −(un−1 , vd ) − τ (f n , vd ), d 0 ud (x, y) = 0, (x, y) ∈ Ω.
∀q d ∈ H d ,
∀vd ∈ M d ,
(3.14)
By using the same approaches as proofs of Theorem 3 in [4], we may prove that Problem (VI) exists a unique solution sequence (und , pnd ) ∈ M d × H d (n = 1, 2, · · · , N ) satisfying kund k20
+τ
n X i=1
kpid k20
n X 0 2 i 2 6 kud k0 + 2τ kf k0 exp(2nτ ),
(3.15)
i=1
which has shown that the solution sequence of Problem (VI) is stable and continuously dependent on source term f (x, y, t) and initial conditions u0 (x, y). Remark 2 If ℑh is a uniformly regular triangulation and Mh and Hh are respectively the spaces of piecewise second order polynomials and piecewise linear polynomial vectors, the number of total degrees of freedom for Problem (VI), i.e., the number of unknown quantities is 3Np +Ns ≈ 5Np (where Np and Ns are respectively the number of vertices and sides of triangles in ℑh ), while the number of total degrees of freedom for Problem (VI) is 3d (d ≪ l 6 L ≪ N ). For scientific engineering problem, the number Np of vertices of triangles in ℑh is more than ten of thousands, even more than a hundred million, while d is only the number of the first few maximal eigenvalues which is chosen the first L snapshots from the N snapshots so that it is very small (for example, in Section 5, d = 6, while Np = 2000 × 2000 = 4 × 106 ). Therefore, Problem (VI) is a reduced-order MFE formulation based on POD method for Problem (IV).
4 Error Analysis of Solutions and Implementation of Algorithm for Problem (VI) 4.1
Error Estimates of Solutions for Problem (VI) In this subsection, we recur to the classical MFE method to derive the errors between the solutions for Problem (IV) and the reduced-order MFE solutions for Problem (VI). To this end, it is necessary to introduce the following discrete Gronwall lemma (see [2, 23]). Lemma 4 (Discrete Gronwall Lemma) If {an }, {bn }, and {cn } are three positive sequences, and {cn } is monotone increase progressively, that satisfy ¯ a n + b n 6 cn + λ
n−1 X
¯ > 0, ai , λ
a0 + b 0 6 c 0 ,
i=0
then ¯ an + bn 6 cn exp(nλ), n > 0. We have the following main results for Problem (VI).
1478
ACTA MATHEMATICA SCIENTIA
Vol.33 Ser.B
Theorem 5 Under hypotheses of Theorem 1, if L2 = O(N ) and τ = O(h), then the following error estimates hold kunh
−
und k0
+k
1/2
X n i=1
kpih
−
pid )k20
1/2
1/2 l X 1/2 6C τ λj , 1 6 n 6 L.
(4.1)
j=d+1
Proof Since M d ⊂ Mh and H d ⊂ Hh , subtracting Problem (VI) from Problem (IV) respectively taking vh = vd ∈ M d and q h = q d ∈ H d yields that (pnh − pnd , q d ) + (q d , ∇(uh − und )) = 0, ∀q d ∈ H d , (uh −
und , vd )
u0h (x, y)
−
−
τ (pnh
u0d (x, y)
=
−
pnd , ∇vd )
u0h (x, y)
−
= (un−1 − un−1 , vd ), ∀vd h d d 0 P uh (x, y), (x, y) ∈ Ω.
(4.2) d
∈M ,
(4.3) (4.4)
By using (3.8), (3.9), (4.2), and (4.3), we have that τ kpnh − pnd k20 = τ (pnh − pnd , pnh − ρd pnh ) + τ (pnh − pnd , ρd pnh − pnd )
= τ (pnh − ρd pnh , pnh − ρd pnh ) − τ (ρd pnh − pnd , ∇(uh − und ))
= τ (pnh − ρd pnh , pnh − ρd pnh ) − τ (ρd pnh − pnh , ∇(uh − und ))
−τ (pnh − pnd , ∇(uh − P d unh )) − τ (pnh − pnd , ∇(P d uh − und ))
= τ (pnh − ρd pnh , pnh − ρd pnh ) − τ (ρd pnh − pnh , ∇(uh − P d unh )) −τ (pnh − ρd pnh , ∇(uh − P d unh )) − (unh − und , P d uh − und )
+(un−1 − un−1 , P d uh − und ) h d
= τ (pnh − ρd pnh , pnh − ρd pnh ) − (unh − und , P d uh − unh )
−(unh − und , uh − und ) + (un−1 − un−1 , P d uh − unh ) h d
+(un−1 − un−1 , uh − und ). h d
(4.5)
By transporting terms for (4.5) and using H¨older inequality, Cauchy inequality, and (3.10), if τ = O(h) we have that τ kpnh − pnd k20 + kunh − und k20
= τ (pnh − ρd pnh , pnh − ρd pnh ) − (unh − und , P d uh − unh )
+(un−1 − un−1 , P d uh − unh ) + (un−1 − un−1 , uh − und ) h d h d
6 τ kpnh − ρd pnh k20 + Chkunh − und k0 k∇(P d uh − unh )k0
+Chkun−1 − un−1 k0 k∇(P d uh − unh )k0 + kun−1 − un−1 k0 kuh − und k0 h d h d
6 Cτ [kpnh − ρd pnh k20 + k∇(P d uh − unh )k20 ] + Cτ kunh − und k20 1 − un−1 k20 + kuh − und k20 ]. +Cτ kun−1 − un−1 k20 + [kun−1 d h d 2 h
(4.6)
That is that τ kpnh − pnd k20 + kunh − und k20
6 2τ kpnh − pnd k20 + kunh − und k20
6 Cτ [kpnh − ρd pnh k20 + k∇(P d uh − unh )k20 ] + Cτ kunh − und k20 +Cτ kun−1 − un−1 k20 + kun−1 − un−1 k20 . h d h d
(4.7)
No.5
Z.D. Luo et al: A REDUCED-ORDER MFE FORMULATION BASED ON POD METHOD
1479
Summing from 1 to n for (4.7), noting that L2 = O(N ) (i.e., τ 1/2 = O(1/L)), and using Lemma 3 could yield that kunh − und k20 + τ
n X i=1
kpih − pid k20 6 Cτ 1/2
l X
λj + Cτ
n X i=1
j=d+1
kuih − uid k20 .
(4.8)
If τ is sufficiently small such that Cτ 6 1/2, we obtain from (4.8) that kunh − und k20 + τ
n X i=1
kpih − pid k20 6 kunh − und k20 + 2τ 6 Cτ 1/2
l X
n X i=1
λj + Cτ
kpih − pid k20 n−1 X i=0
j=d+1
kuih − uid k20 .
(4.9)
Applying Lemma 4 to (4.9) yields that kunh − und k20 + τ
n X i=1
kpih − pid k20 6 Cτ 1/2
l X
λj exp(Cτ n),
(4.10)
j=d+1
which yields (4.1) and completes the proof of Theorem 5. 2 Combining Theorems 1 and 5 yields the following result. Theorem 6 Under hypotheses of Theorem 5, the error estimates between the solution for Problem (II) and the solutions for the reduced-order Problem (VI) are X 21 n 1 ku(tn ) − und k0 + τ 2 kp(ti ) − pid k20 i=1
1
6 C(τ 1+ 2
1/2 l X m+1 1/2 +h )+C τ λj , 1 6 n 6 L. j=d+1
Remark 3 The condition N = O(L2 ) in Theorem 5 shows the relation between the number L of snapshots and the number N of all time instances. Therefore, it is unnecessary to take total transient solutions at all time instances tn as snapshots in [15, 16]. Theorems 5 and 6 have respectively provided the error estimates between the solutions of the reduced-order MFE Problem (VI) and the solutions of the classical MFE Problem (IV) and Problem (II), which have supplied a guide choosing number d of POD basis. Our method here employs the first L MFE solutions (unh , pnh ) (n = 1, 2, · · · , L) for Problem (IV) as assistant analysis. However, when one solves actual problems, he may obtain the ensemble of snapshots from physical system trajectories by drawing samples from experiments with interpolation (or data assimilation) or previous results. Therefore, the assistants (unh , pnh ) (n = 1, 2, · · · , L) could be replaced with the interpolation functions of experimental or previous results, thus rendering it unnecessary to solve Problem (IV) and requiring only to solve directly Problem (VI) such that Theorem 5 is satisfied. And then, time instances are continuously extrapolated forward and POD basis is ceaselessly renewed, the rules of future development and change of natural phenomenon would be very well simulated. 4.2 Implementation of Algorithm for Problem (VI) In the following, we give the implementation of algorithm for solving Problem (VI) which consists of seven steps.
1480
ACTA MATHEMATICA SCIENTIA
Vol.33 Ser.B
Step 1 Generate the snapshots ensemble U i (x, y) = (uih , pih ), i = 1, 2, · · · , L = O(N 1/2 ), which may be the first L solutions for Problem (IV), physical system trajectories by drawing samples from experiments with interpolation (or data assimilation), or previous results. Step 2 Generate the correlation matrix A = (Aik )L×L , where Aik = L1 (U i , U k )H×M , (U i , U k )H×M = (∇uih , ∇ukh ) +(pih , pkh ), and (·, ·) is L2 -inner product. Step 3 Solving the eigenvalue problem Av = λv, v = (a1 , a2 , · · · , aL )T
(4.11)
obtains eigenvalues λ1 > λ2 > · · · > λl > 0 (l = dim{U 1 , U 2 , · · · , U L }) and corresponding eigenvectors v k = (ak1 , ak2 , · · · , akL ) (k = 1, 2, · · · , l). Step 4 For error δ needed, decide on the triangulation parameter h, the time step increment τ , the degree m of interpolation polynomial, and the amount d of POD basis such that 1/2 l P 1 τ 1+ 2 + hm+1 + k 1/2 6 δ. λj j=d+1
Step 5
Generate POD basis ψ k (x, y) = (ψuk (x, y), ψ pk (x, y)) (k = 1, 2, · · · , d): (ψuk (x, y), ψ pk (x, y)) =
L L 1 X k i 1 X k i √ ai u h , √ ai ph . Lλk i=1 Lλk i=1
Step 6 Taking M d = span {ψu1 (x, y), ψu2 (x, y), · · · , ψud (x, y)} and H d =span {ψp1 (x, y), ψ p2 (x, y), · · · , ψ pd (x, y)} and solving Problem (VI) which only includes 3d degrees of freedom yield the solutions (und , pnd ) (n = 1, 2, · · · , L, L + 1, · · · , N ). Step 7 If kun−1 − und k0 > kund − un+1 k0 and kpn−1 − pnd k0 > kpnd − pn+1 k0 (n = L, L + d d d d n n 1, · · · , N −1), then (ud , pd ) (n = 1, 2, · · · , N ) are the solution for Problem (VI), whose errors are 1/2 l P not more than k 1+1/2 + hm+1 + k 1/2 λj . Else, i.e., if kun−1 − und k0 < kund − un+1 k0 d d j=d+1
or kpn−1 − pnd k0 < kpnd − pn+1 k0 (n = L, L + 1, · · · , N − 1), let U i = (uid , pid ) (i = n − L, n − L − d d 1, · · · , n − 1), repeat Step 1 to Step 6 until the solutions (und , pnd ) (n = 1, 2, · · · , L, L + 1, · · · , N ) are obtained.
5
Some Numerical Experiments
In this section, some numerical examples are used to validate the feasibility and efficiency of the reduced-order MFE formulation, i.e., Problem (VI) based on POD method for twodimensional parabolic equations. For the sake of convenience, without loss of generality, herein we consider a two-dimensional parabolic equations by Choosing computing field Ω = [0, 2]×[0, 2], total time T = 0.4, u0 (x, y) = 2 sin πx sin πy, and f (x, y, t) = − sin πx sin πye−2π t in Problem (I) as an example, whose ideas and approaches could directly apply to numerical computation for two-dimensional parabolic equations with real practical applied background. We first divide the field Ω into 2000×2000 small squares with side length △x = △y = 0.001, and then link the diagonal of the square to divide each square into two triangles in the same
No.5
Z.D. Luo et al: A REDUCED-ORDER MFE FORMULATION BASED ON POD METHOD
1481
√ direction which consists of triangulation ℑh . Thus h = 2 × 0.001. In order to make the condition τ = O(h) satisfied in Theorem 5, we take time step size as τ = 0.002. We first find a sequence of numerical solutions (unh , pnh ) ∈ Mh × Hh (n = 1, 2, · · · , 200, i.e., at time t = 1τ, 2τ, · · · , 200τ ) of the classical MFE formulation, i.e., Problem (IV) where Mh is the subspace of piecewise 2-th degree polynomials and Hh the subspace of piecewise 1-th degree polynomial vectors. And then, we chose the first 20 numerical solutions unh (n = 1, 2, · · · , 20) as snapshots U i = (uih , pih ) (i = 1, 2, · · · , 20). Finally, we find 20 eigenvalues which is arranged in a non-increasing order and 20 eigenvectors correspondent to the 20 eigenvalues and we construct 20 POD bases ψ j (j = 1, 2, · · · , 20) with Step 2, 3, and 5 of Subsection 4.2. Take the first 6 POD bases ψ j (j = 1, 2, · · · , 6) from 20 POD bases ψ j = (ψuk (x, y), ψ pk (x, y)) (j = 1, 2, · · · , 20) to span into subspaces M d and H d and compute a numerical solution (und , pnd ) (n = 200, i.e., at t = 200τ ) with the reduced-order MFE formulation, i.e., the reduced-order Problem (VI) according seven steps in Subsection 4.2, without renewing POD basis it is still convergent and depicted graphically on the left charts in Figures 1–3. While the classical MFE solution (unh , pnh ) ∈ Mh × Hh (where Mh is the subspace of piecewise 2-th degree polynomial and Hh the subspace piecewise 1-th degree polynomial vector) is obtained by solving Problem (IV) when n = 200, i.e., at t = 0.4 and is depicted graphically on the right charts in Figures 1–3. Every two charts in Figures 1–3 are exhibiting quasi-identical similarity, respectively, but the reduced-order MFE solutions are better than the classical MFE solutions (it is due to use more initial numerical solutions, i.e., six POD bases which is an open question).
Fig.1
Left chart is the reduced-order MFE solution pn 1d when d = 6 at t = 0.4; right chart is
the classical MFE solution pn 1h at t = 0.4
Fig.2
Left chart is the reduced-order MFE solution pn 2d when d = 6 at t = 0.4; right chart is
the classical MFE solution pn 2h at t = 0.4
1482
ACTA MATHEMATICA SCIENTIA
Fig.3
Vol.33 Ser.B
Left chart is the reduced-order MFE solution un d when d = 6 at t = 0.4; right chart is
the classical MFE solution un h at t = 0.4
Fig.4
When t = 2, the errors (log10 ) between solutions of Problem (VI) with different
number of POD bases for a group of 20 snapshots and the classical MFE formulation Problem (IV)
√ When we take 6 POD bases and τ = 0.002 and h = 2 × 0.001, by computing we h i1/2 20 P obtain that τ 1/2 λj + τ 1+1/2 + h2 6 3 × 10−4 of theoretical errors in Theorem 6. j=7
Figure 4 computationally shows the errors (log10 ) between the reduced-order MFE solutions (und , pnd ) with 20 different POD bases and the classical MFE solution (unh , pnh ) at t = 200τ (i.e., n = 200), respectively. Comparing the classical MFE formulation with the reduced-order MFE formulation containing 6 POD bases implementing the numerical simulation computations when total time t = 200τ , we find that for the classical MFE formulation with piecewise 2-th degree polynomials for unh and piecewise linear polynomials vector for pnh , which has 5 × 2000 × 2000 = 2 × 107 degrees of freedom, the required computing time is 16 minutes, while for the reduced-order MFE formulation with 6 POD bases, which has only 18 degrees of freedom, the corresponding time is only 8 seconds, i.e., the required computing time to solve the classical MFE formulation is about as 120 times as that to do the reduced-order MFE formulation with 6 POD bases, while the errors between their respective solutions do not exceed 3 × 10−4 . It is also shown that finding the approximate solutions for two-dimensional parabolic equations with the reduced-order MFE formulation is computationally very effective. And the results for numerical examples are consistent with those obtained for the theoretical case. Though our examples are in sense recomputing what we have already computed by classical MFE
No.5
Z.D. Luo et al: A REDUCED-ORDER MFE FORMULATION BASED ON POD METHOD
1483
formulation (but only use the first 20 numerical solutions (unh , pnh ) (n = 1, 2, · · · , 20) at the first 20 steps), when we compute actual problems, we may also construct the snapshots and POD basis with interpolation or data assimilation by drawing samples from experiments, then solve directly the reduced-order MFE formulation, while it is unnecessary to solve classical MFE formulation such that the computational load could be alleviated and time-consuming of calculations in the computational process saved.
6
Conclusions
In this paper, we have employed the POD method to establish a new reduced-order MFE formulation for two-dimensional parabolic equations, analyzed the errors between the solutions of its classical MFE formulation and the solutions of the reduced-order MFE formulation based on POD method, and provide the steps of the implementation of algorithm for solving the reduced-order MFE formulation, i.e., Problem (VI), which have shown that our present method has improved and innovated the existing methods (for instance, the methods in [1, 2, 4, 13– 18]). Comparing with the theoretical error estimates with those of numerical examples shows that the results for numerical examples are consistent with those obtained for the theoretical case, further validates both the feasibility and efficiency of the reduced-order MFE formulation. Though snapshots and POD basis of our numerical examples are constructed with the first fewer solutions of the classical MFE formulation, when one computes actual two-dimensional parabolic equations, this process can be omitted in actual applications and one may construct the snapshots and POD basis with interpolation or data assimilation by drawing samples from experiments, then directly solves the reduced-order MFE formulation or Problem (VI), while it is unnecessary to solve Problem (IV) such that the computational load could be alleviated and time-consuming of calculations in the computational process saved. Therefore, the method in this paper holds a good prospect of extensive applications. Future research work in this area will aim at extending the reduced-order MFE formulation, applying it to a set of more complicated PDEs such as the atmosphere quality forecast system and the ocean fluid forecast system. References [1] [2] [3] [4] [5] [6] [7] [8] [9] [10]
Thom´ ee V. Galerkin Finite Element Methods for Parabolic Problems. Berlin: Springer, 1997 Luo Z D. Mixed Finite Element Methods and Applications. Beijing: Science Press, 2006 Brezzi F, Fortin M. Mixed and Hybrid Finite Element Methods. New York: Springer-Verlag, 1991 Li L, Sun P, Luo Z D. A new mixed finite element formulation and error estimates for parabolic equations. Acta Mathematica Scientia, 2012, 32A(6): 1158–1165 Fukunaga K. Introduction to Statistical Recognition. New York: Academic Press, 1990 Jolliffe I T. Principal Component Analysis. Berlin: Springer-Verlag, 2002 Holmes P, Lumley J L, Berkooz G. Turbulence, Coherent Structures, Dynamical Systems and Symmetry. Cambridge, UK: Cambridge University Press, 1996 Sirovich L. Turbulence and the dynamics of coherent structures: Part I-III. Quarterly of Applied Mathematics, 1987, 45(3): 561–590 Lumley J L. Coherent structures in turbulence//Meyer R E, ed. Transition and Turbulence. Academic Press, 1981 Aubry N, Holmes P, Lumley J L, Stone E. The dynamics of coherent structures in the wall region of a turbulent boundary layer. J Fluid Dynamics, 1988, 192: 115–173
1484
ACTA MATHEMATICA SCIENTIA
Vol.33 Ser.B
[11] Moin P, Moser R D. Characteristic-eddy decomposition of turbulence in channel. J Fluid Mech, 1989, 200: 417–509 [12] Rajaee M, Karlsson S K F, Sirovich L. Low dimensional description of free shear flow coherent structures and their dynamical behavior. J Fluid Mech, 1994, 258: 1–29 [13] Joslin R D, Gunzburger M D, Nicolaides R A, Erlebacher G, Hussaini M Y. A self-contained automated methodology for optimal flow control validated for transition delay. AIAA J, 1997, 35(5): 816–824 [14] Kunisch K, Volkwein S. Galerkin proper orthogonal decomposition methods for parabolic problems. Numerische Mathematik, 2001, 90(1): 117–148 [15] Kunisch K, Volkwein S. Galerkin proper orthogonal decomposition methods for a general equation in fluid dynamics. SIAM J Num Anal, 2002, 40(2): 492–515 [16] Sun P, Luo Z D, Zhou Y J. Some reduced finite difference schemes based on a proper orthogonal decomposition technique for parabolic equations. Appl Num Math, 2010, 60(1/2): 15–164 [17] Luo Z D, Chen J, Sun P, et al. Finite element formulation based on proper orthogonal decomposition for parabolic equations. Science in China Series A: Mathematics, 2009, 52(3): 587–596 [18] Luo Z D, Chen J, Xie Z H, et al. A reduced second-order time accurate finite element formulation based on POD for parabolic equations (in Chinese). Science in China Series A: Mathematics, 2011, 41(5): 447–460 [19] Cao Y H, Zhu J, Navon I M, et al. A reduced-order approach to four-dimensional variational data assimilation using proper orthogonal decomposition. Int J Num Methods Fluids, 2007, 53(10): 1571–1583 [20] Tian X J, Xie Z H, Dai A G. An ensemble-based explicit four-dimensional variational assimilation method. J Geophys Res, 2008, 113(D21124), doi: 10.1029/2008JD010358 [21] Tian X J, Xie Z H, Sun Q. A POD-based ensemble four-dimensional variational assimilation method. Tellus, 2011, 63A: 805–816 [22] Adams R A. Sobolev Spaces. New York: Academic Press, 1975 [23] Girault V, Raviart P A. Finite Element Approximations of the Navier-Stokes Equations//Theorem and Algorithms. New York: Springer-Verlag, 1986 [24] Ciarlet P G. The Finite Element Method for Elliptic Problems. Amsterdam: North-Holland, 1978 [25] Rudin W. Functional and Analysis. 2nd ed. McGraw-Hill Companies, Inc, 1973