Main Content

本页对应的英文页面已更新,但尚未翻译。 若要查看最新内容,请点击此处访问英文页面。

equilibrate

缩放矩阵以改善条件

说明

示例

[P,R,C] = equilibrate(A) 置换和重新缩放矩阵 A,使新矩阵 B = R*P*A*C 的对角线上元素的模为 1,非对角线上元素的模不大于 1。

示例

全部折叠

平衡具有较大的条件数的矩阵,以提高用迭代求解器 gmres 求解线性方程组的效率和稳定性。

加载 west0479 矩阵,这是一个实数值 479×479 稀疏矩阵。使用 condest 计算矩阵的估计条件数。

load west0479
A = west0479;
c1 = condest(A)
c1 = 1.4244e+12

尝试使用 gmres 求解线性方程组 Ax=b,迭代 450 次,容差为 1e-11。指定五个输出,以便 gmres 在每次迭代中返回解的残差范数(使用 ~ 隐藏不需要的输出)。在半对数图中绘制残差范数。该图显示,gmres 无法实现合理的残差范数,因此对 x 计算得出的解不可靠。

b = ones(size(A,1),1);
tol = 1e-11;
maxit = 450;
[x,flx,~,~,rvx] = gmres(A,b,[],tol,maxit);
semilogy(rvx)
title('Residual Norm at Each Iteration')

使用 equilibrate 置换和重新缩放 A。创建新矩阵 B = R*P*A*C,它具有更好的条件数且对角线上元素只有 1 和 -1。

[P,R,C] = equilibrate(A);
B = R*P*A*C;
c2 = condest(B)
c2 = 5.1036e+04

使用 equilibrate 返回的输出,您可以将问题 Ax=b 重新表示为 By=d,其中 B=RPACd=RPb。在此形式中,您可以使用 x=Cy 还原原始方程组的解。

使用 gmres 求解 By=d 以得出 y,然后重新绘制每次迭代的残差范数。该图显示,对矩阵进行平衡处理提高了问题的稳定性,gmres 在不到 200 次迭代中收敛于期望容差 1e-11

d = R*P*b;
[y,fly,~,~,rvy] = gmres(B,d,[],tol,maxit);
hold on
semilogy(rvy)
legend('Original', 'Equilibrated', 'Location', 'southeast')
title('Relative Residual Norms (No Preconditioner)')
hold off

使用预条件子改进解

获得矩阵 B 后,要进一步提高问题的稳定性,您可以计算预条件子并将其与 gmres 结合使用。B 的数值属性优于原始矩阵 A 的数值属性,因此您应该使用经平衡处理的矩阵来计算预条件子。

ilu 算出两个不同的预条件子,并将它们用作 gmres 的输入来再次求解问题。在绘制了平衡后范数的同一张图上,叠加绘制应用预条件子后每次迭代的残差范数,以进行比较。该图显示,基于经平衡处理的矩阵算出的预条件子极大地提高了问题的稳定性,gmres 在不到 30 次迭代中达到了所需的容差。

semilogy(rvy)
hold on

[L1,U1] = ilu(B,struct('type','ilutp','droptol',1e-1,'thresh',0));
[yp1,flyp1,~,~,rvyp1] = gmres(B,d,[],tol,maxit,L1,U1);
semilogy(rvyp1)

[L2,U2] = ilu(B,struct('type','ilutp','droptol',1e-2,'thresh',0));
[yp2,flyp2,~,~,rvyp2] = gmres(B,d,[],tol,maxit,L2,U2);
semilogy(rvyp2)

legend('No preconditioner', 'ILUTP(1e-1)', 'ILUTP(1e-2)')
title('Relative Residual Norms with ILU Preconditioner (Equilibrated)')
hold off

输入参数

全部折叠

输入矩阵,指定为方阵。A 可以是稠密矩阵,也可以是稀疏矩阵,但在结构上必须为非奇异矩阵,这由 sprank 决定。

A 为稀疏矩阵时,输出 PRC 也为稀疏矩阵。

数据类型: single | double
复数支持:

输出参数

全部折叠

置换矩阵,以矩阵形式返回。P*AA 的置换,它使对角线元素乘积的绝对值最大化。

行缩放,以对角矩阵形式返回。RC 中对角线上的元素是正实数。

列缩放,以对角矩阵形式返回。RC 中对角线上的元素是正实数。

另请参阅

| | | | |

在 R2019a 中推出