广义Sylvester矩阵方程AX+YA=C一般解的正交投影迭代解法
An Orthogonal Projection Iteration Method for the General Real Solution of Generalized Sylvester Matrix Equations AX+YA=C
DOI: 10.12677/AAM.2023.126283, PDF, HTML, XML, 下载: 154  浏览: 230  科研立项经费支持
作者: 田时宇*, 刘 明:湖南信息学院通识教育学院,湖南 长沙
关键词: 约束矩阵方程正交投影迭代法最佳逼近解极小范数解Constrained Matrix Equations Orthogonal Projection Iteration Method Optimal Approximation Minimum Norm Solution
摘要: 本文讨论了广义Sylvester矩阵方程AX+YA=C的一般实数解及其最佳逼近的正交投影迭代解法,首先利用正交投影及奇异值分解,构造迭代算法,证明了算法的收敛性,得出了收敛速率的估计式;其次给出数值实例,验证了算法的有效性。
Abstract: The general real solution of generalized Sylvester matrix equations AX+YA=C and the orthog-onal projection iteration method to optimal approximation are studied. Firstly, the iterative method is constructed and its convergence is proved by using the theory of orthogonal projection and the singular value decomposition, and the estimation of its convergence rate is obtained; secondly, nu-merical examples are given to verify the validity of the algorithm.
文章引用:田时宇, 刘明. 广义Sylvester矩阵方程AX+YA=C一般解的正交投影迭代解法[J]. 应用数学进展, 2023, 12(6): 2819-2826. https://doi.org/10.12677/AAM.2023.126283

1. 引言

R m × n 表示 m × n 实矩阵的集合, A B 表示矩阵A和B的Kronecker积, v e c ( A ) 表示将矩阵A按行拉直构成的列向量。对 A , B R m × n ,A与B的内积定义为 A , B = t r ( B T A ) ,则由此内积导出的范数 A = t r ( A T A ) 是矩阵A的Frobenius范数。

本文讨论如下问题:

问题I 给定 A R m × n , C R m × n X S 1 R n × n , Y S 2 R m × m ,使得

A X + Y A = C (1)

问题II 设问题I相容,且其解集合为 S E ,给定 X 0 R n × n , Y 0 R m × m ,求解组 [ X ^ , Y ^ ] S E ,使得

X ^ X 0 2 + Y ^ Y 0 2 = min [ X , Y ] S E [ X X 0 2 + Y Y 0 2 ] (2)

广义Sylvester方程在控制论、信号处理、神经网络、模型降阶、图像恢复等领域有着广泛的应用。例如,控制理论及应用领域,在极点配置、观测器设计及构造Sylvester函数中都需要求矩阵方程 A X + Y A = C 的解。正因为矩阵方程的数值解在众多应用中的重要性,近几年国内外众多学者对矩阵方程的数值求解进行了研究,得到了一些有效的数值解法。如文献 [1] [2] [3] [4] 在实矩阵类上不同的方法讨论了 A X + Y A = C 具有一般解、对称解、反对称解的相容性条件和通解表达式。

迄今为止,对矩阵方程 A X + Y A = C 的解及其最佳逼近的迭代算法进行了一些研究,例如顾传青 [5] 给出一种改进的梯度方法,邵新慧 [6] 给出了松弛迭代解法,段复建 [7] 给出了一类Sylvester矩阵方程异类约束解的迭代算法,邓勇 [8] 给出了广义Sylvester矩阵方程反自反解的有限迭代算法,康靖 [9] 给出了改进IO迭代算法,但是对问题I的正交投影迭代算法 [10] 的研究结果尚未见到,并且当前给出的算法都没有有效地进行收敛速率估计。本文利用正交投影的思想构造广义Sylvester矩阵方程一般实数解的迭代算法,并讨论迭代算法的收敛性问题。

2. 问题1的迭代算法

算法1

1) 任取初始值 X 0 , Y 0 ,令 C 0 = C A X 0 Y 0 A

2) 计算 α k = A T C k 2 + C k A T 2 A A T C k + C k A T A 2

3) 令 Δ X k = α k A T C k , Δ Y k = α k C k A T

4) 如果 Δ X k = 0 Δ Y k = 0 ,迭代结束,否则, X k + 1 = X k + Δ X k Y k + 1 = Y k + Δ Y k

5) 令 C k + 1 = C k Δ C k = C k A Δ X k Δ Y k A = C 0 A X k + 1 Y k + 1 A ,则转步骤(2)。

引理1 在迭代算法1中, α k 的选择使得 C k + 1 极小并且使得 C k + 1 Δ C k 相互正交。

证明:根据算法1,我们有

C k + 1 2 = C k Δ C k 2 = C k Δ C k , C k Δ C k = C k 2 2 α k C k , A A T C k + C k A T A + α k 2 A A T C k + C k A T A 2

从上式可知,使得 C k + 1 达到极小的充要条件是

α k = C k , A A T C k + C k A T A A A T C k + C k A T A 2 = A T C k 2 + C k A T 2 A A T C k + C k A T A 2

另一方面,令 C k + 1 , Δ C k = 0 ,我们同样得到 α k = A T C k 2 + C k A T 2 A A T C k + C k A T A 2

引理2 在迭代算法1中,有 C k + 1 2 = C k 2 Δ C k 2

定理1 算法1必定收敛,在问题I中,设矩阵A的最大、最小非零奇异值分别是 σ 1 , σ r ,那么

算法1的收敛速率不小于 0.5 ln ( 1 σ r 2 2 σ 1 2 )

证明:设矩阵A的秩为 r a n k ( A ) = r 。设矩阵A的奇异值分解为 A = U D V T 其中 U , V 都是正交阵, D = [ Σ 0 0 0 ] Σ = d i a g ( σ 1 , σ 2 , , σ r ) σ 1 > σ 2 > > σ r

若问题I相容,则有 C 0 = A X + Y A 成立,而由算法1得 C k = C 0 A X k Y k A ,则

C k = A X + Y A A X k Y k A = U D V T ( X X k ) V V T + U U T ( Y Y k ) U D V T = U ( D V T ( X X k ) V + U T ( Y Y k ) U D ) V T

G = D V T ( X X k ) V + U T ( Y Y k ) U D ,那么 C k = U G V T 其中 G 1 = ( g i j ) g i j R , ( i = 1 , 2 , , m ; j = 1 , 2 , , n ) g i j = 0 , i > r j > r

θ 是矩阵 R k 和矩阵 Δ R k 的夹角,则有

cos ( θ ) = C k , Δ C k C k Δ C k = A T C k 2 + C k A T 2 C k A A T C k + C k A T A = D T G 2 + G D T 2 G D D T G + G D T D

其中

D D T G + G D T D 2 = D D T G 2 + G D T D 2 + 2 t r ( ( D D T G ) T G D T D ) = i = 1 r j = 1 m σ i 4 g i j 2 + i = 1 m j = 1 r σ j 4 g i j 2 + 2 i = 1 r j = 1 r σ i 2 g i j 2 σ j 2 σ 1 2 ( i = 1 r j = 1 m σ i 2 g i j 2 + i = 1 m j = 1 r σ j 2 g i j 2 + i = 1 r j = 1 r σ i 2 g i j 2 + i = 1 r j = 1 r σ j 2 g i j 2 ) σ 1 2 ( i = 1 r j = 1 m σ i 2 g i j 2 + i = 1 m j = 1 r σ j 2 g i j 2 + i = 1 r j = 1 m σ i 2 g i j 2 + i = 1 m j = 1 r σ j 2 g i j 2 ) 2 σ 1 2 ( i = 1 r j = 1 m σ i 2 g i j 2 + i = 1 m j = 1 r σ j 2 g i j 2 )

D T G 2 + G D T 2 = i = 1 r j = 1 m σ i 2 g i j 2 + i = 1 m j = 1 r σ j 2 g i j 2

G 2 = i = 1 r j = 1 r g i j 2 + i = 1 r j = r + 1 m g i j 2 + i = r + 1 m j = 1 r g i j 2

cos ( θ ) i = 1 r j = 1 m σ i 2 g i j 2 + i = 1 m j = 1 r σ j 2 g i j 2 ( i = 1 r j = 1 r g i j 2 + i = 1 r j = r + 1 m g i j 2 + i = r + 1 m j = 1 r g i j 2 ) 1 / 2 ( 2 σ 1 2 ( i = 1 r j = 1 m σ i 2 g i j 2 + i = 1 m j = 1 r σ j 2 g i j 2 ) ) 1 / 2 = ( i = 1 r j = 1 m σ i 2 g i j 2 + i = 1 m j = 1 r σ j 2 g i j 2 ) 1 / 2 ( i = 1 r j = 1 r g i j 2 + i = 1 r j = r + 1 m g i j 2 + i = r + 1 m j = 1 r g i j 2 ) 1 / 2 2 σ 1 σ r ( i = 1 r j = 1 m g i j 2 + i = 1 m j = 1 r g i j 2 ) 1 / 2 ( i = 1 r j = 1 r g i j 2 + i = 1 r j = r + 1 m g i j 2 + i = r + 1 m j = 1 r g i j 2 ) 1 / 2 2 σ 1

= σ r ( 2 i = 1 r j = 1 r g i j 2 + i = 1 r j = r + 1 m g i j 2 + i = r + 1 m j = 1 r g i j 2 ) 1 / 2 ( i = 1 r j = 1 r g i j 2 + i = 1 r j = r + 1 m g i j 2 + i = r + 1 m j = 1 r g i j 2 ) 1 / 2 2 σ 1 σ r ( i = 1 r j = 1 r g i j 2 + i = 1 r j = r + 1 m g i j 2 + i = r + 1 m j = 1 r g i j 2 ) 1 / 2 ( i = 1 r j = 1 r g i j 2 + i = 1 r j = r + 1 m g i j 2 + i = r + 1 m j = 1 r g i j 2 ) 1 / 2 2 σ 1 = σ r 2 σ 1

而由引理2, C k 2 = C k + 1 2 + Δ C k 2 ,我们有 C k + 1 = C k sin ( θ ) C k 1 σ r 2 2 σ 1 2

从上式可知算法1必定收敛,且收敛速率不小于 0.5 ln ( 1 σ r 2 2 σ 1 2 )

引理3 设线性方程组 M y = b 存在解 y 0 R ( M T ) ,则 y 0 必为该线性方程组的唯一的极小范数解。

定理2 若问题I相容,算法1收敛到该问题的极小范数解。

证明:根据算法1,若问题I相容,那么由算法1可以得到问题1的一个迭代解组 [ X , Y ] ,且 X Y 可分别表示为

X = A T H , Y = H A T

下面证明 [ X , Y ] 即为问题I的极小范数解。

将问题I中的矩阵方程(1)两边进行拉直映射运算,并记 v e c ( X ) = x w e c ( Y ) = y v e c ( X ) = x v e c ( Y ) = y v e c ( H ) = h v e c ( C ) = c 则矩阵方程组(1)等价于线性方程组

[ A I I A T ] [ x y ] = c (3)

[ x y ] = [ v e c ( X ) v e c ( Y ) ] = [ v e c ( A T H ) v e c ( H A T ) ] = [ ( A T I ) h ( I A ) h ] = [ A T I I A ] h R ( [ A T I I A ] T )

可知 [ x y ] T 是线性方程组(3)的极小范数解,由于拉直映射是同构的,故 [ X Y ] 也是问题I的极小范数解组。

3. 问题2的解

若问题1相容,则解集 S E 为一非空集,当 [ X , Y ] S E 时,由 A X + Y A = C 可得

A ( X X 0 ) + ( Y Y 0 ) A = C A X 0 Y 0 A

X ˜ = X X 0 , Y ˜ = Y Y 0 , C ˜ = C A X 0 Y 0 A ,则问题2等价于求相容矩阵方程 A X ˜ + Y ˜ A = C ˜ 的极小范数解组 [ X ˜ Y ˜ ]

利用迭代算法1,取特殊初值 X ˜ = A T H , Y ˜ = H A T 时,可得矩阵方程 A X ˜ + Y ˜ A = C ˜ 的唯一的极小范数解组 [ X ˜ Y ˜ ] ,从而得到问题2的解 X ^ = X ˜ + X 0 Y ^ = Y ˜ + Y 0

4. 数值实例

用本文构造的迭代算法求矩阵方程 A X + Y A = C 的一般解及矩阵 X 0 的最佳逼近矩阵。

例设

A = [ 0.8147 0.0975 0.1576 0.1419 0.9058 0.2785 0.9706 0.4218 0.1270 0.5469 0.9572 0.9157 0.9134 0.9575 0.4854 0.7922 0.6324 0.9649 0.8003 0.9595 ]

C = [ 2.2028 2.3979 2.2546 2.0807 3.5648 2.9730 2.6473 2.9713 2.5278 2.2763 1.8380 2.5673 3.6031 3.6455 2.7191 2.7225 4.0793 3.8372 2.8939 3.3182 ]

求问题I的解 X , Y

若问题I相容,给定

X 0 = [ 0.8308 0.2858 0.5678 0.7792 0.5853 0.7572 0.0759 0.9340 0.5497 0.7537 0.0540 0.1299 0.9172 0.3804 0.5308 0.5688 ]

Y 0 = [ 0.4694 0.3112 0.6541 0.2290 0.9961 0.0119 0.5285 0.6892 0.9133 0.0782 0.3371 0.1656 0.7482 0.1524 0.4427 0.1622 0.6020 0.4505 0.8258 0.1067 0.7943 0.2630 0.0838 0.5383 0.9619 ]

求问题II的解

① 由算法1,取初值 X 0 = 0 , Y 0 = 0 ;解得问题I的解为:

X = [ 0 .72470 .87860 .58190 .2200 0 .59310 .85270 .35090 .1406 0 .77060 .22790 .23260 .9298 0 .63660 .60670 .30610 .5422 ]

Y = [ 0 .09240 .40820 .55870 .48490 .6070 0 .38980 .55350 .37980 .57090 .5349 0 .08040 .20850 .44620 .26990 .3700 0 .16010 .45540 .27990 .59650 .5881 0 .33870 .44650 .34570 .74110 .6636 ]

② 由算法1,若初值取特殊矩阵 X 0 = A T H , Y 0 = H A T ,其中

H = [ 0 .25510 .54720 .25430 .1966 0 .50600 .13860 .81430 .2511 0 .69910 .14930 .24350 .6160 0 .89090 .25750 .92930 .4733 0 .95930 .84070 .35000 .3517 ]

则解得问题I的极小范数解为:

X = [ 0 .72470 .87860 .58190 .2200 0 .59310 .85270 .35090 .1406 0 .77060 .22790 .23260 .9298 0 .63660 .60670 .30610 .5422 ]

Y = [ 0 .09240 .40820 .55870 .48490 .6070 0 .38980 .55350 .37980 .57090 .5349 0 .08040 .20850 .44620 .26990 .3700 0 .16010 .45540 .27990 .59650 .5881 0 .33870 .44650 .34570 .74110 .6636 ]

由上述的①,②可知,若问题I相容,则算法1收敛到该问题的唯一极小范数解。

X ˜ = X X 0 , Y ˜ = Y Y 0 , C ˜ = C A X 0 Y 0 A ;取初值 X ˜ 0 = 0 , Y ˜ 0 = 0 ,由算法1可得

X ˜ = [ 0.0662 0.5219 0.1215 0.5786 0.1728 0.3681 0.2329 0.4430 0.1344 0.5747 0.2183 0.6051 0.1063 0.1151 0.0994 0.1558 ]

Y ˜ = [ 0.1995 0.1064 0.1934 0.0145 0.0013 0.2672 0.0680 0.0600 0.0314 0.0378 0.0228 0.0982 0.0527 0.0197 0.0560 0.0641 0.0803 0.1573 0.1054 0.0436 0.1421 0.1111 0.1843 0.1488 0.0066 ]

则由此可得问题II的解为:

X = X ˜ + X 0 = [ 0.7646 0.8077 0.6893 0.2006 0.4125 1.1253 0.3088 0.4910 0.6841 0.1790 0.1643 0.7350 0.8109 0.4955 0.6302 0.4130 ]

Y = Y ˜ + Y 0 = [ 0 .2699 0 .2048 0 .4607 0 .2435 0 .9974 0 .2791 0 .5965 0 .7492 0 .9447 0 .0404 0 .3599 0 .0674 0 .6955 0 .1327 0 .3867 0 .0981 0 .5217 0 .2932 0 .9312 0 .1503 0 .6522 0 .3741 0 .2681 0 .3895 0 .9685 ]

5. 结论

本文对正交投影迭代算法进行改进,构造了求解广义Sylvester矩阵方程 A X + Y A = C 一般实数解的正交投影迭代算法,利用奇异值分解,证明了该算法的收敛性,得出了收敛速率的估计式 0.5 ln ( 1 σ r 2 2 σ 1 2 ) ,为了验证该算法的有效性,进行数值实例,数值实验的结果表明该算法是有效的。本文创新点如下:

1) 对正交投影迭代算法进行改进得到求解广义Sylvester矩阵方程 A X + Y A = C 一般实数解的正交投影迭代算法;

2) 利用奇异值分解、不等式放缩性质证明了该算法的收敛性,得出了收敛速率的估计式为 0.5 ln ( 1 σ r 2 2 σ 1 2 )

3) 进行了数值验证,验证了该算法的有效性。

基金项目

湖南省2022年普通高等学校教学改革研究重点项目(编号:HNJG-2022-0381),湖南信息学院2022年校级科研课题(编号:XXY022YB09)。

参考文献

[1] Hoskins, W.D., Meek, D.S. and Walton, D.J. (1977) The Numerical Solution of the Matrix Equation XA + AY = F. BIT Numerical Mathematics, 17, 184-190.
https://doi.org/10.1007/BF01932289
[2] Hoskins, W.D., Meek, D.S. and Walton, D.J. (1979) High Order Iterative Methods of the Solution of the Matrix Equation XA + AY = F. Linear Algebra and Its Applications, 23, 121-139.
https://doi.org/10.1016/0024-3795(79)90097-1
[3] Chang, X.W. and Wang, J.S. (1993) The Symmetric Solution of the Matrix Equation , and . Linear Algebra and Its Applications, 179, 171-189.
https://doi.org/10.1016/0024-3795(93)90328-L
[4] 邓远北, 胡锡炎. 一类广义Sylvester方程的反对称最小二乘解及其最佳逼近[J]. 系统科学与数学, 2004, 24(3): 382-388.
[5] 顾传青, 蒋祥龙. 求解Sylvester矩阵方程的一种改进的梯度方法[J]. 应用数学与计算数学学报, 2014, 28(4): 432-439.
[6] 邵新慧, 彭程. 一类Sylvester矩阵方程的迭代解法[J]. 东北大学学报(自然科学版), 2017, 38(6): 909-912.
[7] 段复建, 原腾. 一类Sylvester矩阵方程异类约束解的迭代算法[J]. 重庆理工大学学报(自然科学), 2021, 35(6): 247-255.
[8] 邓勇. 关于广义Sylvester矩阵方程反自反解的有限迭代算法[J]. 东北师大学报(自然科学版), 2022, 54(1): 34-43.
[9] 康靖, 马昌凤. Sylvester矩阵方程的改进IO迭代算法[J]. 井冈山大学学报(自然科学版), 2023, 44(2): 1-8.
[10] 郭孔华. 求解约束矩阵方程的正交投影迭代法研究[D]: [博士学位论文]. 长沙: 湖南大学, 2007.