我们考虑如下形式的线性方程

$ \begin{equation} Kx = y, \end{equation} $

(1.1)

其中$ K $是希尔伯特空间$ X $和$ Y $之间的线性紧算子, 且$ R(K) $是$ Y $的无限维子空间, 这是非常典型的线性不适定反问题^[1]. 为方便起见, 设紧算子$ K $是单射的, 那么对于线性算子$ K $的值域$ R\left(K\right) $范围内的每一个$ y $, 算子方程都有唯一解$ x^{\dagger}{\in}X $.

在实际问题中, 通常无法得到方程$ (1.1) $精确的右端项$ y=Kx^\dagger $, 只可观测到满足$ \lVert {y^\delta-y}\rVert \leq\delta $的扰动数据$ y^\delta\in Y $, 其中$ \delta>0 $为误差水平, $ y^\delta $为扰动右端. 相应的方程$ (1.1) $变为

$ \begin{equation} Kx^\delta=y^\delta. \end{equation} $

(1.2)

我们希望找到一个$ x^\dagger $的近似解$ x^{\alpha, \delta}\in X $, 也就是说构造无界逆算子$ K^{-1}:R\left(K\right)\to X $的一族合适的近似算子$ R_{\alpha}:Y\to X $, 将$ x^{\alpha, \delta}:=R_{\alpha}y^{\delta} $视为$ x^{\dagger} $的一种近似, 我们通过正则化方法实现该目标.

对于线性不适定问题, 正则化方法主要分为两类: 变分正则化方法和迭代正则化方法. Tikhonov迭代正则化方法^[2-3]是经典的变分正则化方法, 常见的迭代正则化方法有Landweber迭代正则化方法^[3-5], 隐式迭代正则化方法^[6], 非稳态Tikhonov迭代正则化方法^[7]. 有另外一类正则化方法: 连续正则化方法, 又称渐进正则化方法, Showalter方法^{[1, 8]}, 该方法考虑一阶演化方程

$ \begin{equation} \dot{x}\left(t\right)+K^*Kx\left(t\right)=K^*y, \; \; \; x\left(0\right)=0, \end{equation} $

(1.3)

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

许多迭代正则化方法可以看作对(1.3)进行不同格式数值离散得到的. Landweber迭代正则化方法可以看作对其进行向前Euler离散^[9]; 隐式迭代正则化方法是对其进行向后Euler离散^[9]得到的; 非稳态Tikhonov正则化方法是变步长的隐式迭代正则化方法. 论文[10]证明了对(1.3)使用Runge-Kutta离散得到的迭代方法是正则化方法同时可以达到最优收敛阶.

对于上述迭代正则化方法, 为保证正则化性质, 必须选择合适的有限终止步数$ m=m(\delta) $（先验选择）或$ m=m(\delta, y^\delta) $（后验选择）. 已知在光滑性假设$ x^\dagger=(K^*K)^{v/2}z, \lVert z\rVert\leq E $下, 若迭代步数按照先验方式选取$ m=m(\delta)=O(\delta^{-\frac{2}{v+1}}) $, 则对所有的$ v>0 $^[7], 这些方法都能达到最优收敛阶

$ \begin{equation*} \lVert x^{m\left(\delta\right), \delta}-x^\dagger\rVert \leq cE^\frac{1}{v+1}\delta^\frac{v}{v+1}. \end{equation*} $

若终止迭代步数$ m=m(\delta, y^{\delta}) $是依据Morozov残差准则^[11-12]选取的, 则当$ \delta \to 0 $时, 有如下关系式

$ \begin{equation*} \lVert x^{m\left(\delta\right), \delta}-x^\dagger\rVert \leq cE^\frac{1}{v+1}\delta^\frac{v}{v+1}, \; \; m=O(\delta^{-\frac{2}{v+1}}). \end{equation*} $

在常微分方程的数值离散方法中, 最常见的几类方法分别是向前Euler法、向后Euler法、Crank-Nicolson方法. 一个自然的问题是如果我们尝试对常微分方程(1.3)使用Crank-Nicolson方法离散, 会不会导出一个迭代正则化方法, 该方法能否达到最优收敛阶, 以及变步长的时候是否有可能使得迭代步数在一定条件下类似变步长的向后Euler方法那样达到$ O(|log\delta|) $？本文的主要贡献如下: 首先我们将证明对(1.3)进行Crank-Nicolson离散得到的迭代方法是一种迭代正则化方法, 并且该方法在常数步长下可以达到最优收敛阶; 其次我们发现在非常数步长下, 因为高频扰动的原因导致其不能够像非稳态Tikhonov迭代正则化方法达到$ O(|log\delta|) $的迭代步数, 但是在较大噪声设置下, 其计算效率和非稳态Tikhonov迭代正则化相似.

下面简单描述下本文的结构. 第二节中给出所需要的预备知识, 第三节对(1.3)采用Crank-Nicolson方法离散, 证明该方法是正则化方法, 给出收敛阶并与说明与非稳态Tikhonov迭代正则化方法相比的区别, 第四节用两个数值算例验证理论, 第五节做出总结.

定义2.1^[3]  正则化策略是指一族线性有界算子

$ \begin{equation*} R_\alpha:Y\to X, \; \; \alpha>0, \end{equation*} $

满足对所有的$ x\in X $

$ \begin{equation*} \lim\limits_{\alpha\to 0}R_\alpha Kx=x. \end{equation*} $

即算子$ R_\alpha K $逐点收敛到恒等算子.

由该定义和$ K $的紧性, 可以得到下面的定理

定理2.1^[3]  设$ R_{\alpha} $是紧算子$ K:X\to Y(dimX=\infty) $的正则化策略, 则

1. 算子$ R_{\alpha} $非一致有界: 存在序列$ \alpha_j $, 当$ j\to\infty $时, $ \lVert R_{\alpha_j}\rVert \to\infty $.

2. 序列（$ R_{\alpha}Kx $）在$ X $的有界子集上不一致收敛; 即$ R_\alpha K $不会在算子范数下收敛到恒等算子$ I $.

当用$ x^{\alpha, \delta} $做为$ Kx=y $的解$ x^\dagger $的近似后, 将误差分为两部分:

$ \begin{equation*} \begin{split} \lVert x^{\alpha, \delta}-x^\dagger\rVert &\leq\lVert R_{\alpha}y^\delta-R_{\alpha} y\rVert +\lVert R_\alpha y-x^\dagger\rVert \\ &\leq\lVert R_{\alpha}\rVert \lVert y^\delta-y\rVert +\lVert R_\alpha Kx^\dagger-x^\dagger\rVert, \end{split} \end{equation*} $

因此

$ \begin{equation} \lVert x^{\alpha, \delta}-x^\dagger\rVert \leq\delta\lVert R_\alpha\rVert +\lVert R_\alpha Kx^\dagger-x^\dagger\rVert. \end{equation} $

(2.1)

由定理2.1知, 条件数$ \lVert R_\alpha\rVert $随着$ \alpha $趋于0而趋于无穷大, 根据正则化策略的定义, 第二项随着$ \alpha $趋于0而趋于0. 我们需要选择依赖$ \delta $的参数$ \alpha=\alpha(\delta) $使总误差最小, 即最小化

$ \begin{equation*} \delta\lVert R_{\alpha}\rVert+\lVert R_\alpha Kx^\dagger-x^\dagger\rVert. \end{equation*} $

定义2.2^[3]  若正则化策略$ R_{\alpha} $的参数选择$ \alpha=\alpha\left(\delta\right) $满足$ \lim_{\delta\to 0}\alpha(\delta)=0 $, 且对所有的$ x\in X $有

$ \begin{equation*} \sup\left\{\lVert R_{\alpha\left(\delta\right)}y^\delta-x\rVert :y^\delta\in Y, \lVert Kx-y^\delta\rVert \leq\delta\right\}\to0, \; \; \delta\to0, \end{equation*} $

则称该参数选择是可以接受的.

由估计(2.1)可知: 若$ \alpha(\delta)\to0 $, $ \delta\lVert R_\alpha\rVert \to0(\delta\to0) $且$ \lVert R_\alpha Kx^\dagger-x^\dagger\rVert \to 0(\alpha\to0) $, 则参数选择$ \alpha=\alpha(\delta) $是可接受的. 在迭代正则化方法中, 正则化参数$ \alpha=\frac{1}{m} $.

为了方便构造一族可达的正则化策略, 引入奇异值系统. 令$ K:X\to Y $是一个线性紧算子, 并且令$ \left(\mu_j, x_j, y_j\right) $是$ K $的奇异值系统. 该级数收敛时, 方程$ Kx = y $的解$ x $由Picard定理^[3]给出, 即

$ \begin{equation*} x=\sum\limits_{j=1}^\infty\frac{1}{\mu_j}\left(y_j, y\right)x_j. \end{equation*} $

这个结果说明了$ y $的误差的影响, 小的奇异值会将误差放大, 我们通过阻尼化因子$ 1/\mu_j $来构建正则化策略.

定理2.2  设$ K:X\to Y $是紧算子, 且$ \left(\mu_j, x_j, y_j\right) $为其奇异值系统, 若

$ \begin{equation*} q:\left(0, \infty \right)\times \left(0, \lVert K\rVert \right]\to\mathbb{R} \end{equation*} $

满足以下性质:

1. 对所有的$ m>0 $和满足$ 0<\mu\leq\lVert K\rVert $的$ \mu $, 均有$ |q\left(m, \mu\right)|\leq1 $.

2. 对所有的$ m>0 $, 存在$ c\left(m\right) $, 使得对所有满足$ 0<\mu\leq\lVert K\rVert $的$ \mu $, 都有

$ \begin{equation*} \left|q\left(m, \mu\right)\right|\leq c\left(m\right)\mu. \end{equation*} $

3. 对所有的$ 0<\mu\leq\lVert K\rVert $, 有$ \lim_{m\to\infty}q\left(m, \mu\right)=1 $.

则由

$ \begin{equation} R_{m}y:=\sum\limits_{j=1}^\infty\frac{q\left(m, \mu_j\right)}{\mu_j}\left(y, y_j\right)x_j, \quad y\in Y \end{equation} $

(2.2)

定义的算子$ R_m:Y\to X, m>0 $是一个正则化策略, 且满足$ \lVert R_m\rVert \leq c\left(m\right) $与$ \lVert KR_m\rVert \leq 1 $. 如果参数选择$ m=m(\delta) $满足当$ \delta \to 0 $时, $ m(\delta) \to \infty $且$ \delta c(m(\delta)) \to 0 $, 则该参数选择是可接受的. 函数$ q $被称为算子$ K $的正则化滤子.

下面给出几种迭代正则化方法的正则化策略和正则化滤子:

例1  Landweber迭代正则化方法

$ \begin{equation*} x^{m, \delta} = x^{m-1, \delta}+aK^*\left(y^\delta-Kx^{m-1, \delta}\right), \; \; x^{0}:=0, \; \; a\in(0, 1/\lVert K\rVert^2], \end{equation*} $

$ \begin{equation*} R_m:=a\sum\limits_{k=0}^{m-1}\left(I-aK^*K\right)^kK^*, \; \; m=1, 2, . . ., \end{equation*} $

$ \begin{equation*} q\left(m, \mu\right)=1-\left(1-a\mu^2\right)^m. \end{equation*} $

例2   隐式迭代正则化方法

$ \begin{equation*} x^{m, \delta} = x^{m-1, \delta}+aK^*\left(y^\delta-Kx^{m, \delta}\right), \; \; x^{0}:=0, \end{equation*} $

$ \begin{equation*} R_m:= a\sum\limits_{k=1}^{m}\left(I+aK^*K\right)^{-k}K^*\; \; m=1, 2, . . ., \end{equation*} $

$ \begin{equation*} q\left(m, \mu\right)=1-\frac{1}{\left(1+a\mu^2\right)^m}. \end{equation*} $

例3  中点法离散

$ \begin{equation*} \begin{cases} x^{m, \delta} = x^{m-1, \delta} + a K^* \left(y - K x^{m-\frac{1}{2}, \delta}\right), \\ x^{m-\frac{1}{2}, \delta} = x^{m-1, \delta} + \frac{a}{2} K^* \left(y - K x^{m-1, \delta}\right), \quad x^{0} := 0, \; \; a\in(0, 2/\lVert K\rVert^2], \end{cases} \end{equation*} $

$ \begin{equation*} R_m=\sum\limits_{k=0}^{m-1}\left[I-aK^*K+\frac{a^2}{2}\left(K^*K\right)^2\right]^k\left(a-\frac{a^2}{2}K^*K\right)K^*, \; \; m=1, 2, . . ., \end{equation*} $

$ \begin{equation*} q\left(m, \mu\right)=1-\left(1-a\mu^2+\frac{a^2}{2}\mu^4\right)^m. \end{equation*} $

我们对常微分方程(2.1)使用Crank-Nicolson离散格式^[13]:

$ \begin{equation*} x^{m, \delta}=x^{m-1, \delta}+\frac{a}{2}\left(K^*\left(y^\delta-Kx^{m, \delta}\right)+K^*\left(y^\delta-Kx^{m-1, \delta}\right)\right), \; \; x^0:=0, \; \; m=1, 2. . ., \end{equation*} $

其中$ a $为步长, 整理上式可以得到

$ \begin{equation} x^{m, \delta}=\left(I+\frac{a}{2}K^*K\right)^{-1}\left(I-\frac{a}{2}K^*K\right)x^{m-1, \delta}+a\left(I+\frac{a}{2}K^*K\right)^{-1}K^*y^\delta. \end{equation} $

(3.1)

迭代格式(3.1)有显示格式$ x^{m, \delta}=R_my^\delta $, 其正则化策略$ R_m:Y\to X $由下式定义:

$ \begin{equation*} R_m=a\sum\limits_{k=0}^{m-1}\left(\left(I-\frac{a}{2}K^*K\right)\left(I+\frac{a}{2}K^*K\right)^{-1}\right)^k\left(I+\frac{a}{2}K^*K\right)^{-1}K^*, \; \; m=1, 2, . . .\; . \end{equation*} $

引入奇异值系统$ \left(\mu_j, x_j, y_j\right) $, 则对所有的$ y\in Y $有

$ \begin{equation} \begin{split} R_my&=a\sum\limits_{j=1}^{\infty}\mu_j\sum\limits_{k=0}^{m-1}\left(\frac{1-\frac{a}{2}\mu_j^2}{1+\frac{a}{2}\mu_j^2}\right)^m\left(y, y_j\right)x_j\\ &=\sum\limits_{j=1}^{\infty}\frac{1}{\mu_j}\left[1-\left(\frac{1-\frac{a}{2}\mu_j^2}{1+\frac{a}{2}\mu_j^2}\right)^m\right]\left(y, y_j\right)x_j\\ &=\sum\limits_{j=1}^\infty\frac{q\left(m, \mu_j\right)}{\mu_j}\left(y, y_j\right)x_j, \; \; \; y\in Y, \end{split} \end{equation} $

(3.2)

其中$ q\left(m, \mu\right)=1-\left(\frac{1-\frac{a}{2}\mu^2}{1+\frac{a}{2}\mu^2}\right)^m $.

下面证明当$ a\in(0, 2/\lVert K\rVert^2] $时, $ q(m, \mu) $是正则化滤子.

定理3.1   $ q\left(m, \mu\right)=1-\left(\frac{1-\frac{a}{2}\mu^2}{1+\frac{a}{2}\mu^2}\right)^{m} $满足定理2.2中的性质$ 1, 2, 3 $.

证  性质$ 1, 3 $都是显然的.

性质$ 2 $可以直接用伯努利不等式得到

$ \begin{equation*} 1-\left(\frac{1-\frac{a}{2}\mu^2}{1+\frac{a}{2}\mu^2}\right)^{m}\leq 1-\left(1-\frac{a}{2}\mu^2\right)^{m}\leq 1-\left(1-\frac{am\mu^2}{2}\right)=\frac{am\mu^2}{2}. \end{equation*} $

因此$ \left|q\left(m, \mu\right)\right|\leq\sqrt{\left|q\left(m, \mu\right)\right|}\leq\sqrt{am/2}\mu $.

引理3.1  对所有的$ m\geq1 $和$ 0\leq\mu\leq\sqrt{2/a} $, 有

$ \begin{equation*} m^{v}\mu^{2v}\left(\frac{1-\frac{a}{2}\mu^2}{1+\frac{a}{2}\mu^2}\right)^{2m-2} \leq m^{v}\mu^{2v}\left(1-\frac{a}{2}\mu^2\right)^{2m-2}\\ \leq \left(\frac{2v}{a}\right)^v. \end{equation*} $

证  只需证右边的不等式即可.

$ \begin{equation*} \begin{split} m^{v}\mu^{2v}\left(1-\frac{a}{2}\mu^2\right)^{2m-2} &=\left(\frac{mv}{a\left(m-1\right)}\right)^{v}\left(\frac{m-1}{v}a\mu^2\right)^{v}\left(1-\frac{a}{2}\mu^2\right)^{2m-2}\\ &\leq\left(\frac{mv}{a\left(m-1\right)}\right)^{v}\left(\frac{2m-2}{2m-2+v}\right)^{2m-2+v}\\ &=\left(\frac{v}{a}\right)^{v}\left(\frac{m}{m-1}\right)^v\left(\frac{m-1}{m-1+\frac{v}{2}}\right)^v\left(\frac{m-1}{m-1+\frac{v}{2}}\right)^{2m-2}\\ &\leq\left(\frac{v}{a}\right)^{v}\left(\frac{m}{m-1+\frac{v}{2}}\right)^v\\ &\leq \left(\frac{2v}{a}\right)^v, \end{split} \end{equation*} $

其中第一个小于等于号处使用了均值不等式:

当$ a_1, a_2, ... , a_n>0 $时有

$ \begin{equation*} a_1a_2... a_n\leq\left(\frac{a_1+a_2+...+a_n}{n}\right)^n. \end{equation*} $

在求解$ Kx = y $的每种迭代算法中, 采用Morozov残差准则做为停止条件, 设$ r > 1 $是一个固定数, 当首次出现自然数$ m \in \mathbb{N}_{0} $满足$ \left\|K x^{m, \delta}-y^{\delta}\right\| \leq r \delta $时, 终止算法. 下面证明对(1.3)用Crank-Nicolson离散得到的迭代方法, 按Morozov残差准则给出的停止规则是良定义的; 并且该方法在适当的光滑性假设下构成可接受的正则化策略并能达到最优收敛阶.

我们给出以下前提条件, 本节后面的定理中均使用该前提条件.

设$ K:X\to Y $是一个线性, 紧的, 具有稠密值域的单射算子. 设$ r>1 $, 对所有的$ \delta\in\left(0, \delta_0\right) $, $ y^\delta\in Y $是满足$ \lVert y-y^\delta\rVert \leq\delta $和$ \lVert y^\delta\rVert \geq r\delta $的扰动数据. 设序列$ x^{m, \delta}, m=0, 1, 2, ... $是由下面迭代格式得到的

$ \begin{equation*} \begin{split} x^{m+1, \delta}&=x^{m, \delta}+\frac{a}{2}K^*\left[\left(y-Kx^{m, \delta}\right)+\left(y-Kx^{m+1, \delta}\right)\right], \\ x^{0, \delta}&=0, \; \; m=0, 1, 2, ... , \end{split} \end{equation*} $

其中$ 0<a\leq2/\lVert K\rVert ^2 $.

定理3.2  对任意的$ \delta>0 $均有$ \lim_{m\to\infty}\lVert Kx^{m, \delta}-y^{\delta}\rVert =0 $.

证  在(3.2)中, 对所有的$ y\in Y $, 有表达式

$ \begin{equation*} R_my=\sum\limits_{j=1}^{\infty}\frac{1}{\mu_j}\left[1-\left(\frac{1-\frac{a}{2}\mu_j^2}{1+\frac{a}{2}\mu_j^2}\right)^m\right]\left(y, y_j\right)x_j, \end{equation*} $

因此

$ \begin{equation*} \lVert KR_my-y\rVert ^2=\sum\limits_{j=1}^{\infty}\left(\frac{1-\frac{a}{2}\mu_j^2}{1+\frac{a}{2}\mu_j^2}\right)^{2m}\left|\left(y, y_j\right)\right|^2. \end{equation*} $

由$ \left|\frac{1-\frac{a}{2}\mu_j^2}{1+\frac{a}{2}\mu_j^2}\right|<1 $, 有$ \lVert KR_m-I\rVert \leq1 $, 用$ y^{\delta} $代替$ y $可得

$ \begin{equation*} \lVert Kx^{m, \delta}-y^\delta\rVert ^2=\sum\limits_{j=1}^{\infty}\left(\frac{1-\frac{a}{2}\mu_j^2}{1+\frac{a}{2}\mu_j^2}\right)^{2m}\left|\left(y^\delta, y_j\right)\right|^2. \end{equation*} $

设$ \epsilon>0 $给定, 选择$ J\in\mathbb{N} $满足

$ \begin{equation*} \sum\limits_{j=J+1}^{\infty}\left|\left(y^\delta, y_j\right)\right|^2<\frac{\epsilon^2}{2}. \end{equation*} $

当$ m\to\infty $时, 由于对$ j=1, ..., J $一致的有$ \left|\frac{1-\frac{a}{2}\mu_j^2}{1+\frac{a}{2}\mu_j^2}\right|^{2m}\to0 $, 因此可以找到一个$ m_0\in\mathbb{N} $对所有的$ m\geq m_0 $满足

$ \begin{equation*} \sum\limits_{j=1}^J\left(\frac{1-\frac{a}{2}\mu_j^2}{1+\frac{a}{2}\mu_j^2}\right)^{2m}\left|\left(y^\delta, y_j\right)\right|^2\leq\max\limits_{j=1, ... , J}\left(\frac{1-\frac{a}{2}\mu_j^2}{1+\frac{a}{2}\mu_j^2}\right)^{2m}\sum\limits_{j=1}^J\left|\left(y^\delta, y_j\right)\right|^2\leq\frac{\epsilon^2}{2}, \; \; m\geq m_0. \end{equation*} $

这可以推出对所有的$ m\geq m_0 $, 有$ \lVert Kx^{m, \delta}-y^\delta\rVert ^2\leq\epsilon^2 $; 即该方法是可接受的.

定理3.3  当$ \delta\to0 $时, $ \delta^2m\to0 $, 即: $ m\left(\delta\right) $的选择是可接受的, 序列$ x^{m\left(\delta\right), \delta} $收敛到$ x^\dagger $.

证  只需考虑$ m\left(\delta\right)\to\infty $的情况, 记$ m:=m\left(\delta\right) $. 通过$ m\left(\delta\right) $的选择, 得到对$ y=Kx^\dagger $, 有

$ \begin{equation*} \begin{split} \lVert KR_{m-1}y-y\rVert &\geq\lVert KR_{m-1}y^\delta-y^\delta\rVert -\lVert \left(KR_{m-1}-I\right)\left(y-y^{\delta}\right)\rVert \\ &\geq r\delta-\lVert KR_{m-1}-I\rVert \delta\geq\left(r-1\right)\delta, \end{split} \end{equation*} $

因此

$ \begin{equation} \begin{split} m\left(r-1\right)^2\delta^2&\leq m\sum\limits_{j=1}^\infty\left(\frac{1-\frac{a}{2}\mu_j^2}{1+\frac{a}{2}\mu_j^2}\right)^{2m-2}\left|\left(y, y_j\right)\right|^2\\ &=\sum\limits_{j=1}^\infty m\left(\frac{1-\frac{a}{2}\mu_j^2}{1+\frac{a}{2}\mu_j^2}\right)^{2m-2}\mu_j^2\left|\left(x^\dagger, x_j\right)\right|^2. \end{split} \end{equation} $

(3.3)

注意到该级数在$ \delta\to0 $时收敛到0($ \delta $的依赖性隐藏在m中). 首先由引理3.1, 取$ v=0 $, 对所有的$ m\geq1 $和$ 0\leq\mu\leq\sqrt{2/a} $, 有$ m\mu^2\left(\frac{1-\frac{a}{2}\mu^2}{1+\frac{a}{2}\mu^2}\right)^{2m-2}\leq2/a $. 将级数分成一个有限和与一个剩余级数: 在剩余级数中, 用$ 2/a $代替$ m\left(\frac{1-\frac{a}{2}\mu_j^2}{1+\frac{a}{2}\mu_j^2}\right)^{2m-2}\mu_j^2 $, 并且注意到当$ m\to\infty $时, 对$ j\in 1, ..., J $, 都有$ m\left(\frac{1-\frac{a}{2}\mu_j^2}{1+\frac{a}{2}\mu_j^2}\right)^{2m-2}\to0 $, 这就证明了收敛性.

定理3.4  若存在某个$ v>0 $和某个$ z $满足$ \lVert z\rVert \leq E $, 使得$ x^\dagger=K^*z\in R\left(K^*\right) $或$ x^\dagger=\left(K^*K\right)^{v/2}z\in R\left(\left(K^*K\right)^{v/2}\right) $, 则有如下收敛阶

$ \begin{equation} \lVert x^{m\left(\delta\right), \delta}-x^\dagger\rVert \leq c\sqrt{E\delta} \end{equation} $

(3.4)

或

$ \begin{equation} \lVert x^{m\left(\delta\right), \delta}-x^\dagger\rVert \leq cE^\frac{1}{v+1}\delta^\frac{v}{v+1}, \end{equation} $

(3.5)

其中$ c>0 $为一常数. 这表明对所有$ v>0 $, $ m\left(\delta\right) $的选择是最优的.

证  由定理2.2和定理3.1可得$ \lVert R_m\rVert \leq\sqrt{am/2} $, 代入基本不等式(2.1)

$ \begin{equation} \lVert x^{m, \delta}-x^\dagger\rVert \leq\delta\sqrt{am/2}+\lVert R_mKx^\dagger-x^\dagger\rVert. \end{equation} $

(3.6)

下面严格证明$ x^\dagger=(K^*K)^{v/2}z $的情况: 记$ m=m\left(\delta\right) $, 由不等式(3.3), 有

$ \begin{equation*} \begin{split} \left(r-1\right)^2\delta^2\leq\lVert KR_{m-1}y-y\rVert ^2 &=\sum\limits_{j=1}^\infty\left(\frac{1-\frac{a}{2}\mu_j^2}{1+\frac{a}{2}\mu_j^2}\right)^{2m-2}\left|\left(y, y_j\right)\right|^2\\ &=\sum\limits_{j=1}^\infty\left(\frac{1-\frac{a}{2}\mu_j^2}{1+\frac{a}{2}\mu_j^2}\right)^{2m-2}\mu_j^2\left|\left(x^\dagger, x_j\right)\right|^2\\ &=\sum\limits_{j=1}^\infty\left(\frac{1-\frac{a}{2}\mu_j^2}{1+\frac{a}{2}\mu_j^2}\right)^{2m-2}\mu_j^{2v+2}\left|\left(z, x_j\right)\right|^2. \end{split} \end{equation*} $

由引理3.1, 对所有的$ m\geq1 $, $ v>0 $和$ 0\leq\mu\leq\sqrt{2/a} $, 有

$ \begin{equation*} m^{v+1}\mu^{2v+2}\left(\frac{1-\frac{a}{2}\mu^2}{1+\frac{a}{2}\mu^2}\right)^{2m-2} \leq m^{v+1}\mu^{2v+2}\left(1-\frac{a}{2}\mu^2\right)^{2m-2} \leq\left(\frac{2v+2}{a}\right)^{v+1}. \end{equation*} $

从而可以得到

$ \begin{equation*} \left(r-1\right)^2\delta^2\leq\left(\frac{2v+2}{am}\right)^{\left(v+1\right)}\lVert z\rVert ^2, \end{equation*} $

即

$ \begin{equation*} m(\delta)\leq cE^{2/{v+1}}\delta^{-2/{v+1}}. \end{equation*} $

其中$ c>0 $为某一常数. 为了证明不等式(3.6)的第二项, 用Hölder不等式

$ \begin{equation*} \sum\limits_{j=1}^\infty\left|a_jb_j\right|\leq\left[\sum\limits_{j=1}^\infty\left|a_j\right|^p\right]^{1/p}\left[\sum\limits_{j=1}^\infty\left|b_j\right|^q\right]^{1/q}, \end{equation*} $

其中$ p, q>1 $且$ 1/p+1/q=1 $. 取$ p=\frac{v+1}{v} $, $ q=v+1 $, 可以得出

$ \begin{equation*} \begin{split} \lVert(I-R_mK)x^\dagger\rVert^2&=\sum\limits_{j=1}^\infty\mu_j^{2v}\left(\frac{1-\frac{a}{2}\mu_j^2}{1+\frac{a}{2}\mu_j^2}\right)^{2m}|(z, x_j)|^2\\ &=\sum\limits_{j=1}^\infty\left[\mu_j^{2v}\left(\frac{1-\frac{a}{2}\mu_j^2}{1+\frac{a}{2}\mu_j^2}\right)^\frac{2vm}{v+1}\left|\left(z, x_j\right)\right|^\frac{2v}{v+1}\right]\left[\left(\frac{1-\frac{a}{2}\mu_j^2}{1+\frac{a}{2}\mu_j^2}\right)^\frac{2m}{v+1}\left|\left(z, x_j\right)\right|^\frac{2}{v+1}\right]\\ &\leq\left[\sum\limits_{j=1}^\infty\mu_j^{2v+2}\left(\frac{1-\frac{a}{2}\mu_j^2}{1+\frac{a}{2}\mu_j^2}\right)^{2m}\left|\left(z, x_j\right)\right|^{2}\right]^\frac{v}{v+1}\left[\sum\limits_{j=1}^\infty\left(\frac{1-\frac{a}{2}\mu_j^2}{1+\frac{a}{2}\mu_j^2}\right)^{2m}\left|\left(z, x_j\right)\right|^{2}\right]^\frac{1}{v+1}\\ &\leq\lVert KR_my-y\rVert ^{2v/v+1}\lVert z\rVert ^{2/v+1}\\ &\leq E^{2/v+1}\left[\lVert \left(I-KR_m\right)\left(y-y^{\delta}\right)\rVert +\lVert \left(I-KR_m\right)y^\delta\rVert \right]^{2v/v+1}. \end{split} \end{equation*} $

对第一项使用不等式$ \lVert KR_m-I\rVert \leq1 $, 对第二项使用停机准则可得

$ \begin{equation*} \lVert \left(I-R_mK\right)x^\dagger\rVert \leq E^{1/v+1}\left(1+r\right)^{v/v+1}\delta^{v/v+1}. \end{equation*} $

代入(3.6)得

$ \begin{equation*} \lVert x^{m\left(\delta\right), \delta}-x^\dagger\rVert \leq\delta\sqrt{am(\delta)/2}+\lVert R_{m(\delta)}Kx^\dagger-x^\dagger\rVert \leq cE^{1/v+1}\delta^{v/v+1}. \end{equation*} $

由文献[7]知, 当对隐式迭代正则化方法使用变步长并且取步长为公比大于1的等比数列时, 非稳态Tikhonov迭代正则化方法的迭代步数可以达到$ O(|log\delta|) $级别, 下面说明对Crank-Nicolson离散格式使用同样的步长并不能使迭代步数达到$ O(|log\delta|) $级别.

对非稳态Tikhonov迭代正则化方法, 有

$ \begin{equation} x^\dagger-x^m = r_m(K^*K)x^\dagger, \end{equation} $

(3.7)

其中

$ \begin{equation*} r_m(\lambda) = \prod\limits_{i=1}^m\left(\frac{1}{1+a_i\lambda}\right). \end{equation*} $

文献[7]中说明(3.7)收敛当且仅当

$ \begin{equation*} \sum\limits_{i=1}^ma_i=\infty. \end{equation*} $

类似的, 对于Crank-Nicolson离散格式, 有

$ \begin{equation*} x^\dagger-x^m = r_m(K^*K)x^\dagger, \end{equation*} $

其中

$ \begin{equation*} r_m(\lambda) = \prod\limits_{i=1}^{m}\left(\frac{2-a_i\lambda}{2+a_i\lambda}\right). \end{equation*} $

令$ x_i=a_i\lambda $, 由$ \lambda\in(0, \lVert K^*K\rVert ] $, 级数$ \sum_{i=1}^\infty x_i $与级数$ \sum_{i=1}^\infty a_i\lambda $同敛散, 为简化证明, 考虑级数$ \sum_{i=1}^\infty x_i $敛散性与$ r_m=\prod_{i=1}^{m}\left(\frac{2-x_i}{2+x_i}\right) $收敛性之间的关系.

当$ x_i\to0 $且$ x_0<2 $时,

$ \begin{equation*} \lim\limits_{i\to\infty}\prod\limits_{i=1}^\infty\frac{2-x_i}{2+x_i}=0\Leftrightarrow \lim\limits_{i\to\infty}\sum\limits_{i=1}^\infty ln(1-\frac{2x_i}{2+x_i})=-\infty\Leftrightarrow \sum\limits_{i=1}^\infty x_i = \infty. \end{equation*} $

当$ x_i\to\infty $且$ x_0>2 $时, 只能得到

$ \begin{equation*} \lim\limits_{i\to\infty}\prod\limits_{i=1}^\infty\frac{x_i-2}{2+x_i}=0\Leftrightarrow \lim\limits_{i\to\infty}\sum\limits_{i=1}^\infty ln(1-\frac{4}{2+x_i})=-\infty\Rightarrow \sum\limits_{i=1}^\infty x_i = \infty, \end{equation*} $

还可以得到

$ \begin{equation} \lim\limits_{i\to\infty}\prod\limits_{i=1}^\infty\frac{x_i-2}{2+x_i}=0\Leftrightarrow \lim\limits_{i\to\infty}\sum\limits_{i=1}^\infty ln(1-\frac{4}{2+x_i})=-\infty\Leftrightarrow \sum\limits_{i=1}^\infty \frac{1}{x_i} = \infty. \end{equation} $

(3.8)

当选取步长为公比大于1的等比数列$ \{a_i\} $时, $ \sum_{i=1}^\infty \frac{1}{a_i}<\infty $, 级数和不满足(3.8). 这说明Crank-Nicolson离散格式的收敛条件比非稳态Tikhonov迭代正则化要更严格, 对于此类等比数列, Crank-Nicolson离散格式不能达到对数阶迭代步数, 但在数值实验中发现当误差水平$ \delta $较大时, 该离散格式的迭代步数可以像非稳态Tikhonov迭代正则化方法一样达到对数阶.

求解如下线性积分方程:

$ \begin{equation} Kx(s) = \int_0^1k(s, t)x(t)dt=y(s), \end{equation} $

(4.1)

其中

$ \begin{equation} k(s, t) = \begin{cases} s(1 - t), & s \le t, \\[6pt] t(1 - s), & s > t. \end{cases} \end{equation} $

(4.2)

当$ X=Y=L^2[0, 1] $时, 算子$ K $是紧的, 自伴的, 单射算子. 若$ y\in H^2[0, 1]\cap H_0^1[0, 1] $, 则方程(4.1)有解$ x=-y'' $.

考虑线性积分方程(4.1)的两个不同的精确右端项:

例1  $ y(s)=s^4(1-s)^5 $, $ x^\dagger(t)=-4t^2(1-t)^3(18t^2-16t+3) $, $ K $由梯形法离散, 取离散尺度为100.

例2  $ y(s)=\frac{sin(2\pi s)}{\pi^2} $, $ x^\dagger(t) = 4sin(2\pi t) $, $ K $由梯形法离散, 取离散尺度为100.

使用带噪声的数据$ y^\delta = y+\delta\frac{\epsilon}{\rVert\epsilon\lVert} $, 其中$ y=Kx^\dagger $, $ \epsilon $为服从标准正态分布$ N(0, I) $的向量. 采用Morozov残差准则作为停机条件, 即$ \rVert Kx^{m, \delta}-y^\delta\lVert\leq r\delta $, 取$ r=1.01 $. 迭代算法中的步长$ a $在条件范围内取尽量大的值, 这里取$ a=\frac{1}{\rVert K\lVert^2} $. 记相对均方误差$ err $为:

$ \begin{equation*} err = \frac{\|x^{m, \delta} - x^\dagger\|}{\|x^\dagger\|}. \end{equation*} $

记Landweber迭代正则化方法为LW, 隐式迭代正则化方法为IM, Crank-Nicolson离散格式为C-N.

表 1和表 2分别对例1和例2使用了三种固定步长的迭代正则化方法, 实验表明Crank-Nicolson离散格式结果和我们的理论分析相符, 是正则化方法并且可以达到最优收敛阶, 但是其用时相比于前两种方法较长.

对于非稳态Tikhonov迭代正则化方法和变步长的Crank-Nicolson离散格式, 我们取步长$ a_i=aq^{i-1} $, 其中$ a=1 $, $ q=1.1 $. 记非稳态Tikhonov迭代正则化方法为NT, 变步长的Crank-Nicolson离散格式为NC-N.

表 3和表 4中, 当误差$ \delta $较大时, NC-N方法和NT方法的实验结果相近. 对于例1, NC-N在相对误差和计算用时上优于NT方法, 但当$ \delta\leq10^{-4} $时, NC-N方法在数值实验中不再收敛, 和我们的理论推导相符.

本文对连续正则化方法进行Crank-Nicolson离散, 证明了该迭代格式为正则化方法并且可以达到最优收敛阶, 说明其正则化性能和Landweber迭代正则化方法一致. 进一步的理论分析表明, 与非稳态Tikhonov迭代正则化方法不同, Crank-Nicolson格式在变步长策略下具有更严格的收敛条件, 因此一般的等比步长不能使迭代步数达到$ O(|log\delta|) $级别, 限制了其在变步长情形下的应用. 数值实验表明固定步长的Crank-Nicolson方法具有良好的正则化效果并实现最优误差阶, 但计算耗时相对更高; 变步长情况下在高噪声水平下表现良好, 而在小噪声时可能出现不收敛, 与理论分析一致. 综上, Crank-Nicolson离散格式是一种可行且具有最优误差性质的迭代正则化方法, 但在效率和变步长策略方面仍存在一定局限.