拟谱法 (Pseudo-spectral methods)[ 1] 也称为离散变数表示 (discrete variable representation、DVR)法,是在应用数学 及计算科学 中求解偏微分方程 用的数值分析 方法。拟谱法和谱方法 有密切关系,但在谱方法中基底 函数中使用了拟谱的基底函数,也就是可以在分割网格上表示的函数。此作法简化一些运算子的计算,在使用快速演算法(例如快速傅里叶变换 )时可以加速计算速度。
说明
考虑以下的初值问题
i
∂
∂
t
ψ
(
x
,
t
)
=
[
−
∂
2
∂
x
2
+
V
(
x
)
]
ψ
(
x
,
t
)
,
ψ
(
t
0
)
=
ψ
0
{\displaystyle i{\frac {\partial }{\partial t}}\psi (x,t)={\Bigl [}-{\frac {\partial ^{2}}{\partial x^{2}}}+V(x){\Bigr ]}\psi (x,t),\qquad \qquad \psi (t_{0})=\psi _{0}}
有周期性条件
ψ
(
x
+
1
,
t
)
=
ψ
(
x
,
t
)
{\displaystyle \psi (x+1,t)=\psi (x,t)}
。这个例子是势能为
V
(
x
)
{\displaystyle V(x)}
之粒子的薛定谔方程 ,不过其结构可适用到其他应用。在许多实务上的偏微分方程中,其中有一项是和导数(例如是动能相关的量)有关,另一项则是和另一个函数(此处为势能)的乘积。
在谱方法中,其解
ψ
{\displaystyle \psi }
会展开为一组适合基底函数的组合,例如平面波。
ψ
(
x
,
t
)
=
1
2
π
∑
n
c
n
(
t
)
e
2
π
i
n
x
.
{\displaystyle \psi (x,t)={\frac {1}{\sqrt {2\pi }}}\sum _{n}c_{n}(t)e^{2\pi inx}.}
将解代入,并且计算系数的方程,可以得到系数的常微分方程
i
d
d
t
c
n
(
t
)
=
(
2
π
n
)
2
c
n
+
∑
k
V
n
−
k
c
k
,
{\displaystyle i{\frac {d}{dt}}c_{n}(t)=(2\pi n)^{2}c_{n}+\sum _{k}V_{n-k}c_{k},}
而其元素
V
n
k
{\displaystyle V_{nk}}
可以透过显式的傅立叶转换求得
V
n
−
k
=
∫
0
1
V
(
x
)
e
2
π
i
(
k
−
n
)
x
d
x
.
{\displaystyle V_{n-k}=\int _{0}^{1}V(x)\ e^{2\pi i(k-n)x}dx.}
若将
N
{\displaystyle N}
基底函数的展开到一定项次后截断,并且找
c
n
(
t
)
{\displaystyle c_{n}(t)}
的解,就可以得到偏微分方程的解。一般而言会用数值方法 进行,例如龙格-库塔法 。为了数值解,常微分方程的右侧需重复计算在许多不同时间间隔下的值。此时,谱方法有个和势能项
V
(
x
)
{\displaystyle V(x)}
有关的主要问题。
在谱方法下,和势能函数
V
(
x
)
{\displaystyle V(x)}
的相乘会转换为向量和矩阵的乘法,其复杂度是
N
2
{\displaystyle N^{2}}
,而且在求解系数的微分方程时,需要另外去计算矩阵元素
V
n
k
{\displaystyle V_{nk}}
,这也需要时间。
在拟谱法中,会用不同的方式来计算。给定系数
c
n
(
t
)
{\displaystyle c_{n}(t)}
,会用反离散傅立叶转换来计算函数
ψ
{\displaystyle \psi }
在离散格点
x
j
=
2
π
j
/
N
{\displaystyle x_{j}=2\pi j/N}
下的值。在格点上,计算函数的乘积
ψ
′
(
x
i
,
t
)
=
V
(
x
i
)
ψ
(
x
i
,
t
)
{\displaystyle \psi '(x_{i},t)=V(x_{i})\psi (x_{i},t)}
,再用傅立叶转换转换回来,可以得到一组新的系数
c
n
′
(
t
)
{\displaystyle c'_{n}(t)}
,来代替矩阵乘积运算
∑
k
V
n
−
k
c
k
(
t
)
{\displaystyle \sum _{k}V_{n-k}c_{k}(t)}
。
可以证明二个方法有类似的精准度,而且拟谱法可以使用快速傅立叶转换,其时间复杂度为
O
(
N
ln
N
)
{\displaystyle O(N\ln N)}
,理论上比矩阵乘法要快很多。而且,可以直接计算函数
V
(
x
)
{\displaystyle V(x)}
,不用再经过额外的积分运算。
技术讨论
若用比较抽象的方式来描述,拟谱法是处理在偏微分方程中二个函数
V
(
x
)
{\displaystyle V(x)}
和
f
(
x
)
{\displaystyle f(x)}
的乘积。为了简化表示式,省略掉函数中时间的自变数。在概念上,拟谱法包括三个步骤:
将
f
(
x
)
{\displaystyle f(x)}
以及
f
~
(
x
)
=
V
(
x
)
f
(
x
)
{\displaystyle {\tilde {f}}(x)=V(x)f(x)}
扩展为由有限多个基底函数形成的组合(此为谱方法 )。
针对一组给定的基底函数,找到一组分割方式可以将基底函数的纯量积转换为在格点上的加权和。
在每一个格点将
V
{\displaystyle V}
和
f
{\displaystyle f}
相乘,来计算函数的乘积。
基底的展开
函数
f
{\displaystyle f}
和
f
~
{\displaystyle {\tilde {f}}}
可以用一组有限基底的函数来扩展
{
ϕ
n
}
n
=
0
,
…
,
N
{\displaystyle \{\phi _{n}\}_{n=0,\ldots ,N}}
:
f
(
x
)
=
∑
n
=
0
N
c
n
ϕ
n
(
x
)
{\displaystyle f(x)=\sum _{n=0}^{N}c_{n}\phi _{n}(x)}
f
~
(
x
)
=
∑
n
=
0
N
c
~
n
ϕ
n
(
x
)
{\displaystyle {\tilde {f}}(x)=\sum _{n=0}^{N}{\tilde {c}}_{n}\phi _{n}(x)}
为了简化起见,令基底是正交且正规化的,
⟨
ϕ
n
,
ϕ
m
⟩
=
δ
n
m
{\displaystyle \langle \phi _{n},\phi _{m}\rangle =\delta _{nm}}
,利用内积
⟨
f
,
g
⟩
=
∫
a
b
f
(
x
)
g
(
x
)
¯
d
x
{\displaystyle \langle f,g\rangle =\int _{a}^{b}f(x){\overline {g(x)}}dx}
配合适当的边界
a
,
b
{\displaystyle a,b}
,其系数为
c
n
=
⟨
f
,
ϕ
n
⟩
{\displaystyle c_{n}=\langle f,\phi _{n}\rangle }
c
~
n
=
⟨
f
~
,
ϕ
n
⟩
{\displaystyle {\tilde {c}}_{n}=\langle {\tilde {f}},\phi _{n}\rangle }
配合一些计算可得
c
~
n
=
∑
m
=
0
N
V
n
−
m
c
m
{\displaystyle {\tilde {c}}_{n}=\sum _{m=0}^{N}V_{n-m}c_{m}}
而
V
n
−
m
=
⟨
V
ϕ
m
,
ϕ
n
⟩
{\displaystyle V_{n-m}=\langle V\phi _{m},\phi _{n}\rangle }
。这是谱方法的基础。为了区分
ϕ
n
{\displaystyle \phi _{n}}
的基底以及正交的基底,有时会将上述扩展称为有限基底表示(Finite Basis Representation、FBR)。
分割
针对给定基底
{
ϕ
n
}
{\displaystyle \{\phi _{n}\}}
以及
N
+
1
{\displaystyle N+1}
个基底函数,可以设法找到分割方式,也就是
N
+
1
{\displaystyle N+1}
个点以及加权,使得
⟨
ϕ
n
,
ϕ
m
⟩
=
∑
i
=
0
N
w
i
ϕ
n
(
x
i
)
ϕ
m
(
x
i
)
¯
n
,
m
=
0
,
…
,
N
{\displaystyle \langle \phi _{n},\phi _{m}\rangle =\sum _{i=0}^{N}w_{i}\phi _{n}(x_{i}){\overline {\phi _{m}(x_{i})}}\qquad \qquad n,m=0,\ldots ,N}
特别的例子包括多项式的高斯求积 以及平面波的离散傅里叶变换 ,特别需注意的是格点及加权
x
i
,
w
i
{\displaystyle x_{i},w_{i}}
都是基底及数量
N
{\displaystyle N}
的函数。
利用分割方式,可以透过格点上的值,以另一种方式来表示函数
f
(
x
)
,
f
~
(
x
)
{\displaystyle f(x),{\tilde {f}}(x)}
的数值。此表示法有时称为离散变数表示法(Discrete Variable Representation、DVR),完全等效于基底的展开。
f
(
x
i
)
=
∑
n
=
0
N
c
n
ϕ
n
(
x
i
)
{\displaystyle f(x_{i})=\sum _{n=0}^{N}c_{n}\phi _{n}(x_{i})}
c
n
=
⟨
f
,
ϕ
n
⟩
=
∑
n
=
0
N
w
i
f
(
x
i
)
ϕ
n
(
x
i
)
¯
{\displaystyle c_{n}=\langle f,\phi _{n}\rangle =\sum _{n=0}^{N}w_{i}f(x_{i}){\overline {\phi _{n}(x_{i})}}}
相乘
和函数
V
(
x
)
{\displaystyle V(x)}
的相乘会在格点上进行
f
~
(
x
i
)
=
V
(
x
i
)
f
(
x
i
)
.
{\displaystyle {\tilde {f}}(x_{i})=V(x_{i})f(x_{i}).}
一般来说这里会有一些近似,可以计算其中一个系数
c
~
n
{\displaystyle {\tilde {c}}_{n}}
:
c
~
n
=
⟨
f
~
,
ϕ
n
⟩
=
∑
i
w
i
f
~
(
x
i
)
ϕ
n
(
x
i
)
¯
=
∑
i
w
i
V
(
x
i
)
f
(
x
i
)
ϕ
n
(
x
i
)
¯
{\displaystyle {\tilde {c}}_{n}=\langle {\tilde {f}},\phi _{n}\rangle =\sum _{i}w_{i}{\tilde {f}}(x_{i}){\overline {\phi _{n}(x_{i})}}=\sum _{i}w_{i}V(x_{i})f(x_{i}){\overline {\phi _{n}(x_{i})}}}
利用谱方法,对应的系数会是
c
~
n
=
⟨
V
f
,
ϕ
n
⟩
{\displaystyle {\tilde {c}}_{n}=\langle Vf,\phi _{n}\rangle }
。拟谱法则会用以下皂的近似来处理
⟨
V
f
,
ϕ
n
⟩
≈
∑
i
w
i
V
(
x
i
)
f
(
x
i
)
ϕ
n
(
x
i
)
¯
.
{\displaystyle \langle Vf,\phi _{n}\rangle \approx \sum _{i}w_{i}V(x_{i})f(x_{i}){\overline {\phi _{n}(x_{i})}}.}
若乘积可以用
V
f
{\displaystyle Vf}
给定的有限基底函数组合来表现,则上式在给定分割方式上会完全正确。
特殊的拟谱架构
傅立叶法
若问题中有周期性的边界条件,其周期为
[
0
,
L
]
{\displaystyle [0,L]}
,基底函数可以用平面波来产生,
ϕ
n
(
x
)
=
1
L
e
−
ı
k
n
x
{\displaystyle \phi _{n}(x)={\frac {1}{\sqrt {L}}}e^{-\imath k_{n}x}}
其中
k
n
=
(
−
1
)
n
⌈
n
/
2
⌉
2
π
/
L
{\displaystyle k_{n}=(-1)^{n}\lceil n/2\rceil 2\pi /L}
,而
⌈
⋅
⌉
{\displaystyle \lceil \cdot \rceil }
是取整函数 。
在
n
max
=
N
{\displaystyle n_{\text{max}}=N}
处截断的分割为离散傅立叶转换 ,格点会平均分布
x
i
=
i
Δ
x
{\displaystyle x_{i}=i\Delta x}
,间隔为
Δ
x
=
L
/
(
N
+
1
)
{\displaystyle \Delta x=L/(N+1)}
,各点的加权会是相同旳定值
w
i
=
Δ
x
{\displaystyle w_{i}=\Delta x}
。
在讨论误差时,需注意到平面波的乘积也是平面波,
ϕ
a
+
ϕ
b
=
ϕ
c
{\displaystyle \phi _{a}+\phi _{b}=\phi _{c}}
,
c
≤
a
+
b
{\displaystyle c\leq a+b}
。因此,定量来说,若函数
f
(
x
)
,
V
(
x
)
{\displaystyle f(x),V(x)}
可以用
N
f
,
N
V
{\displaystyle N_{f},N_{V}}
基底函数足够准确的呈现,只要用
N
f
+
N
V
{\displaystyle N_{f}+N_{V}}
个基底函数,即可用拟谱法得到足够准确的结果。
平面波的扩展一般效果较差,需要许多的基底函数才能收敛。不过。基底展开和格点表示的转换可以用快速傅立叶转换 进行,其时间复杂度较低,为
N
ln
N
{\displaystyle N\ln N}
。因此,平面波是拟谱法中常用的一种基底函数。
多项式
另一种常见的展开方式是多项式,此处会使用高斯求积 (Gaussian quadrature),其中提到可以找到加权系数
w
i
{\displaystyle w_{i}}
及格点
x
i
{\displaystyle x_{i}}
使得
∫
a
b
w
(
x
)
p
(
x
)
d
x
=
∑
i
=
0
N
w
i
p
(
x
i
)
{\displaystyle \int _{a}^{b}w(x)p(x)dx=\sum _{i=0}^{N}w_{i}p(x_{i})}
对任意的
2
N
+
1
{\displaystyle 2N+1}
次或是更低次的多项式
p
(
x
)
{\displaystyle p(x)}
都成立。一般而言,加权函数
w
(
x
)
{\displaystyle w(x)}
及范围
a
,
b
{\displaystyle a,b}
都是根据特定问题所选定的,因此会选择几种分割方式中的一种。若要用在拟谱法,需选择基底函数为
ϕ
n
(
x
)
=
w
(
x
)
P
n
(
x
)
{\displaystyle \phi _{n}(x)={\sqrt {w(x)}}P_{n}(x)}
,其中
P
n
{\displaystyle P_{n}}
为
n
{\displaystyle n}
阶多项式,有以下特性
∫
a
b
w
(
x
)
P
n
(
x
)
P
m
(
x
)
d
x
=
δ
m
n
.
{\displaystyle \int _{a}^{b}w(x)P_{n}(x)P_{m}(x)dx=\delta _{mn}.}
在上述条件下,the
ϕ
n
{\displaystyle \phi _{n}}
会是正交基底,其内积为
⟨
f
,
g
⟩
=
∫
a
b
f
(
x
)
g
(
x
)
¯
d
x
{\displaystyle \langle f,g\rangle =\int _{a}^{b}f(x){\overline {g(x)}}dx}
。此基底以及分割点可以用在拟谱法中。
有关其误差,若
f
{\displaystyle f}
可以用
N
f
{\displaystyle N_{f}}
个基底函数很好的呈现,而
V
{\displaystyle V}
可以用
N
V
{\displaystyle N_{V}}
阶多项式很好的呈现,则其积可以用前
N
f
+
N
V
{\displaystyle N_{f}+N_{V}}
个基底函数很好的呈现。且拟谱法用该数量的基底函数,会有足够准确的结果。
在一些标准问题中,会出现这些多项式。例如量子简谐振动子可以扩展为埃尔米特多项式 ,而在转动问题中,会用雅可比多项式 来定义相关的勒壤得多项式 。
参考资料
Orszag, Steven A. Numerical Methods for the Simulation of Turbulence. Physics of Fluids. 1969, 12 (12): II–250. doi:10.1063/1.1692445 .
Gottlieb, David; Orszag, Steven A. Numerical analysis of spectral methods : theory and applications 5. print. Philadelphia, Pa.: Society for Industrial and Applied Mathematics. 1989. ISBN 978-0898710236 .
Hesthaven, Jan S.; Gottlieb, Sigal; Gottlieb, David. Spectral methods for time-dependent problems 1. publ. Cambridge [u.a.]: Cambridge Univ. Press. 2007. ISBN 9780521792110 .
Trefethen, Lloyd N. Spectral methods in MATLAB 3rd. repr. Philadelphia, Pa: SIAM. 2000. ISBN 978-0-89871-465-4 .
Fornberg, Bengt. A Practical Guide to Pseudospectral Methods. Cambridge: Cambridge University Press. 1996. ISBN 9780511626357 .
Boyd, John P. Chebyshev and Fourier spectral methods 2nd ed., rev. Mineola, N.Y.: Dover Publications. 2001 [2019-07-17 ] . ISBN 978-0486411835 . (原始内容存档 于2019-06-24).
Funaro, Daniele. Polynomial approximation of differential equations . Berlin: Springer-Verlag. 1992 [2019-07-17 ] . ISBN 978-3-540-46783-0 . (原始内容存档 于2011-07-22).
de Frutos, Javier; Novo, Julia. A Spectral Element Method for the Navier--Stokes Equations with Improved Accuracy. SIAM Journal on Numerical Analysis. January 2000, 38 (3): 799–819. doi:10.1137/S0036142999351984 .
Claudio, Canuto; M. Yousuff, Hussaini; Alfio, Quarteroni; Thomas A., Zang. Spectral methods fundamentals in single domains. Berlin: Springer-Verlag. 2006. ISBN 978-3-540-30726-6 .
Press, WH; Teukolsky, SA; Vetterling, WT; Flannery, BP. Section 20.7. Spectral Methods . Numerical Recipes: The Art of Scientific Computing 3rd. New York: Cambridge University Press. 2007 [2019-07-17 ] . ISBN 978-0-521-88068-8 . (原始内容存档 于2011-08-11).