52．gmres函数

在MATLAB中，提供了gmres函数实现广义最小残差法。函数的语法格式为：

x＝gmres（A，b）：针对x对线性方程组A∗x＝b求解。n×n系数矩阵A必须是方阵，并且应为大型稀疏矩阵。列向量b的长度必须为n。A可以是函数句柄afun，这样，afun（x）返回A∗x。对于此语法，gmres不会重新启动；最大迭代次数为min（n，10）。

如果gmres收敛，则会显示一条有关该结果的消息。如果gmres无法在达到最大迭代次数后收敛或出于任何原因暂停，则会输出一条显示相对残差norm（b－A∗x）/norm（b）以及该方法停止或失败时所达到的迭代数的警告消息。

gmres（A，b，restart）：每restart次内迭代重新启动该方法一次。最大外迭代次数为min（n/restart，10）。最大总迭代次数为restart∗min（n/restart，10）。如果restart为n或[]，则gmres不重新启动，最大总迭代次数为min（n，10）。

gmres（A，b，restart，tol）：指定该方法的容差。如果tol为[]，gmres使用默认值1e－6。

gmres（A，b，restart，tol，maxit）：指定最大外迭代次数，即，总迭代次数不超过restart∗maxit。如果maxit为[]，gmres使用默认值min（n/restart，10）。如果restart为n或[]，则最大总迭代次数为maxit（而非restart∗maxit）。

gmres（A，b，restart，tol，maxit，M）和gmres（A，b，restart，tol，maxit，M1，M2）：使用预设子条件M或M＝M1∗M2，并高效求解关于x的方程组inv（M）∗A∗x＝inv（M）∗b。如果M为[]，gmres不会应用预设子条件。M可以是函数句柄mfun，这样，mfun（x）返回M\x。

gmres（A，b，restart，tol，maxit，M1，M2，x0）：指定第一个初始估计值。如果x0为[]，gmres使用默认值（即全部为零的向量）。

[x，flag]＝gmres（A，b，…）：也返回一个收敛标志，其取值及含义见表1-5。

如果flag不为0，返回的解x具有在所有迭代中最小的范数残差。如果指定flag输出，则不会显示消息。

[x，flag，relres]＝gmres（A，b，…）：还返回相对残差norm（b－A∗x）/norm（b）。如果flag＝0，则relres≤tol。第三个输出relres是预调节方程组的相对残差。

[x，flag，relres，iter]＝gmres（A，b，…）：还返回计算x时所达到的外和内迭代数，其中0≤iter（1）≤maxit且0≤iter（2）≤restart。

[x，flag，relres，iter，resvec]＝gmres（A，b，…）：还返回每次内迭代中的残差范数的向量。这些是预调节方程组的残差范数。

【例1-53】此示例演示如何使用预设子条件并重新启动gmres。

     %加载west0479,它是一个非对称的479×479实稀疏矩阵
     load west0479;
     A = west0479;
     %定义b以使实际解是全为1的向量
     b = full(sum(A,2));
     %构造一个不完全LU预设子条件
     [L,U] = ilu(A,struct('type','ilutp','droptol',1e-6));

使用重新启动的gmres的好处是限制执行该方法所需的内存量。如果不重新启动，gmres需要存储maxit向量来保存基本Krylov子空间。此外，gmres必须在每一步与以前的所有向量正交。重新启动可限制使用的工作区间以及每次外迭代执行的工作量。请注意，即使预调节的gmres在以上6次迭代中收敛，该算法也允许多达20个基向量，因此应事先分配所有这些空间。

     %请执行gmres(3)、gmres(4)和gmres(5)
     tol = 1e-12;
     maxit = 20;
     re3 = 3;
     [x3,fl3,rr3,it3,rv3] = gmres(A,b,re3,tol,maxit,L,U);
     re4 = 4;
     [x4,fl4,rr4,it4,rv4] = gmres(A,b,re4,tol,maxit,L,U);
     re5 = 5;
     [x5,fl5,rr5,it5,rv5] = gmres(A,b,re5,tol,maxit,L,U);

fl3、fl4和fl5均为0，因为在每种情况下重新启动的gmres使相对残差趋向小于1e－12的规定公差。

     %下列绘图显示了每个重新启动的gmres方法的收敛历史记录。gmres(3)收敛于外迭代5、
     %内迭代3(it3 = [5,3]),这将与外迭代6、内迭代0相同,因此最终刻度线上的标记是6
     subplot(3,1,1);semilogy(1:1/3:6,rv3/norm(b),'-o');
     h1 = gca;
     h1.XTick = [1:1/3:6];
     h1.XTickLabel = ['1';'';'';'2';'';'';'3';'';'';'4';'';'';'5';'';'';'6';];
     title('gmres(3)')
     xlabel('迭代次数');
     ylabel('相对残差');
     subplot(3,1,2);semilogy(1:1/4:3,rv4/norm(b),'-o');
     h2 = gca;
     h2.XTick = [1:1/4:3];
     h2.XTickLabel = ['1';'';'';'';'2';'';'';'';'3'];
     title('gmres(4)')
     xlabel('迭代次数');
     ylabel('相对残差');
     subplot(3,1,3);semilogy(1:1/5:2.8,rv5/norm(b),'-o');
     h3 = gca;
     h3.XTick = [1:1/5:2.8];
     h3.XTickLabel = ['1';'';'';'';'';'2';'';'';'';''];
     title('gmres(5)')
     xlabel('迭代次数');
     ylabel('相对残差');

运行程序，效果如图1-21所示。