Advances in Applied Mathematics
Vol. 10  No. 11 ( 2021 ), Article ID: 46244 , 9 pages
10.12677/AAM.2021.1011383

求解复对称线性系统的改进 EPGS迭代法

隆雪琴,张乃敏

温州大学数学系,浙江 温州

收稿日期:2021年10月2日;录用日期:2021年10月23日;发布日期:2021年11月4日

摘要

本文主要讨论求解复对称线性系统的EPGS迭代法。基于矩阵块三角分裂思想,我们对该方法引入一个加速参数后,得到了改进EPGS迭代法(IEPGS)。我们对IEPGS迭代法进行收敛性分析,给出了收敛性条件,并得到了最优迭代参数及相应的最优收敛因子。数值实验的结果进一步表明了该方法的有效性。

关键词

复对称线性系统,IEPGS迭代法,参数,收敛

Improved EPGS Iterative Method for Solving Complex Symmetric Linear Systems

Xueqin Long, Naimin Zhang

Department of Mathematics, Wenzhou University, Wenzhou Zhejiang

Received: Oct. 2nd, 2021; accepted: Oct. 23rd, 2021; published: Nov. 4th, 2021

ABSTRACT

In this paper, we mainly discuss the EPGS iterative method for solving complex symmetric linear systems. Based on the idea of triangular splitting of matrix blocks, we introduce an acceleration parameter to this method and obtain the improved EPGS iterative method (IEPGS). We analyze the convergence of IEPGS iteration method, give the condition of convergence, and obtain the optimal iteration parameters and the corresponding optimal convergence factors. Numerical experiments show that the method is effective.

Keywords:Complex Symmetric Linear System, IEPGS Iterative Method, Parameter, Convergence

Copyright © 2021 by author(s) and Hans Publishers Inc.

This work is licensed under the Creative Commons Attribution International License (CC BY 4.0).

http://creativecommons.org/licenses/by/4.0/

1. 引言

本文考虑求解如下复对称线性系统:

A u = ( W + i T ) u = b , A C n × n , W , T R n × n , u , b C n (1)

其中 W , T R n × n 分别为对称正定矩阵与对称半正定矩阵, i = 1 为虚数单位。众所周知复对称线性系统(1)可以等价地表示为以下2 × 2块实值线性系统:

B z [ W T T W ] [ x y ] = [ f g ] p , B R 2 n × 2 n , z , p R 2 n , x , y , f , g R n (2)

上述复对称线性系统广泛应用于科学计算与工程应用中的许多实际问题,例如漫散光学层析成像 [1]、结构动力学 [2]、电力建模 [3]、分子散射 [4] 等问题。复对称系统的求解近年来也得到了很大的发展。为了有效求解复对称线性系统(1)和(2),最近,Li和Lu [5] 提出等距参数化的Guass-Seidel (EPGS)迭代法。本文继续讨论EPGS迭代法,我们通过一个矩阵块三角分裂,并添加一个加速参数后,得到改进的EPGS (IEPGS)迭代法,以期提高该方法的收敛速度。本文使用的一些基本符号如下: ( ) T 表示矩阵或向量的转置, σ ( ) ρ ( ) 分别表示矩阵的谱集和谱半径。

2. IEPGS迭代法

在本节中,我们将给出IEPGS迭代法的迭代格式。

基于预处理思想,首先引入一个Givens正交矩阵 G θ

G θ = [ cos ( θ ) I sin ( θ ) I sin ( θ ) I cos ( θ ) I ] ,

其中,I表示n维单位矩阵, θ ( 0 , π 2 ) G θ 1 = G θ T 。对(2)式两边同时左乘矩阵 G θ T 得到

B ˜ z = G θ T B z = [ W ˜ θ T ˜ θ T ˜ θ W ˜ θ ] [ x y ] = G θ T [ f g ] = [ f ˜ g ˜ ] = p ˜ (3)

其中, W ˜ θ = cos ( θ ) W + sin ( θ ) T T ˜ θ = cos ( θ ) T sin ( θ ) W f ˜ = cos ( θ ) f + sin ( θ ) g g ˜ = cos ( θ ) g sin ( θ ) f

显然地, W ˜ θ T ˜ θ 分别为对称正定矩阵和对称矩阵。

算法1 (EPGS迭代法):给定初始向量 x ( 0 ) R n y ( 0 ) R n ,按如下迭代过程计算直至迭代序列 { ( x ( k ) T , y ( k ) T ) } k = 0 收敛到(1)的精确解:

{ W ˜ θ x k + 1 = T ˜ θ y k + f ˜ , W ˜ θ y k + 1 = T ˜ θ x k + 1 + g ˜ (4)

受到文章 [6] [7] 的启发,为了提高EPGS迭代方法的收敛速度,基于矩阵块三角分裂思想,我们引入一个参数来获得其改进版本。对(3)式中的系数矩阵 B ˜ 做如下分裂:

B ˜ = [ α W ˜ θ 0 T ˜ θ W ˜ θ ] [ ( α 1 ) W ˜ θ T ˜ θ 0 0 ] = M ( θ , α ) N ( θ , α )

其中, α 是一个正实常数。由该分裂我们可得如下迭代格式,也即IEPGS迭代:

[ x ( k + 1 ) y ( k + 1 ) ] = Γ ( θ , α ) [ x ( k ) y ( k ) ] + ϕ ( θ , α ) [ f ˜ g ˜ ] , (5)

其中, ϕ ( θ , α ) = M ( θ , α ) 1 且迭代矩阵 Γ ( θ , α ) = M 1 ( θ , α ) N ( θ , α )

Γ ( θ , α ) = [ 1 α W ˜ θ 1 0 1 α W ˜ θ 1 T ˜ θ W ˜ θ 1 W ˜ θ 1 ] [ ( α 1 ) W ˜ θ T ˜ θ 0 0 ] = [ ( 1 1 α ) I 1 α W ˜ θ 1 T ˜ θ ( 1 1 α ) W ˜ θ 1 T ˜ θ 1 α W ˜ θ 1 T ˜ θ W ˜ θ 1 T ˜ θ ] (6)

上述迭代格式(5)也可以表示为如下算法2:

算法2 (IEPGS迭代法):给定任意初始向量 x ( 0 ) R n y ( 0 ) R n ,按如下迭代过程计算直至迭代序列 { ( x ( k ) T , y ( k ) T ) } k = 0 收敛到(1)的精确解:

{ α W ˜ θ x k + 1 = ( α 1 ) W ˜ θ x k + T ˜ θ y k + f ˜ , W ˜ θ y k + 1 = T ˜ θ x k + 1 + g ˜ . (7)

注意到,当 α = 1 时,上述IEPGS迭代方法就退化为EPGS迭代方法,从而IEPGS是EPGS的一种推广。

3. IEPGS迭代法的收敛性分析

定理 1. W , T R n × n 分别是对称正定阵与对称半正定阵, α 是一个正常数且 θ ( 0 , π 2 ) ,则

1) 迭代矩阵 Γ ( θ , α ) 的特征值 λ 取值如下:

λ = 0 , λ = 1 1 α , λ = 1 1 + η i 2 α (8)

其中 η i ( 1 i r ) 是矩阵 S ˜ = W ˜ θ 1 T ˜ θ 的特征值, r = r a n k ( W ˜ θ 1 T ˜ θ ) n

2) 谱半径 ρ ( Γ ( θ , α ) ) < 1 当且仅当 α > 1 + η max 2 2

3) ρ ( Γ ( θ , α ) ) 可表示为:

ρ ( Γ ( θ , α ) ) = max { | 1 1 + η max 2 α | , | 1 1 + η min 2 α | , | 1 1 α | } (9)

此处以及后文中, η min = min η i σ ( S ˜ ) { | η i | } η max = max η i σ ( S ˜ ) { | η i | }

证明:因为 W ˜ θ 是对称正定矩阵,我们可以定义

P 1 = [ W ˜ θ 1 2 0 0 W ˜ θ 1 2 ] (10)

那么(6)式中的EPGS迭代矩阵 Γ ( θ , α ) 相似于矩阵 Γ ^ ( θ , α ) = P 1 1 Γ ( θ , α ) P 1 ,即

Γ ^ ( θ , α ) = [ W ˜ θ 1 2 0 0 W ˜ θ 1 2 ] [ ( 1 1 α ) I 1 α W ˜ θ 1 T ˜ θ ( 1 1 α ) W ˜ θ 1 T ˜ θ 1 α W ˜ θ 1 T ˜ θ W ˜ θ 1 T ˜ θ ] [ W ˜ θ 1 2 0 0 W ˜ θ 1 2 ] = [ ( 1 1 α ) I 1 α Q ( 1 1 α ) Q 1 α Q 2 ]

其中,矩阵 Q = W ˜ θ 1 2 T ˜ θ W ˜ θ 1 2 。设实对称矩阵Q的特征值分解为

Q = U T [ Σ r 0 0 0 ] U , (11)

其中, U R n × n 是一个正交矩阵且 Σ r = d i a g ( η 1 , η 2 , , η r ) 是一个以矩阵Q的非零特征值 η 1 η 2 η r 为对角元的对角矩阵。此处, r = r a n k ( Q ) n 。令

P 2 = [ U T 0 0 U T ] (12)

那么,容易得到 Γ ^ ( θ , α ) 相似于矩阵 Γ ˜ ( θ , α ) = P 2 1 Γ ^ ( θ , α ) P 2 ,即

Γ ˜ ( θ , α ) = [ U 0 0 U ] [ ( 1 1 α ) I 1 α Q ( 1 1 α ) Q 1 α Q 2 ] [ U T 0 0 U T ] = [ ( 1 1 α ) I 1 α U Q U T ( 1 1 α ) U Q U T 1 α U Q U T U Q U T ]

再利用(11)式可得

Γ ˜ ( θ , α ) = [ ( 1 1 α ) I r 0 1 α Σ r 0 0 ( 1 1 α ) I n r 0 0 ( 1 1 α ) Σ r 0 1 α Σ r 2 0 0 0 0 0 ]

注意到,矩阵 Γ ˜ ( θ , α ) 中的非零块矩阵都为对角矩阵。下求矩阵 Γ ˜ ( θ , α ) 的特征值:

| λ I Γ ˜ ( θ , α ) | = | ( λ 1 + 1 α ) I r 0 1 α Σ r 0 0 ( λ 1 + 1 α ) I n r 0 0 ( 1 1 α ) Σ r 0 λ I r 1 α Σ r 2 0 0 0 0 λ I n r | = ( 1 ) r λ r [ ( λ 1 + 1 α ) α η 1 + η 1 ] [ ( λ 1 + 1 α ) α η r + η r ] ( λ 1 + 1 α ) n r ( 1 ) r η 1 α η r α λ n r = λ n [ ( λ 1 + 1 α ) α η 1 + η 1 ] [ ( λ 1 + 1 α ) α η r + η r ] ( λ 1 + 1 α ) n r η 1 α η r α

容易得到矩阵 Γ ˜ ( θ , α ) 的特征值为如下三种形式:

λ i 1 = 0 , λ i 2 = 1 1 + η i 2 α , λ i 3 = 1 1 α , (13)

其中 i 1 = 1 , 2 , , n i 2 = 1 , 2 , , r i 3 = 1 , 2 , , n r 。故结论(1)得证。

下证结论(2),由(13)式知:

ρ ( Γ ( θ , α ) ) < 1 { | λ i 2 | < 1 | λ i 3 | < 1 { | 1 1 + η i 2 α | < 1 | 1 1 α | < 1 { 1 < 1 1 + η i 2 α < 1 1 < 1 1 α < 1 { α > 1 + η i 2 2 α > 1 2 α > 1 + η max 2 2

故结论(2)得证。下证结论(3)。

x 0 时, S ( x ) = 1 1 + x α 关于x单调递减知

max η i | 1 1 + η i 2 α | = max { | 1 1 + η max 2 α | , | 1 1 + η min 2 α | } ,

故迭代矩阵 Γ ( θ , α ) 的谱半径为:

ρ ( Γ ( θ , α ) ) = max { | 1 1 + η max 2 α | , | 1 1 + η min 2 α | , | 1 1 α | }

定理1得证。

上述定理给出了IEPGS迭代法的收敛条件。下面我们讨论使 ρ ( Γ ( θ , α ) ) 极小化的拟最优参数 α * θ * ,以及相应的最优收敛因子。在此之前,我们先给出以下引理。

引理2 [5]. W , T R n × n 分别是对称正定阵与对称半正定阵, θ ( 0 , π 2 ) ,则EPGS迭代法的迭代矩阵 Η ( θ ) = [ 0 W ˜ θ 1 T ˜ θ 0 W ˜ θ 1 T ˜ θ W ˜ θ 1 T ˜ θ ] 有n重0特征值,且 Η ( θ ) 的其余n个特征值 λ 满足以下方程:

λ + ( μ i cos θ sin θ ) 2 ( cos θ + μ i sin θ ) 2 = 0 , ( 1 i n ) ,

其中, μ i ( 1 i n ) 是矩阵 S = W 1 T 的特征值。

μ max , μ min 分别为矩阵 S = W 1 T 的最大和最小特征值,则有如下定理成立。

定理3. W , T R n × n 分别是对称正定阵与对称半正定阵,则对给定的 θ ( 0 , π 2 ) ,极小化 ρ ( Γ ( θ , α ) ) 的最优参数 α * 取值为 α * = 2 + η max 2 2 ;进一步极小化 ρ ( Γ ( θ , α * ) ) 的最优参数 θ * 取值为

θ = arctan μ min μ max 1 + ( 1 + μ min 2 ) ( 1 + μ max 2 ) μ min + μ max ,

且相应的最优收敛因子为:

ρ ( Γ ( θ * , α * ) ) = η max 2 2 + η max 2 (14)

证明:对给定的 θ ,首先考虑最优参数 α * 的选取,即

α * = arg min α ρ ( Γ ( θ , α ) ) = arg min α max { | 1 1 + η max 2 α | , | 1 1 + η min 2 α | , | 1 1 α | }

f η ( α ) = 1 1 + η 2 α ,其中 η 分别取 0 , η min , η max

则有

α * = arg min α max η { | f η ( α ) | } .

由于 d d α f η ( α ) = 1 + η 2 α 2 > 0 ,所以 f η ( α ) 关于 α 单调递增, f η ( α ) 关于 α 单调递减。从而函数 | f η ( α ) | ( η = 0 , η min , η max )关于 α 先递减后递增,它们各自的图像如图1所示。

Figure 1. The graph of | f η ( α ) | , where η = 0 , η min , η max

图1. | f η ( α ) | 的图像,其中 η = 0 , η min , η max

易知图1交点纵坐标取得最大值时相应的横坐标即为最优参数 α * ,经过简单计算得到 α * 满足 1 1 α * = 1 + η max 2 α * 1 α * = 2 + η max 2 2 ,此时, min α ρ ( Γ ( θ , α ) ) = ρ ( Γ ( θ , α * ) ) = η max 2 2 + η max 2

下面我们考虑最优参数 θ * 极小化 ρ ( Γ ( θ , α * ) )

g ( η max 2 ) = η max 2 2 + η max 2 ,则 g ( η max 2 ) 关于 η max 2 单调递增,故要极小化 ρ ( Γ ( θ , α * ) ) ,我们选取最优参数 θ * 极小化 η max 2 即可。由 η i = μ i cos θ sin θ cos θ + μ i sin θ

θ * = arg min θ η max 2 ( θ ) = arg min θ max μ i σ ( W 1 T ) ( cos θ μ i cos θ ) 2 ( sin θ + μ i sin θ ) 2 .

由文 [5] 中定理3.3得

θ * = arctan μ min μ max 1 + ( 1 + μ min 2 ) ( 1 + μ max 2 ) μ min + μ max ,

且取此最优参数 θ * 后,依然有 ρ ( Γ ( θ * , α * ) ) = η max 2 2 + η max 2 。定理3得证。

注:由引理2知EPGS迭代法的迭代矩阵为 Η ( θ ) ,由文 [5] 中定理3.3知其最优收敛因子为 ρ ( Η ( θ * ) ) = ( μ max cos θ * sin θ * ) 2 ( cos θ * + μ max sin θ * ) 2 = η max 2 ,由于不等式 η max 2 2 + η max 2 < η max 2 显然成立。故有 ρ ( Γ ( θ * , α * ) ) < ρ ( H ( θ * ) ) 成立,这表明IEPGS迭代法的收敛速度要快于EPGS迭代法。

4. 数值实验

这一节,我们通过数值试验来说明IEPGS迭代法求解复对称线性系统的有效性。本文所有实验都在内存为4G,操作系统为Windows10的电脑上完成,计算软件为MATLAB (R2018b)。在下面的实验中,取零向量 ( ( x 0 ) T , ( y 0 ) T ) T 为初始迭代向量,记 ( ( x ( k ) ) T , ( y ( k ) ) T ) T 为当前的迭代向量, I T 表示当前迭代步数, R E S 表示残差比,其停机准则为:

R E S = b A u ( k ) 2 b 2 10 9

算例1 [6] 复对称线性系统(1)的形式如下:

[ ( ω 2 M + K ) + i ( ω C V + C H ) ] u = b (15)

其中M和K分别是惯性矩阵和刚度矩阵; C V C H 分别是粘滞阻尼矩阵和迟滞阻尼矩阵; ω 是一个驱动圆频率常数。在本例子中,我们采用 M = I C V = 10 I C H = μ K μ = 0.02 为阻尼系数。

K = I V m + V m I ,其中 V m = h 2 t r i d i a g ( 1 , 2 , 1 ) R m × m n = m 2 。此外,设置 ω = π 并选择右侧向量 b = ( 1 + i ) A 1 ,1表示分量全为1的n维列向量,并通过两边同时乘以 h 2 来正规化(15)式。

通过对m的不同取值,能得到不同大小的系数矩阵。表1是关于MHSS [8]、E-HS [9]、EPGS和IEPGS迭代法的实验参数选取。图2图3对应表1中的实验参数,给出了各迭代法的相对残差曲线。

通过比较图2图3中各迭代方法的递减速度,我们可以看出IEPGS的迭代收敛速度更快,对求解复对称线性系统更有效。

Table 1. Selection of experimental parameters for MHSS, E-HS, EPGS and IEPGS iterative methods

表1. MHSS、E-HS、EPGS和IEPGS迭代法的实验参数选取

Figure 2. Relative residual curves with m = 16 (left) and m = 32 (right)

图2. m = 16 (左)和m = 32 (右)时,MHSS、E-HS、EPGS和IEPGS迭代法的相对残差曲线

Figure 3. Relative residual curves with m = 64 (left) and m = 96 (right)

图3. m = 64 (左)和m = 96 (右)时,MHSS、E-HS、EPGS和IEPGS迭代法的相对残差曲线

文章引用

隆雪琴,张乃敏. 求解复对称线性系统的改进EPGS迭代法
Improved EPGS Iterative Method for Solving Complex Symmetric Linear Systems[J]. 应用数学进展, 2021, 10(11): 3623-3631. https://doi.org/10.12677/AAM.2021.1011383

参考文献

  1. 1. Arridge, S.R. (1999) Optical Tomography in Medical Imaging. Inverse Problems, 15, 41-93. https://doi.org/10.1088/0266-5611/15/2/022

  2. 2. Feriani, A., Perotti, F. and Simoncini, V. (2000) Iterative System Solvers for the Frequency Analysis of Linear Mechanical Systems. Computer Methods in Applied Mechanics and Engineering, 190, 1719-1739. https://doi.org/10.1016/S0045-7825(00)00187-0

  3. 3. Howle, V.E., and Vavasis, S.A. (2005) An Iterative Method for Solving Complex-Symmetric Systems Arising in Electrical Power Modeling. SIAM Journal on Matrix Analysis and Applications, 26, 1150-1178. https://doi.org/10.1137/S0895479800370871

  4. 4. Poirier, B.(2000) Efficient Preconditioning Scheme for Block Partitioned Matrices with Structured Sparsity. Numerical Linear Algebra with Applications, 7, 715-726. https://doi.org/10.1002/1099-1506(200010/12)7:7/8<715::AID-NLA220>3.0.CO;2-R

  5. 5. Li, X.A. and Lu, J. (2020) An Equidistant Parameterized Gauss-Seidel Iteration Method for a Class of Block Two-by-Two Linear Systems. Computational and Applied Mathematics, 39, 292. https://doi.org/10.1007/s40314-020-01341-1

  6. 6. Li, X.A., Zhang, W.H. and Wu, Y.J. (2018) On Symmetric Block Triangular Splitting Iteration Method for a Class of Complex Symmetric System of Linear Equations. Applied Mathematics Letters, 79, 131-137. https://doi.org/10.1016/j.aml.2017.12.008

  7. 7. Salkuyeh, D.K., Hezari, D. and Edalatpour, V. (2014) Generalized SOR Iterative Method for a Class of Complex Symmetric Linear System of Equations. International Journal of Computer Mathematics, 92, 802-815. https://doi.org/10.1080/00207160.2014.912753

  8. 8. Bai, Z.Z., Benzi, M. and Chen, F. (2010) Modified HSS Iteration Methods for a Class of Complex Symmetric Linear Systems. Computing, 87, 93-111. https://doi.org/10.1007/s00607-010-0077-0

  9. 9. Li, C.L. and Ma, C.F. (2018) On Euler-Extrapolated Hermitian/Skew-Hermitian Splitting Method for Complex Symmetric Linear Systems. Applied Mathematics Letters, 86, 42-48. https://doi.org/10.1016/j.aml.2018.06.016

期刊菜单