本文考虑如下形式的第一类Volterra积分方程

其中核函数$k(x, s), (x, s)\in \Delta: = \{(x, s)| -1 < s < x < 1\}$,以及源函数$g(x), x\in I: = [-1, 1]$是给定的,而$u(s)$是一个要求解的未知函数.

令

通过对核函数$k(x, s)$在$s > x$上补充定义(即$k(x, s) = 0, s > x$),我们可以将核函数$k(x, s)$延拓为$[-1, 1]\times [-1, 1]$区域上的函数.则方程(1.1)转化为第一类Fredholm方程

众所周知,具有弱奇性核或连续核的第一类积分算子是紧的.因此,方程(1.1)是不适定的,即:解(如果存在)不连续依赖于数据$g(x)$.但是,因为在使用数值积分公式处理积分项时,不可避免的要产生离散误差.因此,不能直接简单通过离散化直接处理第一类积分方程.随着离散的加细,数值解变得"更糟"不收敛于精确解.在实际问题中,右侧$g(x)$是不能准确知道的,只知道带有一定误差的测量数据,即假设已知$g^\delta$为

其中$\delta > 0$是已知的误差水平.我们的目标是解扰动的方程

在这种情况下,必须考虑一种正则化策略来近似求解.至今已有许多求解Volterra型积分方程的数值方法,见文献[1-4, 7-8, 10-12, 15-16]及其参考文献.然而这些文献中很少有工作是求解带扰动数据的第一类Volterra型积分方程.本文的主要目的是通过同时使用Legendre配置方法和正则化策略来求解方程(1.5).我们对所提出的方法给出了严格的收敛性分析,通过一些数值实验验证了所提方法的有效性.

本文结构如下:第2节中我们介绍第一类Volterra积分方程的谱方法;第3节中提供了收敛性分析中所需要的引理;第4节中分别给出了$L^{\infty}$和$L^{2}$空间中的收敛性分析;在第5节中进行了数值实验;最后,在第6节中给出了简要结论.

如引言中所示,方程(1.1)是不适定的.在这种情况下,解决问题的自然想法是考虑以下正则化策略:

其中$\alpha > 0$是正则化参数.

引理2.1  令$k(x, s)\in C(\bar{\Delta})$,则方程(2.1)存在唯一解$u^{\alpha, \delta}(x)\in C(I)$.

证  首先,众所周知,如果核函数$k(x, s)\in L^2(\Delta)$,那么由(1.2)式定义的算子$K$是从$L^2(I)$到$L^2(I)$的紧算子,并且若$k(x, s)\in C(\bar{\Delta})$,则由(1.2)式定义的算子$K$也是由$C(I)$到$C(I)$的紧算子.下面,我们将证明以下齐次方程只有平凡解$u(x) = 0$,

事实上,通过迭代过程,很容易看出来

其中$\mu = \max_{-1\leq s\leq x\leq1}\mid k(x, s)\mid.$

令$n\to\infty$,我们得到$u(x) = 0, x\in I$.最后,利用Riesz定理,非齐次方程(2.1)有唯一解,且$(\alpha I+K)^{-1}$是有界的.证毕.

方程(2.1)的解$u^{\alpha, \delta}(x)$可以写成

注2.1  求解问题(1.1)的另一种方法是考虑Tikhonov正则化策略,即

其中$K^*$表示$K$的伴随算子.

下面,我们使用Legendre配置方法求解方程(2.1).配置点取成$(N+1)$次Legendre-Gauss或Gauss-Radau或Gauss-Lobatto节点$\{x_{i}\}_{i = 0}^{N}$.假设方程(2.1)在$x_{i}$处成立,即

对于较小的$x_{i}$, $u^{\alpha, \delta}(s)$的可用信息较少.为了克服这个困难,我们将积分区间$[-1, x_{i}]$变换成固定区间$[-1, 1]$,然后利用适当的求积公式处理.准确地说,我们首先做一个简单的线性变换$\tau = s(x_{i}, \theta), \theta\in [-1, 1]$,其中

则(2.4)式变为

使用$(N+1)$个节点的Gauss求积公式可得出

其中$\omega_{j}, 0\leq j\leq N$是Legendre权重.

我们将使用$u_i^{\alpha, \delta}, 0\leq i\leq N$来表示$u^{\alpha, \delta}(s(x_i, \theta_j))$,其中$u_i^{\alpha, \delta} = u^{\alpha, \delta}(x_i)$.为此,我们使用Lagrange插值多项式逼近$u^{\alpha, \delta}(x)$,即

其中$L_p$是$p$次Lagrange基函数.联立(2.7)和(2.8)式得

交换求和顺序得

下面讨论谱配置算法如何实现.记$U_N^{\alpha, \delta} = (u_{0}^{\alpha, \delta}, u_{1}^{\alpha, \delta}, \cdots u_{N}^{\alpha, \delta})^{T}$和$g^\delta_{N} = (g^{\delta}(x_{0})$, $g^\delta(x_{1}),$$\cdots,$$g^\delta(x_{N}))^{T}$,我们可以得到一个矩阵形式的方程

其中矩阵$A$的项由$a_{i, j}$给出:

我们定义第一类Volterra积分方程的正则解$u_N^{\alpha, \delta}(x)$如下:

注2.2  此处我们给出计算$L_j$的一种有效算法.这个想法是用Legendre函数$\{P_p(x)\}_{p = 0}^N$来表示$L_j$:

其中$\beta_{jp}$称为$L_j$的离散多项式系数,它可表示为(见文献[5])

其中$r_p = \sum\limits_{k = 0}^{N}P_{p}^{2}(x_{k})\omega_{k} = (p+\frac{1}{2})^{-1}$, $p < N$;对于Gauss和Gauss-Radau公式$r_N = (N+1/2)^{-1}$;对于Gauss-Lobatto公式$r_N = 2/N$.从(2.9)和(2.10)式可以得出

结合$L_p(s)$递推公式可以有效计算$L_j$,进而

注2.3  对于另一个正则化策略(2.3),我们可以使用类似的过程获得方程

其中, $A^*$表示$A$的伴随.

本节给出一些引理,这些引理将应用于以后的收敛性分析.

引理3.1^[5]  设$u\in H^{m}(I)$, $I: = [-1, 1]$, $m\geq 1$, $\phi \in{\mathcal P}_{N}$.使用带有Legendre权的$(N+1)$点Gauss,或Gauss-Radau,或Gauss-Lobatto求积公式计算$u\phi$的积分,则存在一个与$N$无关的常数$C$,满足

其中

引理3.2^[5]  设$u\in H^{m}(I)$, $I_{N}u$表示$u$的$(N+1)$点Gauss,或Gauss-Radau,或Gauss-Lobatto点插值多项式,即

则以下估计式成立

引理3.3 (Gronwall)不等式  如果可积函数$u(x)$满足

其中$f(x)$是可微函数, $\lambda > 0$是一个给定的常数.那么

并且

其中$\mu = \max_{x, s\in I}\mid k(x, s)\mid$并且${\cal C}_{\lambda} = \lambda(2-2\lambda\mu-\frac{1}{2}+\frac{e^{4\lambda\mu}}{2})^{\frac{1}{2}}.$

证  直接计算可得

令

我们可得

即

在区间$[-1, x]$上,对上面不等式求积分可得

进而

于是

即

简单的计算可得

因此

引理3.3证毕.

在本节中,我们将在$L^{2}$和$L^{\infty}$空间中分别进行收敛性分析.为此,我们假设无扰动方程(1.1)存在唯一解.通过应用三角形不等式,精确解$u(x)$与正则解$u_{N}^{\alpha, \delta}(x)$之间的误差分成两部分

接下来,我们将分别估计$\|u(x)-u^{\alpha, \delta}(x)\|$和$\|u^{\alpha, \delta}(x)-u_{N}^{\alpha, \delta}(x)\|$.

对任意常数$\alpha > 0$,令$u_{i}^{\alpha, \delta}$为方程(2.10)的解,即

按照(3.2)式的表示,令

那么(4.2)式可以写成

于是有

其中

根据引理3.1可得

从(4.5), (2.7)和(2.9)式可得

在(4.5)式两边同时乘以$L_{i}(x)$并从$0$到$N$求和得到

其中$I_N$是引理3.2中定义的插值算子, $e^{\alpha, \delta}$表示误差函数,即

由(4.7)和(2.7)式得出

其中

应用Gronwall不等式得出

其中

下面我们分别估计(4.9)式的最后三项.

根据引理3.2得

应用以上估计及(4.9)式可得

令$N$充分大,我们有

定理4.1  令$u(x)$是Volterra方程(1.1)的精确解,并假设

其中$u_{p}^{\alpha}$由(2.9)式给出,而$L_{p}(x)$是与Gauss点$\{x_{j}\}_{j = 0}^{N}$相关的$p$次Lagrange基函数.如果$u\in H^{m}(I)$,则对于$m\geq 1$, $N$充分大时有

其中$s(x_{i}, \theta)$由(2.5)式定义,并且$C$是一个与$N$无关的常数.

注4.1  按照与定理4.1相同的步骤,我们还可以得到

令$(\mu_j, \phi_j, \psi_j)^{+\infty}_{j = 1}$为算子$K$的奇异系统.由Picard定理方程(1.1)的解$u(x)$可写为如下收敛级数的形式

利用$K$的奇异系统$(\mu_j, \phi_j, \psi_j)$的性质,方程(2.1)的解$u^{\alpha, \delta}(x)$具有如下形式

利用三角不等式简单计算得

显然,当$\alpha\to 0$时第二项趋于零,而$\|R_\alpha\|$趋向于$\infty$,因此,我们需要选择依赖于$\delta$的$\alpha = \alpha(\delta)$以使总误差尽可能小,即

$\alpha = \alpha(\delta)$的选取要根据核函数$k(x, s)$的具体表达式去确定.对于取定的$\alpha = \alpha(\delta)$,从误差表达式(4.14)可知,当$N$趋于无穷大时$\parallel u^{\alpha, \delta}(x)-u_{N}^{\alpha, \delta}(x)\parallel_{L^{2}(I)}$趋于零.因此,当$\delta$趋于零时,精确$u(x)$与正则解$u_{N}^{\alpha, \delta}(x)$之间的总误差趋于零,即我们的方法是收敛的.

注4.2  对于另一种正则化策略(2.3),即$u^{\alpha, \delta}(x) = (\alpha I+K^*K)^{-1}K^*g^\delta(x): = R_\alpha g^\delta$, $\|u(x)-u^{\alpha, \delta}(x)\|$的估计如下

本节的目的是给出一些简单的示例说明所提方法的有效性.数值实验是使用Matlab语言实现的.不失一般性,我们只使用Legendre-Gauss节点(即$L_{N+1}(x)$的零点作为配置点).对于Legendre-Gauss点,相应的权重为$\omega_{j} = \frac{2}{(1-x_{j}^{2})[p'_{N+1}(x_{j})]^{2}}, \; 0\leq j\leq N$.离散噪声数据$g^{\delta}$是通过在精确数据$g$中加入随机噪声而得到的,即

$\delta$是在均方根意义下的误差水平,

在这里,函数${\rm randn}(size(g))$返回一个随机数组(其元素是由均值$0$,方差$\sigma^{2} = 1$的正态分布生成),其大小与$g$相同.

RMSE表示

其中$N$是测试点的总数.

例5.1  第一个例子是一个数值微分问题.更准确地说,在方程(2.1)中核函数$k(x, s)$和$g(x)$取为

相应的精确解由$u(x) = 3x^{2}$给出.

从表 1可以看出,如果数据无扰动,通过使用谱方法求解第一类Volterra方程可以获得谱精度.我们还可以从表 1中看到,误差不随$N$的增加而减小,因为第一类Volterra积分算子是紧的,即,我们解决的问题是不适定的.

从表 2可以看出,如果数据扰动时,使用谱方法(未使用正则化时)计算误差显著增加了.

表 3给出的是误差级数$\epsilon$给定条件下,取不同正则化参数时的正则化解的误差.由表 3可以看出正则化参数要适当的选取,既不能太大,也不能太小,要与$N$有关.

例5.2  第二个例子中,考虑方程(2.1)中

$u(x) = e^{4x}$是相应的精确解.表 4-5显示了具有不同参数的正则化解的误差.

本文利用Legendre配点法和正则化策略,给出了第一类Volterra积分方程的数值解法.这项工作最重要的贡献在于,我们对所提出的方法进行了严格的收敛分析.通过数值试验验证了该方法的有效性.然而,本文尚未给出收敛速度.这一点将在我们未来的工作中加以考虑.