求解大型稀疏线性系统的贪婪双子空间随机Kaczmarz方法

摘要

基于一种有效的从系数矩阵中选取两个工作行的贪婪概率准则， 提出一类求解大型稀疏线性系统的贪婪双子空间随机Kaczmarz方法。理论证明该方法收敛到相容线性系统的最小范数解， 而且该方法的理论收敛因子小于原始双子空间随机Kaczmarz方法的收敛因子。数值实验表明，该方法在求解性能方面较原始双子空间随机Kaczmarz方法更具优势。

关键词

双子空间投影法; 随机Kaczmarz方法; 大型稀疏线性系统; 收敛性分析

考虑求解具有如下形式的相容线性系统：

$$\mathit{A}\mathit{x}=\mathit{b}$$

（1）

其中，系数矩阵$\mathit{A}\in {\mathbf{R}}^{m\times n}(m>n)$可为满秩或秩亏矩阵， 且$\mathit{b}\in {\mathbf{R}}^m$。通常考虑求式（1）的最小范数解

$${\mathit{x}}_{\mathrm{*}}:=\mathrm{a}\mathrm{r}\mathrm{g}\underset{\mathit{x}}{\mathrm{m}\mathrm{i}\mathrm{n}}{‖\mathit{x}‖}_{\mathrm{2}}\hspace{1em}\mathrm{s}.\mathrm{t}.\mathit{A}\mathit{x}=\mathit{b}$$

当系数矩阵$\mathit{A}$为列满秩时，${\mathit{x}}_{\mathrm{*}}$为式（1）的惟一解；当线性系统（1）有无穷多组解时，${\mathit{x}}_{\mathrm{*}}$为式（1）的最小范数解（这里的范数指欧式范数）。

为求解线性系统（1）， 许多迭代方法^［1-5］已被开发和进一步研究。其中， Kaczmarz^［5］于1937年首次提出的Kaczmarz方法， 因其成本低廉、易于操作而迅速得到数值计算领域的专家学者的认可和广泛关注，进而使其获得巨大的理论发展^［6-9］。因其是一种具有代表性的行处理迭代方法且在计算机上易于实现和并行化，因此被广泛应用于计算机断层扫描^［10-13］，图像重建^［14-16］，分布式计算^［17-18］和信号处理^［19-21］等领域。

若以随机顺序而不是以给定顺序选用系数矩阵的行可以极大地提高Kaczmarz方法的收敛速度^［13，16，21］。尽管这些随机选行的Kaczmarz方法在应用中颇具吸引力，但尚不能保证其收敛速率。Strohmer等^［22］首次采用选取工作行的概率与其所对向量的欧式范数的平方成正比这样的选行准则，提出能保证收敛速率的随机Kaczmarz方法，并证明其误差的期望具有指数收敛速率。Dai等^［23］通过求解最小化收敛速度的上限这个凸优化问题来获得选择行的最佳概率分布，提出最优的随机Kaczmarz方法；Bai等^［24］结合贪婪和随机的数学思想引入一种有效的行选择准则提出贪婪随机Kaczmarz方法。此后， 松弛贪婪随机Kaczmarz方法^［25］和贪婪距离随机Kaczmarz方法^［26］也被提出并深入研究。随机Kaczmarz方法的各种变种^［27-39］和收敛理论^［30-35］获得了大量研究和发展。另一方面， Leventhal等^［36］采用同样的随机思想， 提出了求解最小二乘问题的随机坐标下降方法。因随机Kaczmarz方法和随机坐标下降方法存在一定的共性，后期有学者将两者同时比对研究^［37-38］。

尽管随机Kaczmarz方法适用于任何相容线性系统的求解，但当系数矩阵有许多相关行时（常出现于地球物理学中），收敛可能停滞。为克服这一缺点，Needell和Ward采用等可能地随机选择两个不同工作行的选择准则，提出了随机Kaczmarz方法的双子空间拓展‒双子空间随机Kaczmarz方法^［27］，并从理论和数值上说明了该方法至少与随机Kaczmarz方法有相同的指数收敛速率。此外，当系数矩阵高度相干时，它能极大程度地提高随机Kaczmarz方法的收敛速率。此后同时选用多个工作行的随机块方法^［39-40］也被陆续提出并加以深入研究。

Nutini等^［41］研究表明，采用非等可能概率准则选取工作行的Kaczmnarz方法收敛速度至少与采用等可能概率准则选取工作行的Kaczmarz方法收敛速度一样快。在原始双子空间随机Kaczamrz方法^［27］中， 等可能地随机选择两个不同的工作行，然后将当前解投影到由这两行所确定的解的超平面上获得下一个近似解。基于上一次迭代所生成的残差或基于当前解离系数矩阵各行所形成的超平面的距离来选取工作行都能较大程度地提高Kaczmarz方法的收敛速率^［41］。受该思想启发，笔者通过引入控制参数$\theta$建立了一种基于残差的准则来选择工作行， 导出一类贪婪双子空间随机Kaczmarz方法。理论证明新方法收敛到相容线性方程组的最小范数解， 而且新方法的理论收敛因子小于原始双子空间随机Kaczmarz方法的收敛因子。数值实验验证了新方法的有效性， 其在迭代步数和计算时间上均优于原始双子空间随机Kaczmarz方法。

1 贪婪双子空间随机Kaczmarz方法

在经典的随机Kaczmarz方法（RK）中，Strohmer和Vershynin将各行所对向量的欧式范数与${‖\mathit{A}‖}_{\mathrm{F}}$比值的平方作为概率， 然后根据此概率分布随机选取工作行。具体地说， 如果定义${\mathit{A}}^{\left(i\right)}$代表系数矩阵$\mathit{A}$的第$i$行， ${\mathit{b}}^{\left(i\right)}$代表向量$\mathit{b}$的第$i$个分量， 初始向量为${\mathit{x}}_{\mathrm{0}}$， 则随机Kaczmarz方法的具体过程见算法1。

算法1  ^［22］     随机Kaczmarz方法。①置$k:=\mathrm{0}$。②根据概率$P(r=i_k)=\frac{{‖{\mathit{A}}^{(i_k)}‖}_{\mathrm{2}}^{\mathrm{2}}}{{‖\mathit{A}‖}_{\mathrm{F}}^{\mathrm{2}}}$ 选取指标$i_k\in \left\{\mathrm{1,2},\dots ,m\right\}$。

③计算${\mathit{x}}_{k+\mathrm{1}}={\mathit{x}}_k+\frac{{\mathit{b}}^{(i_k)}-{\mathit{A}}^{(i_k)}{\mathit{x}}_k}{{‖{\mathit{A}}^{(i_k)}‖}_{\mathrm{2}}^{\mathrm{2}}}{\left({\mathit{A}}^{(i_k)}\right)}^{\mathrm{T}}$。 ④置$k=k+\mathrm{1}$， 转步骤②。

算法1中， $P(r=i_k)=\frac{{‖{\mathit{A}}^{(i_k)}‖}_{\mathrm{2}}^{\mathrm{2}}}{{‖\mathit{A}‖}_{\mathrm{F}}^{\mathrm{2}}}$代表选取矩阵$\mathit{A}$的第$i_k$行作为本次迭代工作行的概率为$\frac{{‖{\mathit{A}}^{(i_k)}‖}_{\mathrm{2}}^{\mathrm{2}}}{{‖\mathit{A}‖}_{\mathrm{F}}^{\mathrm{2}}}$。

关于随机Kaczmarz方法， 有如下收敛性定理。

定理1^［22，24］ 若线性系统（1）相容， 其中系数矩阵$\mathit{A}\in {\mathbf{R}}^{m\times n}$且右端项$\mathit{b}\in {\mathbf{R}}^m$。初始向量${\mathit{x}}_{\mathrm{0}}\in \mathrm{R}\mathrm{a}\mathrm{n}({\mathit{A}}^{\mathrm{T}})$，其中$\mathrm{R}\mathrm{a}\mathrm{n}({\mathit{A}}^{\mathrm{T}})$表示${\mathit{A}}^{\mathrm{T}}$的列空间， 令${\mathit{x}}_k$为通过随机Kaczmarz方法生成的第$k$个迭代值，则

$$E{‖{\mathit{x}}_{\mathrm{*}}-{\mathit{x}}_k‖}_{\mathrm{2}}^{\mathrm{2}}\le {\left(\mathrm{1}-\frac{{\sigma }_{{}_{\mathrm{m}\mathrm{i}\mathrm{n}}}^{\mathrm{2}}(\mathit{A})}{{‖\mathit{A}‖}_{\mathrm{F}}^{\mathrm{2}}}\right)}^k{‖{\mathit{x}}_{\mathrm{*}}-{\mathit{x}}_{\mathrm{0}}‖}_{\mathrm{2}}^{\mathrm{2}}$$

其中${\sigma }_{\mathrm{m}\mathrm{i}\mathrm{n}}(\cdot )$表示矩阵的最小非零奇异值。

Needell和Ward发现当线性系统是超定且高度相干时，随机Kaczmarz方法的求解效率低下甚至无效，因此提出同时启用两个工作行的双子空间随机Kaczmarz方法（2S‒RK）去求解这类特殊系统。

为了简便，Needell和Ward假设系数矩阵是标准化的，这意味着其每行都具有单位欧几里得范数，在接下来的部分，仍沿用该假设。在原始双子空间随机Kaczmarz方法中， 等可能地随机选取两个不同工作行， 再将当前解投影到由两个工作行所确定的超平面上获得下一个估计值， 具体过程见算法2。

算法2  ^［31］     双子空间随机Kaczmarz方法。①置$k:=\mathrm{0}$。②随机均匀地选择行指标$r_k$和$s_k$。③计算${\mu }_{r_k,s_k}=<{\mathit{A}}^{(r_k)},{\mathit{A}}^{(s_k)}>$。④计算${\mathit{y}}_k={\mathit{x}}_k+\left({\mathit{b}}^{(s_k)}-{\mathit{A}}^{(s_k)}{\mathit{x}}_k\right){\left({\mathit{A}}^{(s_k)}\right)}^{\mathrm{T}}$。 ⑤计算${\mathit{\nu }}_k=\frac{{\mathit{A}}^{(r_k)}-{\mu }_{r_k,s_k}{\mathit{A}}^{(s_k)}}{\sqrt[]{\mathrm{1}-{\left|{\mu }_{r_k,s_k}\right|}^{\mathrm{2}}}}$。⑥计算${\beta }_k=\frac{{\mathit{b}}^{(r_k)}-{\mu }_{r_k,s_k}{\mathit{b}}^{(s_k)}}{\sqrt[]{\mathrm{1}-{\left|{\mu }_{r_k,s_k}\right|}^{\mathrm{2}}}}$。⑦计算${\mathit{x}}_{k+\mathrm{1}}={\mathit{y}}_k+\left({\beta }_k-{\mathit{\nu }}_k{\mathit{y}}_k\right){\left({\mathit{\nu }}_k\right)}^{\mathrm{T}}\mathrm{。}$置$k=k+\mathrm{1}$， 转步②。

关于双子空间随机Kaczmarz方法， 有如下收敛性定理。

定理2  ^［31］     若线性系统（1）相容， 其中系数矩阵$\mathit{A}\in {\mathbf{R}}^{m\times n}$且右端项$\mathit{b}\in {\mathbf{R}}^m$。初始向量${\mathit{x}}_{\mathrm{0}}\in \mathrm{R}\mathrm{a}\mathrm{n}({\mathit{A}}^{\mathrm{T}})$， 其中$\mathrm{R}\mathrm{a}\mathrm{n}({\mathit{A}}^{\mathrm{T}})$表示${\mathit{A}}^{\mathrm{T}}$的列空间， 令${\mathit{x}}_{k+\mathrm{1}}$为通过双子空间随机Kaczmarz方法生成的第$(k+\mathrm{1})$个迭代值，则

$$\begin{array}{l}E{‖{\mathit{x}}_{\mathrm{*}}-{\mathit{x}}_{k+\mathrm{1}}‖}_{\mathrm{2}}^{\mathrm{2}}\le {\left(\mathrm{1}-\frac{{\sigma }_{{}_{\mathrm{m}\mathrm{i}\mathrm{n}}}^{\mathrm{2}}(A)}{m}\right)}^{\mathrm{2}}{‖{\mathit{x}}_{\mathrm{*}}-{\mathit{x}}_k‖}_{\mathrm{2}}^{\mathrm{2}}-\\ \frac{\mathrm{1}}{m(m-\mathrm{1})}\sum _{r\ne s}\left|\frac{{\mu }_{r,s}^{\mathrm{2}}}{\sqrt[]{\mathrm{1}-{\mu }_{r,s}^{\mathrm{2}}}}<{\mathit{e}}_k^{\mathrm{T}},{\mathit{A}}^{(r)}>-\right.\\ {\left.\frac{{\mu }_{r,s}}{\sqrt[]{\mathrm{1}-{\mu }_{r,s}^{\mathrm{2}}}}<{\mathit{e}}_k^{\mathrm{T}},{\mathit{A}}^{(s)}>\right|}^{\mathrm{2}}\end{array}$$

（2）

这里，${\mathit{e}}_k={\mathit{x}}_{\mathrm{*}}-{\mathit{x}}_k$，${\mu }_{r,s}=<{\mathit{A}}^{(r)},{\mathit{A}}^{(s)}>$。

通过算法2，可知双子空间随机Kaczmarz方法的相应误差向量的2‒范数的平方满足

$$\begin{array}{l}{‖{\mathit{x}}_{\mathrm{*}}-{\mathit{x}}_{k+\mathrm{1}}‖}_{\mathrm{2}}^{\mathrm{2}}={‖{\mathit{x}}_{\mathrm{*}}-{\mathit{x}}_k‖}_{\mathrm{2}}^{\mathrm{2}}-\\ {\left|<{\mathit{e}}_k^{\mathrm{T}},{\mathit{A}}^{(s_k)}>\right|}^{\mathrm{2}}-{\left|<({\mathit{x}}_{\mathrm{*}}-{\mathit{y}}_k{)}^{\mathrm{T}},{\mathit{\nu }}_k>\right|}^{\mathrm{2}}=\\ {‖{\mathit{x}}_{\mathrm{*}}-{\mathit{x}}_k‖}_{\mathrm{2}}^{\mathrm{2}}-{\left|{\mathit{b}}^{(s_k)}-{\mathit{A}}^{(s_k)}{\mathit{x}}_k\right|}^{\mathrm{2}}-\\ \frac{{\left|{\mathit{b}}^{(r_k)}-{\mathit{A}}^{(r_k)}{\mathit{y}}_k\right|}^{\mathrm{2}}}{\mathrm{1}-{\left|{\mu }_{r_k,s_k}\right|}^{\mathrm{2}}}\end{array}$$

（3）

若式（3）的最后两项的和取到最大值，则${\mathit{x}}_{\mathrm{*}}$和${\mathit{x}}_{k+\mathrm{1}}$之间的距离可以最小化。基于这个想法，设计最优的贪婪距离选行准则为

$$(r_k,s_k)=\underset{\mathrm{1}\le s,r\le m}{\mathrm{m}\mathrm{a}\mathrm{x}}\left\{{\left|{\mathit{b}}^{(s)}-{\mathit{A}}^{(s)}{\mathit{x}}_k\right|}^{\mathrm{2}}+\frac{{\left|{\mathit{b}}^{(r)}-{\mathit{A}}^{(r)}{\mathit{y}}_k\right|}^{\mathrm{2}}}{\mathrm{1}-{\left|<{\mathit{A}}^{(r)},{\mathit{A}}^{(s)}>\right|}^{\mathrm{2}}}\right\}$$

（4）

若双子空间随机Kaczmarz方法采用式（4）的选行准则，则其定能以一个极快的速率收敛到线性系统（1）的解。因此，可采用式（4）的规则来建立双子空间随机Kaczmarz方法的最佳版本。但在实践中计算满足式（4）行指标$r_k$和$s_k$的成本非常昂贵。为克服其耗时的缺点，笔者将通过依次构造两个具有控制参数$\theta$的行索引集来构造次优版本的贪婪双子空间随机Kaczmarz方法（2S‒GRK（θ）），具体过程见算法3。

首先，构造与$x_k$有关的指标集$U_k$

$$U_k=\left\{s\left||\right.{\mathit{b}}^{(s)}-{\mathit{A}}^{(s)}{\mathit{x}}_k{|}^{\mathrm{2}}\ge \left(\mathrm{1}-\theta \right)\underset{\mathrm{1}\le i\le m}{\mathrm{m}\mathrm{a}\mathrm{x}}\left\{{\left|{\mathit{b}}^{(i)}-{\mathit{A}}^{(i)}{\mathit{x}}_k\right|}^{\mathrm{2}}\right\}\right\},\theta \in [\mathrm{0,1}]$$

（5）

并以概率$P(r=s_k)=\frac{{\left|{\widehat{\mathit{r}}}_{x_k}^{{}^{(s_k)}}\right|}^{\mathrm{2}}}{{‖{\widehat{\mathit{r}}}_{x_k}‖}_{\mathrm{2}}^{\mathrm{2}}}$选取指标$s_k\in U_k$，其中向量${\widehat{\mathit{r}}}_{x_k}$的定义见算法3的步3。将当前解${\mathit{x}}_k$投影到解空间$\left\{\mathit{z}|{\mathit{A}}^{(s_k)}\mathit{z}={\mathit{b}}^{(s_k)}\right\}$得到中间解向量${\mathit{y}}_k$

$${\mathit{y}}_k={\mathit{x}}_k+\left({\mathit{b}}^{(s_k)}-{\mathit{A}}^{(s_k)}{\mathit{x}}_k\right){\left({\mathit{A}}^{(s_k)}\right)}^{\mathrm{T}}$$

由此可得，对于$k=\mathrm{1,2},\mathrm{3}\cdots$成立

$${\mathit{r}}_{y_k}^{(s_k)}={\mathit{b}}^{(s_k)}-{\mathit{A}}^{(s_k)}{\mathit{y}}_k={\mathit{b}}^{(s_k)}-{\mathit{A}}^{(s_k)}{\mathit{x}}_k-\left({\mathit{b}}^{(s_k)}-{\mathit{A}}^{(s_k)}{\mathit{x}}_k\right){\mathit{A}}^{(s_k)}{\left({\mathit{A}}^{(s_k)}\right)}^{\mathrm{T}}=\mathrm{0}$$

即是${‖\mathit{b}-\mathit{A}{\mathit{y}}_k‖}_{\mathrm{2}}^{\mathrm{2}}={\sum _{j\ne s_k}\left|{\mathit{b}}^{(j)}-{\mathit{A}}^{(j)}{\mathit{y}}_k\right|}^{\mathrm{2}}$。

其次，构造与${\mathit{y}}_k$有关的指标集${\widehat{U}}_k$

$${\widehat{U}}_k=\left\{r\left||\right.{\mathit{b}}^{(r)}-{\mathit{A}}^{(r)}{\mathit{y}}_k{|}^{\mathrm{2}}\ge \left(\mathrm{1}-\theta \right)\underset{\mathrm{1}\le j\le m}{\mathrm{m}\mathrm{a}\mathrm{x}}\left\{{\left|{\mathit{b}}^{(j)}-{\mathit{A}}^{(j)}{\mathit{y}}_k\right|}^{\mathrm{2}}\right\}\right\},\theta \in [\mathrm{0,1}]$$

并以概率$P(r=r_k)=\frac{{\left|{\widehat{\mathit{r}}}_{y_k}^{{}^{(r_k)}}\right|}^{\mathrm{2}}}{{‖{\widehat{\mathit{r}}}_{y_k}‖}_{\mathrm{2}}^{\mathrm{2}}}$选取指标${\mathit{r}}_k\in {\widehat{U}}_k$，其中向量${\widehat{\mathit{r}}}_{y_k}$的定义见算法3的步8。将当前解${\mathit{y}}_k$投影到解空间$\left\{\mathit{z}|{\mathit{A}}^{(r_k)}\mathit{z}={\mathit{b}}^{(r_k)}\right\}$得到满足${\mathit{A}}^{(s_k)}{\mathit{x}}_{k+\mathrm{1}}={\mathit{b}}^{(s_k)}$的解向量${\mathit{x}}_{k+\mathrm{1}}$

$$\begin{array}{l}{\mathit{x}}_{k+\mathrm{1}}={\mathit{y}}_k+\left(\frac{{\mathit{b}}^{(r_k)}-<{\mathit{A}}^{(r_k)},{\mathit{A}}^{(s_k)}>{\mathit{b}}^{(s_k)}}{\sqrt[]{\mathrm{1}-{\left|<{\mathit{A}}^{(r_k)},{\mathit{A}}^{(s_k)}>\right|}^{\mathrm{2}}}}-\right.\\ \left.\frac{{\mathit{A}}^{(r_k)}-<{\mathit{A}}^{(r_k)},{\mathit{A}}^{(s_k)}>{\mathit{A}}^{(s_k)}}{\sqrt[]{\mathrm{1}-{\left|<{\mathit{A}}^{(r_k)},{\mathit{A}}^{(s_k)}>\right|}^{\mathrm{2}}}}{\mathit{y}}_k\right)\times \\ {\left(\frac{{\mathit{A}}^{(r_k)}-<{\mathit{A}}^{(r_k)},{\mathit{A}}^{(s_k)}>{\mathit{A}}^{(s_k)}}{\sqrt[]{\mathrm{1}-{\left|<{\mathit{A}}^{(r_k)},{\mathit{A}}^{(s_k)}>\right|}^{\mathrm{2}}}}\right)}^{\mathrm{T}}\end{array}$$

由指标集${\widehat{U}}_k$的定义，对于$\forall r\in {\widehat{U}}_k$，可得$\begin{array}{l}{\left|{\mathit{r}}_{y_k}^{(r)}\right|}^{\mathrm{2}}={\left|{\mathit{b}}^{(r)}-{\mathit{A}}^{(r)}{\mathit{y}}_k\right|}^{\mathrm{2}}\ge \\ (\mathrm{1}-\theta )\underset{\mathrm{1}\le j\le m}{\mathrm{m}\mathrm{a}\mathrm{x}}\left\{{\left|{\mathit{b}}^{(j)}-{\mathit{A}}^{(j)}{\mathit{y}}_k\right|}^{\mathrm{2}}\right\}=\end{array}$

$$\begin{array}{l}(\mathrm{1}-\theta ){‖\mathit{b}-\mathit{A}{\mathit{y}}_k‖}_{\mathrm{2}}^{\mathrm{2}}\frac{\underset{\mathrm{1}\le j\le m}{\mathrm{m}\mathrm{a}\mathrm{x}}\left\{{\left|{\mathit{b}}^{(j)}-{\mathit{A}}^{(j)}{\mathit{y}}_k\right|}^{\mathrm{2}}\right\}}{{‖\mathit{b}-\mathit{A}{\mathit{y}}_k‖}_{\mathrm{2}}^{\mathrm{2}}}=\\ (\mathrm{1}-\theta ){‖\mathit{b}-\mathit{A}{\mathit{y}}_k‖}_{\mathrm{2}}^{\mathrm{2}}\frac{\underset{\mathrm{1}\le j\le m}{\mathrm{m}\mathrm{a}\mathrm{x}}\left\{{\left|{\mathit{b}}^{(j)}-{\mathit{A}}^{(j)}{\mathit{y}}_k\right|}^{\mathrm{2}}\right\}}{\sum _{j\ne s_k}{\left|{\mathit{b}}^{(j)}-{\mathit{A}}^{(j)}{\mathit{y}}_k\right|}^{\mathrm{2}}}\ge \\ (\mathrm{1}-\theta ){‖\mathit{b}-\mathit{A}{\mathit{y}}_k‖}_{\mathrm{2}}^{\mathrm{2}}\frac{\underset{\mathrm{1}\le j\le m}{\mathrm{m}\mathrm{a}\mathrm{x}}\left\{{\left|{\mathit{b}}^{(j)}-{\mathit{A}}^{(j)}{\mathit{y}}_k\right|}^{\mathrm{2}}\right\}}{(m-\mathrm{1})\underset{\mathrm{1}\le j\le m}{\mathrm{m}\mathrm{a}\mathrm{x}}\left\{{\left|{\mathit{b}}^{(j)}-{\mathit{A}}^{(j)}{\mathit{y}}_k\right|}^{\mathrm{2}}\right\}}=\\ \frac{\mathrm{1}-\theta }{m-\mathrm{1}}{‖\mathit{b}-\mathit{A}{\mathit{y}}_k‖}_{\mathrm{2}}^{\mathrm{2}}\end{array}$$

（6）

对于$k=\mathrm{1,2},\mathrm{3}\cdots$，同样可得${\mathit{r}}_{x_{k+\mathrm{1}}}^{(s_k)}={\mathit{b}}^{(s_k)}-{\mathit{A}}^{(s_k)}{\mathit{x}}_{k+\mathrm{1}}=\mathrm{0}$和${\mathit{r}}_{x_{k+\mathrm{1}}}^{(r_k)}={\mathit{b}}^{(r_k)}-{\mathit{A}}^{(r_k)}{\mathit{x}}_{k+\mathrm{1}}=\mathrm{0}$，即是${‖\mathit{b}-\mathit{A}{\mathit{x}}_{k+\mathrm{1}}‖}_{\mathrm{2}}^{\mathrm{2}}={\sum _{i\ne s_k,r_k}\left|{\mathit{b}}^{(i)}-{\mathit{A}}^{(i)}{\mathit{x}}_{k+\mathrm{1}}\right|}^{\mathrm{2}}$。同理， 由指标集$U_k$的定义，对于$\forall s\in U_k$，可得

$$\begin{array}{l}{\left|{\mathit{r}}_{x_k}^{(s)}\right|}^{\mathrm{2}}={\left|{\mathit{b}}^{(s)}-{\mathit{A}}^{(i)}{\mathit{x}}_k\right|}^{\mathrm{2}}\ge \\ (\mathrm{1}-\theta )\underset{\mathrm{1}\le i\le m}{\mathrm{m}\mathrm{a}\mathrm{x}}\left\{{\left|{\mathit{b}}^{(i)}-{\mathit{A}}^{(i)}{\mathit{x}}_k\right|}^{\mathrm{2}}\right\}=\\ (\mathrm{1}-\theta ){‖\mathit{b}-\mathit{A}{\mathit{x}}_k‖}_{\mathrm{2}}^{\mathrm{2}}\frac{\underset{\mathrm{1}\le i\le m}{\mathrm{m}\mathrm{a}\mathrm{x}}\left\{{\left|{\mathit{b}}^{(i)}-{\mathit{A}}^{(i)}{\mathit{x}}_k\right|}^{\mathrm{2}}\right\}}{{‖\mathit{b}-\mathit{A}{\mathit{x}}_k‖}_{\mathrm{2}}^{\mathrm{2}}}=\\ (\mathrm{1}-\theta ){‖\mathit{b}-\mathit{A}{\mathit{x}}_k‖}_{\mathrm{2}}^{\mathrm{2}}\frac{\underset{\mathrm{1}\le i\le m}{\mathrm{m}\mathrm{a}\mathrm{x}}\left\{{\left|{\mathit{b}}^{(i)}-{\mathit{A}}^{(i)}{\mathit{x}}_k\right|}^{\mathrm{2}}\right\}}{\sum _{i\ne s_{k-\mathrm{1}},r_{k-\mathrm{1}}}{\left|{\mathit{b}}^{(i)}-{\mathit{A}}^{(i)}{\mathit{x}}_k\right|}^{\mathrm{2}}}\ge \\ (\mathrm{1}-\theta ){‖\mathit{b}-\mathit{A}{\mathit{x}}_k‖}_{\mathrm{2}}^{\mathrm{2}}\frac{\underset{\mathrm{1}\le i\le m}{\mathrm{m}\mathrm{a}\mathrm{x}}\left\{{\left|{\mathit{b}}^{(i)}-{\mathit{A}}^{(i)}{\mathit{x}}_k\right|}^{\mathrm{2}}\right\}}{(m-\mathrm{2})\underset{\mathrm{1}\le i\le m}{\mathrm{m}\mathrm{a}\mathrm{x}}\left\{{\left|{\mathit{b}}^{(i)}-{\mathit{A}}^{(i)}{\mathit{x}}_k\right|}^{\mathrm{2}}\right\}}=\\ \frac{\mathrm{1}-\theta }{m-\mathrm{2}}{‖\mathit{b}-\mathit{A}{\mathit{x}}_k‖}_{\mathrm{2}}^{\mathrm{2}}\end{array}$$

（7）

算法3  　贪婪双子空间随机Kaczmarz方法。

步1 置$k:=\mathrm{0}$。计算${\epsilon }_k=\left(\mathrm{1}-\theta \right)\underset{\mathrm{1}\le i\le m}{\mathrm{m}\mathrm{a}\mathrm{x}}\left\{{\left|{\mathit{b}}^{(i)}-{\mathit{A}}^{(i)}{\mathit{x}}_k\right|}^{\mathrm{2}}\right\}\mathrm{。}$定义正整数指标集$U_k=\left\{s\left||\right.{\mathit{b}}^{(s)}-{\mathit{A}}^{(s)}{\mathit{x}}_k{|}^{\mathrm{2}}\ge {\epsilon }_k\right\}\mathrm{。}$ 计算 ${\widehat{r}}_{x_k}$ 的第$s$个分量${\widehat{r}}_{x_k}^{(s)}=\left\{\begin{array}{ll}{\mathit{b}}^{(s)}-{\mathit{A}}^{(s)}{\mathit{x}}_k,& \hspace{1em}s\in U_k\\ \mathrm{0},& \hspace{1em}s\notin U_k\end{array}\right.$

步2 根据概率$P(r=s_k)=\frac{{\left|{\widehat{\mathit{r}}}_{x_k}^{(s_k)}\right|}^{\mathrm{2}}}{{‖{\widehat{\mathit{r}}}_{x_k}‖}_{\mathrm{2}}^{\mathrm{2}}}$选取指标$s_k\in U_k$。

步3 计算${\mathit{y}}_k={\mathit{x}}_k+\left({\mathit{b}}^{(s_k)}-{\mathit{A}}^{(s_k)}{\mathit{x}}_k\right){\left({\mathit{A}}^{(s_k)}\right)}^{\mathrm{T}}$， ${\widehat{\epsilon }}_k=\left(\mathrm{1}-\theta \right)\underset{\mathrm{1}\le j\le m}{\mathrm{m}\mathrm{a}\mathrm{x}}\left\{{\left|{\mathit{b}}^{(j)}-{\mathit{A}}^{(j)}{\mathit{y}}_k\right|}^{\mathrm{2}}\right\}\mathrm{。}$

步4 定义正整数指标集${\widehat{U}}_k=\left\{\mathit{r}||{\mathit{b}}^{(r)}-{\mathit{A}}^{(r)}{\mathit{y}}_k{|}^{\mathrm{2}}\ge {\widehat{\epsilon }}_k\right\}\mathrm{。}$计算 ${\widehat{\mathit{r}}}_{y_k}$ 的第$r$个分量${\widehat{r}}_{y_k}^{(r)}=\left\{\begin{array}{ll}{\mathit{b}}^{(r)}-{\mathit{A}}^{(r)}{\mathit{y}}_k,& \hspace{1em}r\in {\widehat{U}}_k\\ \mathrm{0},& \hspace{1em}r\notin {\widehat{U}}_k\end{array}\right.$

步5 根据概率$P(r=r_k)=\frac{{\left|{\widehat{\mathit{r}}}_{y_k}^{(r_k)}\right|}^{\mathrm{2}}}{{‖{\widehat{\mathit{r}}}_{y_k}‖}_{\mathrm{2}}^{\mathrm{2}}}$选取指标${\mathit{r}}_{\mathit{k}}\in {\widehat{U}}_k$。

步6 计算${\mu }_{r_k,s_k}=<{\mathit{A}}^{(r_k)},{\mathit{A}}^{(s_k)}>;$ ${\mathit{\nu }}_k=\frac{{\mathit{A}}^{(r_k)}-{\mu }_{r_k,s_k}{\mathit{A}}^{(s_k)}}{\sqrt[]{\mathrm{1}-{\left|{\mu }_{r_k,s_k}\right|}^{\mathrm{2}}}}$； ${\beta }_k=\frac{{\mathit{b}}^{(r_k)}-{\mu }_{r_k,s_k}{\mathit{b}}^{(s_k)}}{\sqrt[]{\mathrm{1}-{\left|{\mu }_{r_k,s_k}\right|}^{\mathrm{2}}}}\mathrm{。}$

步7 计算${\mathit{x}}_{k+\mathrm{1}}={\mathit{y}}_k+\left({\beta }_k-{\mathit{\nu }}_k{\mathit{y}}_k\right){\left({\mathit{\nu }}_k\right)}^{\mathrm{T}}\mathrm{。}$

步8 置$k=k+\mathrm{1}$， 转步1。

此外，由${\mathit{y}}_k$和${\mathit{x}}_{k+\mathrm{1}}$的递推式可以推得相对应的残量的递推公式如下：

$$\begin{array}{l}{\mathit{r}}_{y_k}=\mathit{b}-\mathit{A}{\mathit{y}}_k=\\ \mathit{b}-\mathit{A}\left({\mathit{x}}_k+\left({\mathit{b}}^{(s_k)}-{\mathit{A}}^{(s_k)}{\mathit{x}}_k\right){\left({\mathit{A}}^{(s_k)}\right)}^{\mathrm{T}}\right)=\\ \mathit{b}-\mathit{A}{\mathit{x}}_k-\left({\mathit{b}}^{(s_k)}-{\mathit{A}}^{(s_k)}{\mathit{x}}_k\right)\mathit{A}{\left({\mathit{A}}^{(s_k)}\right)}^{\mathrm{T}}=\\ {\mathit{r}}_{x_k}-{\mathit{r}}_{x_k}^{(s_k)}\mathit{A}{\left({\mathit{A}}^{(s_k)}\right)}^{\mathrm{T}}=\\ {\mathit{r}}_k-{\mathit{r}}_{x_k}^{(s_k)}{\mathit{B}}_{(s_k)}\end{array}$$

$$\begin{array}{l}{\mathit{r}}_{x_{k+\mathrm{1}}}=\mathit{b}-\mathit{A}{x}_{k+\mathrm{1}}=\\ \mathit{b}-\mathit{A}\left({\mathit{y}}_k+\left({\beta }_k-{\mathit{\nu }}_k{\mathit{y}}_k\right){\left({\mathit{\nu }}_k\right)}^{\mathrm{T}}\right)=\\ \mathit{b}-\mathit{A}{\mathit{y}}_k-\left({\beta }_k-{\mathit{\nu }}_k{\mathit{y}}_k\right)\mathit{A}{\left({\mathit{\nu }}_k\right)}^{\mathrm{T}}=\\ {\mathit{r}}_{y_k}-\frac{{\mathit{r}}_{{}^{y_k}}^{(r_k)}-{\mu }_k{\mathit{r}}_{{}^{y_k}}^{(s_k)}}{\sqrt[]{\mathrm{1}-{\left|{\mu }_k\right|}^{\mathrm{2}}}}\frac{\mathit{A}{\mathit{A}}_{(r_k)}^T-{\mu }_k\mathit{A}{\mathit{A}}_{(s_k)}^{\mathrm{T}}}{\sqrt[]{\mathrm{1}-{\left|{\mu }_k\right|}^{\mathrm{2}}}}=\\ {\mathit{r}}_{y_k}-\frac{{\mathit{r}}_{{}^{y_k}}^{(r_k)}-{\mu }_k{\mathit{r}}_{{}^{y_k}}^{(s_k)}}{\mathrm{1}-{\left|{\mu }_k\right|}^{\mathrm{2}}}\left({\mathit{B}}_{(r_k)}-{\mu }_k{\mathit{B}}_{(s_k)}\right)\end{array}$$

其中$\mathit{B}=\mathit{A}{\mathit{A}}^{\mathrm{T}}$且${\mathit{B}}_{(s_k)}$、${\mathit{B}}_{(r_k)}$表示矩阵$\mathit{B}$的第$s_k$、$r_k$列。若迭代前已知$\mathit{A}{\mathit{A}}^{\mathrm{T}}$， 则可进一步提高贪婪双子空间随机Kaczmarz方法的求解性能， 参见文献［24］。关于采用随机策略近似计算矩阵$\mathit{A}{\mathit{A}}^{\mathrm{T}}$， 参见文献［42］。

若每次迭代仅选择一个工作行，则算法3成为一类广义贪婪随机Kaczmarz（GRK（θ））方法。通过改变参数$\theta$，可得到一系列常用算法作为特殊情况。取$\theta =\mathrm{0}$，$\theta =\frac{\mathrm{1}}{\mathrm{2}}$，$\theta =\frac{\mathrm{1}}{\mathrm{2}}\left(\mathrm{1}-\frac{‖\mathit{b}-\mathit{A}{\mathit{x}}_k‖}{m\underset{\mathrm{1}\le i\le m}{\mathrm{m}\mathrm{a}\mathrm{x}}\left\{{\left|{\mathit{b}}^{(i)}-{\mathit{A}}^{(i)}{\mathit{x}}_k\right|}^{\mathrm{2}}\right\}}\right)$，$\theta =(\mathrm{1}-\vartheta )\left(\mathrm{1}-\frac{‖\mathit{b}-\mathit{A}{\mathit{x}}_k‖}{m\underset{\mathrm{1}\le i\le m}{\mathrm{m}\mathrm{a}\mathrm{x}}\left\{{\left|{\mathit{b}}^{(i)}-{\mathit{A}}^{(i)}{\mathit{x}}_k\right|}^{\mathrm{2}}\right\}}\right)$ $\left(\vartheta \in [\mathrm{0,1}]\right)$时，GRK（θ）方法分别简化为贪婪Kaczmarz方法^［41］，贪婪距离随机Kaczmarz方法^［25］，贪婪随机Kaczmarz方法^［24］和松弛贪婪随机Kaczmarz方法^［33］。

2 贪婪双子空间随机Kaczmarz方法的收敛性分析

相容线性系统（1）的最小范数解${\mathit{x}}_{\mathrm{*}}$形如

$${\mathit{x}}_{\mathrm{*}}:=\mathrm{a}\mathrm{r}\mathrm{g}\underset{\mathbf{A}\mathbf{x}=\mathit{b}}{\mathrm{m}\mathrm{i}\mathrm{n}}{‖\mathit{x}‖}_{\mathrm{2}}\triangleq {\mathit{A}}^{†}\mathit{b}\in \mathrm{R}\mathrm{a}\mathrm{n}({\mathit{A}}^{\mathrm{T}})$$

其中 $\mathrm{R}\mathrm{a}\mathrm{n}({\mathit{A}}^{\mathrm{T}})$表示${\mathit{A}}^{\mathrm{T}}$的列空间。若初始向量${\mathit{x}}_{\mathrm{0}}\in \mathrm{R}\mathrm{a}\mathrm{n}({\mathit{A}}^{\mathrm{T}})$， 由算法3的步5和步13可知， 算法3生成的迭代序列${\left\{{\mathit{x}}_k\right\}}_{k=\mathrm{0}}^{\mathrm{\infty }}$一定在$\mathrm{R}\mathrm{a}\mathrm{n}({\mathit{A}}^{\mathrm{T}})$中， 故若算法3收敛， 则一定收敛到相容线性系统（1）的最小范数解。

关于算法3特殊情况（每次迭代仅选一个工作行）——广义贪婪随机Kaczmarz方法， 有如下收敛性定理。

定理 3  　若线性系统（1）相容， 其中系数矩阵$\mathit{A}\in {\mathbf{R}}^{m\times n}$且右端项$\mathit{b}\in {\mathbf{R}}^m$。初始向量${\mathit{x}}_{\mathrm{0}}\in \mathrm{R}\mathrm{a}\mathrm{n}({\mathit{A}}^{\mathrm{T}})$， 令${\mathit{x}}_{k+\mathrm{1}}$为通过广义贪婪随机Kaczmarz方法生成的第$(k+\mathrm{1})$个迭代值，则有

$$E{‖{\mathit{x}}_{\mathrm{*}}-{\mathit{x}}_{k+\mathrm{1}}‖}_{\mathrm{2}}^{\mathrm{2}}\le \left(\mathrm{1}-\frac{\mathrm{1}-\theta }{m-\mathrm{1}}{\sigma }_{{}_{\mathrm{m}\mathrm{i}\mathrm{n}}}^{\mathrm{2}}(\mathit{A})\right){‖{\mathit{x}}_{\mathrm{*}}-{\mathit{x}}_k‖}_{\mathrm{2}}^{\mathrm{2}}\hspace{1em}$$

（8）

和

$$E{‖{\mathit{x}}_{\mathrm{*}}-{\mathit{x}}_{k+\mathrm{1}}‖}_{\mathrm{2}}^{\mathrm{2}}\le {\left(\mathrm{1}-\frac{\mathrm{1}-\theta }{m-\mathrm{1}}{\sigma }_{{}_{\mathrm{m}\mathrm{i}\mathrm{n}}}^{\mathrm{2}}(\mathit{A})\right)}^{k+\mathrm{1}}{‖{\mathit{x}}_{\mathrm{*}}-{\mathit{x}}_{\mathrm{0}}‖}_{\mathrm{2}}^{\mathrm{2}},\hspace{1em}k=\mathrm{0,1},\mathrm{2}\cdots$$

证明：　由广义贪婪随机Kaczmarz方法（算法3的前5步）可得

$${\mathit{x}}_{k+\mathrm{1}}={\mathit{x}}_k+\left({\mathit{b}}^{(s_k)}-{\mathit{A}}^{(s_k)}{\mathit{x}}_k\right){\left({\mathit{A}}^{(s_k)}\right)}^{\mathrm{T}}$$

根据概率$P(r=s_k)=\frac{{\left|{\widehat{\mathit{r}}}_{x_k}^{(s_k)}\right|}^{\mathrm{2}}}{{‖{\widehat{\mathit{r}}}_{x_k}‖}_{\mathrm{2}}^{\mathrm{2}}}$选指标$s_k\in U_k$。由${\mathit{A}}^{(s_k)}\left({\mathit{x}}_{\mathrm{*}}-{\mathit{x}}_{k+\mathrm{1}}\right)=\mathrm{0}$，即${\mathit{x}}_{\mathrm{*}}-{\mathit{x}}_k$与${\mathit{x}}_{k+\mathrm{1}}-{\mathit{x}}_k$正交，可得

$${‖{\mathit{x}}_{\mathrm{*}}-{\mathit{x}}_{k+\mathrm{1}}‖}_{\mathrm{2}}^{\mathrm{2}}={‖{\mathit{x}}_{\mathrm{*}}-{\mathit{x}}_k‖}_{\mathrm{2}}^{\mathrm{2}}-{‖{\mathit{x}}_k-{\mathit{x}}_{k+\mathrm{1}}‖}_{\mathrm{2}}^{\mathrm{2}}$$

（9）

对等式（9）两侧取期望， 有

$$\begin{array}{l}E{‖{\mathit{x}}_{\mathrm{*}}-{\mathit{x}}_{k+\mathrm{1}}‖}_{\mathrm{2}}^{\mathrm{2}}=\\ {‖{\mathit{x}}_{\mathrm{*}}-{\mathit{x}}_k‖}_{\mathrm{2}}^{\mathrm{2}}-E{‖{\mathit{x}}_k-{\mathit{x}}_{k+\mathrm{1}}‖}_{\mathrm{2}}^{\mathrm{2}}=\\ {‖{\mathit{x}}_{\mathrm{*}}-{\mathit{x}}_k‖}_{\mathrm{2}}^{\mathrm{2}}-\sum _{s\in U_k}\frac{{\left|{\widehat{\mathit{r}}}_{x_k}^{(s_k)}\right|}^{\mathrm{2}}}{{‖{\widehat{\mathit{r}}}_{x_k}‖}_{\mathrm{2}}^{\mathrm{2}}}{‖\left({\mathit{b}}^{(s_k)}-{\mathit{A}}^{(s_k)}{\mathit{x}}_k\right){\left({\mathit{A}}^{(s_k)}\right)}^{\mathrm{T}}‖}_{\mathrm{2}}^{\mathrm{2}}=\\ {‖{\mathit{x}}_{\mathrm{*}}-{\mathit{x}}_k‖}_{\mathrm{2}}^{\mathrm{2}}-\sum _{s\in U_k}\frac{{\left|{\widehat{\mathit{r}}}_{x_k}^{(s_k)}\right|}^{\mathrm{2}}}{{‖{\widehat{\mathit{r}}}_{x_k}‖}_{\mathrm{2}}^{\mathrm{2}}}{\left|{\mathit{b}}^{(s_k)}-{\mathit{A}}^{(s_k)}{\mathit{x}}_k\right|}^{\mathrm{2}}\end{array}$$

（10）

由${\mathit{r}}_{x_{k+\mathrm{1}}}^{(s_k)}={\mathit{b}}^{(s_k)}-{\mathit{A}}^{(s_k)}{\mathit{x}}_{k+\mathrm{1}}=\mathrm{0},k=\mathrm{0,1},\mathrm{2}\cdots$，则${‖\mathit{b}-\mathit{A}{\mathit{x}}_{k+\mathrm{1}}‖}_{\mathrm{2}}^{\mathrm{2}}={\sum _{i\ne s_k}\left|{\mathit{b}}^{(i)}-{\mathit{A}}^{(i)}{\mathit{x}}_{k+\mathrm{1}}\right|}^{\mathrm{2}}$。类似不等式（6），对于$\forall r\in {\widehat{U}}_k$，可以推得

$${\left|{\mathit{r}}_{x_k}^{(s)}\right|}^{\mathrm{2}}={\left|{\mathit{b}}^{(s)}-{\mathit{A}}^{(s)}{\mathit{x}}_k\right|}^{\mathrm{2}}\ge \frac{\mathrm{1}-\theta }{m-\mathrm{1}}{\sigma }_{{}_{\mathrm{m}\mathrm{i}\mathrm{n}}}^{\mathrm{2}}(\mathit{A}){‖{\mathit{x}}_{\mathrm{*}}-{\mathit{x}}_k‖}_{\mathrm{2}}^{\mathrm{2}}$$

（11）

其中式（11）最后一个不等式需要用到不等式（12），即对于任意$\mathit{u}\in \mathrm{R}\mathrm{a}\mathrm{n}({\mathit{A}}^{\mathrm{T}})$， 成立

$${‖\mathit{A}\mathit{u}‖}_{\mathrm{2}}^{\mathrm{2}}\ge {\sigma }_{{}_{\mathrm{m}\mathrm{i}\mathrm{n}}}^{\mathrm{2}}(\mathit{A}){‖\mathit{u}‖}_{\mathrm{2}}^{\mathrm{2}}$$

（12）

联立等式（11）和不等式（12）可得

$$\begin{array}{l}E{‖{\mathit{x}}_{\mathrm{*}}-{\mathit{x}}_{k+\mathrm{1}}‖}_{\mathrm{2}}^{\mathrm{2}}={‖{\mathit{x}}_{\mathrm{*}}-{\mathit{x}}_k‖}_{\mathrm{2}}^{\mathrm{2}}-\\ \sum _{s\in U_k}\frac{{\left|{\widehat{\mathit{r}}}_{x_k}^{(s)}\right|}^{\mathrm{2}}}{{‖{\widehat{\mathit{r}}}_{x_k}‖}_{\mathrm{2}}^{\mathrm{2}}}{\left|{\mathit{b}}^{(s)}-{\mathit{A}}^{(s)}{\mathit{x}}_k\right|}^{\mathrm{2}}\le \\ {‖{\mathit{x}}_{\mathrm{*}}-{\mathit{x}}_k‖}_{\mathrm{2}}^{\mathrm{2}}-\frac{\mathrm{1}-\theta }{m-\mathrm{1}}{\sigma }_{{}_{\mathrm{m}\mathrm{i}\mathrm{n}}}^{\mathrm{2}}(\mathit{A})‖{\mathit{x}}_{\mathrm{*}}-\\ {{\mathit{x}}_k‖}_{\mathrm{2}}^{\mathrm{2}}\sum _{s\in U_k}\frac{{\left|{\mathit{b}}^{(s)}-{\mathit{A}}^{(s)}{\mathit{x}}_k\right|}^{\mathrm{2}}}{{‖{\widehat{\mathit{r}}}_{x_k}‖}_{\mathrm{2}}^{\mathrm{2}}}=\\ {‖{\mathit{x}}_{\mathrm{*}}-{\mathit{x}}_k‖}_{\mathrm{2}}^{\mathrm{2}}-\frac{\mathrm{1}-\theta }{m-\mathrm{1}}{\sigma }_{{}_{\mathrm{m}\mathrm{i}\mathrm{n}}}^{\mathrm{2}}(\mathit{A}){‖{\mathit{x}}_{\mathrm{*}}-{\mathit{x}}_k‖}_{\mathrm{2}}^{\mathrm{2}}=\\ \left(\mathrm{1}-\frac{\mathrm{1}-\theta }{m-\mathrm{1}}{\sigma }_{{}_{\mathrm{m}\mathrm{i}\mathrm{n}}}^{\mathrm{2}}(\mathit{A})\right){‖{\mathit{x}}_{\mathrm{*}}-{\mathit{x}}_k‖}_{\mathrm{2}}^{\mathrm{2}}\end{array}$$

通过归纳法，可得

$$\begin{array}{l}E{‖{\mathit{x}}_{\mathrm{*}}-{\mathit{x}}_{k+\mathrm{1}}‖}_{\mathrm{2}}^{\mathrm{2}}\le {\left(\mathrm{1}-\frac{\mathrm{1}-\theta }{m-\mathrm{1}}{\sigma }_{{}_{\mathrm{m}\mathrm{i}\mathrm{n}}}^{\mathrm{2}}(\mathit{A})\right)}^{k+\mathrm{1}}{‖{\mathit{x}}_{\mathrm{*}}-{\mathit{x}}_{\mathrm{0}}‖}_{\mathrm{2}}^{\mathrm{2}}\\ K=\mathrm{0,1},\mathrm{2}\cdots \end{array}$$

关于贪婪双子空间随机Kaczmarz方法， 有如下的收敛性定理。

定理4  　若线性系统（1）相容， 其中系数矩阵$\mathit{A}\in {\mathbf{R}}^{m\times n}$且右端项$\mathit{b}\in {\mathbf{R}}^m$。初始向量${\mathit{x}}_{\mathrm{0}}\in \mathrm{R}\mathrm{a}\mathrm{n}({\mathit{A}}^{\mathrm{T}})$， 令${\mathit{x}}_{k+\mathrm{1}}$为通过贪婪双子空间随机Kaczmarz方法生成的第$(k+\mathrm{1})$个迭代值，则

$$\begin{array}{l}E{‖{\mathit{x}}_{\mathrm{*}}-{\mathit{x}}_{k+\mathrm{1}}‖}_{\mathrm{2}}^{\mathrm{2}}\le \\ M(\theta ){‖{\mathit{x}}_{\mathrm{*}}-{\mathit{x}}_k‖}_{\mathrm{2}}^{\mathrm{2}}-\\ {\displaystyle \sum _{r\in {\widehat{U}}_k}^{}}\frac{{\left|{\widehat{\mathit{r}}}_{y_k}^{(r)}\right|}^{\mathrm{2}}}{{‖{\widehat{\mathit{r}}}_{y_k}‖}_{\mathrm{2}}^{\mathrm{2}}}{\displaystyle \sum _{s\in U_k}^{}}\frac{{\left|{\widehat{\mathit{r}}}_{x_k}^{(s)}\right|}^{\mathrm{2}}}{{‖{\widehat{\mathit{r}}}_{x_k}‖}_{\mathrm{2}}^{\mathrm{2}}}{\left|{\mu }_{r,s}\right|}^{\mathrm{2}}{\left|<{\mathit{e}}_k^{\mathrm{T}},{\mathit{\nu }}_{r,s}>\right|}^{\mathrm{2}}\end{array}$$

（13）

这里，$M(\theta )=\left(\mathrm{1}-\frac{\mathrm{1}-\theta }{m-\mathrm{1}}{\sigma }_{{}_{\mathrm{m}\mathrm{i}\mathrm{n}}}^{\mathrm{2}}(\mathit{A})\right)\left(\mathrm{1}-\frac{\mathrm{1}-\theta }{m-\mathrm{2}}{\sigma }_{{}_{\mathrm{m}\mathrm{i}\mathrm{n}}}^{\mathrm{2}}(\mathit{A})\right)$，${\mathit{e}}_k={\mathit{x}}_{\mathrm{*}}-{\mathit{x}}_k$，${\mu }_{r,s}=<{\mathit{A}}^{(r)},{\mathit{A}}^{(s)}>$，${\mathit{\nu }}_{r,s}=\frac{{\mathit{A}}^{(r)}-{\mu }_k{\mathit{A}}^{(s)}}{\sqrt[]{\mathrm{1}-{\left|{\mu }_{r,s}\right|}^{\mathrm{2}}}}$。进一步放缩可得一个松弛但更简单的收敛界

$$E{‖{\mathit{x}}_{\mathrm{*}}-{\mathit{x}}_{k+\mathrm{1}}‖}_{\mathrm{2}}^{\mathrm{2}}\le M_k(\theta ){‖{\mathit{x}}_{\mathrm{*}}-{\mathit{x}}_k‖}_{\mathrm{2}}^{\mathrm{2}}$$

其中$M_k(\theta )=\left(\mathrm{1}-\frac{(\mathrm{1}+D_k)(\mathrm{1}-\theta )}{m-\mathrm{1}}{\sigma }_{{}_{\mathrm{m}\mathrm{i}\mathrm{n}}}^{\mathrm{2}}(\mathit{A})\right)\left(\mathrm{1}-\frac{\mathrm{1}-\theta }{m-\mathrm{2}}{\sigma }_{{}_{\mathrm{m}\mathrm{i}\mathrm{n}}}^{\mathrm{2}}(\mathit{A})\right)$，$D_k=\underset{r\in {\widehat{U}}_k,s\in U_k}{\mathrm{m}\mathrm{i}\mathrm{n}}\left\{\frac{{\left|{\mu }_{r,s}\right|}^{\mathrm{2}}}{\mathrm{1}-{\left|{\mu }_{r,s}\right|}^{\mathrm{2}}}\right\}$。

证明：　由算法3的步5和步13可得

$$\begin{array}{l}{\mathit{y}}_k={\mathit{x}}_k+\left({\mathit{b}}^{(s_k)}-{\mathit{A}}^{(s_k)}{\mathit{x}}_k\right){\left({\mathit{A}}^{(s_k)}\right)}^{\mathrm{T}}=\\ {\mathit{x}}_k+<{\mathit{e}}_k^{\mathrm{T}},{\mathit{A}}^{(s_k)}>{\left({\mathit{A}}^{(s_k)}\right)}^{\mathrm{T}}\\ {\mathit{x}}_{k+\mathrm{1}}={\mathit{y}}_k+\left({\beta }_k-{\mathit{\nu }}_k{\mathit{y}}_k\right){\left({\mathit{\nu }}_k\right)}^{\mathrm{T}}=\\ {\mathit{x}}_k+<{\mathit{e}}_k^{\mathrm{T}},{\mathit{A}}^{(s_k)}>{\left({\mathit{A}}^{(s_k)}\right)}^{\mathrm{T}}+<{\mathit{e}}_k^{\mathrm{T}},{\mathit{\nu }}_k>{\left({\mathit{\nu }}_k\right)}^{\mathrm{T}}\end{array}$$

记${\xi }_k=<{\mathit{A}}^{(r_k)},{\mathit{\nu }}_k>$。为了便利，在后续的证明中，${\mathit{A}}^{(r_k)}$、${\mathit{A}}^{(s_k)}$、${\mu }_k$、${\mathit{\nu }}_k$、${\xi }_k$分别简记为${\mathit{A}}^{(r)}$、${\mathit{A}}^{(s)}$、$\mu$、$\mathit{\nu }$、$\xi$。由于$<{\mathit{A}}^{(s)},\mathit{\nu }>=\mathrm{0}$，即向量$\mathit{\nu }$正交于向量${\mathit{A}}^{(s)}$，展开${‖{\mathit{x}}_{\mathrm{*}}-{\mathit{x}}_{k+\mathrm{1}}‖}_{\mathrm{2}}^{\mathrm{2}}$可得

$${‖{\mathit{x}}_{\mathrm{*}}-{\mathit{x}}_{k+\mathrm{1}}‖}_{\mathrm{2}}^{\mathrm{2}}={‖{\mathit{e}}_k‖}_{\mathrm{2}}^{\mathrm{2}}-{\left|<{\mathit{e}}_k^{\mathrm{T}},{\mathit{A}}^{(s)}>\right|}^{\mathrm{2}}-{\left|<{\mathit{e}}_k^{\mathrm{T}},\mathit{\nu }>\right|}^{\mathrm{2}}$$

（14）

在算法3中，一次迭代需执行两次投影。因此，将式（14）中的误差与广义贪婪随机Kaczmarz 方法两次迭代后获得的误差进行比较。设$\mathit{z}$是广义贪婪随机Kaczmarz 方法基于${\mathit{y}}_k$的下一次迭代值， 则有

$$\begin{array}{l}\mathit{z}={\mathit{y}}_k+\left({\mathit{b}}^{(r)}-{\mathit{A}}^{(r)}{\mathit{y}}_k\right){\left({\mathit{A}}^{(r)}\right)}^{\mathrm{T}}=\\ {\mathit{x}}_k+<{\mathit{e}}_k^{\mathrm{T}},{\mathit{A}}^{(s)}>{\left({\mathit{A}}^{(s)}\right)}^{\mathrm{T}}+\\ <{\mathit{x}}_{\mathrm{*}}^{\mathrm{T}}-{\mathit{y}}_k^{\mathrm{T}},{\mathit{A}}^{(r)}>{\left({\mathit{A}}^{(r)}\right)}^{\mathrm{T}}\end{array}$$

（15）

由${\mathit{y}}_k$的递推式可得

$$\begin{array}{l}<{\mathit{x}}_{\mathrm{*}}^{\mathrm{T}}-{\mathit{y}}_k^{\mathrm{T}},{\mathit{A}}^{(r)}>=\\ <{\mathit{e}}_k^{\mathrm{T}},{\mathit{A}}^{(r)}>-<{\mathit{e}}_k^{\mathrm{T}},{\mathit{A}}^{(s)}><{\mathit{A}}^{(s)},{\mathit{A}}^{(r)}>=\\ <{\mathit{e}}_k^{\mathrm{T}},{\mathit{A}}^{(r)}>-\mu <{\mathit{e}}_k^{\mathrm{T}},{\mathit{A}}^{(s)}>=\\ <{\mathit{e}}_k^{\mathrm{T}},{\mathit{A}}^{(r)}-\mu {\mathit{A}}^{(s)}>\end{array}$$

（16）

由$\mu$、$\mathit{\nu }$、$\xi$的定义可得

$${\mathit{A}}^{(r)}=\mu {\mathit{A}}^{(r)}+\mathit{\nu },{\xi }^{\mathrm{2}}+{\mu }^{\mathrm{2}}=\mathrm{1}$$

（17）

联立式（15）~（17）可得

$$\begin{array}{l}\mathit{z}={\mathit{x}}_k+<{\mathit{e}}_k^{\mathrm{T}},{\mathit{A}}^{(s)}>{\left({\mathit{A}}^{(s)}\right)}^{\mathrm{T}}+\\ <{\mathit{e}}_k^{\mathrm{T}},{\mathit{A}}^{(r)}-\mu {\mathit{A}}^{(s)}>{\left({\mathit{A}}^{(r)}\right)}^{\mathrm{T}}=\\ {\mathit{x}}_k+<{\mathit{e}}_k^{\mathrm{T}},{\mathit{A}}^{(s)}>{\left({\mathit{A}}^{(s)}\right)}^{\mathrm{T}}+\\ <{\mathit{e}}_k^{\mathrm{T}},\xi \mathit{\nu }>{\left(\mu {\mathit{A}}^{(s)}+\xi \mathit{\nu }\right)}^{\mathrm{T}}=\\ {\mathit{x}}_k+<{\mathit{e}}_k^{\mathrm{T}},{\mathit{A}}^{(s)}>{\left({\mathit{A}}^{(s)}\right)}^{\mathrm{T}}+\\ \xi \mu <{\mathit{e}}_k^{\mathrm{T}},\mathit{\nu }>{\left({\mathit{A}}^{(s)}\right)}^{\mathrm{T}}+\\ {\xi }^{\mathrm{2}}<{\mathit{e}}_k^{\mathrm{T}},\mathit{\nu }>{\left(\mathit{\nu }\right)}^{\mathrm{T}}=\\ {\mathit{x}}_k+<{\mathit{e}}_k^{\mathrm{T}},{\mathit{A}}^{(s)}>{\left({\mathit{A}}^{(s)}\right)}^{\mathrm{T}}+\\ \xi \mu <{\mathit{e}}_k^{\mathrm{T}},\mathit{\nu }>{\left({\mathit{A}}^{(s)}\right)}^{\mathrm{T}}+<{\mathit{e}}_k^{\mathrm{T}},\mathit{\nu }>{\left(\mathit{\nu }\right)}^{\mathrm{T}}-\\ {\mu }^{\mathrm{2}}<{\mathit{e}}_k^{\mathrm{T}},\mathit{\nu }>{\left(\mathit{\nu }\right)}^{\mathrm{T}}\end{array}$$

利用向量$\mathit{\nu }$与向量${\mathit{A}}^{(s)}$的正交性，展开${‖{\mathit{x}}_{\mathrm{*}}-\mathit{z}‖}_{\mathrm{2}}^{\mathrm{2}}$可得

$$\begin{array}{l}{‖{\mathit{x}}_{\mathrm{*}}-\mathit{z}‖}_{\mathrm{2}}^{\mathrm{2}}=‖{\mathit{e}}_k-<{\mathit{e}}_k^{\mathrm{T}},{\mathit{A}}^{(s)}>{\left({\mathit{A}}^{(s)}\right)}^{\mathrm{T}}-\\ \xi \mu <{\mathit{e}}_k^{\mathrm{T}},\mathit{\nu }>{\left({\mathit{A}}^{(s)}\right)}^{\mathrm{T}}-\\ {{\left|<{\mathit{e}}_k^{\mathrm{T}},\mathit{\nu }>\right|}^{\mathrm{2}}{\mathit{\nu }}^{\mathrm{T}}+{\mu }^{\mathrm{2}}<{\mathit{e}}_k^{\mathrm{T}},\mathit{\nu }>{\mathit{\nu }}^{\mathrm{T}}‖}_{\mathrm{2}}^{\mathrm{2}}=\\ {‖{\mathit{e}}_k‖}_{\mathrm{2}}^{\mathrm{2}}-{\left|<{\mathit{e}}_k^{\mathrm{T}},{\mathit{A}}^{(s)}>\right|}^{\mathrm{2}}-\\ {\left|<{\mathit{e}}_k^{\mathrm{T}},\mathit{\nu }>\right|}^{\mathrm{2}}+{\mu }^{\mathrm{2}}{\left|<{\mathit{e}}_k^{\mathrm{T}},\mathit{\nu }>\right|}^{\mathrm{2}}\end{array}$$

（18）

联立式（14）和式（18）可得

$${‖{\mathit{x}}_{\mathrm{*}}-{\mathit{x}}_{k+\mathrm{1}}‖}_{\mathrm{2}}^{\mathrm{2}}={‖{\mathit{x}}_{\mathrm{*}}-\mathit{z}‖}_{\mathrm{2}}^{\mathrm{2}}-{\left|\mu <{\mathit{e}}_k^{\mathrm{T}},\mathit{\nu }>\right|}^{\mathrm{2}}$$

（19）

对式（19）两端取期望可得

$$E{‖{\mathit{x}}_{\mathrm{*}}-{\mathit{x}}_{k+\mathrm{1}}‖}_{\mathrm{2}}^{\mathrm{2}}=E{‖{\mathit{x}}_{\mathrm{*}}-\mathit{z}‖}_{\mathrm{2}}^{\mathrm{2}}-E{\left|\mu <{\mathit{e}}_k^{\mathrm{T}},\mathit{\nu }>\right|}^{\mathrm{2}}$$

（20）

由定理3可得

$$E{‖{\mathit{x}}_{\mathrm{*}}-\mathit{z}‖}_{\mathrm{2}}^{\mathrm{2}}\le \left(\mathrm{1}-\frac{\mathrm{1}-\theta }{m-\mathrm{1}}{\sigma }_{{}_{\mathrm{m}\mathrm{i}\mathrm{n}}}^{\mathrm{2}}(\mathit{A})\right)E{‖{\mathit{x}}_{\mathrm{*}}-{\mathit{y}}_k‖}_{\mathrm{2}}^{\mathrm{2}}$$

（21）

再由

$$\begin{array}{l}E{‖{\mathit{x}}_{\mathrm{*}}-{\mathit{y}}_k‖}_{\mathrm{2}}^{\mathrm{2}}=\\ {‖{\mathit{x}}_{\mathrm{*}}-{\mathit{x}}_k‖}_{\mathrm{2}}^{\mathrm{2}}-E{\left|{\mathit{b}}^{(s_k)}-{\mathit{A}}^{(s_k)}{\mathit{x}}_k\right|}^{\mathrm{2}}=\\ {‖{\mathit{x}}_{\mathrm{*}}-{\mathit{x}}_k‖}_{\mathrm{2}}^{\mathrm{2}}-\sum _{s\in U_k}\frac{{\left|{\widehat{\mathit{r}}}_{x_k}^{(s)}\right|}^{\mathrm{2}}}{{‖{\widehat{\mathit{r}}}_{x_k}‖}_{\mathrm{2}}^{\mathrm{2}}}{\left|{\mathit{b}}^{(s)}-{\mathit{A}}^{(s)}{\mathit{x}}_k\right|}^{\mathrm{2}}\le \\ {‖{\mathit{x}}_{\mathrm{*}}-{\mathit{x}}_k‖}_{\mathrm{2}}^{\mathrm{2}}-\frac{\mathrm{1}-\theta }{m-\mathrm{2}}{‖\mathit{b}-\mathit{A}{\mathit{x}}_k‖}_{\mathrm{2}}^{\mathrm{2}}(\mathrm{由}\mathrm{等}\mathrm{式}(\mathrm{7}))\le \\ \left(\mathrm{1}-\frac{\mathrm{1}-\theta }{m-\mathrm{2}}{\sigma }_{{}_{\mathrm{m}\mathrm{i}\mathrm{n}}}^{\mathrm{2}}(\mathit{A})\right){‖{\mathit{x}}_{\mathrm{*}}-{\mathit{x}}_k‖}_{\mathrm{2}}^{\mathrm{2}}\end{array}$$

因此，可得

$$E{‖{\mathit{x}}_{\mathrm{*}}-\mathit{z}‖}_{\mathrm{2}}^{\mathrm{2}}\le \left(\mathrm{1}-\frac{\mathrm{1}-\theta }{m-\mathrm{1}}{\sigma }_{{}_{\mathrm{m}\mathrm{i}\mathrm{n}}}^{\mathrm{2}}(\mathit{A})\right)\left(\mathrm{1}-\frac{\mathrm{1}-\theta }{m-\mathrm{2}}{\sigma }_{{}_{\mathrm{m}\mathrm{i}\mathrm{n}}}^{\mathrm{2}}(\mathit{A})\right){‖{\mathit{x}}_{\mathrm{*}}-{\mathit{x}}_k‖}_{\mathrm{2}}^{\mathrm{2}}$$

（22）

再由期望的定义可得

$$E{\left|\mu <{\mathit{e}}_k^{\mathrm{T}},\mathit{\nu }>\right|}^{\mathrm{2}}={\displaystyle \sum _{r\in {\widehat{U}}_k}^{}}\frac{{\left|{\widehat{\mathit{r}}}_{y_k}^{(r)}\right|}^{\mathrm{2}}}{{‖{\widehat{\mathit{r}}}_{y_k}‖}_{\mathrm{2}}^{\mathrm{2}}}{\displaystyle \sum _{s\in U_k}^{}}\frac{{\left|{\widehat{\mathit{r}}}_{x_k}^{(s)}\right|}^{\mathrm{2}}}{{‖{\widehat{\mathit{r}}}_{x_k}‖}_{\mathrm{2}}^{\mathrm{2}}}{\left|{\mu }_{r,s}\right|}^{\mathrm{2}}{\left|<{\mathit{e}}_k^{\mathrm{T}},{\mathit{\nu }}_{r,s}>\right|}^{\mathrm{2}}$$

（23）

把式（22）和式（23）带入式（20）可得式（13）。

把式（16）和式（17）带入式（23），同时再结合式（6）可得

$$\begin{array}{l}E{\left|\mu <{\mathit{e}}_k^{\mathrm{T}},\mathit{\nu }>\right|}^{\mathrm{2}}=\\ {\displaystyle \sum _{r\in {\widehat{U}}_k}^{}}\frac{{\left|{\widehat{\mathit{r}}}_{y_k}^{(r)}\right|}^{\mathrm{2}}}{{‖{\widehat{\mathit{r}}}_{y_k}‖}_{\mathrm{2}}^{\mathrm{2}}}{\displaystyle \sum _{s\in U_k}^{}}\frac{{\left|{\widehat{\mathit{r}}}_{x_k}^{(s)}\right|}^{\mathrm{2}}}{{‖{\widehat{\mathit{r}}}_{x_k}‖}_{\mathrm{2}}^{\mathrm{2}}}{\left|{\mu }_{r,s}\right|}^{\mathrm{2}}{\left|<{\mathit{e}}_k^{\mathrm{T}},{\mathit{\nu }}_{r,s}>\right|}^{\mathrm{2}}=\\ {\displaystyle \sum _{r\in {\widehat{U}}_k}^{}}\frac{{\left|{\widehat{\mathit{r}}}_{y_k}^{(r)}\right|}^{\mathrm{2}}}{{‖{\widehat{\mathit{r}}}_{y_k}‖}_{\mathrm{2}}^{\mathrm{2}}}\frac{{\left|{\mu }_{r,s}\right|}^{\mathrm{2}}}{\mathrm{1}-{\left|{\mu }_{r,s}\right|}^{\mathrm{2}}}{\left|<{\mathit{x}}_{\mathrm{*}}^{\mathrm{T}}-{\mathit{y}}_k^{\mathrm{T}},{\mathit{A}}^{(r)}>\right|}^{\mathrm{2}}\ge \\ D_k{\displaystyle \sum _{r\in {\widehat{U}}_k}^{}}\frac{{\left|{\widehat{\mathit{r}}}_{y_k}^{(r)}\right|}^{\mathrm{2}}}{{‖{\widehat{\mathit{r}}}_{y_k}‖}_{\mathrm{2}}^{\mathrm{2}}}{\left|<{\mathit{x}}_{\mathrm{*}}^{\mathrm{T}}-{\mathit{y}}_k^{\mathrm{T}},{\mathit{A}}^{(r)}>\right|}^{\mathrm{2}}\ge \\ D_k\frac{\mathrm{1}-\theta }{m-\mathrm{1}}{\sigma }_{{}_{\mathrm{m}\mathrm{i}\mathrm{n}}}^{\mathrm{2}}(\mathit{A}){‖{\mathit{x}}_{\mathrm{*}}-{\mathit{y}}_k‖}_{\mathrm{2}}^{\mathrm{2}}\end{array}$$

（24）

其中$D_k=\underset{r\in {\widehat{U}}_k,s\in U_k}{\mathrm{m}\mathrm{i}\mathrm{n}}\left\{\frac{{\left|{\mu }_{r,s}\right|}^{\mathrm{2}}}{\mathrm{1}-{\left|{\mu }_{r,s}\right|}^{\mathrm{2}}}\right\}$。把式（21）和式（24）带入式（20）可得$\begin{array}{l}E{‖{\mathit{x}}_{\mathrm{*}}-{\mathit{x}}_{k+\mathrm{1}}‖}_{\mathrm{2}}^{\mathrm{2}}=\\ \left(\mathrm{1}-\frac{\mathrm{1}-\theta }{m-\mathrm{1}}{\sigma }_{{}_{\mathrm{m}\mathrm{i}\mathrm{n}}}^{\mathrm{2}}(\mathit{A})\right)‖{\mathit{x}}_{\mathrm{*}}-\end{array}$

$${{\mathit{y}}_k‖}_{\mathrm{2}}^{\mathrm{2}}-D_k\frac{\mathrm{1}-\theta }{m-\mathrm{1}}{\sigma }_{{}_{\mathrm{m}\mathrm{i}\mathrm{n}}}^{\mathrm{2}}(\mathit{A}){‖{\mathit{x}}_{\mathrm{*}}-{\mathit{y}}_k‖}_{\mathrm{2}}^{\mathrm{2}}=$$

$$\begin{array}{l}\left(\mathrm{1}-\frac{(\mathrm{1}+D_k)(\mathrm{1}-\theta )}{m-\mathrm{1}}{\sigma }_{{}_{\mathrm{m}\mathrm{i}\mathrm{n}}}^{\mathrm{2}}(\mathit{A})\right){‖{\mathit{x}}_{\mathrm{*}}-{\mathit{y}}_k‖}_{\mathrm{2}}^{\mathrm{2}}=\\ \left(\mathrm{1}-\frac{(\mathrm{1}+D_k)(\mathrm{1}-\theta )}{m-\mathrm{1}}{\sigma }_{{}_{\mathrm{m}\mathrm{i}\mathrm{n}}}^{\mathrm{2}}(\mathit{A})\right)E{‖{\mathit{x}}_{\mathrm{*}}-{\mathit{y}}_k‖}_{\mathrm{2}}^{\mathrm{2}}\le \\ \left(\mathrm{1}-\frac{(\mathrm{1}+D_k)(\mathrm{1}-\theta )}{m-\mathrm{1}}{\sigma }_{{}_{\mathrm{m}\mathrm{i}\mathrm{n}}}^{\mathrm{2}}(\mathit{A})\right)\\ \left(\mathrm{1}-\frac{\mathrm{1}-\theta }{m-\mathrm{2}}{\sigma }_{{}_{\mathrm{m}\mathrm{i}\mathrm{n}}}^{\mathrm{2}}(\mathit{A})\right){‖{\mathit{x}}_{\mathrm{*}}-x_k‖}_{\mathrm{2}}^{\mathrm{2}}\end{array}$$

收敛因子越小， 方法收敛得越快。通过定理4收敛界可知，参数$\theta$越小，收敛因子$M_k(\theta )$越小，即贪婪双子空间随机Kaczmarz方法收敛速度越快。为便于比较，用$\theta$取值特殊的贪婪双子空间随机Kaczmarz方法与原双子空间随机Kaczmarz方法作比较， 但对实际问题， 因贪婪双子空间随机Kaczmarz方法选行概率更具优势性， 对于$\forall \theta \in [\mathrm{0,1}]$，贪婪双子空间随机Kaczmarz方法都比原双子空间随机Kaczmarz方法收敛更快。

事实 ①令$t_{\mathrm{1}},t_{\mathrm{2}},\cdots ,t_l$是在概率分别为$p_{\mathrm{1}},p_{\mathrm{2}},\cdots ,p_l$的某个概率空间上定义的一组随机变量。若越大的$\left|t_i\right|$所对应的概率$p_i$越大，则值${\displaystyle \sum _{i=\mathrm{1}}^l}\left|t_i\right|p_i$越大；②对于任意数组$\Gamma =\left\{t_{\mathrm{1}},t_{\mathrm{2}},\cdots ,t_l\right\}$，$t_i\in \mathbf{R}$，$\mathrm{a}\mathrm{v}\mathrm{e}(\Gamma )$表示数组$\Gamma$的平均值。令$\widehat{\Gamma }=\left\{t_i\ge \mathrm{a}\mathrm{v}\mathrm{e}(\Gamma ),t_i\in \Gamma \right\}$，有$\mathrm{a}\mathrm{v}\mathrm{e}(\widehat{\Gamma })\ge \mathrm{a}\mathrm{v}\mathrm{e}(\Gamma )$。

当$\theta \le \frac{\mathrm{1}}{m}$时，显然成立

$$\begin{array}{l}\mathrm{1}-\frac{\mathrm{1}-\theta }{m-\mathrm{2}}{\sigma }_{{}_{\mathrm{m}\mathrm{i}\mathrm{n}}}^{\mathrm{2}}(\mathit{A})\le \mathrm{1}-\frac{\mathrm{1}-\theta }{m-\mathrm{1}}{\sigma }_{{}_{\mathrm{m}\mathrm{i}\mathrm{n}}}^{\mathrm{2}}(\mathit{A})\le \mathrm{1}-\frac{{\sigma }_{{}_{\mathrm{m}\mathrm{i}\mathrm{n}}}^{\mathrm{2}}(\mathit{A})}{m}\\ \mathrm{即}M(\theta )\le {\left(\mathrm{1}-\frac{{\sigma }_{{}_{\mathrm{m}\mathrm{i}\mathrm{n}}}^{\mathrm{2}}(\mathit{A})}{m}\right)}^{\mathrm{2}}\mathrm{。}\end{array}$$

同时，当$\theta \le \frac{\mathrm{1}}{m}$时，对于任意$r\in {\widehat{U}}_k$，显然成立

$$\begin{array}{l}{\left|{\mathit{r}}_{y_k}^{(r)}\right|}^{\mathrm{2}}\ge \\ (\mathrm{1}-\theta )\underset{\mathrm{1}\le j\le m}{\mathrm{m}\mathrm{a}\mathrm{x}}\left\{{\left|{\mathit{b}}^{(j)}-{\mathit{A}}^{(j)}{\mathit{y}}_k\right|}^{\mathrm{2}}\right\}\ge \\ \frac{m-\mathrm{1}}{m}\underset{\mathrm{1}\le j\le m}{\mathrm{m}\mathrm{a}\mathrm{x}}\left\{{\left|{\mathit{b}}^{(j)}-{\mathit{A}}^{(j)}{\mathit{y}}_k\right|}^{\mathrm{2}}\right\}\ge \\ \frac{{‖\mathit{b}-\mathit{A}{\mathit{y}}_k‖}_{\mathrm{2}}^{\mathrm{2}}}{m}\mathrm{。}\end{array}$$

基于事实①和②可得$\begin{array}{l}E{\left|\mu <{\mathit{e}}_k^{\mathrm{T}},\mathit{\nu }>\right|}^{\mathrm{2}}=\\ {\displaystyle \sum _{r\in {\widehat{U}}_k}^{}}\frac{{\left|{\widehat{\mathit{r}}}_{y_k}^{(r)}\right|}^{\mathrm{2}}}{{‖{\widehat{\mathit{r}}}_{y_k}‖}_{\mathrm{2}}^{\mathrm{2}}}{\displaystyle \sum _{s\in U_k}^{}}\frac{{\left|{\widehat{\mathit{r}}}_{x_k}^{(s)}\right|}^{\mathrm{2}}}{{‖{\widehat{\mathit{r}}}_{x_k}‖}_{\mathrm{2}}^{\mathrm{2}}}{\left|{\mu }_{r,s}\right|}^{\mathrm{2}}{\left|<{\mathit{e}}_k^{\mathrm{T}},{\mathit{\nu }}_{r,s}>\right|}^{\mathrm{2}}=\end{array}$

$$\begin{array}{l}{\displaystyle \sum _{r\in {\widehat{U}}_k}^{}}\frac{{\left|{\widehat{\mathit{r}}}_{y_k}^{(r)}\right|}^{\mathrm{2}}}{{‖{\widehat{\mathit{r}}}_{y_k}‖}_{\mathrm{2}}^{\mathrm{2}}}{\displaystyle \sum _{s\in U_k}^{}}\frac{{\left|{\widehat{\mathit{r}}}_{x_k}^{(s)}\right|}^{\mathrm{2}}}{{‖{\widehat{\mathit{r}}}_{x_k}‖}_{\mathrm{2}}^{\mathrm{2}}}\frac{{\left|{\mu }_{r,s}\right|}^{\mathrm{2}}}{\mathrm{1}-{\left|{\mu }_{r,s}\right|}^{\mathrm{2}}}{\left|<{\mathit{x}}_{\mathrm{*}}^{\mathrm{T}}-{\mathit{y}}_k^{\mathrm{T}},{\mathit{A}}^{(r)}>\right|}^{\mathrm{2}}\ge \\ \frac{\mathrm{1}}{\left|{\widehat{U}}_k\right|}\frac{\mathrm{1}}{\left|U_k\right|}\sum _{r\in {\widehat{U}}_k}{\displaystyle \sum _{s\in U_k}^{}}\frac{{\left|{\mu }_{r,s}\right|}^{\mathrm{2}}}{\mathrm{1}-{\left|{\mu }_{r,s}\right|}^{\mathrm{2}}}{\left|<{\mathit{x}}_{\mathrm{*}}^{\mathrm{T}}-{\mathit{y}}_k^{\mathrm{T}},{\mathit{A}}^{(r)}>\right|}^{\mathrm{2}}\ge \\ \frac{\mathrm{1}}{m-\mathrm{1}}\frac{\mathrm{1}}{m}\sum _{r\ne s}{\displaystyle \sum _{s\in m}^{}}\frac{{\left|{\mu }_{r,s}\right|}^{\mathrm{2}}}{\mathrm{1}-{\left|{\mu }_{r,s}\right|}^{\mathrm{2}}}{\left|<{\mathit{x}}_{\mathrm{*}}^{\mathrm{T}}-{\mathit{y}}_k^{\mathrm{T}},{\mathit{A}}^{(r)}>\right|}^{\mathrm{2}}=\\ \frac{\mathrm{1}}{m(m-\mathrm{1})}\sum _{r\ne s}\left|\frac{{\mu }_{r,s}^{\mathrm{2}}}{\sqrt[]{\mathrm{1}-{\mu }_{r,s}^{\mathrm{2}}}}<{\mathit{e}}_k^{\mathrm{T}},{\mathit{A}}^{(r)}>-\right.\\ {\left.\frac{{\mu }_{r,s}^{}}{\sqrt[]{\mathrm{1}-{\mu }_{r,s}^{\mathrm{2}}}}<{\mathit{e}}_k^{\mathrm{T}},{\mathit{A}}^{(s)}>\right|}^{\mathrm{2}}\end{array}$$

因此，贪婪双子空间随机Kaczmarz方法收敛因子小于原双子空间随机Kaczmarz方法的收敛因子。

浮点运算量方面， 原双子空间随机Kaczmarz方法每次迭代需要$\mathrm{12}n+\mathrm{4}$个浮点运算量。若残差由递推公式计算， 则贪婪双子空间随机Kaczmarz方法每次迭代需要$\mathrm{16}m+\mathrm{6}n+\mathrm{7}$个浮点运算量。因此， 当$m<\frac{\mathrm{6}n-\mathrm{3}}{\mathrm{16}}$时， 贪婪双子空间随机Kaczmarz方法计算量小于原双子空间随机Kaczmarz方法。

3 数值实验

在本节中， 分别通过数值实验比较随机Kaczmarz方法（RK）、双子空间随机Kaczmarz方法（2S‒RK）和贪婪双子空间随机Kaczmarz方法（2S‒GRK（θ））， 并显示后者无论在迭代步数（IT）还是计算时间（CPU）上均优于前两者。这里，迭代步数（IT）和计算时间（CPU）取30次计算结果的平均值。

在实验中， 通过Matlab函数$\mathbf{r}\mathbf{a}\mathbf{n}\mathbf{d}\mathbf{n}$随机生成解$x_{\mathrm{*}}\in {\mathbf{R}}^n$， 右端项$b\in {\mathbf{R}}^m$由$\mathit{A}{\mathit{x}}_{\mathrm{*}}$给出。此外，需要指出2S‒GRK（θ）方法按照算法3中定义的过程精确执行， 并没有显式地计算矩阵$\mathit{B}=\mathit{A}{\mathit{A}}^{\mathrm{T}}$来计算残量。本节的所有数值实验， 均取初始向量为${\mathit{x}}_{\mathrm{0}}=\mathrm{0}$， 停机准则为近似解的相对误差$R_{\mathrm{S}\mathrm{E}}$满足

$$R_{\mathrm{S}\mathrm{E}}=\frac{{‖x_{\mathrm{*}}-x_k‖}_{\mathrm{2}}^{\mathrm{2}}}{{‖x_{\mathrm{*}}‖}_{\mathrm{2}}^{\mathrm{2}}}<{\mathrm{10}}^{-\mathrm{6}}$$

或者迭代步数超过30万。在实验表格中， 若在30万步内未达到指定精度， 即用“$--$”表示。

为了测试这些方法，选取了两类矩阵进行实验。一类是通过Matlab函数$\mathbf{u}\mathbf{n}\mathbf{i}\mathbf{f}\mathbf{r}\mathbf{n}\mathbf{d}$生成的500×100的相干矩阵。矩阵元素是在区间$[d,\mathrm{1}]$上独立且均匀分布的随机变量。通过改变$d$的值，可以得到具有不同相干系数对的矩阵，相干系数对$(\sigma ,\Delta )$定义如下

$$\sigma =\sigma (A)=\underset{i\ne j}{\mathrm{m}\mathrm{i}\mathrm{n}}\left|<A^{(i)},A^{(j)}>\right|,\Delta =\Delta (A)=\underset{i\ne j}{\mathrm{m}\mathrm{a}\mathrm{x}}\left|<A^{(i)},A^{(j)}>\right|$$

表1列出了测试的相干矩阵的相关信息。包含不同相干程度的500×100矩阵。

表1 相干矩阵信息

Tab.1 Information of coherent matrices

矩阵名称	阶数	取值区间	相干系数对
A1	$500\times 100$	$[-\mathrm{0.4,1}]$	$(0.008\mathrm{4,1})$
A2	$500\times 100$	$[-\mathrm{0.1,1}]$	$(0.477\mathrm{9,1})$
A3	$500\times 100$	$[\mathrm{0.2,1}]$	$(0.799\mathrm{3,1})$
A4	$500\times 100$	$[\mathrm{0.5,1}]$	$(0.944\mathrm{8,1})$
A5	$500\times 100$	$[\mathrm{0.8,1}]$	$(0.993\mathrm{5,1})$

对于测试的矩阵， 表2给出了RK， 2S‒RK， 2S‒GRK（θ）3种方法的迭代步数（IT）与计算时间（CPU）， 以及2S‒GRK（θ）方法对于2S‒RK方法的加速比。

表2 相干矩阵数值结果

Tab.2 Numerical results of coherent matrices

算法	A1		A2		A3		A4		A5
	IT	CPU	IT	CPU	IT	CPU	IT	CPU	IT	CPU
RK	2 530.0	0.133 8	5 020.0	0.255 2	12 120.0	0.618 7	43 870.0	2.236 4	$--$	$--$
2S‒RK	1 130.6	0.034 5	1 474.7	0.036 0	1 644.8	0.039 1	1 752.6	0.041 5	1 745.8	0.041 1
2S‒GRK（θ）	124.0	0.013 7	136.0	0.011 0	142.0	0.011 3	150	0.011 7	141.0	0.011 2
加速比	2.48		3.27		3.46		3.55		3.64

由表2可知， 2S‒RK与2S‒GRK（θ）方法总能计算出符合精度的解， 而RK方法有时在迭代步数达到30万步后仍不能计算出符合精度的解。而且即便RK方法收敛， 2S‒RK与2S‒GRK（θ）方法在迭代步数与计算时间上也均少于RK方法。此外，2S‒GRK（θ）方法进一步优化了2S‒RK方法， 其关于2S‒RK方法迭代时间的加速比最大可以达到3.64， 最小可以达到2.48。值得注意的是， 无论矩阵的相干程度如何，2S‒GRK（θ）方法均收敛， 且其迭代步数与计算时间均优于原始的2S‒RK方法。

图1描绘了相干矩阵A2（图1a）和A4（图1b）的近似解的相对误差$R_{\mathrm{S}\mathrm{E}}$以10为底的对数随着迭代步数变化的曲线， 进一步验证了2S‒GRK（θ）方法比经典的2S‒RK方法收敛更快。

图1 以A2和A4为系数矩阵的线性系统的近似解的相对误差随着迭代步数变化的曲线

Fig.1 $\mathbf{l}\mathbf{g}({\mathit{R}}_{\mathbf{S}\mathbf{E}})$ versus IT for linear systems with coefficient matrices:A2 and A4

表3列出了测试的稀疏矩阵的相关信息。

表3 稀疏矩阵信息

Tab.3 Information of sparse matrices

矩阵名称	阶数	非零元个数	相干系数对
Abtaha1	$14596\times 209$	51 307	$(0.816\mathrm{5,1})$
Ash219	$219\times 85$	438	$(\mathrm{0,1})$
Ash331	$331\times 104$	662	$(\mathrm{0,1})$
Ash608	$608\times 188$	1 216	$(\mathrm{0,1})$
Ash958	$958\times 292$	1 916	$(\mathrm{0,1})$

对于测试的稀疏矩阵， 表4给出了RK、2S‒RK、2S‒GRK（θ）3种方法的迭代步数（IT）与计算时间（CPU），以及2S‒GRK（θ）方法对于2S‒RK方法的加速比。

表4 稀疏矩阵数值结果

Tab.4 Numerical results of sparse matrices

算法	Abtaha1		Ash219		Ash331		Ash608		Ash958
	IT	CPU	IT	CPU	IT	CPU	IT	CPU	IT	CPU
RK	81 883	26.572 9	1 896	0.123 1	2 183	0.145 1	4 119	0.302 7	5.885	0.512 8
2S‒RK	24 931	12.346 3	901	0.065 0	1 050	0.078 4	2 056	0.184 3	3.005	0.316 3
2S‒GRK（θ）	254	7.059 6	127	0.010 6	129	0.012 7	237	0.030 2	353	0.055 7
加速比	1.75		6.13		6.17		6.10		5.68

由表4可得与表2类似结论， 即 2S‒GRK（θ）与2S‒RK方法在迭代步数与计算时间上也均少于RK方法。此外， 2S‒GRK（θ）方法进一步优化了2S‒RK方法， 其关于2S‒RK方法迭代时间的加速比最大可以达到6.17， 最小可以达到1.75。因此，无论对于相干矩阵还是稀疏（不相干）矩阵， 2S‒GRK（θ）方法均收敛， 且其迭代步数与计算时间均优于传统的 2S‒RK方法。

图2描绘了矩阵Ash219（图2a）和Ash958（图2b）的近似解的相对误差以10为底的对数随着迭代步数变化的曲线， 进一步验证了 2S‒GRK（θ）方法方法比经典的2S‒RK方法收敛更快。

图2 以Ash219和Ash958为系数矩阵的线性系统的近似解的相对误差随着迭代步数变化的曲线

Fig.2 $\mathbf{l}\mathbf{g}({\mathit{R}}_{\mathbf{S}\mathbf{E}})$ versus IT for linear systems with coefficient matrices:Ash219 and Ash958

对于4个测试矩阵A2、A4、Ash219和Ash958，图3描绘了2S‒GRK（θ）方法的计算时间（图3a）和迭代步数（图3b）随控制参数$\theta$变化的曲线。由图3可知，只要在计算时采用合适的控制参数$\theta$， 2S‒GRK（θ）方法的性能就会得到很大提升。结合表2和表4可知，无论$\theta \in [\mathrm{0,1}]$取何值， 2S‒GRK（θ）方法方法总比经典的2S‒RK方法收敛更快。

图3 2S‒GRK（θ）方法的计算时间和迭代步数随控制参数θ变化的曲线

Fig.3 CPU and IT versus θ of 2S‒GRK(θ)

4 结论

提出一类贪婪双子空间随机Kaczmarz方法， 理论分析证明新方法的收敛性， 还表明收敛因子小于传统的双子空间随机Kaczmarz方法的收敛因子。数值实验结果也表明所提出的新方法在迭代步数和计算时间上均优于传统的双子空间随机Kaczmarz方法。贪婪双子空间随机Kaczmarz方法能够快速求解大规模稀疏相容线性方程组， 主要原因在于新方法定义的随机选取工作行的概率更具优势。如何随机选取工作行以加速Kaczmarz方法及其应用^［43］仍然是一个值得研究的问题。

作者贡献声明

荆燕飞：核心思想提炼、论文修改。

李彩霞：论文撰写。

胡少亮：论文修改。