基于几何平滑动量的大规模线性方程组自适应确定性块Kaczmarz方法

An Adaptive Deterministic Block Kaczmarz Algorithm Based on Geometrically Smoothed Momentum for Solving Large-Scale Linear Systems

摘要: 为高效求解大规模相容线性方程组，本文在已有自适应确定性块Kaczmarz (ADBK)方法的基础上，引入几何平滑动量加速机制，提出了基于几何平滑动量的自适应确定性块Kaczmarz方法。该方法继承了ADBK通过残差范数自适应选取块索引、避免显式计算子矩阵伪逆的计算优势，并在不增加额外计算复杂度的前提下进一步提升了收敛速度。数值实验结果表明，与原始ADBK方法及多种典型算法相比，所提出的方法在迭代次数、计算时间和收敛速度方面均具有明显优势。

Abstract: To efficiently solve large-scale consistent systems of linear equations, this paper proposes an adaptive deterministic block Kaczmarz method based on geometric smooth momentum, building upon the existing Adaptive Deterministic Block Kaczmarz (ADBK) method. The proposed method inherits the computational advantages of ADBK in adaptively selecting block indices through residual norms and avoiding explicit computation of submatrix pseudoinverses. Furthermore, by introducing a geometric smooth momentum acceleration mechanism, the convergence speed is enhanced without increasing additional computational complexity. Numerical experimental results demonstrate that compared with the original ADBK method and various typical algorithms, the proposed method exhibits significant advantages in terms of iteration count, computational time, and convergence rate.

1. 引言

线性方程组的高效求解是大数据分析、机器学习和图像处理等领域中高性能数值算法和软件开发的核心环节，常常受到问题规模、不相容等因素制约。因此，构造格式简单、内存耗费低和收敛速度快的迭代方法在大规模线性方程组、线性反问题和优化等问题中尤为重要。

波兰数学家Stefan Kaczmarz于1937年提出Kaczmarz方法[1]，用于求解线性方程组和线性不适定问题，如电子计算机断层扫描和信号处理[2]等。然而，Kaczmarz方法的收敛性质依赖于行指标的选择顺序，收敛速度慢且很难做出收敛分析结果。在此方法基础上，2009年Strohmer和Vershynin [3]提出了基于工作行指标随机选取的随机Kaczmarz (RK)方法，并证明了其依期望指数收敛，但是RK方法易受大型矩阵稀疏度的影响。

块Kaczmarz (BK)方法最早可以追溯到Elfving [4]的工作，是Kaczmarz方法的自然推广。理论分析和数值实验表明构造块Kaczmarz类方法是提升Kaczmarz方法和RK方法数值性能的有效途径。目前国内外学者主要从两方面来设计新型高效的块方法：一种是基于投影思想建立的BK方法[5]-[7]，这类方法在每一次迭代过程中都需要计算伪逆或求解一个最小二乘问题，且难以实施并行计算。另一种是基于免伪逆思想建立的BK方法，该类方法能避免投影类方法的主要缺陷。比较典型的方法有求解相容线性方程组的随机平均块Kaczmarz (RaBK)方法[8]和快速确定块Kaczmarz (FDBK)方法[9]。

本文研究基于免伪逆思想建立的BK方法，在自适应确定性块Kaczmarz方法的基础上将几何平滑动量运用到该方法中得到基于动量的块Kaczmarz方法。利用数值实验验证带几何平滑动量的自适应确定性块Kaczmarz方法的有效性和高效性。

2. 自适应确定性块Kaczmarz法

考虑如下的大规模线性方程组：

$Ax=b$ ， (1)

其中$A\in {ℝ}^{m\times n}$ ，未知向量$x\in {ℝ}^n$ 和右端项$b\in {ℝ}^m$ 。对于相容的线性方程组，若解不唯一，往往关注方程组(1)的最小欧式范数解

$x_{LN}=\underset{x\in {ℝ}^n}{\arg \min }{\Vert x\Vert }_2$ ，使得$Ax=b$ ，

其中${\Vert \cdot \Vert }_2$ 表示欧几里得范数。对于不相容的线性方程组，一般考虑最小二乘解

$x_{LS}=\underset{x\in {ℝ}^n}{\arg \min }{\Vert b-Ax\Vert }_2$ ，

重点研究的是最小二乘最小范数解：$x_{\ast }=A^{†}b$ 。对于相容线性方程组，其最小二乘最小范数解就是它的最小范数解[10]。

2.1. 块Kaczmarz型方法

设$\left[m\right]$ 表示集合$\left\{1,2,\cdots ,m\right\}$ ，如果一组索引集$J_1,J_2,\cdots ,J_p$ 满足：对于任意$i\ne j\left(i,j=1,2,\cdots ,p\right)$ ，有$J_i\cap J_j=\varnothing$ ，并且${\displaystyle {\cup }_{i=1}^p{J}_i}=\left[m\right]$ ，那么可以称这组集合$J_1,J_2,\cdots ,J_p$ 是对集合$\left[m\right]$ 的一个划分。给定一个正整数$p$ ，可以将矩阵$A$ 的行索引划分为$J=\left\{J_1,J_2,\cdots ,J_p\right\}$ ，即将$A$ 表示为$A={\left[A_{J_1}^{\text{T}},A_{J_2}^{\text{T}},\cdots ,A_{J_p}^{\text{T}}\right]}^{\text{T}}$ ，同时将向量$b$ 表示为$b={\left[b_{J_1}^{\text{T}},b_{J_2}^{\text{T}},\cdots ,b_{J_p}^{\text{T}}\right]}^{\text{T}}$ 。在第$k+1$ 次迭代中，选择一个子矩阵$A_{J_{i_k}}$ 和对应的子向量$b_{J_{i_k}},$ 其中$i_k=\text{mod}\left(k,p\right)+1$ 。在该步中，当前迭代点$x_k$ 被正交投影到解空间$\left\{x\in {ℝ}^n:A_{J_{i_k}}x=b_{J_{i_k}}\right\}$ 上。因此从初始向量$x_0$ 出发，块Kaczmarz方法的迭代格式可以表示为：

$x_{k+1}=x_k+A_{J_{i_k}}^{†}\left(b_{J_{i_k}}-A_{J_{i_k}}{x}_k\right)$ ，$k=0,1,2,\cdots$ (2) \

其中$A_{J_{i_k}}^{†}$ 表示矩阵$A_{J_{i_k}}$ 的Moore-Penrose伪逆。

在给定集合$\left[m\right]$ 的一个预定划分$J$ 之后，Needell和Tropp [5]提出随机块Kaczmarz (RBK)方法。为降低计算复杂度，Necoara [8]提出了一个统一的随机平均块Kaczmarz (RaBK)方法，该方法不需要计算伪逆，其迭代格式如下：

$x_{k+1}=x_k+{\alpha }_k\left({\displaystyle \sum _{i\in J_{i_k}}{\omega }_i^k}\frac{b^{\left(i\right)}-A^{\left(i\right)}{x}_k}{{\Vert A^{\left(i\right)}\Vert }_2^2}{\left(A^{\left(i\right)}\right)}^{\text{T}}\right)$ ，(3)

其中${\omega }_i^k\in \left[0,1\right]$ 为权重，满足${\displaystyle \sum _{i\in J_{i_k}}{\omega }_i^k}=1$ ，${\alpha }_k\in \left(0,2\right)$ 为步长。为了进一步提高收敛速度，Chen和Huang提出了一个快速确定性块Kaczmarz (FDBK)方法[9]，数值实验表明，FDBK方法在性能上优于RaBK方法。

2.2. 自适应确定性块Kaczmarz方法

首先提出一种基于残差向量的欧几里得范数的新指标集序列

$U_k=\left\{i_k|{\left|b^{\left(i_k\right)}-A^{\left(i_k\right)}{x}_k\right|}^2\ge \frac{{\Vert b-A{x}_k\Vert }_2^2}{m}\right\}$ ，(4)

上述指标序列$U_k$ 非空，具体内容见算法1。

算法1. ADBK方法[10]

1：输入：$A,b,l,x_0\in range\left(A^{\text{T}}\right)$ ；

2：输出：$x_l$

3：for $k=0,1,\cdots ,l-1$ do

4：确定指标集序列

$U_k=\left\{i_k|{\left|b^{\left(i_k\right)}-A^{\left(i_k\right)}{x}_k\right|}^2\ge \frac{{\Vert b-A{x}_k\Vert }_2^2}{m}\right\}$ ；

5：计算

${\eta }_k={\displaystyle \sum }_{i_k\in U_k}\left(b^{\left(i_k\right)}-A^{\left(i_k\right)}{x}_k\right)e_{i_k}$ ；(5)

6：计算

$x_{k+1}=x_k+\frac{{\eta }_k^{\text{T}}\left(b-A{x}_k\right)}{{\Vert A^{\text{T}}{\eta }_k\Vert }_2^2}{A}^{\text{T}}{\eta }_k.$ (6)

7：end for

定理2.2.1. 对$A^{\text{T}}$ 列空间中的任意初始猜测向量$x_0$ ，由ADBK方法生成的迭代解序列${\left\{x_k\right\}}_{k=0}^{\infty }$ 收敛到相容线性方程组(1)的最小范数解$x_{\ast }=A^{†}b$ ，且迭代解序列${\left\{x_k\right\}}_{k=0}^{\infty }$ 的误差满足

${\Vert x_{k+1}-x_{\ast }\Vert }_2^2\le \left(1-\frac{\mu {\lambda }_{\text{min}}\left(A^{\text{T}}A\right)}{m{\lambda }_{\text{max}}^{block}}\right){\Vert x_k-x_{\ast }\Vert }_2^2$ ， (7)

其中$\mu =\underset{k}{\text{inf}}\left\{N\left(U_k\right)\right\}$ ，$N\left(U_k\right)$ 表示指标集$U_k$ 中所含元素的个数以及${\lambda }_{\text{max}}^{block}={\text{max}}_{U_k}\left\{{\lambda }_{\text{max}}\left(A_{U_k}{A}_{U_k}^{\text{T}}\right)\right\}$ 。

3. 基于几何平滑动量的自适应确定性块Kaczmarz方法

考虑使用几何平滑动量来加快ADBK方法的收敛速度。类似于文献[11]中的工作，对于参数$\beta \in \left[0,1\right)$ 和$M\in \left[0,1\right]$ ，定义几何平滑动量Kaczmarz方法的算法，其迭代形式如下：

$\{\begin{array}{l}{x}_{k+1}=x_k+\frac{b^{\left(i_k\right)}-A^{\left(i_k\right)}{x}_k}{{\Vert A^{\left(i_k\right)}\Vert }_2^2}{\left(A^{\left(i_k\right)}\right)}^{\text{T}}+M{y}_k\hfill \\ y_{k+1}=\beta y_k+\left(1-\beta \right)\left(x_{k+1}-x_k\right)\hfill \end{array}$ ，(8)

下面给出几何平滑动量Kaczmarz (KGSM)方法的收敛理论。

定理3.1. 固定参数$\beta \in \left[0,1\right)$ ，$M\in \left[0,1\right]$ ，以及$l\in \left\{1,\cdots ,n\right\}$ 。设$v_l$ 是矩阵$A$ 的第$l$ 大奇异值${\sigma }_l$ 所对应的右奇异向量。假设$x_k$ 是由KGSM方法所定义的。那么对于所有$k\ge 0,$ 有

$\mathbb{E}\langle x_{k+1}-x,v_l\rangle ={\left[\begin{array}{c}r\\ \varsigma \end{array}\right]}^{\text{T}}{\left[\begin{array}{cc}r& \varsigma \\ -1& \beta \end{array}\right]}^k\left[\begin{array}{c}1\\ -1/\left(1-\beta \right)\end{array}\right]\langle x_0-x,v_l\rangle ,$

其中

$r:=1-\frac{{\sigma }_l^2}{{\Vert A\Vert }_F^2}+M\left(1-\beta \right),\varsigma :=M{\left(1-\beta \right)}^2$ 。

当$\beta$ 接近于1时，动量项会变得高度平滑。所以将几何平滑动量更新公式嵌入ADBK方法的更新迭代公式中有

$\{\begin{array}{l}{x}_{k+1}=x_k+\frac{{\eta }_k^{\text{T}}\left(b-A{x}_k\right)}{{\Vert A^{\text{T}}{\eta }_k\Vert }_2^2}{A}^{\text{T}}{\eta }_k+M{y}_k\hfill \\ y_{k+1}=\beta y_k+\left(1-\beta \right)\left(x_{k+1}-x_k\right)\hfill \end{array}$ (9)

综上，可得带几何平滑动量的ADBK (gsmADBK)方法，具体内容见如下算法2。

算法2. gsmADBK方法

1：输入：$A,b,l,0\le \beta \le 1$ ，并且$x_0\in range\left(A^{\text{T}}\right)$ ，$y_0=x_0$ ；

2：输出：$x_l$

3：for $k=0,1,\cdots ,l-1$ do

4：确定指标集序列

$U_k=\left\{i_k|{\left|b^{\left(i_k\right)}-A^{\left(i_k\right)}{x}_k\right|}^2\ge \frac{{\Vert b-A{x}_k\Vert }_2^2}{m}\right\}$ ；

5：计算

${\eta }_k={\displaystyle \sum }_{i_k\in U_k}\left(b^{\left(i_k\right)}-A^{\left(i_k\right)}{x}_k\right)e_{i_k}$ ；

6：计算

$x_{k+1}=x_k+\frac{{\eta }_k^{\text{T}}\left(b-A{x}_k\right)}{{\Vert A^{\text{T}}{\eta }_k\Vert }_2^2}{A}^{\text{T}}{\eta }_k+M{y}_k$

7：计算

$y_{k+1}=\beta y_k+\left(1-\beta \right)\left(x_{k+1}-x_k\right).$

8：end for

注记3.1. 当参数$M=0$ 时，gsmADBK方法退化为ADBK方法。

gsmADBK方法的收敛性分析由于${\eta }_k$ 的选择不确定，因此给出一个准确的收敛速度的上界是困难的，有待今后的进一步研究。

4. 数值实验

下面通过几组数值算例来验证gsmADBK方法求解大规模相容线性方程组的有效性以及高效性，实验均通过MATLAB (版本R2020b)编程实现，初始向量均设为$x_0=zeros\left(n,1\right)$ 。停机准则是相对解误差满足

$RSE=\frac{{\Vert x_k-x_{\ast }\Vert }_2}{{\Vert x_{\ast }\Vert }_2}\le {10}^{-3}$

或者迭代次数超过100,000次，或者计算时间超过3600秒。在后两种情况下，相应的IT迭代次数和CPU时间记为“--”。IT和CPU是50次重复计算所需迭代步数和CPU时间的平均值(单位：秒)。令右端向量$b=A{x}_{\ast },x_{\ast }\in {ℝ}^n$ 由MATLAB函数randn随机生成。

4.1. 超定和欠定线性方程组

首先考虑第一类系数矩阵高斯矩阵，由MATLAB函数randn随机生成，即$A=randn\left(m,n\right)$ 。使用RaBK、FDBK、ADBK和gsmADBK方法求解大规模相容线性方程组$Ax=b$ ，其中RaBK方法的采样方式为分块采样，块大小为$\left[m/{\Vert A\Vert }_2^2\right]$ ，步长${\alpha }_k=1.95L_k$ ，$L_k=\frac{{\displaystyle {\sum }_{i\in J_k}{\overline{\omega }}_i^k{\left(A^{\left(i\right)}{x}_k-b^{{}^{\left(i\right)}}\right)}^2}}{{\Vert {\displaystyle {\sum }_{i\in J_k}{\overline{\omega }}_i^k{\left(A^{\left(i\right)}{x}_k-b^{{}^{\left(i\right)}}\right)}^2{\left(A^{\left(i\right)}\right)}^{\text{T}}}\Vert }_2^2}$，${\overline{\omega }}_i^k=\frac{{\omega }_i^k}{{\Vert A^{\left(i\right)}\Vert }_2^2}$ 。数值结果由表1和表2给出，其中gsmADBK方法采用的参数是数值实验过程中确定的最优参数，即通过CPU值最小获取。为了显示收敛程度，gsmADBK方法对其他方法的加速定义为

$\text{speed-up\_ADBK}=\frac{\text{CPU of Method}}{\text{CPU of gsmADBK}}$ 。

由表1和表2的数值结果可以看出，RaBK、FDBK、ADBK和gsmADBK方法求解超定及欠定线性方程组(1)都是有效的。从迭代步数和计算时间来看，gsmADBK方法优于RaBK、FDBK和ADBK方法，speed-up_RaBK、speed-up_FDBK和speed-up_ADBK的取值范围分别为895.56~5419.41、8.69~103.16、1.05~3.39，即增加动量项给ADBK方法的收敛特性带来有效的改进并且优于传统方法。

Table 1. Numerical results of four iterative methods when $A=randn\left(m,n\right)$ , $m>n$

表1. $A=randn\left(m,n\right)$ ，$m>n$ 时四种迭代方法的数值结果

m		1000	3000	5000	7000	9000	12,000	15,000
n		500	1000	2000	3500	5000	6500	8000
RaBK	IT	12,085	12,606	32,123	86,189	>105	>105	>105
	CPU	2.1848	16.6574	173.0970	1243.2086	--	--	--
FDBK	IT	278	128	230	446	689	636	650
	CPU	0.0289	0.1617	1.2585	5.2504	15.0544	23.5612	36.3798
ADBK	IT	70	32	42	70	100	87	88
	CPU	0.0062	0.0378	0.2238	0.8485	2.1926	3.2254	4.9407
gsmADBK	IT	23	15	17	23	29	28	28
	CPU	0.0022	0.0186	0.0915	0.2919	0.6470	1.0641	1.6170
	$M_{opt}$	0.5	0.3	0.4	0.5	0.6	0.6	0.5
	${\beta }_{opt}$	0.2	0.2	0.1	0.2	0.1	0.1	0.2
speed-up_RaBK		993.09	895.56	1891.77	4259.02	--	--	--
speed-up_FDBK		13.14	8.69	13.75	17.99	23.27	22.14	22.50
speed-up_ADBK		2.82	2.03	2.45	2.91	3.39	3.03	3.06

Table 2. Numerical results of four iterative methods when $A=randn\left(m,n\right)$ , $m<n$

表2. $A=randn\left(m,n\right)$ ，$m<n$ 时四种迭代方法的数值结果

m		500	1000	2000	3500	5000	6500	8000
n		10,000	3000	5000	7000	9000	12,000	15,000
RaBK	IT	11,740	12,467	32,609	86,022	>105	>105	>105
	CPU	4.3706	39.2750	319.8674	1764.5594	--	--	--
FDBK	IT	320	174	269	506	753	709	721
	CPU	0.2579	0.3220	1.7888	6.4280	17.2074	26.1081	40.1144
ADBK	IT	70	32	42	70	100	87	88
	CPU	0.0065	0.0372	0.2247	0.8590	2.1213	3.3347	4.7675
gsmADBK	IT	23	15	17	23	29	28	28
	CPU	0.0025	0.0226	0.1096	0.3256	0.7025	1.1920	1.7674
	$M_{opt}$	0.5	0.4	0.4	0.5	0.6	0.6	0.5
	${\beta }_{opt}$	0.3	0.2	0.2	0.2	0.2	0.1	0.3
speed-up_RaBK		1748.24	1737.83	2918.50	5419.41	--	--	--
speed-up_FDBK		103.16	14.25	16.32	19.74	24.49	21.90	22.70
speed-up_ADBK		2.60	1.65	1.05	2.64	3.02	2.80	2.70

图1展示了两个不同参数$M$ 和$\beta$ 对gsmADBK方法计算效率的影响。确定gsmADBK方法的最优参数主要是通过计算不同参数对应的CPU值进行选择的。

(a) $m=7500,n=3000$ (b) $m=3500,n=7000$

Figure 1. CPU values of solving overdetermined and underdetermined linear systems using gsmADBK method with different $\left(M,\beta \right)$

图1. 不同$\left(M,\beta \right)$ 的gsmADBK方法求解超定和欠定线性方程组的CPU值

在图2和图3中，对于超定和欠定线性方程组，有RSE与迭代步数(左)和计算时间(右)的关系曲线。从中反映出gsmADBK方法收敛最快，综上所述，gsmADBK方法的数值性能表现优秀。

4.2. 稀疏线性方程组

考虑第二类系数矩阵稀疏矩阵，均是从SuiteSparse矩阵集合[12]中选取的。表3和表4分别列出了六个细长满秩矩阵和六个扁平满秩矩阵的相关信息，表5列出了六个秩亏稀疏矩阵的相关信息。使用RaBK、FDBK、ADBK、和gsmADBK方法求解初始向量为$x_0$ 的大规模相容线性方程组$Ax=b$ 。数值结果由表6，表7以及表8给出，其中不同方法中采用的参数是数值实验过程中确定的最优参数，即通过CPU值最小获取。

Figure 2. The relationship of RSE with IT (left) and CPU (right) for RaBK, FDBK, ADBK and gsmADBK methods of$A=randn\left(7000,3500\right)$

图2. $A=randn\left(7000,3500\right)$ 的RaBK、FDBK、ADBK和gsmADBK方法的RSE与IT (左)和CPU (右)的关系

Figure 3. The relationship of RSE with IT (left) and CPU (right) for RaBK, FDBK, ADBK and gsmADBK methods of$A=randn\left(3500,7000\right)$

图3. $A=randn\left(3500,7000\right)$ 的RaBK、FDBK、ADBK和gsmADBK方法的RSE与IT (左)和CPU (右)的关系

Table 3. Related properties of full-rank sparse matrices, $m\ge n$

表3. 满秩稀疏矩阵的相关性质，$m\ge n$

Name	Cities	ash219	ash331	ash608	cage9	cage10
$m\times n$	55 × 46	219 × 85	331 × 104	608 × 188	3534 × 3534	11,397 × 11,397
density	53.04%	2.35%	1.92%	1.06%	0.33%	0.16%
cond(A)	207.15	3.02	3.10	3.37	12.60	11.02

Table 4. Related properties of full-rank sparse matrices, $m<n$

表4. 满秩稀疏矩阵的相关性质，$m<n$

Name	Trec8	crew1	bibd_17_8	model1	nemsafm	flower_5_4
$m\times n$	23 × 84	135 × 6469	136 × 24,310	362 × 798	334 × 2348	5226 × 14,721
density	39.44%	5.38%	20.59%	1.05%	0.36%	0.057%
cond(A)	26.89	18.20	9.04	17.57	4.77	14.9

Table 5. Relevant properties of low-rank sparse matrices

表5. 秩亏稀疏矩阵的相关性质

Name	relat5	relat6	rel5	rel6	Sandi_sandi	flower_4_4
$m\times n$	340 × 35	2340 × 157	340 × 35	2340 × 157	314 × 360	1837 × 5529
density	8.89%	2.21%	5.51%	1.39%	0.54%	0.20%
cond(A)	Inf	Inf	Inf	Inf	1.47e+17	1.05e+19

Table 6. Numerical results for full-rank sparse matrices $A$ , $m\ge n$

表6. 满秩稀疏矩阵$A$ ，$m\ge n$ 时的数值结果

Name		Cities	ash219	ash331	ash608	cage9	cage10
RaBK	IT	112488	846	814	1970	31,359	>105
	CPU	0.1753	0.0111	0.0138	0.0810	438.8697	--
FDBK	IT	45593	46	30	50	356	290
	CPU	0.4229	4.7216e−04	3.6335e−04	7.6797e−04	0.0698	0.6507
ADBK	IT	67770	21	19	23	252	107
	CPU	0.1366	5.5380e−05	7.2638e−05	8.5800e−05	0.0140	0.0187
gsmADBK	IT	2084	12	11	13	52	37
	CPU	0.0045	3.2180e−05	3.3107e−05	5.5590e−05	0.0032	0.0063
	$M_{opt}$	0.9	0.2	0.2	0.3	0.7	0.7
	${\beta }_{opt}$	0.6	0.1	0.1	0.2	0.4	0.4
speed-up_RaBK		38.96	344.93	416.83	1457.10	137146.78	--
speed-up_FDBK		93.98	14.67	10.98	13.81	21.81	103.29
speed-up_ADBK		30.36	1.72	2.19	1.54	4.38	2.97

Table 7. Numerical results for full-rank sparse matrices $A$ , $m<n$

表7. 满秩稀疏矩阵$A$ ，$m<n$ 时的数值结果

Name		Trec8	crew1	bibd_17_8	model1	nemsafm	flower_5_4
RaBK	IT	2303	7478	1293	11,357	2614	>105
	CPU	0.0046	2.0872	2.4240	1.5826	1.1650	--
FDBK	IT	1333	808	257	646	136	3501
	CPU	0.0111	0.1113	2.5492	0.0136	0.0036	0.9751
ADBK	IT	1240	319	119	375	56	463
	CPU	0.0019	0.0182	0.0492	0.0021	0.0007	0.0515
gsmADBK	IT	134	92	38	69	24	57
	CPU	0.0004	0.0055	0.0174	0.0004	0.0003	0.0063
	$M_{opt}$	0.8	0.8	0.6	0.8	0.5	0.8
	${\beta }_{opt}$	0.5	0.2	0.1	0.2	0.1	0.2
speed-up_RaBK		11.50	379.49	139.31	3956.50	3883.33	--
speed-up_FDBK		27.75	20.24	146.51	34.00	12.00	154.78
speed-up_ADBK		4.75	3.31	2.83	5.25	2.33	8.17

由表6和表7的数值结果可以得出以下结论：

1) RaBK、FDBK、ADBK和gsmADBK方法求解满秩稀疏线性方程组$Ax=b$ 都是有效的。

2) RaBK、FDBK、ADBK和gsmADBK方法呈现出与表1和表2中类似的数值结果。在所有测试算例求解中，gsmADBK方法消耗的迭代步数和计算时间优于RaBK、FDBK和ADBK方法，speed-up_RaBK、speed-up_FDBK、speed-up_ADBK的取值范围分别为11.50~137146.78、10.98~154.78、1.54~30.96，表明添加动量后加速效果显著且优于传统方法。但搜索不同稀疏矩阵的最优参数是较为耗时的。

Table 8. Numerical results of low-rank sparse matrices

表8. 秩亏稀疏矩阵的数值结果

Name		relat5	relat6	rel5	rel6	Sandi_sandi	flower_4_4
FDBK	IT	75	252	79	322	505	561
	CPU	0.0009	0.0122	0.0008	0.0117	0.0059	0.0481
ADBK	IT	58	84	70	192	353	147
	CPU	0.0002	0.0010	0.0002	0.0019	0.0011	0.0048
gsmADBK	IT	22	28	34	38	88	39
	CPU	6.6070e−05	0.0004	8.4310e−05	0.0004	0.0003	0.0014
	$M_{opt}$	0.4	0.6	0.4	0.7	0.8	0.6
	${\beta }_{opt}$	0.3	0.2	0.2	0.1	0.3	0.3
speed-up_FDBK		13.62	30.50	9.49	29.25	19.67	34.36
speed-up_ADBK		3.03	2.50	2.37	4.75	3.67	3.43

由表8的数值结果可以得出以下结论：

1) RaBK方法无法有效求解病态秩亏稀疏线性方程组，问题在于秩亏系数矩阵使RaBK方法步长的选择失效。其他FDBK、ADBK和gsmADBK方法都能够有效求解病态秩亏稀疏线性方程组。

2) 在所有测试算例求解中，gsmADBK方法消耗的迭代步数和计算时间优于FDBK和ADBK方法，speed-up_FDBK、speed-up_ADBK的取值范围分别为9.49~34.36和2.37~4.75，表明添加动量加速效果显著且优于传统方法。但搜索不同稀疏矩阵的最优参数是较为耗时的。

5. 结论

本文基于数值代数与随机优化等领域的理论知识，在自适应确定性块Kaczmarz (ADBK)方法的基础上结合动量思想提出了基于几何平滑动量的自适应确定性块Kaczmarz (gsmADBK)方法。将该方法应用于超定、欠定及稀疏线性方程组的求解问题，数值实验结果表明该方法能够有效求解大规模相容线性方程组并且引入动量项能够显著改善 ADBK 方法的收敛特性，在迭代次数、计算时间以及收敛速度方面均取得提升。

基金项目

国家自然科学基金(12172186: 11772166)。

NOTES

^*第一作者。

^#通讯作者。