ABSTRACT

The paper uses the Lagrange’s form of one-dimensional radial basis function (abbr. RBF) interpolation with zero-degree algebraic precision to give the finite difference formulas (abbr. RBF-FD) for approximation of first derivative and second derivative of the interpolated function at the central node of the three equidistant nodes. In particular, the approximation error of second derivative of the RBF difference formula is analyzed, and the optimal parameter of the radial basis function, which makes the approximation order of RBF-FD reach the highest level, is obtained. Then, using these RBF-FD formulas, the RBF-FD scheme of the non-homogeneous two-point boundary value problem with small viscosity coefficient is given. The convergence order of the RBF-FD scheme constructed in this paper is two times of the polynomial finite difference scheme under the same node stencil. While, the calculating time of the RBF-FD scheme is only slightly increased. Numerical experiments indicate that the RBF-FD scheme can maintain the fourth-order convergence order under viscosity coefficient greater than or equal to 10⁻³.

Keywords:The Finite Difference Based on Radial Basis Function, Best Parameter Value, Convergence Order, Non-Homogeneous Two-Point Boundary Value Problems, Small Viscosity Coefficient

基于径向基有限差分法求解带小粘性系数 非齐次两点边值问题

段俊娜，郭子滔，冯仁忠

北京航空航天大学数学与系统科学学院，数学、信息与行为教育部重点实验室，北京

收稿日期：2017年12月15日；录用日期：2018年1月11日；发布日期：2018年1月18日

摘 要

本文利用具有零次代数精度的一元径向基函数(Radial Basis Function，简记RBF)插值的Lagrange形式，给出在三等距节点的中心节点处逼近被插函数的一阶导数和二阶导数的有限差分(简记RBF-FD)公式。特别地，对二阶导函数径向基差分的逼近误差进行分析，给出了使其逼近阶达到最高的径向基函数的最佳参数值。然后，利用这些RBF-FD公式给出了求解带小粘性系数的非齐次两点边值问题的差分格式，使所构造格式的收敛阶是同节点模板的多项式差分格式收敛阶的2倍，而计算时间略有增加。数值实验表明所构造的RBF-FD格式在小粘性系数大于等于10⁻³的情况下能够保持四阶收敛速度。

关键词 :径向基函数有限差分，最佳参数，收敛阶，非齐两点边值问题，小粘性系数

Copyright © 2018 by authors and Hans Publishers Inc.

This work is licensed under the Creative Commons Attribution International License (CC BY).

http://creativecommons.org/licenses/by/4.0/

1. 引言

有限差分(finite difference，简记FD)是光滑函数导数逼近的最经典方法 [1] ，它是用函数在某点附近的函数值的线性组合来逼近函数在该点的导数值，一般是通过递推方法来定义，广泛应用于数值求解微分方程。常用的逼近一元函数导数的差分公式为基于三等距节点 $\left\{x_{i-1},x_i,x_{i+1}:x_{i+1}-x_i=x_i-x_{i-1}=h\right\}$ 的一阶中心差分和二阶中心差分，其截断误差分别为

$\frac{f\left(x_{i+1}\right)-f\left(x_{i-1}\right)}{2h}=f^{\prime }\left(x_i\right)+\frac{h^2}{6}{f}^{‴}\left(x_i\right)+O\left(h^3\right),$ (1.1)

$\frac{f\left(x_{i+1}\right)-2f\left(x_i\right)+f\left(x_{i-1}\right)}{h^2}=f^{″}\left(x_i\right)+\frac{h^2}{12}{f}^{\left(4\right)}\left(x_i\right)+O\left(h^4\right).$ (1.2)

显然，中心差分公式(1.1)，(1.2)的截断误差仅达到二阶精度，且由于差分公式中的系数项均是固定常数，因此其逼近阶不可能再提高。

径向基函数(Radial Basis Function，简记RBF)插值是针对散乱数据插值提出的 [2] [3] [4] 。表1列举了一些常用的无限光滑径向基函数，它们都含有一个形状参数 $\epsilon$ ，该形状参数的不同取值会影响径向基插值的精度和稳定性 [5] 。

G. B. Wright和B. Fornberg在文献 [7] 中利用文献 [6] 的结果证明了由径向基函数插值产生的有限差分收敛到经典的多项式差分，当 $\epsilon \to \text{0}$ 时。不同的形状参数 $\epsilon$ 的取值会影响此类差分公式的逼近精度。但到目前为止，几乎没有理论分析使差分公式逼近精度达到最佳时参数 $\epsilon$ 的精确表达式。人们一般通过数值测试找出拟最佳 $\epsilon$ 值，如 [7] [8] 。我们对本文构造的RBF-FD公式的截断误差进行Taylor展开分析，给出使其逼近阶达到最高的参数 $\epsilon$ 的精确表达式，并将此类差分公式用于求解带小粘性系数的非齐两点边值问题。

表1. 几类典型的全局径向基函数

在应用科学与工程技术领域中，许多实际的问题经常涉及到求解含有极小参数的微分方程，当方程中的参数值有很小的变化时，会导致方程的解产生很大变动，如带小粘性系数的非齐次两点边值问题：

$\{\begin{array}{l}{\left(\nu \left(x\right)u^{\prime }\left(x\right)\right)}^{\prime }+p\left(x\right)u^{\prime }\left(x\right)+q\left(x\right)u\left(x\right)=f\left(x\right),a<x<b\\ u\left(a\right)=\alpha ,u\left(b\right)=\beta ，\end{array}$ (1.3)

其中函数 $\nu \left(x\right),p\left(x\right),q\left(x\right),f\left(x\right)\in C^2\left[a,b\right]$ ，且 $0<\nu \left(x\right)\le 1$ ， $q\left(x\right)\ge 0$ 。由于求解带小粘性系数的边值问题不论在理论上还是在数值方法上都要比一般两点边值问题困难 [9] [10] ，目前对于带小粘性系数两点边值问题，主要的数值解法有差分余项反向修正法 [11] 、高精度差分法 [12] 、隐式紧致差分法 [13] 、混合型高精度紧致差分法 [14] 等。本文中，我们利用构造的径向基函数有限差分(简记RBF-FD)，设计求解带小粘性系数非齐次两点边值问题(1.3)的一种四阶精度的径向基有限差分格式，使所构造格式的收敛阶是同节点模板的多项式差分格式收敛阶的2倍，此外，数值实验表明所构造的RBF-FD格式在小粘性系数大于等于 ${10}^{-3}$ 的情况下能够保持四阶收敛速度。

文章余下章节安排如下：第二节构造了在三等距节点下逼近一元函数的一阶和二阶导数的RBF-FD公式，并通过对其截断误差分析得到最佳参数 $\epsilon$ 的表达式，以保证RBF-FD公式达到最佳逼近阶；第三节利用第二节构造的逼近二阶导函数的RBF-FD公式，给出了求解带小粘性系数的非齐两点边值问题的RBF-FD格式，并巧妙地利用微分方程给出了最佳参数 $\epsilon$ 的计算式，避免了直接估值带来的困难；第四节是结论部分。

2. 基于三等距节点的RBF-FD公式的构造及最佳参数e的选择

2.1. 基于三等距节点的RBF-FD公式

本文中，我们使用一元径向基函数插值来构造差分公式。假设给定三节点数据组为 ${\left\{\left(x_i,u\left(x_i\right)\right)\right\}}_{i=1}^3\subset {ℝ}^2$ ，且节点是等距的，即 $h=x_1-x_3=x_2-x_1$ 。这里使用带常数项的径向基函数

$s\left(x\right)={\displaystyle \underset{j=1}{\overset{\text{3}}{\sum }}{\lambda }_j\phi \left(\Vert x-x_j\Vert \right)}+\beta ,$ (2.1)

作为节点处的插值函数，其中 $\phi \left(r\right)$ 是径向基函数， $\Vert \cdot \Vert$ 是 ${ℝ}^{\text{s}}$ 空间中的欧氏距离。方程(2.1)中的系数 $\left\{{\lambda }_1,{\lambda }_2,{\lambda }_3,\beta \right\}$ 由下列插值条件和额外条件

$\{\begin{array}{l}s\left(x_i\right)=u\left(x_i\right),i=1,2,3,\\ {\displaystyle \underset{j=1}{\overset{3}{\sum }}{\lambda }_j}=0,\end{array}$ (2.2)

决定。方程组(2.2)的矩阵向量形式为

$\underset{A}{\underset{︸}{\left(\begin{array}{cccc}\phi \left(0\right)& \phi \left(\Vert x_1-x_2\Vert \right)& \phi \left(\Vert x_1-x_3\Vert \right)& 1\\ \phi \left(\Vert x_2-x_1\Vert \right)& \phi \left(0\right)& \phi \left(\Vert x_2-x_3\Vert \right)& 1\\ \phi \left(\Vert x_3-x_1\Vert \right)& \phi \left(\Vert x_3-x_2\Vert \right)& \phi \left(0\right)& 1\\ 1& 1& 1& 0\end{array}\right)}}\left(\begin{array}{c}{\lambda }_1\\ {\lambda }_2\\ {\lambda }_3\\ \beta \end{array}\right)=\left(\begin{array}{c}u\left(x_1\right)\\ u\left(x_2\right)\\ u\left(x_3\right)\\ 0\end{array}\right)$ (2.3)

插值函数(2.1)在节点 $x_1$ 处的 $p$ 阶导数为

$s^{\left(p\right)}\left(x_1\right)={\displaystyle \underset{j=1}{\overset{3}{\sum }}{\lambda }_j{\phi }^{\left(p\right)}\left(\Vert x_1-x_j\Vert \right)},p=1,2.$

我们将 $s^{\left(p\right)}\left(x_1\right)$ 称为基于径向基函数插值的 $p$ 阶有限差分，简记为 $p$ 阶RBF-FD。一般有 $s^{\left(p\right)}\left(x_1\right)\approx u^{\left(p\right)}\left(x_1\right)$ 。如果要得到 $s^{\left(p\right)}\left(x_1\right)$ 的表达式，则需要符号求解方程(2.3)，这在一般情况下比较困难。下面我们通过插值函数 $s\left(x\right)$ 的Lagrange形式

$s\left(x\right)={\displaystyle \underset{j=1}{\overset{3}{\sum }}{\psi }_j\left(x\right)u\left(x_j\right)}.$ (2.4)

来给出 $s^{\left(p\right)}\left(x_1\right)$ 的精确表达式，其中 ${\psi }_j\left(x\right),j=1,2,3$ 是Lagrange基函数。根据文献 [15] 有

${\psi }_j\left(x\right)=\frac{\det \left(A_j\left(x\right)\right)}{\det \left(A\right)}.$ (2.5)

其中 $A_j\left(x\right),j=1,2,3$ 是(2.3)中矩阵 $A$ 的第 $j$ 个行向量被以下向量

$B\left(x\right)=\left(\phi \left(\Vert x-x_1\Vert \right)\phi \left(\Vert x-x_2\Vert \right)\phi \left(\Vert x-x_3\Vert \right)1\right)$ .

置换后的矩阵。借助(2.4)有

$s^{\left(p\right)}\left(x_1\right)={\displaystyle \underset{j=1}{\overset{3}{\sum }}{\psi }_j^{\left(p\right)}\left(x_1\right)u\left(x_j\right)},p=1,2.$

在文章余下部分，我们取径向基函数为 $\phi \left(r\right)=\sqrt{1+{\left(\epsilon r\right)}^2}$ (Multiquadric，MQ函数)来进行讨论。根据公式(2.5)， $s^{″}\left(x_1\right)$ 中的系数 ${{\psi }^{″}}_j\left(x_1\right),j=1,2,3$ 的表达式分别为

${{\psi }^{″}}_{\text{3}}\left(x_1\right)=-\frac{{\epsilon }^{2}h}{a\left(b-1\right)},{{\psi }^{″}}_{\text{1}}\left(x_1\right)=0,{{\psi }^{″}}_{\text{2}}\left(x_1\right)=\frac{{\epsilon }^{2}h}{a\left(b-1\right)}.$

其中

$a=\sqrt{1+{\epsilon }^2{h}^2},b=\sqrt{1+4{\epsilon }^2{h}^2}.$ (2.6)

于是有

$s^{″}\left(x_1\right)=\frac{{\epsilon }^2\left(a^3-1\right)}{a^3\left(b-4a+3\right)}\left(2u\left(x_1\right)-u\left(x_3\right)-u\left(x_2\right)\right).$ (2.7)

2.2. 最佳参数e的选择

假设被插值函数 $u\left(x\right)$ 在包含节点组 $\left\{x_3,x_1,x_2\right\}$ 的区间 $I$ 上充分光滑，且不妨记 $x_{\epsilon }=\epsilon h$ ，则 ${{\psi }^{″}}_{\text{2}}\left(x_1\right),2u\left(x_1\right)-u\left(x_3\right)-u\left(x_2\right),s^{″}\left(x_1\right)$ 关于 $x_{\epsilon }$ 和 $h$ 的Taylor展开式分别为

$-\frac{{\epsilon }^2\left(a^3-1\right)}{a^3\left(b-4a+3\right)}=-{\epsilon }^2\left(-\frac{1}{x_{\epsilon }^2}-\frac{5}{4}+\frac{95}{48}{x}_{\epsilon }^2-\frac{365}{96}{x}_{\epsilon }^4+O\left(x_{\epsilon }^6\right)\right),$

$2u\left(x_1\right)-u\left(x_3\right)-u\left(x_2\right)=-2\frac{u^{″}\left(x_1\right)}{2!}{h}^2-2\frac{u^{\left(4\right)}\left(x_1\right)}{4!}{h}^4+O\left(h^6\right),$

$s^{″}\left(x_1\right)=u^{″}\left(x_1\right)+\left(\frac{5u^{″}\left(x_1\right)}{4}{\epsilon }^2+\frac{u^{\left(4\right)}\left(x_1\right)}{12}\right)h^2+O\left(h^4\right).$ (2.8)

由式(2.8)可以看出，当 ${\epsilon }^{\text{2}}=-\frac{\text{1}}{\text{15}}\cdot u^{\left(4\right)}\left(x_1\right)/u^{″}\left(x_1\right)$ 时(称为 $\epsilon$ 的最佳值)，有 $s^{″}\left(x_1\right)=u^{″}\left(x_1\right)+O\left(h^4\right)$ 。下一节，我们利用(2.8)构造求解带小粘性系数的非齐两点边值问题的具有四阶逼近精度的差分格式。

3. 带小粘性系数的非齐次两点边值问题的求解

带小粘性系数的非齐次两点边值问题：

$\{\begin{array}{l}v{u}^{″}\left(x\right)+p\left(x\right)u^{\prime }\left(x\right)+q\left(x\right)u\left(x\right)=f\left(x\right),a<x<b,\\ u\left(a\right)=\alpha ,u\left(b\right)=\beta \text{,}\end{array}$ (3.1)

其中 $v>0$ 为粘性常系数，函数 $p\left(x\right),q\left(x\right),f\left(x\right)\in C^2\left[a,b\right]$ ，且 $q\left(x\right)\ge 0$ 。

3.1. 求解两点边值问题的RBF-FD格式

区间 $a\le x\le b$ 被剖成 $n$ 等份，节点 $x_i=a+\left(i-1\right)h,h=\left(b-a\right)/n,1\le i\le n+1$ 。在节点 $x_i$ 处，一方面，我们利用问题(3.1)的解 $u\left(x_i\right)$ 有

$s^{″}\left(x_i\right)=\frac{{\epsilon }^2\left(a^3-1\right)}{a^3\left(b-4a+3\right)}\left(2u\left(x_i\right)-u\left(x_{i-1}\right)-u\left(x_{i+1}\right)\right),$ (3.2)

其中 $a=\sqrt{1+{\epsilon }^2{h}^2},b=\sqrt{1+4{\epsilon }^2{h}^2}$ 。由2.2节的分析可知，当 ${\epsilon }^{\text{2}}=-\frac{\text{1}}{\text{15}}\cdot u^{\left(4\right)}\left(x_i\right)/u^{″}\left(x_i\right)$ (最佳参数)时， $s^{″}\left(x_i\right)=u^{″}\left(x_i\right)+O\left(h^4\right)$ 。这表明为了使 $s^{″}\left(x_i\right)$ 为 $u^{″}\left(x_i\right)$ 的四阶逼近，形状参数 ${\epsilon }^2$ 必须取与 $u^{″}\left(x_i\right),u^{\left(4\right)}\left(x_i\right)$ 有关的值，从而 $a,b$ 也与 $u^{″}\left(x_i\right),u^{\left(4\right)}\left(x_i\right)$ 有关，不妨将这类系数 $\epsilon ,a$ 和 $b$ 分别记为 ${\epsilon }_i,a_i$ 和 $b_i$ ，同时引入记号 ${\gamma }_{\text{i}}=\left({\epsilon }_i^2\left(a_i^3-1\right)\right)/\left(a_i^3\left(b_i-4a_i+3\right)\right)$ 。

另一方面，欲根据中心差分公式

$u^{\prime }\left(x_i\right)=\frac{u\left(x_{i+1}\right)-u\left(x_{i-1}\right)}{2h}-\frac{h^2}{6}{u}^{‴}\left(x_i\right)+O\left(h^4\right),$ (3.3)

给出一个逼近 $u^{\prime }\left(x_i\right)$ 达到 $O\left(h^4\right)$ 精度的近似式。为了达到这个目的，需要将(3.3)右端的三阶导数 $u^{‴}\left(x_i\right)$ 用具有二阶精度的差分逼近公式替代。做到这一点并不难，由两点边值问题(3.1)中的微分方程，可知

$u^{″}\left(x\right)=\frac{1}{v}\left(f\left(x\right)-p\left(x\right)u^{\prime }\left(x\right)-q\left(x\right)u\left(x\right)\right),$ (3.4)

从而

$\begin{array}{c}{u}^{‴}\left(x\right)=\frac{1}{v}{\left(f\left(x\right)-p\left(x\right)u^{\prime }\left(x\right)-q\left(x\right)u\left(x\right)\right)}^{\prime }\\ =\frac{1}{v}\left(f^{\prime }\left(x\right)-p\left(x\right)u^{″}\left(x\right)-\left(p^{\prime }\left(x\right)+q\left(x\right)\right)u^{\prime }\left(x\right)-q^{\prime }\left(x\right)u\left(x\right)\right).\end{array}$ (3.5)

将式(3.5)代入(3.3)右边，并利用中心差分公式，易得到

$\begin{array}{c}{u}^{\prime }\left(x_i\right)=\frac{u\left(x_{i+1}\right)-u\left(x_{i-1}\right)}{2h}-\frac{h^2}{6v}(f^{\prime }\left(x_i\right)-p^{\prime }\left(x_i\right)\frac{u\left(x_{i+1}\right)-u\left(x_{i-1}\right)}{2h}\\ -p\left(x_i\right)\frac{u\left(x_{i+1}\right)-2u\left(x_i\right)+u\left(x_{i-1}\right)}{h^2}-q^{\prime }\left(x_i\right)u\left(x_i\right)\\ -q\left(x_i\right)\frac{u\left(x_{i+1}\right)-u\left(x_{i-1}\right)}{2h}+O\left(h^2\right))+O\left(h^4\right).\end{array}$ (3.6)

利用式(3.2)，(3.6)得到求解两点边值问题(3.1)的具有四阶精度的RBF-FD方程：

$\begin{array}{l}v{\gamma }_i\left(2u_i-u_{i+1}-u_{i-1}\right)+p_i(\frac{u_{i+1}-u_{i-1}}{2h}-\frac{h^2}{6v}({f^{\prime }}_i-{p^{\prime }}_i\frac{u_{i+1}-u_{i-1}}{2h}\\ -p_i\frac{u_{i+1}-2u_i+u_{i-1}}{h^2}-{q^{\prime }}_i{u}_i-q_i\frac{u_{i+1}-u_{i-1}}{2h}))+q_i{u}_i=f_i,2\le i\le n\end{array}$

且满足边值条件 $u_1=\alpha ,u_{n+1}=\beta$ 。进一步整理，得

$\{\begin{array}{l}{\tilde{a}}_i{u}_{i-1}+{\tilde{b}}_i{u}_i+{\tilde{c}}_i{u}_{i+1}={\tilde{d}}_i,2\le i\le n,\\ u_1=\alpha ,u_{n+1}=\beta ,\end{array}$ (3.7)

其中

$\begin{array}{l}{\tilde{a}}_i=-v{\gamma }_i-\frac{p_i}{2h}+\frac{p_i^2}{6v}-\frac{p_{i}h}{12v}\left({p^{\prime }}_i+q_i\right),\\ {\tilde{b}}_i=2v{\gamma }_i-\frac{p_i^2}{3v}+\frac{p_i{q^{\prime }}_i{h}^2}{6v}+q_i,\\ {\tilde{c}}_i=-v{\gamma }_i+\frac{p_i}{2h}+\frac{p_i^2}{6v}+\frac{p_{i}h}{12v}\left({p^{\prime }}_i+q_i\right),\\ {\tilde{d}}_i=f_i+\frac{p_i{h}^2}{6v}{f^{\prime }}_i.\end{array}$

格式(3.7)的矩阵向量形式记为 $AU=F$ ，这里

$A=\left(\begin{array}{cccccccc}{\tilde{b}}_2& {\tilde{c}}_2& 0& 0& \cdots & 0& 0& 0\\ {\tilde{a}}_3& {\tilde{b}}_3& {\tilde{c}}_3& 0& \cdots & 0& 0& 0\\ ⋮& ⋮& ⋮& ⋮& \ddots & ⋮& ⋮& ⋮\\ 0& 0& 0& 0& \cdots & {\tilde{a}}_{n-1}& {\tilde{b}}_{n-1}& {\tilde{c}}_{n-1}\\ 0& 0& 0& 0& \cdots & 0& {\tilde{a}}_n& {\tilde{b}}_n\end{array}\right),U=\left(\begin{array}{c}{u}_2\\ u_3\\ ⋮\\ \begin{array}{c}{u}_{n-1}\\ u_n\end{array}\end{array}\right),F=\left(\begin{array}{c}{\tilde{d}}_2-\alpha {\tilde{a}}_2\\ {\tilde{d}}_3\\ ⋮\\ \begin{array}{c}{\tilde{d}}_{n-1}\\ {\tilde{d}}_n-\beta {\tilde{c}}_n\end{array}\end{array}\right).$

明显地， $A,F$ 都是 $U$ 的函数，即方程组 $AU=F$ 是非线性的。我们构造如下迭代格式

$U^{k+1}=A^{-1}\left(U^k\right)F\left(U^k\right),k=0,1,2,\cdots .$ (3.8)

进行求解，直到满足迭代终止条件 ${\Vert U^{k+!}-U^k\Vert }_2/{\Vert U^k\Vert }_2\le {10}^{-5}$ 时，迭代终止。

3.2. 最佳参数e²的计算

根据3.1节的分析可知最佳参数 ${\epsilon }^{\text{2}}=-\frac{\text{1}}{\text{15}}\cdot u^{\left(4\right)}\left(x_i\right)/u^{″}\left(x_i\right)$ ，而 $u^{″}\left(x_i\right),u^{\left(4\right)}\left(x_i\right)$ 是未知的，如果能给出最佳参数 ${\epsilon }_i^2$ 的一个具有二阶精度的近似值，那么仍有 $\left|s^{″}\left(x_i\right)-u^{″}\left(x_i\right)\right|=O\left(h^4\right)$ 成立。如何给出 ${\epsilon }_i^2$ 的一个具有二阶精度的近似值呢？

首先，对方程(3.4)两边连续求导两次，得

$\begin{array}{l}{u}^{\left(4\right)}\left(x_i\right)=\frac{1}{v}(f^{″}\left(x_i\right)-p\left(x_i\right)u^{‴}\left(x_i\right)-\left(2p^{\prime }\left(x_i\right)+q\left(x_i\right)\right)u^{″}\left(x_i\right)\\ -\left(p^{″}\left(x_i\right)+2q^{\prime }\left(x_i\right)\right)u^{\prime }\left(x_i\right)-q^{″}\left(x_i\right)u\left(x_i\right)).\end{array}$ (3.9)

将式(3.4)，(3.5)代入(3.9)得

$\begin{array}{l}{u}^{\left(4\right)}\left(x_i\right)=\frac{1}{v}(\left(-p^{″}\left(x_i\right)+\frac{3p\left(x_i\right)p^{\prime }\left(x_i\right)}{v}+\frac{2p\left(x_i\right)q\left(x_i\right)}{v}-\frac{p^3\left(x_i\right)}{v^2}-2q^{\prime }\left(x_i\right)\right)u^{\prime }\left(x_i\right)\\ +\left(\frac{2p^{\prime }\left(x_i\right)q\left(x_i\right)}{v}+\frac{p\left(x_i\right)q^{\prime }\left(x_i\right)}{v}-\frac{p^2\left(x_i\right)q\left(x_i\right)}{v^2}+\frac{q^2\left(x_i\right)}{v}-q^{″}\left(x_i\right)\right)u\left(x_i\right)\\ +\left(-\frac{2p\left(x_i\right)}{v}+\frac{p^2\left(x_i\right)}{v^2}-\frac{q\left(x_i\right)}{v}\right)f\left(x_i\right)+f^{″}\left(x_i\right)-f^{\prime }\left(x_i\right)\frac{p\left(x_i\right)}{v}).\end{array}$

不妨记作 $u^{\left(4\right)}\left(x_i\right)=\frac{1}{v}\left(g_1\left(x_i\right)u^{\prime }\left(x_i\right)+g_2\left(x_i\right)u\left(x_i\right)+g_3\left(x_i\right)f\left(x_i\right)+g_4\left(x_i\right)\right)$ ，其中

$\begin{array}{l}{g}_1\left(x_i\right)=-p^{″}\left(x_i\right)+\frac{3p\left(x_i\right)p^{\prime }\left(x_i\right)}{v}+\frac{2p\left(x_i\right)q\left(x_i\right)}{v}-\frac{p^3\left(x_i\right)}{v^2}-2q^{\prime }\left(x_i\right)，\\ g_2\left(x_i\right)=\frac{2p^{\prime }\left(x_i\right)q\left(x_i\right)}{v}+\frac{p\left(x_i\right)q^{\prime }\left(x_i\right)}{v}-\frac{p^2\left(x_i\right)q\left(x_i\right)}{v^2}+\frac{q^2\left(x_i\right)}{v}-q^{″}\left(x_i\right)，\\ g_3\left(x_i\right)=-\frac{2p\left(x_i\right)}{v}+\frac{p^2\left(x_i\right)}{v^2}-\frac{q\left(x_i\right)}{v},\\ g_4\left(x_i\right)=f^{″}\left(x_i\right)-f^{\prime }\left(x_i\right)\frac{p\left(x_i\right)}{v}.\end{array}$

根据RBF-FD公式(2.8)的截断误差，可知

$s^{″}\left(x_i\right)-u^{″}\left(x_i\right)=\left(\frac{5u^{″}\left(x_i\right)}{4}{\epsilon }^2+\frac{u^{\left(4\right)}\left(x_i\right)}{12}\right)h^2+O\left(h^4\right).$

从而，当

$\begin{array}{c}{\epsilon }_i^2=-\frac{1}{15}\frac{u^{\left(4\right)}\left(x_i\right)}{u^{″}\left(x_i\right)}\\ \approx -\frac{1}{15}\frac{\left(u\left(x_{i+1}\right)-u\left(x_{i-1}\right)\right)g_1\left(x_i\right)+2h\left(u\left(x_i\right)g_2\left(x_i\right)+f\left(x_i\right)g_3\left(x_i\right)+g_4\left(x_i\right)\right)}{2h\left(f\left(x_i\right)-q\left(x_i\right)u\left(x_i\right)\right)-p\left(x_i\right)\left(u\left(x_{i+1}\right)-u\left(x_{i-1}\right)\right)}.\end{array}$

时，有 $s^{″}\left(x_i\right)=u^{″}\left(x_i\right)+O\left(h^4\right)$ 成立。

Algorithm I：差分方程(3.7)的迭代求解过程：

STEP 1：迭代格式(3.8)的初始值选择：利用二阶中心差分格式

$\{\begin{array}{l}v\cdot \frac{u_{i+1}-2u_i+u_{i-1}}{h^2}+p_i\cdot \frac{u_{i+1}-u_{i-1}}{2h}+q_i{u}_i=f_i,i=2,\cdots ,n,\\ u_1=\alpha ,u_{n+1}=\beta \end{array}$ (3.10)

计算出一组具有 $O\left(h^2\right)$ 精度的近似解 ${\left\{u_i\right\}}_{i=1}^{n+1}$ 。

STEP 2：以 ${\left\{u_i\right\}}_{i=1}^{n+1}$ 作为格式(3.8)的初始值 $U^0$ 进行迭代计算，直到满足终止条件 ${\Vert U^{k+!}-U^k\Vert }_2/{\Vert U^k\Vert }_2\le {10}^{-5}$ 时，迭代终止。

3.3. 数值实验

例1：利用Algorithm I求解下列两点边值问题：

$\{\begin{array}{l}v{u}^{″}-u^{\prime }=-v{\text{π}}^2\sin \left(\text{π}x\right)-\cos \left(\text{π}x\right),0<x<1\\ u\left(0\right)=0,u\left(1\right)=0.\end{array}$

其精确解为 $u\left(x\right)=\sin \left(\text{π}x\right)$ 。

下表分别给出了，当小粘性系数 $v$ 取不同值时，本文构造的RBF-FD格式(3.7)的解的迭代次数(Iter)、误差( $l_1$ -error)、收敛阶(order)和计算时间(Time(s))，其中

$l_1\text{-error}=\frac{1}{n+1}{\displaystyle \underset{i=1}{\overset{n+1}{\sum }}\left|u_i-u\left(x_i\right)\right|}.$

同时，表2中还列出了二阶中心差分格式(3.10)的解的误差、收敛阶及计算时间。

表2~表5分别给出了 $v=1,0.1,0.01,0.001$ 时，本文构造的RBF-FD格式(3.7)在不同步长下解的误差

表2. v = 1时RBF-FD格式(3.7)与中心差分格式(3.10)在不同步长下解的误差、收敛阶及计算时间

表3. v = 0.1时RBF-FD格式(3.7)与中心差分格式(3.10)在不同步长下解的误差、收敛阶及计算时间

表4. v = 0.01时RBF-FD格式(3.7)与中心差分格式(3.10)在不同步长下解的误差、收敛阶及计算时间

表5. v = 0.001时RBF-FD格式(3.7)与中心差分格式(3.10)在不同步长下解的误差、收敛阶及计算时间

( $l_1$ -error)及收敛阶(order)。由数值结果可以看出，对于例1中的问题，本文构造的格式在 $v=1~0.001$ 时能够达到四阶收敛速度，与传统的基于三节点的中心差分格式(3.10)相比，收敛阶提高了2阶，而计算时间略有增加。但是在 $v$ 小于0.001的粗网格情况下，RBF-FD格式的精度与收敛阶会有所下降。导致这一

问题的本质原因是，因为格式(3.7)的截断误差中含有 $\frac{1}{v}$ 和 $\frac{1}{v^2}$ 项，从而在步长 $h$ 较大时，截断误差会随着粘性系数 $v$ 的减小而增大，这就导致了当 $v$ 取值较小，而网格较粗时，格式(3.7)的精度与收敛阶也会随之下降。

4. 结论

本文利用基于三等距节点的具有零次代数精度的径向基函数插值的Lagrange形式，给出了逼近被插函数在中心节点处一阶导数和二阶导数的RBF-FD公式，通过对其截断误差进行分析，给出了使逼近阶达到最高阶的最佳参数 $\epsilon$ 值。然后，利用这些具有最佳参数值的RBF-FD公式，给出了求解带小粘性系数的非齐两点边值问题的RBF-FD差分格式，所构造的RBF-FD格式的收敛阶是同节点模板的二阶多项式中心差分格式的2倍，而计算时间略有增加，显示出带有最佳参数的RBF-FD格式的优势。此外，本文构造的RBF-FD格式，保证了粘性系数 $v\ge 0.001$ 时，数值解能够达到4阶收敛速度。未来的工作是要解决如何修正现有的差分格式，使得当粘性系数 $v$ 取值更小时，问题的数值解仍能够达到高阶的逼近精度。

基金项目

国家自然科学基金(91630203, 11271041)和民机专项(MJ-F-2012-04)资助项目。

文章引用

段俊娜,郭子滔,冯仁忠. 基于径向基有限差分法求解带小粘性系数非齐次两点边值问题
A Finite Difference Method Based on Radial Basis Function for Non-Homogeneous Two-Point Boundary Value Problems with Small Viscosity Coefficient[J]. 应用数学进展, 2018, 07(01): 20-29. http://dx.doi.org/10.12677/AAM.2018.71004

参考文献 (References)