Reduced-order extrapolation spectral-finite difference scheme based on POD method and error estimation for three-dimensional parabolic equation

Jing AN1, Zhendong LUO2, Hong LI3, Ping SUN1
1 School of Mathematics and Computer Science, Guizhou Normal University,
Guiyang 550001, China
2 School of Mathematics and Physics, North China Electric Power University,
Beijing 102206, China
3 School of Mathematical Sciences, Inner Mongolia University, Hohhot 010021, China
c
Higher Education Press and Springer-Verlag Berlin Heidelberg 2015
Abstract In this study, a classical spectral-finite difference scheme (SFDS)
for the three-dimensional (3D) parabolic equation is reduced by using proper
orthogonal decomposition (POD) and singular value decomposition (SVD).
First, the 3D parabolic equation is discretized in spatial variables by using
spectral collocation method and the discrete scheme is transformed into
matrix formulation by tensor product. Second, the classical SFDS is obtained
by difference discretization in time-direction. The ensemble of data are
comprised with the first few transient solutions of the classical SFDS for the
3D parabolic equation and the POD bases of ensemble of data are generated
by using POD technique and SVD. The unknown quantities of the classical
SFDS are replaced with the linear combination of POD bases and a reducedorder
extrapolation SFDS with lower dimensions and sufficiently high accuracy
for the 3D parabolic equation is established. Third, the error estimates between
the classical SFDS solutions and the reduced-order extrapolation SFDS
solutions and the implementation for solving the reduced-order extrapolation
SFDS are provided. Finally, a numerical example shows that the errors of
numerical computations are consistent with the theoretical results. Moreover,
it is shown that the reduced-order extrapolation SFDS is effective and feasible
to find the numerical solutions for the 3D parabolic equation.
Keywords Singular value decomposition (SVD), proper orthogonal
decomposition (POD) bases, spectral-finite difference scheme (SFDS),
error estimation, parabolic equation
Received August 29, 2012; accepted March 26, 2015
Corresponding author: Zhendong LUO, E-mail: zhdluo@ncepu.edu.cn
1026 Jing AN et al.
MSC 65N12, 65N30, 65M60
1 Introduction
The three-dimensional (3D) parabolic equation has been an extremely vital
significance in the field of physics. For example, it is widely applied to
describe many physical phenomenon such as process of heat conduction and
spread of electromagnetic field. It is not easy to find their exact solutions
for the 3D parabolic equation describing practical engineering problems which
include generally many complex physical nature or irregular computing field.
On the contrary, an efficient approach is to find their numerical solutions.
The spectral-finite difference scheme (SFDS) is one of the most important and
efficient numerical methods to solve 3D parabolic equation due to SFDS with
high precision in spatial directions. However, SFDS includes many degrees of
freedom and is a large-scale system of linear equations at each time level so that
solving such system of linear equations needs a lot of computational load. Thus,
an essential problem is how to alleviate the computational load and save timeconsuming
calculations and resource demands in the computational process in
a way that guarantees sufficiently accurate numerical solutions.
Proper orthogonal decomposition (POD) combined with some numerical
methods can reduce infinite-dimensional partial differential equations (PDEs)
into some adequate approximate models with fewer degrees of freedom and
sufficiently high accuracy so as to alleviate the computational load and memory
requirements savings (see [10]). POD is also known as Karhunen-Lo`eve
expansions in signal analysis and pattern recognition (see [8]), or principal
component analysis in statistics (see [11]), or the method of empirical
orthogonal functions in geophysical fluid dynamics and meteorology (see [34]).
POD method has been successfully applied in many fields. In fluid
dynamics, Lumley [18] first employed the POD technique to capture the large
eddy coherent structures in a turbulent boundary layer. This technique was
further extended in [4], where a link between the turbulent structure and the
dynamics of a chaotic system was investigated. The method of snapshot was
introduced by Sirovich [36] and was widely used in applications to reduce the
order of the POD eigenvalue problem. Examples of these were applications
to optimal flow control problems [12] and turbulence [33]. POD method was
also applied to data compression, signal analysis (see [2]), simulation, and
control in chemistry reaction system (see [9,35]). Until resent ten years, some
researches that POD method was used to reduce some Galerkin methods of
PDEs have not been presented (see [13,14]). More recently, some reduced
finite difference models and finite element formulations have been studied
(see [5,6,17,19–32,37]). Moreover, there are related works available for POD
applications in optimization, for instance, adaptive POD [1], trust-region-POD
[3], optimality systems-POD [15], and POD-a-posteriori error estimates [39].
However, almost all existing reduced-order numerical computational
Reduced-order extrapolation spectral-finite difference scheme 1027
methods (see, e.g., [1–3,5,6,13–15,17,19–32,37]) based on POD technique
utilize numerical solutions obtained from classical numerical methods on the
total time span [0, T] to construct POD bases and reduced-order models, and
then recompute the solutions on the same time span [0, T], which actually
belongs to re-computations on the same time span [0, T]. Especially, to the
best of our knowledge, there are no published results addressing that POD
method is used to reduce SFDS for the 3D parabolic equation or providing error
estimates between the classical SFDS solutions and the reduced SFDS
solutions based on POD method. In this work, we establish a reduced-order
extrapolation SFDS based on POD method and singular value decomposition
(SVD) for the 3D parabolic equation. It is different existing some works. First,
the 3D parabolic equation is discretized in spatial directions by using spectral
collocation method and the discrete scheme is transformed into matrix
formulation by tensor product. Second, a classical SFDS is obtained by
difference discretization in time-direction. The ensemble of data are comprised
with the first few transient solutions of the classical SFDS for the 3D parabolic
equation and the POD bases of ensemble of data are generated by using POD
technique and SVD. The unknown quantities of SFDS are replaced with the
linear combination of POD bases and a reduced-order extrapolation model,
i.e., a reduced-order extrapolation SFDS with lower dimensions and sufficiently
high accuracy for the 3D parabolic equation is established. Third, the error
estimates between the classical SFDS solutions and the reduced-order
extrapolation SFDS solutions and the implementation for solving the reducedorder
extrapolation SFDS are provided. Finally, a numerical example has shown
that the errors of numerical computations are consistent with the theoretical
results. Moreover, it is shown that the reduced-order extrapolation SFDS is
effective and feasible to find the numerical solutions for the 3D parabolic
equation.
2 Classical SFDS for 3D parabolic equation and generation of snapshots
Spectral method is one of the ‘big three’ technologies for the numerical solution
of PDEs (see [38]). Naturally, the origins of each technology can be traced
further back. For spectral method, some of the ideas are as old as interpolation
and expansion, and more specifically algorithmic developments arrived with
Lanczos as early as 1938 (see [16]) and with Clenshaw et al. and others in the
1960s (see [7]). Then, in the 1970s, a transformation of the field took place
initiated by work by Orszag and others on problems in fluid dynamics and
meteorology, and spectral method became famous. In particular, collocation
method is the most popular form of the spectral methods among practitioners.
It is very easy to implement and generally leads to satisfactory results as long as
the problems possess sufficient smoothness. Collocation method is derived by
asking that the approximation solutions satisfy exactly the boundary conditions
and the equation at the interior collocation points. In this work, we mainly take
1028 Jing AN et al.
the following 3D parabolic equation into account by use of collocation method.
⎧⎪⎪⎨
⎪⎪⎩
∂u
∂t
− αΔu = f, (x, y, z, t) ∈ Ω × (0, T),
u|∂Ω = 0, t∈ (0, T),
u(x, y, z, 0) = u0(x, y, z), (x, y, z) ∈ Ω,
(1)
where u = u(x, y, z, t) is unknown function, α > 0 is a given positive parameter,
f(x, y, z, t) and u0(x, y, z) are two given functions, T is the total time, and for
the sake of convenience, without loss of generality, we may as well suppose
Ω = (−1, 1)3.
First, the 3D parabolic equation is discretized in spatial directions by using
spectral collocation method. For the sake of simplicity, we shall use the same
number of points N at the x, y, and z directions, although in practical
applications, one may wish to use different numbers of points at each
direction. Let PN denote the space of all algebraic polynomials of degree less
than or equal to N, let
XN = {p ∈ PN × PN × PN : p|∂Ω = 0},
and let {ξi}N
i=0 be Legendre Gauss-Lobatto points. Then the Legendre
collocation method is to look for uN ∈ C1(0, T;XN) such that, for any t ∈
(0, T),
⎧⎨
⎩
∂uN(ξi, ξj, ξk, t)
∂t
− αΔuN(ξi, ξj, ξk, t) = f(ξi, ξj, ξk, t), 1 i, j, k N − 1,
uN(ξi, ξj, ξk, 0) = u0(ξi, ξj, ξk), 0 i, j, k N.
(2)
Let
U(t) = (uN(ξ1, ξ1, ξ1, t), . . . , uN(ξN−1, ξ1, ξ1, t), ξ1, t), . . . ,
uN(ξN−1, ξN−1, ξN−1, t))T,
U0 = (u0(ξ1, ξ1, ξ1), . . . , u0(ξN−1, ξ1, ξ1), . . . , u0(ξN−1, ξN−1, ξ1), . . . ,
u0(ξN−1, ξN−1, ξN−1))T,
F(t) = (f(ξ1, ξ1, ξ1, t), . . . , f(ξN−1, ξ1, ξ1, t), . . . , f(ξN−1, ξN−1, ξN−1, t))T.
Then (2) can be transformed into following matrix form by tensor product
⎧⎨
⎩
∂U(t)
∂t
− α(I ⊗ I ⊗ D2 +I ⊗ D2 ⊗ I + D2 ⊗ I ⊗ I)U(t) = F(t),
U(0) = U0,
(3)
where D2 is the second-order differentiation matrix of order N−1 forChebyshev
case (see [40]), I is the unit matrix of order N − 1, and ⊗ denotes the tensor
product of matrixes, i.e.,
A ⊗ B = (Abi,j)N−1
i,j=1.
Reduced-order extrapolation spectral-finite difference scheme 1029
Now, a classical SFDS will be obtained by difference discretization in timedirection.
Let τ be the time step increment, tn = nτ (n = 0, 1, . . .,M,M =
[T/τ]). The implicit difference scheme of (3) is as follows:
Un+1 − Un
τ
− α(I ⊗ I ⊗ D2 +I ⊗ D2 ⊗ I + D2 ⊗ I ⊗ I)Un+1 = Fn+1, (4)
where
Un = (uN(ξ1, ξ1, ξ1, tn), . . . , uN(ξN−1, ξ1, ξ1, tn), . . . ,
uN(ξN−1, ξN−1, ξN−1, tn))T,
Fn = (f(ξ1, ξ1, ξ1, tn), . . . , f(ξN−1, ξ1, ξ1, tn), . . . , f(ξN−1, ξN−1, ξ1, tn), . . . ,
f(ξN−1, ξN−1, ξN−1, tn))T.
Then, we get by simplifying (4) that
AUn+1 = Un + τFn+1, (5)
where
A = I1 − ατ(I ⊗ I ⊗ D2 + I ⊗ D2 ⊗ I + D2 ⊗ I ⊗ I),
I1 is the unit matrix of order (N −1)3. Thus, if α > 0, f, and τ are given, then
we can get the set of SFDS solutions, i.e., {U0, U1, . . . , UM} by solving (5). A
subset, i.e., {Ui}L
i=0 (usually, L M), known as snapshots which are useful
and of interest to us, is chosen from first L elements of {U0, U1, . . . , UM}.
It is well known (see [40]) that the errors between the solution U(x, y, z, t)
to equation (1) and the solutions un to equation (4) are as follows:
U(tn) − Un 0 = O(τ,N
−3), (6)
where · 0 is the standard norm of vector.
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 date assimilation). For example, when
one finds numerical solutions to PDEs representing weather forecast, one can
use the previous weather result, then restructure the optimal basis for the
ensemble of snapshots by following POD technique, and finally combine it
with spectral method to derive a reduced order dynamical system. Thus, the
forecast of future weather change can be quickly simulated, which is of major
importance for actual real-life applications.
3 FormulatePODbases
Let Au = (U0, U1, . . . , UL). For matrix Au, there exists the SVD as follows:
Au = Uu

Su O
O O

V T
u ,
1030 Jing AN et al.
where Uu ∈ R(N−1)3
, Vu ∈ R(L+1)×(L+1) are all orthogonal matrices. Su =
diag(σ1, σ2, . . . , σr) ∈ Rr×r is the diagonal matrix correspondent to Au. σi > 0
(i = 1, 2, . . . , r) are the positive singular values. The matrices Uu = (φ1, φ2, . . . ,
φ(N−1)3) and Vu = (ϕ1, ϕ2, . . . , ϕL+1) consist of the orthogonal eigenvectors to
the AuATu
and ATu
Au, respectively. The columns of these eigenvector matrices
are organized so as to correspond to the singular values σi comprised in Su
in a non-increasing order. And the singular values of the decomposition are
connected to the eigenvalues of the matrices AuATu
and ATu
Au in a manner such
that λi = σ2
i (i = 1, 2, . . . , r). Since the number of mesh points is far larger
than that of transient moment points, i.e.,
(N − 1)3 L + 1,
that is also that the order (N − 1)3 for matrix AuATu
is far larger than the
order L+1 for matrix ATu
Au, however, their non-zero eigenvalues are identical,
therefore, we may first solve the eigenequation corresponding to matrix ATuAu
to find the eigenvectors ϕj (j = 1, 2, . . . , L + 1), and then by relationship
φj = Auϕj
σj
, j= 1, 2, . . . , r,
we may obtain r (r L + 1) eigenvectors corresponding to the non-zero eigenvalues
for matrix AuATu
.
Define matrix norm
A a,β = sup
x =0
Ax α
x β
.
Let
AMu =
Mu
i=1
σiφiϕTi
,
where σi (i = 1, 2, . . .,Mu) are the first Mu diagonal elements of matrix Su, φi
and ϕi (i = 1, 2, . . .,Mu) be the first Mu column vectors of matrices Uu and
Vu, respectively. According the relationship properties of spectral radius and
· 2,2 for matrix, if
Mu < r = rank(Au),
then there holds the following formula (see [26,27]):
min
rank(B) Mu
Au − B 2,2 = Au − AMu
2,2 = σ(Mu+1) =

λ(Mu+1) , (7)
which shows that AMu is an optimal representation of Au. Denote the L + 1
column vectors of matrices Au by
al
u = (ul
1, ul
2, . . . , ul
(N−1)2 )T, l= 0, 1, . . . , L,
Reduced-order extrapolation spectral-finite difference scheme 1031
by εl (l = 1, 2, . . . , L+1) the unit column vectors except that a vector component
is 1, while the other components are 0. Then, by the compatibility of the norm
for matrixes and vectors, we obtain
al
u
− PMu(al
u) 2 = (Au − AMu)εl 2
(Au − AMu) 2,2 εl 2
= σ(Mu+1)
=

λ(Mu+1) , (8)
where
PMu(al
u) =
Mu
j=1
(φj, al
u)φj ,
and (φj, al
u) is the canonical inner product for vectors φj and al
u. Inequality (8)
shows that PMu(al
u) is the optimal approximation to al
u. It follows that a group
of optimal base Φu = (φ1, φ2, . . . , φMu) is obtained from (7) and (8), which is
known as a group of POD base.
4 Reduced-order extrapolation SFDS
If we approximate Un with Φuαn
Mu, i.e.,
Un ≈ Φuαn
Mu, (9)
where αn
Mu is a constant vector, then, by inserting (9) into (5), we could obtain
an approximating formula of (5) as follows:
Φuαn+1
Mu
≈ A
−1Φuαn
Mu + τA
−1Fn+1. (10)
Then we may obtain an iterative scheme of least square approximation of (10)
by multiplying ΦTu
as follows:
βn
Mu = ΦTu
Un, 0 n L, (11)
ΦTu
AΦuβn+1
Mu = βn
Mu +ΦTu
τFn+1, L n M − 1, (12)
U
∗n
Mu = Φuβn
Mu, 1 n M. (13)
The system of equations (11)–(13) is known as a reduced-order extrapolation
SFDS based on POD technique, U∗n
Mu is known as the reduced-order solution
of SFDS based on POD technique, it has no repeating computations.
Especially, the method here is totally different from the existing reduced-order
numerical methods based on POD technique, whose snapshots are taken from
total transient solutions at all time instances tn (1 n M) or taken one
from every ten transient solutions at uniform intervals (see, e.g., [1–3,5,6,13–
15,17,19–32,37]), while snapshots here are taken from the first L + 1 transient
1032 Jing AN et al.
solutions (L M (N −1)3), it is used existing data (the first L+1 transient
solutions) to simulate future development and change of natural phenomenon
(all M solutions, L M), which is an important computational method. Since
the classical SFDS (5) includes (N − 1)3 linear equations on each time level,
while the reduced-order extrapolation SFDS (11)–(13) based on POD technique
includes only Mu unknowns (Mu (N − 1)3). Thus, the reduced-order
extrapolation SFDS (11)–(13) could greatly save the computing time in the
computational process.
5 Error estimates of reduced-order SFDS solutions and implementation for
solving reduced-order extrapolation SFDS
5.1 Error estimates of reduced-order SFDS solutions
In this subsection, we devote to deriving the error estimates for the reducedorder
SFDS solutions above section. To this end, we first introduce two lemmas.
Lemma 1 AMu = ΦuΦTu
Au.
Proof Owing to
AMu =
Mu
i=1
σiφiϕTi
= (φ1, φ2, . . . , φMu)diag(σ1, σ2, . . . , σMu)(ϕ1, ϕ2, . . . , ϕMu)T,
Au = (φ1, φ2, . . . , φr)diag(σ1, σ2, . . . , σr)(ϕ1, ϕ2, . . . , ϕr)T,
we have
ΦuΦTu
Au = Φu
⎛
⎜⎝
φT1
...
φT
Mu
⎞
⎟⎠
(φ1, φ2, . . . , φr)diag(σ1, σ2, . . . , σr)(ϕ1, ϕ2, . . . , ϕr)T
= Φu(IMu O)diag(σ1, σ2, . . . , σr)(ϕ1, ϕ2, . . . , ϕr)T
= Φu(diag(σ1, σ2, . . . , σMu) O)(ϕ1, ϕ2, . . . , ϕr)T
= (φ1, φ2, . . . , φMu)diag(σ1, σ2, . . . , σMu)(ϕ1, ϕ2, . . . , ϕMu)T
= AMu,
which completes the proof.
Lemma 2 al
u
− ΦuΦTu
al
u
2 σ(Mu+1) =

λ(Mu+1) (l = 0, 1, . . . , L).
Proof By using Lemma 1, we obtain AMu = ΦuΦTu
Au. Then, we have
ΦuΦTu
al
u = AMuεl = PMu(al
u).
We derive from (8) that
al
u
−ΦuΦTu
al
u
2 = al
u
−PMu(al
u) 2 σ(Mu+1) =

λ(Mu+1), l= 0, 1, . . . , L,
Reduced-order extrapolation spectral-finite difference scheme 1033
which completes the proof.
We have the following main result for the reduced-order extrapolation SFDS
(11)–(13).
Theorem 1 Let Un (n = 0, 1, . . .,M) be classical SFDS solutions, and let
U∗n
Mu (n = 0, 1, . . .,M) be the reduced SFDS solutions. Then we have
Un − U
∗n
Mu
∞ C(n)

λ(Mu+1) , (14)
where
C(n) =

1, 1 n L,
A−1 n−L
∞ , L+ 1 n M.
Proof By using (5), we have
Un+1 = A
−1Un + τA
−1Fn+1.
While, by using (12), we have
U
∗n+1
Mu
= A
−1U
∗n
Mu + τA
−1Fn+1.
Let
en = Un − U
∗n
Mu.
Then we have
en+1 ∞ = A
−1en ∞ A
−1 ∞ en ∞ (L + 1 n M − 1).
Thus, we obtain that
en ∞ A
−1 n−L
∞ eL ∞ (L n M). (15)
Using Lemma 2, we obtain that
en ∞ = Un − U
∗n
Mu
∞
= Un − ΦuΦTu
Un ∞
Un − ΦuΦTu
Un 2


λ(Mu+1) , 1 n L. (16)
Then it follows from (15) and (16) that
en ∞ A
−1 n−L
∞

λ(Mu+1), L+ 1 n M, (17)
which completes the proof.
Combining (6) with Theorem 1 yields the following result.
1034 Jing AN et al.
Theorem 2 Let U be the solution to equation (1), and let U∗n
Mu (n = 0, 1, . . . ,
M) be the reduced SFDS solutions. Then we have
U(tn) − U
∗n
Mu
∞ = O

τ,N
−3,C(n)

λ(Mu+1)

, (18)
where
C(n) =

1, 1 n L,
A−1 n−L
∞ , L+ 1 n M.
Remark 2 The error estimates of Theorems 1 and 2 have given the guides to
determine the number Mu of POD basic satisfying

λ(Mu+1) = O(τ,N
−3)
and whether or not need to renew POD bases in the computational process,
i.e., if
A
−1 n−L
∞

λ(Mu+1) = O(τ,N
−3) (L + 1 n M),
then we do not need to renew POD bases, otherwise, we need.
5.2 Implementation for solving reduced-order extrapolation SFDS
In this subsection, we provide the implementation for solving the reduced-order
extrapolation SFDS, which consists of the following 5 steps.
Step 1 For given positive parameter α > 0, source term f(x, y, z, t), initial
value function u0(x, y, z), the time step increment τ, the space of all algebraic
polynomials of degree no more than N :
XN = {p ∈ PN × PN × PN : p|∂Ω = 0},
and Legendre Gauss-Lobatto points {ξi}N
i=0, solving the classical FD scheme at
the first fewer L steps (in usually, take L = 20):
AUn+1 = Un + τFn+1, 0 n L − 1,
where
A = I1 − ατ(I ⊗ I ⊗ D2 + I ⊗ D2 ⊗ I + D2 ⊗ I ⊗ I),
I1 is the unit matrix of order (N − 1)3, D2 the second-order differentiation
matrix of order N −1 for Chebyshev case (see [40]), I the unit matrix of order
N − 1,
A ⊗ B = (Abi,j)N−1
i,j=1,
Un = (uN(ξ1, ξ1, ξ1, tn), . . . , uN(ξN−1, ξ1, ξ1, tn), . . . , uN(ξN−1, ξN−1, ξ1, tn),
. . . , uN(ξN−1, ξN−1, ξN−1, tn))T,
Fn = (f(ξ1, ξ1, ξ1, tn), . . . , f(ξN−1, ξ1, ξ1, tn), . . . , f(ξN−1, ξN−1, ξ1, tn), . . . ,
f(ξN−1, ξN−1, ξN−1, tn))T,
Reduced-order extrapolation spectral-finite difference scheme 1035
yields the set of classical SFDS solutions, i.e., {Ui}L
i=0 (usually, L M).
Step 2 Formulate the snapshot matrixes
Au = (U0, U1, . . . , UL)
and solve the linear system of equation
(ATu
Au − λuIL+1)ϕu = 0
obtaining, respectively, the eigenvalues
λ1 λ2 ... λ Mu
> 0
( Mu = rankAu) and corresponding eigenvectors ϕuj (j = 1, 2, ..., Mu).
Step 3 For the error μ = O(τ,N−3) needed, determine the numberMu (Mu
Mu) of POD bases such that

λu(Mu+1) μ and formulate the POD bases
Φu = (φu1, φu2, ..., φuMu) (where φuj = Auϕuj/

λuj , j = 1, 2, ...,Mu).
Step 4 Solve the reduced-order extrapolation SFDS
βn
Mu = ΦTu
Un, 0 n L,
ΦTu
AΦuβn+1
Mu
= βnM
u +ΦTu
τFn+1, L n M − 1,
U
∗n
Mu = ΦuβnM
u, 1 n M.
obtaining the reduced-order solution vectors
U
∗n
Mu = (u
∗
N(ξ1, ξ1, ξ1, tn), . . . , u
∗
N(ξN−1, ξ1, ξ1, tn), . . . , u
∗
N(ξN−1, ξN−1, ξ1, tn),
. . . , u
∗
N(ξN−1, ξN−1, ξN−1, tn))T (0 n M).
Step 5 If
C(n)

λu(Mu+1) μ,
then U∗n
Mu (n = 1, 2, ...,M) are just solutions satisfying accuracy needed. Else,
i.e., if
C(n)

λu(Mu+1) > μ,
then put Ul = U∗l
Mu (l = n − L − 1, n − L, ..., n − 1), return Step 2.
6 Numerical example
To verify the efficiency of the method, we take
u0(x, y, z) = sinπx sin πy sin πz,
f(x, y, z, t) = (3π2 − 1) sin πx sin πy sin πz e−t,
(x, y, z) ∈ [−1, 1]3, t∈ [0, 10], α= 1,
1036 Jing AN et al.
in equation (1). It is easy to obtain the accurate solution as follows:
u(x, y, z, t) = sinπx sin πy sin πz e−t.
Let us denote the accurate solution at discrete points by
u = ( U1, U2, . . . , UM),
where
Ui = (u(ξ1, ξ1, ξ1, ti), . . . , u(ξN−1, ξ1, ξ1, ti), . . . , u(ξN−1, ξN−1, ξ1, ti), . . . ,
u(ξN−1, ξN−1, ξN−1, ti))T.
Denote
u1 = (U1, U2, . . . , UM)
the classical SFDS solution at the same discrete points and
uMu
2 = (U
∗1
Mu, U
∗2
Mu, . . . , U
∗M
Mu )
the reduced-order extrapolation SFDS solution based on POD method at the
same discrete points, which is found with the reduced-order extrapolation SFDS
according five steps of implementation of algorithm in Section 5.2 with the first
L = 10 snapshots. Let
e( u, u1) = max
1 i,j,k N−1, 1 n M
|( u − u1)n
i,j,k
|,
e( u, uMu
2 ) = max
1 i,j,k N−1, 1 n M
|( u − uMu
2 )n
i,j,k
|,
e( u1, uMu
2 ) = max
1 i,j,k N−1,1 n M
|( u1 − uMu
2 )n
i,j,k
|.
Let t1 be the computing time of classical SFDS solution, and let t2 be the
computing time of the reduced-order SFDS solution based on POD method.
When τ = 0.0001 and N = 15, the errors between accurate solution u and
the reduced-order SFDS solution u2 at the same discrete points for different
numbers of POD bases are
e( u, uMu
2 ) = 1.48 × 10−4, Mu = 1, 2, 3. (19)
Since the error
e( u, u1) = 1.48 × 10−4
between accurate solution and classical SFDS solution which is obtained by
computing, it follows that the error between accurate solution and classical
SFDS solution is very proximate to the errors between accurate solution and
the reduced-order SFDS solutions based on POD method with Mu = 1, 2, 3.
Table 1 shows the errors between the classical SFDS solution u1 and the
reduced-order SFDS solutions u2 based on POD method with Mu = 1, 2,3 and
their computing time.
Reduced-order extrapolation spectral-finite difference scheme 1037
Table 1 Error between SFDS solutions and reduced SFDS solutions
and computing time
Mu
1 2 3
e(u1, uMu
2 ) 2.96 × 10
−12 2.08 × 10
−12 2.53 × 10
−13
t2/s 3.43 3.62 3.7606
Since the computing time t1 = 553.2 s of classical SFDS solution when
T = 10, it follows from Table 1 that the computing time finding the reducedorder
SFDS solution based on POD method can be saved greatly and since
the reduced-order extrapolation SFDS only includes Mu degrees of freedom at
every time level, the computational load can be alleviated greatly, too.
(19) and Table 1 show that the numerical results are consistent with the
theoretical result (10−4), and the error between the accurate solution and the
classical SFDS solution and the error between the accurate solution and the
reduced-order SFDS solutions based on POD method are all very small, while
the errors between the classical SFDS solution and the reduced-order SFDS
solution based on POD method is smaller, but the computing time of the
reduced-order SFDS solutions based on POD method can be saved greatly.
Therefore, the reduced-order SFDS solution based on POD method is effective
and feasible to solve high dimensional problems.
7 Conclusions and perspective
In this study, we have employed POD method and SVD to study the classical
SFDS for the 3D parabolic equation. First, the ensemble of data is obtained
by using the classical SFDS of 3D parabolic equation. Then the POD bases
of the ensemble of data are generated by using SVD and POD technique.
The unknown quantities of the classical SFDS are replaced with the linear
combination of POD bases and the reduced-order extrapolation model, i.e.,
the reduced-order extrapolation SFDS with lower dimensions and sufficiently
high accuracy for the 3D parabolic equation is established. Finally, the error
estimates of the reduced-order extrapolation SFDS solutions based on POD
method is provided, which can be used to determine the number Mu of POD
basic and whether or not to renew POD bases in the computational process.
The numerical example has shown that the numerical errors between the
classical SFDS solutions and the reduced-order SFDS solutions based on POD
method are consistent with the theoretical results. Moreover, it is also shown
that the reduced-order extrapolation SFDS based on POD method is effective
and feasible to find the numerical solutions for 3D parabolic equation.
Especially, we here thoroughly improve the existing methods, where we do
only use the first fewer given numerical solutions on very short time span [0, T0]
(T0 T) as snapshots to formulate the POD basis and establish the reducedorder
extrapolation SFDS based on POD technique for finding the numerical
1038 Jing AN et al.
solutions on total time span [0, T]. Thus, we sufficiently adopt the advantage
of POD method, namely, sufficiently utilize the given data (on very short time
span [0, T0] and T0 T) to forecast future physic phenomenons (on time span
[T0, T]). For some practical problems, such as in order to simulate numerical
weather forecast corresponding to weather forecast PDEs, it is only necessary
to structure the POD bases by using the previous data information as snapshots,
and then solve the reduced-order extrapolation SFDS without solving
the classical SFDS. Thus, we can forecast quite precisely the future natural
phenomena with existing information and can greatly save the computing time
solving the reduced-order extrapolation SFDS.
Acknowledgements This work was supported in part by the National Natural Science
Foundation of China (Grant Nos. 11271127, 11361035), the Doctoral Foundation of Guizhou
Normal University, and the Science and Technology Fund of Guizhou Province (Grant No.
7052) in 2014.
References
1. Afanasiev K, Hinze M. Adaptive control of a wake flow using proper orthogonal
decomposition. Lect Notes Pure Appl Math, 2001, 216: 317–332
2. Algazi V, Sakrison D. On the optimality of Karhunen-Lo`eve expansion. IEEE Trans
Inform Theory, 1969, 15: 319–321
3. Arian E, Fahl M, Sachs E W. Trust-region proper orthogonal decomposition models by
optimization method. In: Proceedings of the 41st IEEE Conference on Decision and
Control, Las Vegas, Nevada, 2002. 2002, 3300–3305
4. 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
5. Cao Y H, Zhu J, Luo Z H, Navon I M. Reduced order modeling of the upper
tropical Pacific Ocean model using proper orthogonal decomposition. Comput Math
Appl, 2006, 52: 1373–1386
6. Cao Y H, Zhu J, Navon I M, Luo Z D. A reduced order approach to four-dimensional
variational data assimilation using proper orthogonal decomposition. Int J Numer
Meth Fluids, 2007, 53: 1571–1583
7. Fox L, Parker I B. Chebyshev Polynomials in Numerical Analysis. Oxford: Oxford
University Press, 1968
8. Fukunaga K. Introduction to Statistical Recognition. New York: Academic Press, 1990
9. Graham M, Kevrekidis I. Alternative approaches to the Karhunen-Lo`eve decomposition
for model reduction and data analysis. Comput Chem Eng, 1996, 20: 495–506
10. Holmes P, Lumley J L, Berkooz G. Turbulence, Coherent Structures, Dynamical
Systems and Symmetry. Cambridge: Cambridge University Press, 1996
11. Jolliffe I T. Principal Component Analysis. Berlin: Springer-Verlag, 2002
12. Joslin R D, Gunzburger M D, Nicolaides R, Erlebacher G, Hussaini M Y. A selfcontained
automated methodology for optimal flow control validated for transition
delay. AIAA Journal, 1997, 35: 816–824
13. Kunisch K, Volkwein S. Galerkin proper orthogonal decomposition methods for
parabolic problems. Numer Math, 2001, 90: 117–148
14. Kunisch K, Volkwein S. Galerkin proper orthogonal decomposition methods for a
general equation in fluid dynamics. SIAM J Numer Anal, 2002, 40(2): 492–515
15. Kunisch K, Volkwein S. Proper orthogonal decomposition for optimality systems.
ESAIM: Math Model Numer Anal, 2008, 42(1): 1–23
Reduced-order extrapolation spectral-finite difference scheme 1039
16. Lanczos C. Trigonometric interpolation of empirical and analytical functions. J Math
Phys, 1938, 17: 123–199
17. Li H R, Luo Z D, Chen J. Numerical simulation based on proper orthogonal decomposition
for two-dimensional solute transport problems. Appl Math Model, 2011, 35(5):
2489–2498
18. Lumley J L. Coherent structures in turbulence. In: Meyer R E, ed. Transition and
Turbulence. Proceedings of the Symposium on Transition and Turbulence in Fluids,
Madison, WI, October 13–15, 1980. New York: Academic Press, 1981, 215–242
19. Luo Z D, Chen J, Navon I M, Yang X Z. Mixed finite element formulation and error
estimates based on proper orthogonal decomposition for the non- stationary Navier-
Stokes equations. SIAM J Numer Anal, 2008, 47(1): 1–19
20. Luo Z D, Chen J, Navon I M, Zhu J. An optimizing reduced PLSMFE formulation
for non-stationary conduction-convection problems. Int J Numer Meth Fluids, 2009,
60(4): 409–436
21. Luo Z D, Chen J, Sun P, Yang X Z. Finite element formulation based on proper
orthogonal decomposition for parabolic equations. Sci China Ser A: Math, 2009, 52(3):
585–596
22. Luo Z D, Chen J, Zhu J, Wang R W, Navon I M. An optimizing reduced order FDS
for the tropical Pacific Ocean reduced gravity model. Int J Numer Meth Fluids, 2007,
55(2): 143–161
23. Luo Z D, Du J, Xie Z H, Guo Y. A reduced stabilized mixed finite element formulation
based on proper orthogonal decomposition for the no-stationary Navier-Stokes
equations. Int J Numer Meth Eng, 2011, 88(1): 31–46
24. Luo Z D, Li H, Zhou Y J, Huang X M. A reduced FVE formulation based on POD
method and error analysis for two-dimensional viscoelastic problem. J Math Anal Appl,
2012, 385: 310–321
25. Luo Z D, Li H, Zhou Y J, Xie Z H. A reduced finite element formulation and error
estimates based on POD method for two-dimensional solute transport problems.
J Math Anal Appl, 2012, 385: 371–383
26. Luo Z D, Ou Q L, Xie Z X. A reduced finite difference scheme and error estimates based
on POD method for the non-stationary Stokes equation. Appl Math Mech, 2011, 32(7):
847–858
27. Luo Z D, Wang R W, Zhu J. Finite difference scheme based on proper orthogonal
decomposition for the non-stationary Navier-Stokes equations. Sci China Ser A: Math,
2007, 50(8): 1186–1196
28. Luo Z D, Xie Z H, Chen J. A reduced MFE formulation based on POD for the nonstationary
conduction-convection problems. Acta Math Sci Ser B Engl Ed, 2011, 31(5):
1765–1785
29. Luo Z D, Xie Z H, Shang Y Q, Chen J. A reduced finite volume element formulation
and numerical simulations based on POD for parabolic equations. J Comput Appl
Math, 2011, 235(8): 2098–2111
30. Luo Z D, Yang X Z, Zhou Y J. A reduced finite difference scheme based on singular
value decomposition and proper orthogonal decomposition for Burgers equation.
J Comput Appl Math, 2009, 229(1): 97–107
31. Luo Z D, Zhou Y J, Yang X Z. A reduced finite element formulation based on proper
orthogonal decomposition for Burgers equation. Appl Numer Math, 2009, 59(8): 1933–
1946
32. Luo Z D, Zhu J, Wang R W, Navon I M. Proper orthogonal decomposition approach
and error estimation of mixed finite element methods for the tropical Pacific Ocean
reduced gravity model. Comput Meth Appl Mech Eng, 2007, 196(41-44): 4184–4195
33. Rajaee M, Karlsson S K F, Sirovich L. Low dimensional description of free sheer flow
coherent structures and their dynamical behavior. J Fluid Mech, 1994, 258: 1401–1402
34. Selten F. Baroclinic empirical orthogonal functions as basis functions in an atmospheric
model. J Atmospheric Sci, 1997, 54: 2100–2114
1040 Jing AN et al.
35. Shvartsman S, Kevrekisis I. Low-dimensional approximation and control of periodic
solutions in spatially extended systems. Phys Rev E, 1998, 58(3): 361–368
36. Sirovich L. Turbulence and the dynamics of coherent sructures: part I–III. Quart Appl
Math, 1987, 45(3): 561–590
37. Sun P, Luo Z D, Zhou Y J. Some reduced finite difference schemes based on a proper
orthogonal decomposition technique for parabolic equations. Appl Numer Math, 2010,
60(1-2): 154–164
38. Trefethen L N. Spectral Method in MATLAB. Philadephia: SIAM, 2000
39. Tr¨ultzsch F, Volkwein S. POD a-posteriori error estimates for linear-quadratic optimal
control problems. Comput Optim Appl, 2009, 44(1): 83–115
40. Weideman J A C, Reddy S C. A Matlab differentiation matrix suite. ACM Trans Math
Software, 2000, 26: 465–511