Operator Splitting Methods in PDEs

为什么要使用算符分解方法求解偏微分方程
基本思路
方法的优点
方法的缺点
First Order Splitting Lie-Trotter Splitting
First Order Splitting Additive Splitting

求解偏微分方程时，常用的一个方法，今天才开始好好学习一下，做个笔记和大家分享下。
学习的书籍是：Helge Holden et al. (2010):
Splitting Methods for Partial Differential Equations with Rough Solutions
Analysis and MATLAB programs

为什么要使用算符分解方法求解偏微分方程

对于许多问题的模拟，最后往往归结为求解偏微分方程。由于现在计算能力的大大提高，求解的偏微分方程越来越复杂
求解过程可能涉及到不同的物理过程，这些过程反映在PDE中，表现为不同的数学项。这些不同的项数学上的表现也大不相同，需要做不同的处理，使得不论进行解析求解还是数值模拟都很困难
如果求解的方程是一个新类型的PDE，现有的解算方法可能会失效
利用算符分解方法可以将一个复杂的PDE分解为一些简单的PDE，这些简单的PDE一般都有现成的算法求解，最后再将这些简单PDE的解综合起来得到最终结果（类似divide and conquer的思想）

基本思路

将PDE方程分解成一系列的子方程，这些子方程是容易求解者有现成的解法
使用一定的算符分解算法，将这些子方程的解连起来，形成最终的解

以一种抽象的方式举例，比如求解柯西问题：

dUdt+A(U)=0, U(0)=U0

\frac{\mathrm{d}U}{\mathrm{d}t} + \mathcal{A}(U)=0, ~~ U(0)=U_0
其中 A \mathcal{A}代表一些未指明的算符，以后会举例说明（比如在热流传输过程中的对流项、输运项）

使用分离变量法（我自己的理解，有问题大家请指正），可以获得下列形式的解：

U(t)=e−tAU0

U(t) = e^{-t\mathcal{A}}U_0

假设算符能够分解，即 A=A1+A2 \mathcal{A}=\mathcal{A}_1 + \mathcal{A}_2，并且可以很容易的求解下列子问题：

dUdt+Aj(U)=0, U(0)=U0, j=1,2

\frac{\mathrm{d}U}{\mathrm{d}t} + \mathcal{A}_j(U)=0, ~~ U(0)=U_0, ~~~ j=1,2
获得的解形式为：

U(t)=e−tAjU0

U(t) = e^{-t\mathcal{A}_j}U_0

考虑最简单的形式，假设 tn=nΔt, n=1,2,... t_n=n\Delta t, ~~ n=1,2,..., Δt>0 \Delta t >0，且足够小，则希望获得如下近似：

U(tn+1)=e−ΔtA2e−ΔtA1U(tn)

U(t_{n+1})=e^{-\Delta t\mathcal{A}_2}e^{-\Delta t\mathcal{A}_1}U(t_n)
假设 A1 \mathcal{A}_1和 A2 \mathcal{A}_2是对易算符，即满足

∂A2∂UA1−∂A1∂UA2=0

\frac{\partial \mathcal{A}_2}{\partial U}\mathcal{A}_1 - \frac{\partial \mathcal{A}_1}{\partial U}\mathcal{A}_2 = 0
则 e−ΔtA2e−ΔtA1=e−ΔtA e^{-\Delta t\mathcal{A}_2}e^{-\Delta t\mathcal{A}_1} = e^{-\Delta t\mathcal{A}} , 获得的 U(tn+1) U(t_{n+1})是精确解。
但实际上， A1 \mathcal{A}_1和 A2 \mathcal{A}_2往往是不对易的，此时我们希望

U(tn+1)=e−tAjU0=limΔt→0,t=nΔte−ΔtA2e−ΔtA1U(tn)

U(t_{n+1})=e^{-t\mathcal{A}_j}U_0 = \lim\limits_{\Delta t \to 0,t=n\Delta t} e^{-\Delta t\mathcal{A}_2}e^{-\Delta t\mathcal{A}_1}U(t_n)
这个就是著名的Lie–Trotter–Kato公式。将解中的算符 etAj e^{t\mathcal{A}_j}替换为数值近似，就获得了求解PDE的数值方法。所有更精确的方法，都是在这个基础上发展得到的。采用不同的数值近似，可以获得不同的方法。

方法的优点

这个方法看上去非常粗糙，但是具有下列优点

对子方程，可采用不同的数值或分析方法。在某些情况下，这是非常必要的，可以使求解更简单有效
更方便利用一些典型PDE的现有算法，具有更大的灵活性。
减少计算过程中对内存的需求，提高求解的稳定性
对某些问题，可能是唯一可行的求解方法
可以很方便的扩充模拟问题的模型，加入不同的物理过程

方法的缺点

如果子系统耦合度非常高，比如所代表的物理过程演化的时间尺度非常短，则算符分解方法可能会失效，分解过程中步长的选取会受到严重限制
必须考虑分解造成的误差，以防止出现严重的错误

补充：
对方程： Ut=Ux U_t = U_x，使用前向欧拉公式，对所有的 λ \lambda，解是不稳定的；使用leap-frog方法（2阶时间精度，2阶空间精度），对 λ<1 \lambda ，解是稳定的
然后，不幸的是，对方程： Ut=Uxx U_t = U_{xx}，使用前向欧拉公式，对 λ<1/2 \lambda ，解是稳定的；使用leap-frog方法（2x2精度），对所有 λ \lambda，解都是不稳定的

现在我们想解下列方程： Ut=Ux+Uxx U_t = U_x + U_{xx}，这两种方法还可以用吗？！！

使用算符分解方法，可以做如下处理：
（1）求解 Ut=Ux U_t = U_x，初值为 U0 U_0，获得解 U1 U_1（可以使用leap-frog）
（2）求解 Ut=Uxx U_t = U_{xx}，用上一步的 U1 U_1作为初值，获得解 U2 U_2（可以使用欧拉方法）
（3）将 U2 U_2当作（1）中的初值 U0 U_0，循环执行（1），（2），获得下个步长的解，
END

First Order Splitting: Lie-Trotter Splitting

问题仍然如前面所述。
算法如下：

(1)

∂U(t)∂t=A1U(t),t[tn,tn+1],u(tn)=unsp

\frac{\partial U(t)} { \partial t} = \mathcal{A}_1 U(t), t [t^n,t^{n+1}], u(t^n)=u^n_{sp}
(2)
(3)

First Order Splitting: Additive Splitting

————————————-万恶分割线————————————————-
由于时间和水平所限，错误在所难免，大家批判阅读，特此敬告！