邵新慧, 李晨, 王心怡
东北大学 理学院， 辽宁 沈阳 110819
收稿日期: 2015-10-21
基金项目: 国家自然科学基金资助项目 (11371081)。
作者简介: 邵新慧(1970-)，女，山东青岛人，东北大学副教授。

摘要: 鞍点问题广泛出现在众多的工程研究领域，如流体力学、电磁学、最优化问题、最小二乘问题、椭圆偏微分方程问题等.以SOR类方法为基础，结合HS分裂思想，将经典鞍点问题的求解方法推广到特殊鞍点问题的求解上.给出一种具有新型分裂迭代格式的MSOR-Like方法，用以求解一类含有非对称块的鞍点系统，给出了相应的收敛性分析以及最优松弛参数选取方法.数值算例验证了对于不同的预优矩阵，MSOR-Like方法只有收敛速度的分别，没有收敛性能的影响，且在相同计算精度下，该方法解决特殊鞍点问题的迭代效果优于常规方法解决经典鞍点问题.
关键词：鞍点问题迭代法HS分裂SOR方法收敛
Modified SOR-Like Method for Saddle Point Problems
SHAO Xin-hui, LI Chen, WANG Xin-yi
School of Sciences, Northeastern University, Shenyang 110819, China
Corresponding author: SHAO Xin-hui,E-mail:xinhui1002@126.com
Abstract: Saddle point problems exist in many engineering research areas such as fluid mechanics, electromagnetism, optimization problems, the least squares problems, elliptic partial differential equations, and etc. Based on SOR-Like methods in combination of the concept of HS splitting, a new iteration splitting improvement method was presented so as to apply the classic saddle point solutions to special saddle point problems. Then, the MSOR-Like method was proposed to handle the above special saddle point system containing asymmetric blocks, and the convergence analysis as well as the selection of optimal relaxation parameters were also given. Finally, a numerical example was given to verify different optimal matrix of the modified SOR method, and it was found that the only difference is in the convergence rate while there is no difference in the convergence effect. Furthermore, under the same calculation accuracy, the modified SOR method for solving the special saddle point problems is better than the conventional methods in solving the classical saddle point problems.
Key Words: saddle point problemiterative methodHermitian and Skew-Hermitian (HS) splittingSOR methodconvergence
在工程应用的诸多研究中，如涉及流体力学问题、电磁学问题、最优化问题等相关计算时，常常需要求解一类具有特殊结构的鞍点问题^[1-5]：

(1)

其中:A∈R^m×m为非对称的正定矩阵;B∈R^m×n(m≥n)为列满秩矩阵; x, f∈R^m; y, g∈Rⁿ;B^Τ为B的转置矩阵;f, g为已知向量.鞍点系统呈现一定病态，使得一些经典求解方法难以实现，因此如何行之有效地解决这类问题是近年来国内外学者们研究的热点.通常做法是采用定常迭代与非定常迭代双线并行，结合相应的预条件处理技术, 其中包括Uzawa类方法^[6-9]，SOR类方法^[10-18]，HSS类方法^[19-21]，Krylov子空间类方法^[21]，以及相应的预处理方法.
1 MSOR-Like方法的迭代格式和收敛性分析基于HS分裂，令A=H+S, 其中.由于A为正定矩阵，则H也是正定的.作分裂M=D－L－U, 其中D=这里Q∈R^n×n为对称正定矩阵，若设ω>0为松弛参数，建立迭代格式：
其中:

(2)

整理得如下迭代算法：

(3)

将式(3)称之为求解鞍点问题的MSOR-Like方法.
若令λ为MSOR-Like方法迭代矩阵M_ω的特征值，(x y)^T∈R^m+n为相应特征向量，则有

(4)

引理1?设A∈R^m×m为非对称正定矩阵，为对称正定矩阵，为反对称矩阵，B∈R^m×m(m≥n)为列满秩矩阵，M_ω如式(2)，λ为M_ω的任一特征值，则λ≠1.
证明 ?(反证法)若λ=1，由式(4)可得由于系数矩阵是非奇异阵，则有(x y)^T=0，这与其为迭代矩阵M_ω的特征向量矛盾，故λ≠1.
引理2 ?在引理1条件下，设λ为迭代矩阵M_ω的任一特征值，(x y)^T∈R^m+n为M_ω的相应特征向量，必有x≠0.若y=0，则|λ|＜1.
证明 ?① (反证法)若x=0，由式(4)可得By=0.由于B为列满秩矩阵，因而有y=0, 与(x y)^T为迭代矩阵M_ω的特征向量矛盾，故x≠0.② 若x≠0，y=0，由式(4)可得(1－ω－λ)Hx－λωSx=0.用左乘等式两边，则有(1－ω－λ)×设则有(1－ω－λ)a－λωb=0.由于H是对称正定矩阵，则a>0；S是反对称矩阵，则b=0.因而，有λ=1－ω.由于松弛参数需满足0＜ω＜2, 则有|λ|=|1－ω|＜1.
引理3^[1]? 实一元二次方程λ²－bλ+c=0两根之模均小于1，当且仅当方程系数满足|c|＜1，且 |b|＜1+c.
定理1 ?设A∈R^m×m为非对称的正定矩阵，为对称正定矩阵，为反对称矩阵，B∈R^m×m(m≥n)为列满秩矩阵，Q∈R^n×n为适当的对称正定矩阵，则MSOR-Like方法迭代格式(3)收敛，当且仅当松弛参数满足 其中
证明?设λ为迭代矩阵M_ω的任一特征值，为M_ω的相应特征向量.由于Q的非奇异性，由引理1知λ≠1，将式(4)等价变为
整理并同时左乘到等式两边，则有(1－=0.设则有(1－ω－λ)a－λωb－.由于H和BQ^－1B^T是对称正定矩阵，因而a>0，c>0.而S是反对称矩阵，则b=0，于是有λ²+由引理3知，|λ|＜1当且仅当由于a>0，则因而，若使MSOR－Like方法迭代格式收敛，松弛参数ω的范围需满足
推论1?在定理1条件下，若MSOR-Like方法迭代格式收敛，则松弛参数范围可进一步化简为 其中a=
证明 ?由于H和BQ^－1B^Τ是对称正定矩阵，则分别为H和BQ^－1B^Τ的Rayleigh商.因而，
2 数值实验考虑模型Stokes问题，有对流扩散方程：
其中:Ω为(0, 1)×(0, 1);方程满足Dirichlet边界条件;Δ为Laplace算子;u为表示速度的向量函数;p为表示压力的数值函数.
将上述问题用迎风差分格式离散化处理，可得形式如下的鞍点问题：
其中:

(－1, 1, 0)∈R^p×p, 为Kronecker量积符号，h=为离散网格值，且m=2p²，n=p².将模型中的T作适当修改，得到(－1.5, 2, －0.5)∈R^p×p.于是问题模型重新定义为
令A=H+S, 其中(A－A^Τ)，有 其中
取则m=2p²，n=p².取零向量为初始向量，其精确解为(x^Τ, y^Τ)^Τ=(1, 1, …，1)^Τ.设IT为迭代步数，CPU为上机运算满足收敛要求所需时间，使用MATLAB编程语言实现，当ERR＜10^－6时计算中止, 其中ERR表示迭代误差，定义如下：
数值实验结果表 1~表 3表明，对于不同的预优矩阵Q，MSOR-Like方法只有收敛速度的区别，没有收敛性的影响.因而可以得出，在解决非对称鞍点问题上，该方法是切实有效的.
表 1(Table 1)

表 1 h=1/9(p=8, m=128, n=64)时，ωopt附近区域收敛情况(最优松弛参数ωopt=0.25)Table 1 Converging situation around optimal omega with h=1/9(p=8, m=128, n=64)(ωopt=0.25) ω IT CPU ERR 0.10 276 2.184 0 9.501 0e-07 0.15 193 1.638 0 8.882 2e-07 0.20 155 1.216 8 8.849 8e-07 0.25 141 1.060 8 9.162 5e-07 0.30 158 1.263 6 8.837 1e-07

表 1 h=1/9(p=8, m=128, n=64)时，ωopt附近区域收敛情况(最优松弛参数ωopt=0.25) Table 1 Converging situation around optimal omega with h=1/9(p=8, m=128, n=64)(ωopt=0.25)

表 2(Table 2)

表 2 取h=1/9(p=8, m=128, n=64)时，ωopt附近区域收敛情况(最优松弛参数ωopt=1.00)Table 2 Converging situation around optimal omega with h=1/9(p=8, m=128, n=64)(ωopt=1.00) ω IT CPU ERR 0.98 9 0.093 6 5.698 5e-08 0.99 7 0.062 4 8.099 8e-07 1.00 2 0.031 2 9.488 1e-16 1.01 8 0.062 4 1.952 4e-07 1.02 9 0.093 6 6.221 3e-07

表 2 取h=1/9(p=8, m=128, n=64)时，ωopt附近区域收敛情况(最优松弛参数ωopt=1.00) Table 2 Converging situation around optimal omega with h=1/9(p=8, m=128, n=64)(ωopt=1.00)

表 3(Table 3)

表 3 取h=1/9(p=8, m=128, n=64)时，ωopt附近区域收敛情况(最优松弛参数ωopt=1.15)Table 3 Converging situation around optimal omega with h=1/9 (p=8, m=128, n=64) (ωopt=1.15) ω IT CPU ERR 1.00 51 0.468 0 9.946 3e-07 1.05 50 0.405 6 8.140 7e-07 1.10 47 0.327 6 9.673 2e-07 1.15 47 0.421 2 6.614 1e-07 1.20 88 0.655 2 9.756 8e-07

表 3 取h=1/9(p=8, m=128, n=64)时，ωopt附近区域收敛情况(最优松弛参数ωopt=1.15) Table 3 Converging situation around optimal omega with h=1/9 (p=8, m=128, n=64) (ωopt=1.15)

3 结 语本文以SOR-Like方法为基础，融合HS分裂思想，将经典鞍点问题的求解方法推广到特殊鞍点问题的求解，给出了一种具有新型分裂迭代格式的MSOR-Like方法，理论分析和数值算例验证了MSOR-Like方法在解决非对称鞍点问题时的有效性.
参考文献

[1]	Benzi M, Golub G H, Liesen J. Numerical solution of saddle point problems[J].Acta Numerica , 2005, 14(1): 130–137.
[2]	Duff I S, Gould N I M, Reid J K, et al. The factorization of sparse symmetric indefinite matrices[J].IMA Journal of Numerical Analysis , 1991, 11(2): 181–204.DOI:10.1093/imanum/11.2.181
[3]	Varga R S. Matrix iterative analysis[M]. Englewood Cliffs: Prentice-Hall, 1962.
[4]	Hadjidimos A. Accelerated over relaxation method[J].Mathematics of Computation , 1978, 32(141): 149–157.DOI:10.1090/S0025-5718-1978-0483340-6
[5]	Arrow K, Hurwicz L, Uzawa H. Studies in nonlinear programming[M]. Stanford: Stanford University Press, 1958: 5.
[6]	Elman H C, Golub G H. Inexact and preconditioned Uzawa algorithms for saddle point problems[J].SIAM Journal on Numerical Analysis , 1994, 31(6): 1645–1661.DOI:10.1137/0731085
[7]	Bramble J H, Pasciak J E, Vassilev A T. Analysis of the inexact Uzawa algorithm for saddle point problems[J].SIAM Journal on Numerical Analysis , 1997, 34(3): 1072–1092.DOI:10.1137/S0036142994273343
[8]	Bai Z Z, Wang Z Q. On parameterized inexact Uzawa methods for generalized saddle point problems[J].Linear Algebra and Its Applications , 2008, 428(11): 2900–2932.
[9]	Zhang J, Shang J. A class of Uzawa-SOR methods for saddle point problems[J].Applied Mathematics and Computation , 2010, 216(7): 2163–2168.DOI:10.1016/j.amc.2010.03.051
[10]	Yang A L, Wu Y J. The Uzawa-HSS method for saddle point problems[J].Applied Mathematics Letters , 2014, 38: 38–42.DOI:10.1016/j.aml.2014.06.018
[11]	Yun J H. Variants of the Uzawa method for saddle point problem[J].Computers & Mathematics with Applications , 2013, 65(7): 1037–1046.
[12]	Young D M. Iterative solution for large systems[M]. New York: Academic Press, 1971.
[13]	Golub G H, Wu X, Yuan J Y. SOR-Like methods for augmented systems[J].BIT Numerical Mathematics , 2001, 41(1): 71–85.DOI:10.1023/A:1021965717530
[14]	Bai Z Z, Parlett B N, Wang Z Q. On generalized successive overrelaxation methods for augmented linear systems[J].Numerische Mathematik , 2005, 102(1): 1–38.DOI:10.1007/s00211-005-0643-0
[15]	Li C, Li Z, Nie Y Y, et al. Generalized AOR method for the augmented system[J].International Journal of Computer Mathematics , 2004, 81(4): 495–504.DOI:10.1080/00207160410001661663
[16]	Darvishi M T, Hessari P. Symmetric SOR method for augmented systems[J].Applied Mathematics and Computation , 2006, 183(1): 409–415.DOI:10.1016/j.amc.2006.05.094
[17]	Zhang G F, Lu Q. On generalized symmetric SOR method for augmented systems[J].Journal of Computational and Applied Mathematics , 2008, 219(1): 51–58.DOI:10.1016/j.cam.2007.07.001
[18]	Bai Z Z, Golub G H, Ng M K. Hermitian and Skew-Hermitian splitting methods for non-Hermitian positive definite linear systems[J].SIAM Journal on Matrix Analysis and Applications , 2003, 24(3): 603–626.DOI:10.1137/S0895479801395458
[19]	Benzi M, Gander M J, Golub G H. Optimization of the Hermitian and Skew-Hermitian splitting iteration for saddle-point problems[J].BIT Numerical Mathematics , 2003, 43(5): 881–900.DOI:10.1023/B:BITN.0000014548.26616.65
[20]	Bai Z Z, Golub G H, Pan J Y. Preconditioned Hermitian and Skew-Hermitian splitting methods for non-Hermitian positive semidefinite linear systems[J].Numerische Mathematik , 2004, 98(1): 1–32.DOI:10.1007/s00211-004-0521-1
[21]	Bai Z Z, Golub G H. Accelerated Hermitian and Skew-Hermitian splitting iteration methods for saddle-point problems[J].IMA Journal of Numerical Analysis , 2007, 27(1): 1–23.