求解时间分布阶扩散方程的两个高阶有限差分格式
求解变系数对流扩散方程的高阶紧致差分格式
ห้องสมุดไป่ตู้Ｖｏ １ ．２ ７
Ｎｏ． １ １
重 庆 理 工 大 学 学 报（ 自然科 学 ）
Ｊ ｏ ｕ ｎａ ｒ ｌ ｏ ｆ Ｃ ｈ ｏ ｎ ｇ ｑ ｉ ｎ ｇ Ｕ ｎ ｉ ｖ ｅ ｒ ｓ ｉ ｔ ｙ ｏ ｆ Ｔ ｅ ｃ ｈ ｎ ｏ ｌ ｏ ｇ ｙ （ Ｎ ａ ｔ ｕ ｒ ａ ｌ Ｓ ｃ ｉ ｅ ｎ ｃ ｅ ）
对 流 扩散方 程 是一类 基本 的运 动方 程 ， 是描 述 黏性 流体 的非线 性 方 程 的线 性 化 模 型方 程 ， 它可 以用
Ａ Ｈｉ ｇ ｈ － ｏ ｒ ｄ ｅ ｒ Ｃｏ ｍｐ ａ ｃ ｔ Ｆ ｉ ｎ ｉ ｔ ｅ Ｄｉ ｆｅ ｒ ｅ ｎ ｃ ｅ Ｓ ｃ ｈ ｅ ｍｅ ｆ ｏ ｒ Ｓ ｏ ｌ ｖ ｉ ｎ ｇ ｔ ｈ ｅ Ｖａ ｒ ｉ ａ ｂ ｌ ｅ Ｃｏ ｅ ｉｃ ｆ ｉ ｅ ｎ ｔ Ｃｏ ｎ ｖ ｅ ｃ ｔ ｉ ｏ ｎ Ｄｉ ｆｕ ｓ ｉ ｏ ｎ Ｅｑ ｕ ａ ｔ ｉ ｏ ｎ ｓ
种新方法具有更好的健壮性 ， 并且可有效求解对流 占优问题。
关 键 词： 变量 替换 ； 紧致 差分 格式 ； Ｃ ｒ ａ ｎ ｋ Ｎ ｉ ｃ ｏ ｌ ｓ ｏ ｎ格 式 ； 无条 件稳 定 ； 对流 扩散 方程 文献标 识 码 ： Ａ 文章编 号 ： １ ６ ７ ４— ８ ４ ２ ５ （ ２ ０ １ ３ ） ｌ １ — ０ １ ２ ０— ０ ６ 中图分类 号 ： Ｏ ２ ４ １ ． ８ ２
ｔ ｉ ａ ｌ ｄ ｅ ｉ ｒ ｖ ａ ｔ ｉ ｖ ｅ ．Ｐ ｒ ｏ ｏ ￣ｏ ｆ ｕ ｎ ｃ ｏ ｎ ｄ ｉ ｔ ｉ ｏ ｎ ａ ｌ ｓ ｔ ａ ｂ ｉ ｌ ｉ ｔ ｙ ｏ ｆ ｔ ｈ ｅ ｓ ｅ ｎ ｅ ｗ ｓ ｃ ｈ ｅ ｍｅ ｓ ｗ ｅ ｒ ｅ ｇ ｉ ｖ ｅ ｎ ｉ ｎ ｔ ｈ ｅ ａ ｒ ｔ ｉ ｃ ｌ ｅ ．Ｃ ｏ ｍ． ｐ ａ ｒ ｅ ｄ ｗ ｉ ｔ ｈ ｔ ｈ ｅ ｓ ｔ ａ ｎ ｄ ａ ｒ ｄ ｃ ｅ ｎ ｔ ｒ ａ ｌ ｄ ｉ ｆ ｆ ｅ ｒ ｅ ｎ ｃ ｅ ｓ ｃ ｈ ｅ ｍｅ． ｔ ｈ ｅ ｎ ｅ ｗ ｍｅ ｔ ｈ ｏ ｄ ｓ ａ ｒ ｅ ｍｏ ｒ ｅ ｒ ｏ ｂ ｕ ｓ ｔ ｏｒ ｆ ｔ ｈ ｅ ｃ ｏ ｎ ｖ ｅ ｃ ｔ ｉ ｏ ｎ
时间分数阶亚扩散方程的高阶差分方法
考虑 如下 形式 的时 间分 数 阶亚扩 散方 程 川
ｒ ｃ Ｄ
．
（ ， ）＝ ａ ２ ｕ＋ ／ ． （ ， ｔ ） ，
（ ， ｔ ） ∈ ，× （ ０ ， Ｔ ） ， Ｔ ＞０ ， ０ ＜卢 ＜ １ ，
｛ （ ， ０ ）＝ 。 （ ） ，
Ｄ ｕ （ ）
［ Ｄ ｕ （ ） ］
（
（ ）
目前 ， 已经 有一 些方 法数 值模 拟形 如式 （ １ ） 亚扩 散方 程 ． 时间分 数 阶导 数 常用 的离 散 格 式是 Ｌ １ 方法 ， 具
收 稿 日期 ： ２ ０ １ ２—１ １—１ ２； 修 回 日期 ：２ ０ １ ３— ０ ｌ 一２ ６ 基金项 目： 上 海 市 教 委 科 研 创 新 重 点项 目（ １ ２ Ｚ Ｚ ０ ８ ４ ） 资 助 作者简介 ： 曾 凡海 （ １ ９ ８ ２一）， 男， 博士生． 研 究方 向 ： 分 数 阶微 分方 程 数 值 计 算 ， Ｅ － ｍａ ｉ ｌ ： ｆ ａ ｎ ｈ ａ ｉ ｚ ＠ｓ ｈ ｕ ． ｅ ｄ ｕ ． ｃ ｎ
一
ｃ
ｕ （ ， ￡ ） ＝Ｄ ｏ
，
’ ［ Ｏ ｔ ｕ （ ｘ ， ｔ ） ］
ｊ ｌ （ — ｓ ） 一 Ｏ ， Ｕ （ ， ｓ ） ｄ ｓ ，
（ ２ ）
其中Ｄ 是 分数 阶积 分 （ 或者 Ｒ ｉ ｅ ｍ ａ ｎ ｎ ． Ｌ ｉ ｏ ｕ ｖ ｉ ｌ ｌ ｅ积分 ） 算子 ， 其 定义 为
李 常品（ １ ９ ６ ８一）， 男， 博士 ， 教授 ． 研 究方 向 为分 数 阶微 分 方程 数 值 计 算 、 分数阶动力学 ， Ｅ — ｍａ ｉ ｌ ： ｌ ｃ ｐ ＠ｓ ｈ ｕ ． ｅ ｄ ｕ ． ｃ ｎ ｝通 讯 作 者 ：
对流扩散方程有限差分方法
对流扩散方程有限差分方法求解对流扩散方程的差分格式有很多种，在本节中将介绍以下3种有限差分格式：中心差分格式、Samarskii 格式、Crank-Nicolson 型隐式差分格式。
3.1 中心差分格式时间导数用向前差商、空间导数用中心差商来逼近，那么就得到了（1）式的中心差分格式]6[21111122h u u u vhu u au u nj n j n j nj n j n jn j -+-+++-=-+-τ（3）若令 haτλ=，2h vτμ=，则（3）式可改写为)2()(2111111nj n j n j n j n j n j n ju u u u u u u -+-+++-+--=μλ （4）从上式我们看到，在新的时间层1+n 上只包含了一个未知量1+n j u ，它可以由时间层n 上的值n j u 1-，n j u ，n j u 1+直接计算出来。
因此，中心差分格式是求解对流扩散方程的显示格式。
假定),(t x u 是定解问题的充分光滑的解，将1+n j u ，n j u 1+，n j u 1-分别在),(n j t x 处进行Taylor 展开：)(),(),(211ττO t u t x u t x u u njn j n j n j+⎥⎦⎤⎢⎣⎡∂∂+==++)(2),(),(322211h O x u h x u h t x u t x u u nj nj n j n j n j +⎥⎦⎤⎢⎣⎡∂∂+⎥⎦⎤⎢⎣⎡∂∂+==++)(2),(),(322211h O x u h x u h t x u t x u u njnj n j n j n j +⎥⎦⎤⎢⎣⎡∂∂+⎥⎦⎤⎢⎣⎡∂∂-==--代入(4)式，有 21111122),(hu u u vhu u au u t x T nj n j n j nj n j n jn j n j -+-+++---+-=τ)()()(2222h O v x u v h O a x u a O t u nj nj nj ⋅-⎥⎦⎤⎢⎣⎡∂∂-⋅+⎥⎦⎤⎢⎣⎡∂∂++⎥⎦⎤⎢⎣⎡∂∂=τ )()()(222h O v a O x u v x u a t u njnj nj ⋅-++⎥⎦⎤⎢⎣⎡∂∂-⎥⎦⎤⎢⎣⎡∂∂+⎥⎦⎤⎢⎣⎡∂∂=τ)(2h O +=τ显然，当0→τ，0→h 时，0),(→n j t x T ，即中心差分格式与定解问题是相容的。
时间分数阶慢扩散方程的一类有效差分方法
时间分数阶慢扩散方程的一类有效差分方法
在求解时间分数阶慢扩散方程时，可以使用一类有效的差分方法，如下所示：
首先，将时间区间 [0, T] 离散化为 N 个子区间，每个子区间的长度为Δt = T/N。
接下来，使用差分格式逐步逼近时间分数阶慢扩散方程的解。
假设 u(t) 是解的函数，我们可以使用如下的差分格式来近似
u(t) 在每个子区间的值：
u(tn+1) - u(tn) = Δt^α * D^β * u(tn) / Γ(α+1) + Δt^α * f(tn)
其中，tn 和 tn+1 是相邻的离散时间点，Δt = tn+1 - tn，f(tn) 是已知的源项函数，D^β 是时间分数阶微分算子。
根据差分格式，我们可以得到如下的迭代公式：
u(t0) = u0
u(tn+1) = (1 + Δt^α * D^β / Γ(α+1)) * u(tn) + Δt^α * f(tn)
对于初始条件，我们需要给出 u(t0) 的值，可以使用边界条件或者初始值条件来确定。
以上就是时间分数阶慢扩散方程的一类有效差分方法。
需要注
意的是，在实际求解过程中，我们需要根据具体的问题和数值条件来选择合适的差分格式和离散化步长。
求解空间分数阶扩散方程和对流扩散方程的有限差分格式研究
求解空间分数阶扩散方程和对流扩散方程的有限差分格式研究求解空间分数阶扩散方程和对流扩散方程的有限差分格式研究近年来，随着科学技术的不断发展，空间分数阶扩散方程和对流扩散方程的求解方法成为研究的热点之一。
这两类方程在工程和科学中的应用非常广泛，如热传导、空气污染模拟、地下水流动等。
为了更好地理解和应用这两类方程，研究人员通过有限差分格式的方法进行了深入研究。
首先，我们来介绍一下空间分数阶扩散方程。
空间分数阶扩散方程是一种描述扩散现象的方程，它的形式为：∂^αu(x,t)/∂t^α = D∇^2u(x,t)其中，u(x,t)是扩散物质的浓度分布函数，D是扩散系数，∂^α/∂t^α是Riemann-Liouville分数阶导数算子，α是分数阶。
对于空间分数阶扩散方程的求解，常用的方法是有限差分格式。
有限差分格式将连续的空间和时间离散化，将求解问题转化为代数问题。
通过将空间和时间网格分段，可以近似地表示方程中的偏导数。
在求解空间分数阶扩散方程时，有限差分格式需要考虑到方程中的分数阶导数，将其离散化为差分形式。
接下来，我们来介绍对流扩散方程的求解方法。
对流扩散方程是一种描述物质传输过程的方程，它的形式为：∂u(x,t)/∂t = D∇^2u(x,t) - v∇u(x,t)其中，u(x,t)是物质的浓度分布函数，D是扩散系数，v是流速。
对于对流扩散方程的求解，同样可以使用有限差分格式。
有限差分格式将方程中的偏导数转化为差分形式，通过逼近连续的空间和时间变量来得到代数格式。
在求解对流扩散方程时，有限差分格式需要考虑方程中的扩散项和对流项，并将其离散化。
在研究中，我们对空间分数阶扩散方程和对流扩散方程的有限差分格式进行了深入研究。
通过对方程的离散化，并应用差分格式求解代数问题，得到了方程的数值解。
我们通过调整网格大小和时间步长，对结果的精度进行了分析。
此外，我们还对不同边界条件下的方程求解进行了讨论。
研究结果显示，对于空间分数阶扩散方程和对流扩散方程，有限差分格式能够给出合理的数值解，并且能够有效地反映扩散现象和物质传输过程。
解扩散方程的指数时间差分方法
解扩散方程的指数时间差分方法指数时间差分方法（Exponential Time Differencing，简称ETD方法）是一种数值解扩散方程的方法，它通过将时间的离散化与指数函数的特性相结合，提高了计算效率和数值稳定性。
以下将对ETD方法进行详细介绍。
一、基本原理考虑一维扩散方程：∂u/∂t=D∂²u/∂x²其中，u是扩散物质的浓度，D是扩散系数。
二、离散化将时间离散化，令t = nh，其中，n为离散时间步长的索引，h为时间步长。
使用ETD方法后的求解格式如下：u(n+1)=e^(-hDk²)u(n)+[1-e^(-hDk²)]u(n)其中，k为空间离散化步长。
三、指数函数的近似计算ETD方法的关键在于指数函数的近似计算，常用的计算方法有：1. Padé展开：将指数函数在一些点进行泰勒展开，然后用有理函数近似，然后求解所得的微分方程系统。
这种方法具有高精度和高效率的优势。
2. Caley型：将指数函数通过特征多项式进行近似。
这种方法具有高阶精度。
3.向量化方法：将指数函数的计算转化为向量运算，提高计算效率。
四、算法流程使用ETD方法求解扩散方程的基本流程如下：1.确定求解区域和初始条件。
2.选择合适的离散化步长k和时间步长h。
3.将扩散方程的时间部分离散化，并利用指数函数的近似计算方法进行计算。
4.对空间部分进行差分离散化。
5.将时间离散化的方程和空间离散化的方程通过时间推进方法进行求解。
6.循环进行步骤3~5，直到达到所需的时间步数。
五、优缺点ETD方法相对于传统差分方法有以下优点：1.高效性：指数时间差分方法能够对指数函数进行有效计算，提高了计算速度。
2.数值稳定性：该方法具有良好的数值稳定性，可以更准确地求解扩散方程。
3.高精度：ETD方法的数值精度较高，可以减小数值误差。
然而，ETD方法也存在一些缺点：1.适用性有限：ETD方法主要适用于线性扩散方程，对非线性扩散方程的求解效果有限。
有限差分求解扩散的数值解
Jeff
C2
无孔结构
J0 =-D ∂C ∂z
多孔介质
ε ∂C Jeff =- D τ ∂z
其中ε为孔隙率
曲折因子
τ= J0 ε Jeff
3.扩散模拟表征多孔介质物质输运性能（相转换流延制备） SRCT构建孔结构数组
灰度图
二值图
黑色：孔 白色：实体 孔隙率：0.35
三维结构重构图A
3.扩散模拟表征多孔介质物质输运性能
1.有限差分的基本原理（ 一维扩散例子） 计算结果
稳态下
∂C ∂2 C =D 2 = 0 ∂t ∂x
2.边界条件处理（ 二维柱坐标扩散例子）
扩散进口
扩散出口
绝缘 环
圆柱上底面和侧面为高浓度面（红色） 圆柱下底面同心小圆为低浓度面（白色） 圆柱下底面圆环为绝缘面 （蓝色）
2.边界条件处理（ 二维柱坐标扩散例子） 0 Z
r
扩散微分方程： ∂C ∂2 C 1 ∂C ∂2 C =D( 2 + r + 2) ∂t ∂r ∂r ∂Z
圆心对称
柱坐标下取体积微元由质量衡算推导扩散偏微分方程
r+ △r
△Z
r
体积微元质量变化量
∆m=-2π r+∆r ∆Z∆t*J -π r+∆r 2 -r2 ∆t*J
移项化简得
r+∆r+2πr∆Z∆t*J r Z
2 2 + π r+ ∆r -r ∆t*J Z+∆Z
∆m 1 -J = * 2πr∆Z∆r∆z r
r+∆r *
即
r+∆r -J r *r J Z+∆Z -J Z ∆r ∆Z ∂C ，将 J=-D 代入可得 ∂x
有限差分法基本原理课件
离散网格点
有限差分法基本原理
差分和逼近误差
差分概念：
设有x的解析函数 y f(x) ，函数y 对x 的导
数为：
d yli m ylim f(x x )f(x )
dx x 0 x x 0
x
dy dy 、dx 分别是函数及自变量的微分，dx 是函数 对自变量的导数，又称微商。上式中的y 、x 分别称 为函数及其自变量的差分，y 为函数对自变量的差商。
有限差分法基本原理
模型方程
为了抓住问题的实质，同时又不使讨论的问题过于
复杂，常用一些简单的方程来模拟流体力学方程进行讨 论分析，以阐明关于一些离散方法的概念。这些方程就 叫做模型方程。常用的模型方程：
对流方程：
0
t x
对流－扩散方程：
t
x
2x2
热传导方程：
2
t
x2
有限差分法基本原理
Poisson方程：
*n i
为差分方程的近似数
值
解，之间的误差为 。同样，近似数值解也满足同样的方程：
T i* n 1 S* i n T 1 ( 1 2 S )T * in S* i n T 1
in 1 Sn i 1 ( 1 2 S )n i Sn i 1
上式称为误差传播方程。
有限差分法基本原理
x等价定理
2
x2
2
y2
f
2 2
Laplace方程： x2 y2 0
有限差分法基本原理
差分方程的建立过程
以对流方程说明差分方程的建立过程。
0
t x
(x,0) (x)
有限差分法基本原理
差分方程的建立过程
1.划分网格
计算传热学第4讲扩散方程的数值解
w
e
w
( S c S p T ) A ( x ) dx
(6)
( S c A ) P x ( S p A ) P xT P
在上面的积分过程中，我们假定：
–
–
待求变量T在控制容积P上为常数 整个控制容积的A(x)为常数，且等于P点的值。
4.1 一维稳态导热问题的数值解
将（5）和（6）代入方程（4），
其中，
4.1 一维稳态导热问题的数值解
a P T P a W TW a E T E b P (8)
其中，
A aW x w (8a) A aE x e
(8c)
(8b)
a P aW a E ( S p A ) P x
网格是用外节点法划分的
边界上出现半个控制容积 边界节点的差分方程可以用 下述方法推出： qB 一阶精度的Taylor级数展开 法
dT dx
x0
x=0
e 1 e (x)1 (x)2 2 3
T 2 T1 ( x ) 1
(9)
T 2 T1 dT 1 qB ( x ) 1 dx x 0
边界条件的处理
整理后得到：
T1 T 2 ( x ) 1
1
qB
(10)
特点： 最简单的处理方法 只有一阶精度 与控制方程的精度不匹配
边界条件的处理
–
元体能量平衡法：
在研究边界节点所代表的控 制容积（元体）的能量平衡
–
x=0
½(x)1 S
e q1 2 e 3 (x)2
插值法解时间分数阶扩散方程
ｓｍａｌｌ ｉｎｔｅｒｖａｌ Ｅｔ 一１，ｔ ］（２≤ ≤ Ｎ），ｗｈｉｌｅ ｔｈｅ ｌｉｎｅａｒ ｉｎｔｅｒｐｏｌａｔｉｏｎ ａｐｐｒｏｘｉｍａｔｉｏｎ ｉｓ ａｐｐｌｉｅｄ ｏｎ ｔｈｅ ｆｉｒｓｔ ｓｍａｌｌ ｉｎｔｅｒｖａｌ Ｅｔ０，ｔ１］．Ｕｓｉｎｇ ｔｈｅ ｅｎｅｒｇｙ ｎｏｒｍ，ｔｈｅ ｕｎｃｏｎｄｉｔｉｏｎａｌ ｓｔａｂｉｌｉｔｙ ａｎｄ ｃｏｎｖｅｒｇｅｎｃｅ ｏｆ
ＭＲ ｓｕｂｊｅｃｔ ｃｌａｓｓｉｆｉｃａｔｉｏｎ：３５Ｒ１ １
近 年来 ，分数 阶微 积分 理 论和应 用 发展 迅速 ，在 一 定 程度 上 填补 了传 统整 数 阶微积 分 的不 足 。分 数 阶微 积 分 已被用 来 描 述 应用 学 科 里 的许 多 现 象 ，例 如 多孔材 料 中 的液 体 流动 、反 常扩散 运 动 ，力 学里 的 自相 似结 构 、信号 处理 等 ］。
第五章对流扩散问题(假扩散 高阶格式第一个内节点的差分方程)
第五章 对流扩散问题———假扩散
1.0 0.8 0.6 0.4 0.2 0.0 0.0 t=0 AS t=4 UD c=1.0 t=8
0.6 0.8 1.0 x c值进一步增大，我们看到：c=0.010.050.0750.40.81.0 ，而 FD/FDc=0.01=10.960.9340.6060.2020，可见此时UD格式的整体假扩散 已经衰减了100%，使其计算结果和精确解完全一样。
第五章 对流扩散问题———假扩散
以上结果再次证明了我们关于时间项和对流项采用一阶
格式时，对其假扩散的分析。 也证明了我们提出了二阶
迎风格式和QUICK格式确实能够消除由对流项离散而引
起的假扩散。
第五章 对流扩散问题———假扩散
例 设有下列一维非稳态对流问题
u , t x 0 x 1, u 0.1
n i 1 / 2
n i 1 / 2
)
n i
ut c x
第五章 对流扩散问题———假扩散
其中，
n 1 i
c(
n i 1 / 2
n i 1 / 2
)
n i
ut c x
对SUD
n n i 1 / 2 i 1 / 2
n n 3 n 4 i 1 i2 i 2
第五章 对流扩散问题———假扩散
UD格式
x 1 / 100
xi xi 1 x
t cx / u
t=0
i= 1
2
3
4
98
99
100
101
0 i 由初始条件确定
t=t
i= 1
2
3
4
时间分数阶扩散方程的三次样条差分格式
时间分数阶扩散方程的三次样条差分格式时间分数阶扩散方程的三次样条差分格式如下：$$。
\frac{\partial^{\alpha}u(x_i, t^n)}{\partial t^{\alpha}} = D \frac{\partial^2u(x_i, t^n)}{\partial x^2}, \quad 0 < \alpha < 1。
$$。
其中，$\alpha$为时间分数阶，$D$为扩散系数，$u(x_i, t^n)$为在$x_i$处和$t^n$时刻的解。
为方便表示，下面将$t^n$简单表示为$t$。
采用前向差分和后向差分的一阶导数定义，可以得到时间导数的差分格式：$$。
\frac{\partial^{\alpha}u(x_i, t)}{\partial t^{\alpha}}\approx \frac{u(x_i, t) - u(x_i, t - \Delta t)}{\Deltat^{\alpha}}, \quad 0 < \alpha < 1。
$$。
其中，$\Delta t$为时间步长。
由于采用了前向差分，此格式为显式差分格式。
若采用三点中心差分法近似求解二阶导数，则。
$$。
\frac{\partial^2u(x_i, t)}{\partial x^2} \approx\frac{u(x_{i+1},t) - 2u(x_i,t) + u(x_{i-1},t)}{\Delta x^2}。
$$。
其中，$\Delta x$为空间步长。
将上式代入原方程得到：$$。
\frac{u(x_i, t) - u(x_i, t - \Delta t)}{\Delta t^{\alpha}} = D \frac{u(x_{i+1},t) - 2u(x_i,t) + u(x_{i-1},t)}{\Delta x^2}。
$$。
再将时间步长和空间步长合并，得到：$$。
u_i^n - u_i^{n-1} = D\Delta t^{\alpha} \frac{u_{i+1}^n -2u_i^n + u_{i-1}^n}{\Delta x^2}。
第五章对流扩散问题(假扩散 高阶格式差分方程的求解)
第五章 对流扩散问题———假扩散
在P点的控制容积上对上边的 NW 微分方程进行积分，则：
1 6 N 5 w P E NE
u(Te Tw ) v(Tn Ts ) 0
W 2 SW
A 迎风格式
u(Te Tw ) v(Tn Ts ) 0 ue / w 0
n bb a EE n a EE WW WW b
第五章 对流扩散问题———假扩散
A5.3 PDMA算法 对五对角阵，有没有类似三对角阵TDMA那样的直接 求解方法呢？实际上对五对角阵，人们也可以找到相 应的直接求解方法，这个方法就是PDMA算法。下边 以一维为例来介绍这个算法。
第五章 对流扩散问题———假扩散
由此可见，对一维而言，所得到的差分方程不再是可以 用TDMA直接求解的三对角矩阵方程，而变成一个五对 角阵方程。对二维或三维而言，逐线联立求解的方程组 也不再是可以用TDMA直接求解的三对角矩阵方程，也
变成一个五对角阵方程。
那么，针对这样一种五对角方程，通常我们如何来求解 它们呢？
假设经过代入后得到的上三角方程为
i Ai i 2 Bi i 1 Ci
*#
问题的关键就变成为：找出系数Ai, Bi与Ci和系数a, b, c, d
, e及f之间的关系。为此，写出i-2点的上三角方程如下：
i 2 Ai 2i Bi 2i 1 Ci 2
****
(****)*(di+ eiBi-2)，有：
(di ei Bi 2 )i 1 Ai 1 (di ei Bi 2 )i 1 Bi 1 (di ei Bi 2 )i Ci 1 (di ei Bi 2 )
扩散方程的差分解法
扩散方程的差分解法在研究热传导过程、扩散过程、边界层现象时，我们常常遇到抛物型方程，这类方程中最典型、最简单的就是热传导方程。
热传导方程中的自变量中包括时间t ，它是描述一种随时间变化的物理过程，即所谓不定常现象。
这类问题的基本定解问题应是初值问题，即在初始时刻（t=0）时给定定解条件，求解t>0时的解。
本文主要运用有限差分法对一维扩散方程进行求解，并对差分解的适定性、相容性、收敛性及稳定性进行分析，同时与解析解进行对比。
1.扩散方程一维扩散方程为：22u u t xα∂∂=∂∂ （1）式中，u 为因知量，α为扩散系数，x 为坐标，t 为时间。
其定解条件如下： 初始条件： (,0)() 0x u x f x L =≤≤（2）边界条件： 12(0,)() , (,)()u t f t u L t f t ==（3） 一般假定函数()f x ，1()f t ，2()f t 满足连接条件，即1(0)(0) f f =，2()(0) f L f =。
2.有限差分法有限差分法是数值计算解微分方程古老的方法之一，也是系统化地、数值地求解数学物理方法的方程。
其控制方程中的导数用离散点上函数值的差商代替。
差分格式可以分为显格式和隐格式。
所谓显格式是指在任一结点上因变量在新是时间层上的值可以通过之前的时间层上相邻结点变量的值显式解出来。
由于这些层的变量值是已知的，当时间向前推进时，空间点上的新的变量值就只需逐点计算就行了，因此显格式计算起来比较省事。
隐格式则是指任一结点上变量在新的时间层的值，不能通过之前的时间层上相邻结点的值显式解出来，它不仅与之前的时间层上的已知值有关，而且也与新时间层的相邻结点的变量值有关。
因而一个差分方程常常包括几个相邻结点上的未知数，未知数的个数取决于格式的构成形式。
为了解出这些未知数需要联立新的方程，而每引进一个新的方程往往又同时引进了新的未知数。
因此，隐格式总是伴随着求解巨大的代数方程组。
- 1、下载文档前请自行甄别文档内容的完整性，平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益，请联系客服反馈,我们会尽快为您处理(人工客服工作时间：9:00-18:30)。
应用数学和力学，第 ４０ ２０１９ 年 ７ 月 １ 日出版
卷
第
７
期 Leabharlann Ａｐｐｌｉｅｄ Ｍａｔｈｅｍａｔｉｃｓ ａｎｄ Ｍｅｃｈａｎｉｃｓ Ｖｏｌ．４０，Ｎｏ．７，Ｊｕｌ．１，２０１９
文章编号：１０００⁃ ０８８７（２０１９）０７⁃ ０７９１⁃ １０
ⓒ 应用数学和力学编委会，ＩＳＳＮ １０００⁃ ０８８７
［
ｈ（
β
（ ｌ
１）
）
＋
ｈ（
β
（ ２） ｌ
）
］
＋
２５
Δβ４ × １３５
ｈ（４）（ η） ，
η ∈ （０，１） ．
１．２ 时间离散
接下来考虑时间导数的离散．时间分布阶扩散方程的时间导数是 Ｃａｐｕｔｏ 型分数阶导数，为
了使构造的格式在时间方向达到三阶精度，本文采用一种离散左 Ｒｉｅｍａｎｎ⁃Ｌｉｏｕｖｉｌｌｅ 分数阶导
１ 准 备 工 作
１．１ 数值积分
先介绍两个关键的数值积分公式：复化 Ｓｉｍｐｓｏｎ 公式和复化两点 Ｇａｕｓｓ⁃Ｌｅｇｅｎｄｒｅ 公式．设 Ｊ
是正整数，
Δβ ＝ １ ／ Ｊ， βｌ ＝ ｌΔβ， ｌ ＝ ０，１，…，Ｊ；
βｌ＋１／ ２
＝ （βｌ
＋
βｌ＋１）
／ ２，
β（１） ｌ
而在文献13中他们又采用加权?平移的grnwald?let?nikov公式得到了求解分布阶扩散方程的两种高阶差分格式其收敛速率分别为o2h2197应用数学和力学第40卷第7期2019年7月1日出版appliedmathematicsandmechanicsvol40no7jul12019?收稿日期
＝ βｌ＋１／ ２
－ Δβ ／ ２
３，
β（２） ｌ
＝ βｌ＋１／ ２
＋ Δβ ／ ２
３，
ｌ ＝ ０，１，…，Ｊ － １．
引理 １（ 复化 Ｓｉｍｐｓｏｎ 公式） 设 ｈ（ β） ∈ Ｃ４（ ［０，１］ ） ， 有
∫ ∑
１
ｈ（ β） ｄβ
０
＝
Δβ ６
Ｊ－１ ｌ＝０
［ ｈ（ βｌ ）
＋
４ｈ（βｌ＋１／ ２）
＋
ｈ（ βｌ ＋ １ ） ］
－
Δβ４ １８０
ｈ（４）（ η） ，
η ∈ （０，１） ．
引理 ２［１５］（ 复化两点 Ｇａｕｓｓ⁃Ｌｅｇｅｎｄｒｅ 公式） 设 ｈ（ β） ∈ Ｃ４（ ［０，１］ ） ， 有
∫ ∑
１
ｈ（ β） ｄβ
０
＝
Δβ ２
Ｊ－１ ｌ＝０
∗ 收稿日期： ２０１８⁃１２⁃２５； 修订日期： ２０１９⁃０３⁃０９ 基金项目： 国家自然科学基金（１１４７１２６２） 作者简介： 胡嘉卉（１９８０—） ，女，博士生（ Ｅ⁃ｍａｉｌ： ｈｕｊｈ＠ ｍａｉｌ．ｎｗｐｕ．ｅｄｕ．ｃｎ） ； 聂玉峰（１９６８—） ，男，教授，博士生导师（ 通讯作者． Ｅ⁃ｍａｉｌ： ｙｆｎｉｅ＠ ｎｗｐｕ．ｅｄｕ．ｃｎ） ．
求解时间分布阶扩散方程的 两个高阶有限差分格式∗
胡嘉卉１，２， 王俊刚１， 聂玉峰１
（１． 西北工业大学 应用数学系， 西安 ７１０１２９； ２． 河南工业大学 理学院， 郑州 ４５０００１）
摘要： 基于复化 Ｓｉｍｐｓｏｎ 公式和复化两点 Ｇａｕｓｓ⁃Ｌｅｇｅｎｄｒｅ 公式，构造了两个求解时间分布阶扩散 方程的高阶有限差分格式．不同于以往文献中提出的时间一阶或二阶格式，这两种格式在时间方向 都具有三阶精度，而在分布阶和空间方向可达到四 阶 精 度． 数 值 结 果 表 明，两 种 算 法 都 是 稳 定 且 收 敛的，从而是有效的．两种格式的收敛速率也通过数值实验进行了验证，并且通过和文献中的算法 对比可以得出其更为高效．
关 键 词： 时间分布阶扩散方程； 分数阶导数； 高阶差分格式； 收敛速率
中图分类号： Ｏ２４２．２ 文献标志码： Ａ
ＤＯＩ： １０．２１６５６ ／ １０００⁃０８８７．３９０３５８
引 言
近年来，分数阶微分方程（ ｆｒａｃｔｉｏｎａｌ ｄｉｆｆｅｒｅｎｔｉａｌ ｅｑｕａｔｉｏｎ，ＦＤＥ） 因其非局部特性，在研究具 有历史记忆性的过程和具有遗传特性的材料方面得到了广泛的应用［１⁃３］ ． 而分布阶 ＦＤＥ 作为 ＦＤＥ 的一种，最早是由 Ｃａｐｕｔｏ 在 ２０ 世纪 ６０ 年代提出的，用于模拟滞弹性介质中的应力⁃应变 行为［４］ ．不同于常数阶 ＦＤＥ 和多项 ＦＤＥ，分布阶 ＦＤＥ 是通过把 ＦＤＥ 中分数阶导数的阶在某个 给定区间上积分得到的，这类方程可以看作是前面提到的两类 ＦＤＥ 的推广．分布阶 ＦＤＥ 的典 型应用是模拟减速的次扩散过程［５⁃７］ ，在这类过程中，粒子束的均方位移关于时间以对数增长， 最终将形成超慢扩散．近年来，分布阶 ＦＤＥ 还被应用于控制和信号过程［８］ 、电介质感应和扩散 的模拟［９］ 以及识别系统［１０］ 等多个研究领域．
由于分布阶 ＦＤＥ 的广泛应用，方程的求解就成为一个重要问题．在大多数情况下，要求得 解析解非常困难，甚至是不可能的，因此，研究可靠而有效的数值解法就势在必行． 数值求解分 布阶 ＦＤＥ 的相关研究目前还处于初始阶段，其中，有限差分法以其格式便于构造、编程方便等 优点成为采用较多的一种方法．Ｆｏｒｄ 等［１１］ 提出了一种求解时间分布阶扩散方程的隐式差分格 式．Ｇａｏ（ 高广花） 等［１２］ 采用 Ｇｒüｎｗａｌｄ⁃Ｌｅｔｎｉｋｏｖ 公式，并应用外推技术提高近似精度，提出了两 种求解一维分布阶 ＦＤＥ 的差分格式；而在文献［１３］ 中，他们又采用加权⁃平移的 Ｇｒüｎｗａｌｄ⁃Ｌｅｔ⁃ ｎｉｋｏｖ 公式，得到了求解分布阶扩散方程的两种高阶差分格式，其收敛速率分别为 Ｏ（ τ２ ＋ ｈ２ ＋
数的加权⁃平移 Ｇｒüｎｗａｌｄ⁃Ｌｅｔｎｉｋｏｖ 公式［１６］ 离散时间导数．
７９１
７９２
胡 嘉 卉 王 俊 刚 聂 玉 峰
Δβ２） 和 Ｏ（τ２ ＋ ｈ４ ＋ Δβ４）， 这里 τ，ｈ 和 Δβ 分别表示时间、空间和分布阶三个方向的步长．文献 ［１４］ 提出了求解带有非线性源项的二维和三维时间分布阶波动方程的差分格式，并对解进行 了模拟．本文基于复化 Ｓｉｍｐｓｏｎ 公式和复化两点 Ｇａｕｓｓ⁃Ｌｅｇｅｎｄｒｅ 公式，提出了求解一维时间分 布阶扩散方程的两种高阶差分格式．不同于以往文献中提出的时间一阶或二阶格式，本文提出 的格式在时间方向可达到三阶精度，而且在分布阶和空间方向都具有四阶精度．本文在第 ３ 节 中列出的数值结果表明，所提出的两种格式是稳定且收敛的，并且通过与文献［ １３］ 和文献 ［１１］ 中的算法对比，得出其更为高效．