求解泊松方程的高精度紧致差分方法

收稿日期：2021-04-09作者简介：刘昊（1995-），男，宁夏吴忠人，硕士研究生。

泊松方程的有限差分方法及快速实现刘昊，张荣培，霍俊蓉（沈阳师范大学数学与系统科学学院，辽宁沈阳110136）摘要：针对d 维（d =1,2,3）带有Dirichlet 边界的泊松方程，设计一类快速求解方法。

首先采用有限差分方法将方程离散，利用Kronecker 积的性质将离散后的方程进行矩阵分解，进而应用快速离散正弦变换（DST ）方法进行有效求解。

数值实验结果表明，该方法可快速求解d 维泊松方程，并验证了其准确性和有效性。

关键词：泊松方程；有限差分；Crank -Nicolson 方法；离散正弦变换中图分类号：O175文献标识码：A文章编号：1673-1603（2021）04-0091-06DOI ：10.13888/ki.jsie （ns ）.2021.04.018第17卷第4期2021年10月Vol.17No.4Oct.2021沈阳工程学院学报（自然科学版）Journal of Shenyang Institute of Engineering （Natural Science ）泊松方程是常见于静电学、机械工程和理论物理的偏微分方程，在无引力源的情况下，ΔΦ=0（即拉普拉斯方程）；当考虑引力场时，ΔΦ=f （f 为引力场的质量分布）。

该方程通常用格林函数法求解，也可以采用分离变量法和特征线法求解［1］。

带有规则区域Dirichlet 边界条件的泊松方程：ìíîïï-Δu =f (x )u|∂Ω=g (x )（1）式中，Δ表示拉普拉斯算子；x ∈Ω；Ω⊂R d ；当d =1时，Δu =u xx ；当d =2时，Δu =u xx +u yy ；当d =3时，Δu =u xx +u yy +u zz 。

由于带有Dirichlet 边界的泊松方程的解析解不容易求出，一般情况下采用有限差分法、有限元法和有限体积方法的数值方法进行求解［2-4］。

三维泊松方程的七点差分格式
豆桂芳;吴振远;杜艳林
【期刊名称】《工程地球物理学报》
【年(卷),期】2009(006)006
【摘要】泊松方程在地球物理中的应用非常广泛,随着三维勘探的展开,三维泊松方程的快速、精确求解变得越来越重要.本文用有限差分方法推导出长方体域的三维泊松方程的七点差分格式,及其边界条件的差分格式,并对满足第一类边界条件的七点差分格式进行求解,详细推导了七点差分格式的求解过程.
【总页数】4页(P802-805)
【作者】豆桂芳;吴振远;杜艳林
【作者单位】中国地质大学,数学与物理学院,武汉,430074;中国地质大学,数学与物理学院,武汉,430074;中国地质大学,数学与物理学院,武汉,430074
【正文语种】中文
【中图分类】P631
【相关文献】
1.五点差分格式求解泊松方程并行算法的研究 [J], 廖臣;祝大军;刘盛纲
2.二维泊松方程非均匀网格上的高精度紧致差分格式 [J], 郭锐;黄雪芳;葛永斌
3.解抛物型方程的七点半显差分格式 [J], 陈泽龙;张大凯
4.求解薛定谔–泊松方程组的时间分裂紧致差分格式 [J], 姜珊;刘荣华;马秀;王汉权;
5.半导体器件非线性泊松方程的一种新差分格式 [J], 王国彬;严荣良;任迪远
因版权原因，仅展示原文概要，查看原文内容请购买。

三维泊松方程的高精度求解方法三维泊松方程是一种常见的偏微分方程，广泛应用于自然科学与工程领域的建模和计算。

解决三维泊松方程的高精度方法是很有挑战性的，因为它涉及到大规模的计算和内存需求。

本文将介绍几种常见的高精度求解三维泊松方程的方法。

一、直接方法直接方法是通过离散化三维泊松方程，并利用线性代数方法求解线性方程组。

其中最常用的方法是基于有限差分法的离散化方法。

具体而言，我们可以将泊松方程转化为一个线性方程组Ax=b的形式，其中A是一个稀疏矩阵，x是未知量，b是常数向量。

对于高精度求解，可以使用稀疏矩阵的特殊性质进行优化，如共轭梯度法（Conjugate Gradient Method）、预条件共轭梯度法（Preconditioned Conjugate Gradient Method）等。

这些方法可以充分利用矩阵的稀疏性，减少计算复杂度和存储空间。

二、迭代方法迭代方法是通过迭代更新解的数值，直到达到预设的收敛条件。

常见的迭代方法包括雅可比迭代法（Jacobi Method）、高斯-赛德尔迭代法（Gauss-Seidel Method）和多重网格法（Multigrid Method）等。

对于高精度求解，可以选择迭代方法的高阶或者加速版本，如SOR迭代法（Successive Over-Relaxation Method）、MGCG迭代法（Multigrid Conjugate Gradient Method）等。

这些方法可以加快收敛速度，提高数值精度。

三、谱方法谱方法是基于函数的全局性质进行求解的方法。

将输入函数在一组合适的基函数上展开，并选择合适的权重，然后通过将展开系数带入泊松方程，得到一组线性方程，进而求解出解。

对于高精度求解，可以选择具有较高收敛阶的基函数，如拉格朗日插值基函数和切比雪夫基函数等。

此外，合理地选择基函数的数量和位置，可以进一步提高数值精度。

四、混合方法混合方法是将上述不同的方法结合在一起，以充分发挥各自的优点。

五点差分格式求解poisson方程Poisson方程是数学物理中的一个重要方程，广泛应用于电场、热传导、流体力学等领域。

求解Poisson方程的方法有很多，其中五点差分格式是一种常用的数值求解方法。

五点差分格式是一种离散化的方法，将Poisson方程中的二阶导数用差分近似表示。

具体来说，我们将求解区域划分为网格，每个网格点上的解值用一个未知数表示。

然后，根据Poisson方程的离散形式，我们可以得到一个线性方程组，通过求解这个方程组，即可得到Poisson方程的数值解。

在五点差分格式中，我们使用中心差分来近似二阶导数。

对于一个二维的求解区域，我们可以将其划分为若干个网格点，每个网格点的坐标为(i,j)，其中i表示横向的网格编号，j表示纵向的网格编号。

假设网格的步长为h，则(i,j)点的解值为u(i,j)。

根据中心差分的定义，我们可以得到(i,j)点的二阶导数近似为：(u(i+1,j) - 2u(i,j) + u(i-1,j))/h^2 + (u(i,j+1) - 2u(i,j) + u(i,j-1))/h^2将上述近似代入Poisson方程中，我们可以得到：(u(i+1,j) + u(i-1,j) + u(i,j+1) + u(i,j-1) - 4u(i,j))/h^2 = f(i,j)其中f(i,j)为Poisson方程中的源项。

根据上述离散形式，我们可以得到每个网格点上的线性方程，将所有的线性方程组合起来，即可得到一个大的线性方程组。

通过求解这个线性方程组，我们可以得到Poisson方程的数值解。

在实际求解中，我们可以使用迭代法来求解这个线性方程组。

常用的迭代方法有Jacobi迭代法、Gauss-Seidel迭代法等。

这些迭代方法的基本思想是通过不断迭代更新未知数的值，直到收敛为止。

除了迭代法，我们还可以使用直接法来求解线性方程组。

常用的直接法有LU分解法、Cholesky分解法等。

这些直接法的基本思想是通过将线性方程组转化为一个上三角矩阵和一个下三角矩阵的乘积，从而求解出未知数的值。

泊松方程有限差分在这篇文章中，我们将探讨泊松方程的有限差分方法。

有限差分是一种数值解微分方程的方法，它将方程中的微分算子用差分算子来近似表示，从而将连续空间中的问题转化为离散空间中的问题。

首先，我们先来回顾一下泊松方程的一般形式。

泊松方程通常可以写为：∇^2φ = -ρ其中，∇^2表示拉普拉斯算子，φ是待求解的标量场，ρ是源项。

在物理学中，φ通常代表电势、温度或者密度等物理量，而ρ则表示外部的电荷密度、热源或质量密度等。

对于一个给定的区域Ω上的泊松方程问题，我们要求解φ满足泊松方程以及边界条件。

边界条件通常给出了φ在Ω的边界上的数值或者导数信息。

在有限差分方法中，我们首先需要将问题的空间离散化。

设Ω是一个二维区域，我们用一个网格来离散Ω。

假设Ω上有N个网格点，我们用(i,j)来表示第i行第j列的网格点，并假设Ω被水平方向和竖直方向的线段分别分成了M+1和N+1个小区间。

接下来，我们需要定义泊松方程的差分格式。

对于一个给定的网格点(i,j)，我们可以用中心差分来近似拉普拉斯算子∇^2：∇^2φ(i,j) ≈ (φ(i+1,j) - 2φ(i,j) + φ(i-1,j))/Δx^2 + (φ(i,j+1) - 2φ(i,j) + φ(i,j-1))/Δy^2其中Δx和Δy分别是水平和竖直方向上的网格间距。

通过这样的近似，我们可以得到φ在点(i,j)的近似解。

然后，我们将泊松方程中的微分算子用差分算子来替代，得到离散的泊松方程格式：(φ(i+1,j) - 2φ(i,j) + φ(i-1,j))/Δx^2 + (φ(i,j+1) - 2φ(i,j) + φ(i,j-1))/Δy^2 = -ρ(i,j)这就是泊松方程的有限差分格式。

通过对所有网格点应用这样的格式，我们可以得到一个关于φ(i,j)的代数方程组。

通过求解这个方程组，我们就可以得到φ在整个Ω上的近似解。

在实际计算中，我们通常采用迭代方法来求解这个代数方程组。

泊松方程的高效求解方法研究与实现泊松方程是物理学和工程学中非常重要的方程之一。

它描述了许多现象，如流体力学、电势、重力势和热传导等。

在工程领域，泊松方程的解决在很多电子、光学、医学等领域具有广泛的应用。

但是，由于其解析解十分困难，因此研究泊松方程的高效求解方法成为了相当重要的课题。

一、背景泊松方程是一个二阶偏微分方程，由于其数值解比解析解更为适用，因此我们通常使用离散化方法来求解。

离散化方法将连续问题转化为有限数量的离散点问题，然后应用数值计算方法来求解。

在离散化的过程中，我们使用有限差分法或有限元方法来逼近连续的泊松方程，进而可用计算机求解。

二、问题及挑战针对泊松方程，在形式上看，其只涉及二维的空间位置而不涉及时间。

泊松方程不仅解析解困难，而且求解效率也很低。

在高精度计算时，泊松方程的求解难度会明显增加，这就为计算的最终结果质量带来了挑战。

针对这个问题，许多研究者极力探索一些高效求解方法，以了解泊松方程的数值解，以便在各种不同的应用中使用。

三、高效求解方法1.有限差分法有限差分法是一种数值求解方法，它将求解区域离散化为一个网格，并使用近似方式逼近泊松方程。

有限差分法是求解偏微分方程的常用技术，它可以通过发现微分方程的近似解来得到真实解。

然而，有限差分法需要相对更大的计算量来处理问题，因此，该方法通常需要相对更长的时间才能完成计算。

2.快速傅立叶变换（FFT）FFT是一种用于解决线性DFT的优化算法。

由于傅里叶变换是泊松方程的线性算子，因此使用FFT来解决泊松方程是可行的。

但是，FFT的缺点是它需要用到的网格大小通常是2的幂。

因此，在确定网格大小时，需要注意网络的大小和计算所需的时间。

3.并行算法并行算法可提高计算机求解的速度，减少泊松方程的求解时间。

并行算法将问题分成多个部分使得多个处理器可以同时工作。

这为提高计算速度提供了途径。

并行算法可以加速求解大型泊松方程，特别是在处理更大的数据集时，这种算法的优势更加明显。

求解波动方程的2种显式高精度紧致差分格式姜蕴芝;葛永斌【摘要】针对一维波动方程,空间采用四阶Padé逼近,时间采用中心差分离散得到了一种时间二阶、空间四阶精度的显式紧致差分格式,其截断误差为O(τ2+h4).之后采用截断误差余项修正的方法对时间离散进行改进,改进后的格式的截断误差为O(τ4+τ2h2+h4),即格式具有整体四阶精度.然后,通过Fourier方法分析了2种格式的稳定性.最后,通过数值实验验证了本格式的精确性和可靠性.%In this paper,an explicit compact difference scheme is obtained for solving the one dimensional wave equation.The truncation error of the scheme is O(τ2 + h4).It's constructed by applying the fourth-order accurate Padé approximation in space and the second-order accurate central difference in time.Then,the remainder of the truncation error correction method is employed to improve the accuracy of the discretization of time,the truncation error of the improved scheme is O(τ4 + τ2 h2 + h4),which means the scheme has an overall fourth-order accuracy.And then,the stability conditions of the two schemes are obtained by the Fourier method.Finally,the accuracy and the reliability of the present two schemes are verified by numerical experiments.【期刊名称】《四川师范大学学报（自然科学版）》【年(卷),期】2017(040)002【总页数】7页(P177-183)【关键词】波动方程;Padé逼近;紧致格式;显式差分;稳定性【作者】姜蕴芝;葛永斌【作者单位】宁夏大学数学统计学院,宁夏银川750021;宁夏大学数学统计学院,宁夏银川750021【正文语种】中文【中图分类】O242.1波动方程是一类重要的双曲型偏微分方程，在数学、物理、化学等领域内有着广泛的应用[1].对这类方程进行数值求解的方法主要有有限差分法、有限元法和有限体积法等.就有限差分法而言，对该类问题研究的理论成果有：文献[2]采用二阶中心差分格式和非均匀网格离散，提出了一种求解一维波动方程在非均匀时间网格上的三层显式差分格式，该格式具有一阶精度；文献[3]利用泰勒级数展开及待定系数法建立了一种求解一维波动方程的三层隐格式，该格式是条件稳定的，并且具有四阶精度；文献[4]将Runge-Kutta方法应用到哈密顿系统中并与辛格式相结合，提出了求解一维波动方程的一类显式辛方法，该方法具有二阶精度，并且是条件稳定的；文献[5]采用三次样条公式推导出精度分别为O(τ2+h2)、O(τ2+h4)、O(τ4+h2)和O(τ4+h4)的4种差分格式，并采用Fourier方法分析了格式的稳定性；文献[6]对一维二阶波动方程提出了具有二阶精度的精细时程积分法，该方法能够在时间方向上精确计算，空间方向的局部截断误差为O(h2)，并且该方法是无条件稳定的；文献[7]通过加权平均和紧致差分离散的思想提出了2种精度分别为O(τ2+h4)和O(τ4+h4)的隐式紧致差分格式；文献[8]在Crank-Nicolson格式的基础上设计了重叠型区域分解的并行算法，该算法的最优逼近阶为O(τ2+h2)；文献[9]利用三次样条插值，提出了求解一维波动方程的一种三层隐式差分格式，该格式最优能够达到时间二阶，空间四阶精度，并且是无条件稳定的；文献[10]通过四次样条函数与广义梯形算法相结合的方法提出了一维波动方程的一类两层差分格式，其精度为O(τ2+h4)，当选择适当的参数时，其精度可提高到O(τ3+h4)；文献[11]采用四阶紧致差商逼近公式及加权平均思想，提出了2种精度分别为O(τ2+h4)和O(τ4+h4)的交替方向隐式格式，前者是无条件稳定的，后者是条件稳定的；文献[12]提出了一种求解一维波动方程的高精度隐式差分格式，该格式是无条件稳定的，并且具有时间二阶、空间四阶精度.本文将建立2种显式紧致差分格式，为此，考虑如下一维波动方程的初边值问题:其中,u(x,t)是待求未知函数，a为波动系数，f(x,t)为非齐次项，φ(x)、ψ(x)、g0(t)、g1(t)为已知函数，且具有充分的光滑性.设时间步长为τ，空间步长为h，网格节点为(xj,tn)，其中xj=jh，tn=nτ，j=0,1,2,…,N，n>0，用表示u(xj,tn)的近似值.1.1 CTFS格式对初始时间层的离散.利用泰勒展开公式有将(2)式代入(4)式中，进行整理并舍去其截断误差项有计算的值.对空间内部节点采用文献[13]中的四阶紧致差分格式进行逼近，则有对于空间边界节点的处理，由(1)式与(3)式可得时间的推进.对时间方向上的推进，采用中心差分格式有对上式进行整理，并舍去其截断误差项有(10)式即为所构造的显式紧致差分格式，记为CTFS格式.等号右端的项采用(6)～(8)式通过追赶法进行计算，初始时间步由(5)式进行计算.通过差分格式的构造过程不难发现，该格式是显格式，其截断误差为O(τ2+h4).1.2 FTFS格式为了使时间精度与空间精度能够相匹配，使格式整体精度达到四阶，下面对上述差分格式进行改进，进一步提高时间精度，为此，对初始时间层的离散采用文献[7]中的方法，有其中λ=τ/h.对于时间二阶导数项的离散，采用中心差分格式离散后保留其截断误差主项，可得又由(1)式可得将上式代入(12)式中有对(14)式进行整理，并舍去其截断误差项有又由(1)式，对上式进行变形有(15)式即为整体四阶精度的显式紧致差分格式，记为FTFS格式.由格式的构造过程可知，该差分格式的截断误差为O(τ4+τ2h2+h4).下面采用Fourier方法分析本文所提2种格式的稳定性.假设源项f精确无误差，令=ηneiωj，其中ξ、η为振幅，ω为相位角，为虚数单位.引理 1[1] 实系数二次方程μ2-bμ-c=0的根按模不大于1的充分必要条件为|b|≤1-c≤2.对于(6)式有对上式进行化简整理有对于CTFS格式(10)，令，将其写为矩阵的形式有令Uj=(uj,vj)T，并将(17)式代入进行整理有从而可得CTFS格式(10)的误差增长矩阵为令λ=τ/h，r=aλ,则得上述误差增长矩阵的特征方程为其中,，c=-1.因此，可得格式(10)稳定的充要条件为即.进而可得稳定性条件为|r|，即|a|.对于FTFS格式(15)，令，将其写为矩阵的形式有令代入上式，并进行整理有将(17)式代入(24)式进行整理，可得FTFS格式的误差增长矩阵为则误差增长矩阵的特征方程为其中，因此，可得该格式稳定的充要条件为上式等价于如下2个不等式:对于不等式(28)易得由于不等式对任意ω的取值都要成立，所以有即|a|.对于不等式(29)，由求根公式易得其等式解为令，则有在y∈[0,2]时，恒有<0，故有Smin=S(y)|y=2=2，此时cos ω=-1.又由于关于对称，故有另一解为1.进而可得第二个不等式的解为，与第一个不等式取交集即得差分格式(15)的稳定性条件为，即|a|].为了验证本文所提2种格式的精确性和可靠性，现考虑如下2个具有精确解的初边值问题.问题 1[7]:其精确解为u(x,t)=sin(πt)sin(πx).问题 2 :其精确解为u(x,t)=te-πtsin(πx).表1～3给出了问题1的计算结果.表1采用本文CTFS格式与文献[7]中的四阶格式进行了计算.由于这2种格式的精度均为时间二阶、空间四阶，因此取τ=2h2，计算了t=0.5时刻取不同h时(τ也相应不同)的L∞和L2范数误差和收敛阶.eL∞和eL2范数及收敛阶的定义如下：由表1的数据可以看出，2种格式空间均达到了四阶精度，而本文的CTFS格式要比文献[7]中的四阶格式更为精确.表2给出了本文FTFS格式和文献[7]中的四阶格式的计算结果，由于2种格式均具有整体四阶精度，因此取τ=h，计算了t=0.5时刻取不同h时的L∞和L2范数误差和收敛阶.可以看出，2种格式几乎具有相同的精度.表3验证了本文2种格式的稳定性，取h=1/32,给出了不同时间步steps(不同时间t)、不同网格比λ的计算结果.可以看出，当λ=0.8时CTFS格式是稳定的，而当λ≥1时，CTFS格式是发散的，这与本文的理论分析结果|a|是吻合的.当λ=0.8、1.0、1.5、1.6、1.7时，本文FTFS格式是稳定的，这也验证了本文的理论分析结果.而文献[7]中的四阶格式的稳定性条件为因此，当|a|λ大过1之后，计算结果是发散的.表4～6给出了问题2的计算结果.表4采用了本文CTFS格式与文献[7]中的四阶格式进行了计算.取τ=h2，计算了t=1时刻取不同h时的L∞和L2范数误差和收敛阶.同样可以看出，本文的CTFS格式要比文献[7]中的四阶格式更为精确.表5给出了本文FTFS格式和文献[7]中的四阶格式的计算结果，取τ=h，计算了t=1时刻取不同时的L∞和L2范数误差和收敛阶.可以看出，2种格式均达到了四阶精度，但本文FTFS格式比文献[7]中的四阶格式的计算结果更加精确.与问题1齐次问题的结果相比较，说明了本文FTFS格式更加适用于求解非齐次问题.表6验证了本文2种格式的稳定性，取h=1/32，给出了不同时间步、不同网格比λ的计算结果.同样可以看出，当λ=0.8、1.0、1.5、1.6、1.7时，本文FTFS格式是稳定的，这与我们对该格式稳定性的理论分析结果|a|]是吻合的.本文针对一维波动方程提出了2种显式高精度紧致差分格式，2种格式的截断误差分别为O(τ2+h4)和O(τ4+τ2h2+h4)，并通过Fourier分析法分析了2种格式的稳定性，其中前一种CTFS格式的稳定性条件为|a|，后一种FTFS格式的稳定性条件为|a|].然后通过数值实验将本文格式与文献[7]中的格式的计算结果进行对比，可以看出本文格式计算结果更加精确.并且本文格式的精度、稳定性与理论分析相一致，验证了本文格式的精确性与可靠性.此外，文献[14]提出了数值求解二维波动方程的三层全隐式紧致差分格式，其精度分别为O(τ2+h4)和O(τ4+h4)，为了加快迭代收敛速度，采用了多重网格方法进行加速.而本文的显式差分方法可以推广到二维，并且在二维情况下只需进行2次追赶法求解和一次显式递推计算，无需迭代，因此可望本文方法推广到二维后，在保有高精度的情况下会较文献[14]的方法具有更高的求解效率.目前，我们正在进行此方面的研究.另外，文献[15]提出了分数阶波动方程的一种差分方法，发展分数阶波动方程的显式高精度差分方法也是一个有意义的研究方向.致谢宁夏大学自然科学基金项目(ZR1407)和宁夏大学研究生创新项目(GIP2016032)对本文给予了资助,谨致谢意.【相关文献】[1] 戴嘉尊,邱建贤.偏微分方程数值解法[M].南京:东南大学出版社,2002.[2] CHO C H.Stability for the finite difference schemes of the linear wave equation with nonuniform time meshes[J].Num Meth Part Diff Eqns,2013,29(3):1031-1042.[3] 张天德,孙传灼.关于波动方程的差分格式[J].山东工业大学学报,1995,25(3):283-287.[4] 孙耿.波动方程的一类显式辛格式[J].计算数学,1997,19(1):1-10.[5] 齐远节,刘利斌.求解二阶波动方程的三次样条差分方法[J].大学数学,2011,27(1):59-64.[6] 金承日,吕万金.二阶双曲型方程的精细时程积分法[J].计算力学学报,2003,20(1):113-115.[7] 葛永斌,朱琳,田振夫.求解波动方程的高精度紧致隐式差分方法[J].宁夏大学学报(自然科学版),2005,26(4):297-299.[8] 田敏,羊丹平.波动方程的重叠型区域分解并行有限差分算法[J].山东大学学报(理学版),2007,42(2):28-38.[9] RASHIDINIA J,JALILIAN R,KAZEMI V.Spline methods for the solutions of hyperbolic equations[J].Appl Math Comput,2007,190(1):882-886.[10] LIU T,LIU L,XU H,et al.A new two level difference scheme for solving one-dimensional second-order hyperbolic equations[C]//2012 Fifth International Joint Conference on Computational Sciences and Optimization,2010:218-221.[11] 马月珍,李小刚,葛永斌.二维波动方程的高精度交替方向隐式方法[J].四川师范大学学报(自然科学版),2010,33(2):179-183.[12] 崔进,吴宏伟.一类波动方程初边值问题的高阶差分格式[J].应用数学,2014,27(1):166-174.[13] LELE S pact finite difference schemes with spectral-like resolution[J].J Comput Phys,1992,103(1):16-42.[14] 葛永斌,吴文权,田振夫.二维波动方程的高精度隐格式及其多重网格算法[J].厦门大学学报(自然科学版),2003,42(6):691-696.[15] 余跃玉.一种Caputo型时间分数阶波动方程的差分方法[J].四川师范大学学报(自然科学版),2014,37(4):524-528.2010 MSC：35A35; 65M99。