CR和CAB,等级揭示矩阵分解

线性变换的秩是线性代数中的基本概念,矩阵分解是数值线性代数中的基本概念。吉尔·斯特朗2020年线性代数展望试图在线性代数入门课程的早期引入这些概念。

CR矩阵分解提供了一个rref,行简化阶梯形,作为揭示矩阵分解的秩。我讨论了CR在一个一对属于帖子十月。现在我想描述一下驾驶室分解,它使用rref为了以相同的方式处理行和列,两次。

吉尔和我正在写一篇关于这些想法的评论文章,LU和CR消除.这是当前草案的链接,LUCR_43.pdf

内容

驾驶室

这是一个快速的描述驾驶室

[C,W,B]=cab(A)返回C=A(:,cols)、B=A(rows,:)和W=A(rows,cols),以便
一个= C *发票(W) * B
[C,W,B,cols,rows] = cab(A)也返回索引cols和rows。

所有元素CWB从输入矩阵中获取一个它本身C是它的一些列,B它的一些行,和W集合的交点是多少CB. 事实上科尔斯所有这些因数分解都是计算的。列指标来自行简化阶梯形,rref (A),和转置的行指标,rref(的)

排名

一个矩阵的线性无关列数等于线性无关行数这一事实并不平凡,需要证明。这个数字是排名和是非奇异矩阵的大小W

W是r-by-r,其中秩(A)=r=长度(cols)=长度(行)。

如果你试着去找CWR用一个W那就太小了C*inv(W)*R将无法复制一个. 如果W如果太大,它将是单数的。

第一名

如果一个排名第一,,一个= u * v '对于列向量u和行向量v'. 然后CuBv'W是中的第一个非零元素uv. 下面是一个示例,其中第一行故意为零。

A=(0:3)“*(4:6)
A=0 0 4 5 6 8 10 12 15 18
总体安排契约[C,W,B]=驾驶室(A)
C=04812W=4B=456

第二名

在本例中,即使是不经意的检查也会显示等级。注意一个它是长方形的。

A=重塑(1:24,4,6)'
A=12345678910112151617171819212324
[C W B,关口,行]=出租车(A)
C = 1 2 5 6 9 10 13 14 17 18 21 22 W = 1 2 5 6 6 B = 1 2 3 4 5 6 7 8 cols = 1 2 rows = 1 2
排名=长度(峡路)
排名= 2

正式

如果A是正方形且非奇异的,则C=B=W=A。

例如

A=帕斯卡(3)
A=1 1 2 3 1 3 6
[C,W,B]=驾驶室(A)
C=1112336w=1112336b=12112336

富兰克林魔术广场

富兰克林半魔术大约在1752年,该广场归属于本杰明·富兰克林。富兰克林给他的一位同事的一封信包含在下面引用的富兰克林论文中。

富兰克林平方的所有行和列都有所需的魔和,但对角线没有,因此矩阵不是完全魔和。然而,许多其他有趣的子矩阵是神奇的,包括弯曲的对角线和围绕中心对称排列的任何八个元素。

这个8 × 8矩阵的秩是多少?

一个=富兰克林(8)
A=52 61 4 13 20 29 36 45 14 3 62 51 46 35 30 19 53 60 5 12 28 37 44 11 6 59 54 43 27 22 55 58 7 10 23 26 39 42 9 8 57 56 41 40 25 24 50 63 2 15 18 31 34 16 1 64 48 32 17

很难认为富兰克林的平方是线性变换,但是驾驶室仍然可以计算它的等级。

[C,W,B]=驾驶室(A)
C=526141436253605511156555879857063216164 W=526141436253605b=52614132093644362514635305360121283744

因此,排名只有三位。

计算C*inv(W)*B引入了一些舍入误差。

测试=C*inv(W)*B;显示器(“测试”
测试=52.00 61.00 4.00 13.00 20.00 29.00 36.00 45.00 14.00 3.00 62.00 51.00 46.00 35.00 30.00 19.00 53.00 60.00 5.00 12.00 21.00 28.00 37.00 44.00 11.00 6.00 59.00 54.00 43.00 38.00 27.00 22.00 55.00 58.00 7.00 10.00 23.00 26.00 39.00 42.00 9.00 8.00 57.00 56.00 41.00 40.00 25.00 24.00 50.00 20.00 20.00 27.00 31.00 1.0064.00 49.00 48.00 33.00 32.00 17.00

舍入测试精确到最接近的整数可以复制富兰克林的幻方。

反斜杠和正斜杠

通过(C/W)*B或C*(W\B)计算C*inv(W)*B通常更精确。如果A有整数元素,您可以尝试四舍五入(C/W*B)或四舍五入(C*(W\B))。

这个5 × 5矩阵是故意构造的,这样它的秩不足就不明显了。

A=画廊(5)
A=-911-2163-25270-69141-4211684-575575-11493451-138013891-38917782-2334593361024-10242048-614424572
[C,W,B]=驾驶室(A)
C = 9 11 -21 63 70 -69 141 -421 -575 575 -1149 3451 3891 -3891 7782 -23345 1024 -1024 2048 -6144 W = 9 11 -21 63 70 -69 141 -421 -575 575 -1149 3451 3891 -3891 7782 -23345 B = 9 11 -21 63 -252 70 -69 141 -421 1684 -575 575 -1149 3451 -13801 3891 -3891 7782 -23345 93365
秩=大小(W,1)
排名= 4

在这个例子中,W他有轻微的病态。

卡帕=秒(W)
kappa=2.4301e+04

有三种方法可以重新创建一个

invw = C *发票(W) * B;返回= (C/W) * B;forex = C * (W\B);

反演的相对误差受反演条件的影响W,而反斜杠和正斜杠的精确度几乎要高出五个数量级,

总体安排短的erelerr=[norm(A-invw)norm(A-fore)norm(A-back)]/norm(A)
relerr=2.1266e-13 4.7760e-18 3.7217e-17

将这三个数中的任何一个舍入为整数都会产生一个确切地

=轮(C *发票(W) * B)
A=-911-2163-25270-69141-4211684-575575-11493451-138013891-38917782-2334593361024-10242048-614424572

Gauss-Jordan

计算的算法rref是高斯-约当消去法,它在数值分析界的名声不太好。在听我描述了这个算法的恶名之后,吉尔在序言中写道我们的调查高斯·乔丹是“昂贵的,而且在数量上是危险的”。

“昂贵”的部分很容易解释。温度场的计算rref在轴元素上方和下方的列中引入零,而高斯消元仅在轴下方引入零。那是更少的工作。我们可以用一条融合的乘加指令FMA来衡量计算复杂度,FMA将两个浮点值相乘,然后在一次操作中将第三个值添加到结果中。对于一个大的n×n系统,Gauss-Jordan需要大约(1/2)n^3 FMA,而Gaussian消元需要更少,只有(1/3)n^3 FMA。

然而,我们并不打算这样做CR驾驶室用于大型矩阵。就我们的目的而言,即使使用高斯-乔丹两次也不昂贵。

高斯·乔丹名声中“数字危险”的部分更为微妙。它源于求解一个线性系统A*x=b通过计算rref关于增广矩阵,[甲、乙].另一方面,高斯消去法可以找到三角形因子lU然后计算x通过反向替换。我们可以说高斯消去是部分消去,而高斯·乔丹是完成消除。Jim Wilkinson和Gwen Peters在1975年的论文中对这两种方法进行了比较。高斯消去法几乎总是计算附近系统的解,但对于增广矩阵上的高斯约当却不能这样说。

然而,我们并没有使用两种形式的消去法来实际计算驾驶室因素;rref只选择索引。矩阵CWB从…起驾驶室(A)余子式的一个.他们肯定是准确的.有两个潜在的数值困难--W可能是条件不好,或者可能是尺寸不对。

条件

如果…怎么办一个他身体状况很差吗?然后W必须继承这种条件。这里有一个例子。从

n = 12;U = eye(n,n) - triu(ones(n,n),1)
U = 1 1 1 1 1 1 1 1 1 1 1 1 0 1 1 1 1 1 1 1 1 1 1 1 0 0 1 1 1 1 1 1 1 1 1 1 0 0 0 1 1 1 1 1 1 1 1 1 0 0 0 0 1 1 1 1 1 1 1 1 0 0 0 0 0 1 1 1 1 1 1 1 0 0 0 0 0 0 1 1 1 1 1 1 0 0 0 0 0 0 0 1 1 1 1 1 0 0 0 0 0 0 0 0 1 1 1 1 0 0 0 0 0 0 0 0 0 1 1 1 0 0 0 0 0 0 0 0 0 0 1 1 0 0 0 0 0 0 0 00 0 0 1

U他身体状况很差。电导率(U, 1)条件(U,inf)都等于n * 2 ^ (n - 1)

允许

一个= U ' * U
A=1-1-1-1-1-1-1-1-1-1-1-1-1-2-1-1-1-1-1-1-1-1-1-2-2-2-2-2-3-3-3-1-1-2-3-4-1-4-1-2-4-4-5-1-1-2-6-1-1-2-3-6-7-1-1-1-2-3-6-7

一个是对称的,正定的,满秩的。决定因素一个等于一。我看不出一个简洁的公式条件(A),但我知道它确实像4^n

自从一个是正式的,它的一个= C *发票(W) * B因子分解必须是A=A*inv(A)*A

[C W B] =出租车(A);test = [isequal(C,A) isequal(W,A) isequal(B,A)]
测试=1×3逻辑阵列1

W他身体状况很差。投资部(西)有很多元素。

Winv =发票(W)
Winv=第1列到第6列1398102 699051 349526 174764 87384 43696 699051 349526 174763 87382 43692 21848 349526 174763 87382 43691 21846 10924 174764 87382 43691 21846 10923 5462 87384 43692 21846 10923 5462 2731 43696 21848 10924 5462 2731 1366 21856 10928 5464 2732 1366 683 1094 542 43692 2752 13688 344 172 2816 1408 704352 176 88 1536 768 384 192 96 48 1024 512 256 128 64 32列7至12 21856 10944 5504 2816 1536 1024 10928 5472 2752 1408 768 512 5464 2736 1376 704 384 256 2732 1368 688 352 192 128 1366 684 344 176 64 683 342 172 48 342 171 86 24 1617 86 12 86 43 22 11 11 11 44 22 36 12 12 6 2 1 16 8 2 1 1
logWinv=楼层(log2(Winv))
logWinv=2019181717171515141811019181818181818181818171514181818181715141818181717151818181818171717151514141418110101019181717171810191010101010101019181818181818181818181818181818171717181717171717171717171718181818181818171717171717181818101010101019976 5 4 3 2 1 0 0

沿着…的对角线logWinv从…起(n,n)(1,1)

Winv(1,1)>4^(n-2)
ans=逻辑1

一个接近低阶矩阵?答案是肯定的,但低阶近似值来自SVD。

CAB与SVD

SVD是最终秩显示因式分解,但它实际上并不试图计算秩。

[U,S,V]=svd(A)

计算一整套奇异值和向量,就像矩阵具有满秩一样。奇异值,s=诊断(s),按降序索引,

s(1)>=…>=s(r)>=s(r+1)。

您可以将奇异值与公差进行比较,以获得近似的秩和因子分解。

选择一个有效的等级是由Eckart-Young-Mirsky定理:如果基于“增大化现实”技术近似是通过选择秩得到的吗r

基于“增大化现实”技术= U (:, 1: r) * S (1: r, 1: r) * V (:, 1: r) '

误差的范数以省略的第一项为界,

范数(A - Ar,2) <= s(r+1)

如果奇异值迅速减小,则只需几个项即可实现精确的低秩近似。

另一方面驾驶室从矩阵本身的列和行中选择基。这个问题没有埃克特-杨定理CR驾驶室分解。

我们打算CR驾驶室用于课堂上的小矩阵,通常带有整数项。我们将让奇异值分解和随机矩阵分解处理大问题。

工具书类

本杰明·富兰克林,本杰明·富兰克林的论文,第44卷,1750-53页,写给彼得·柯林森。(美国哲学学会和耶鲁大学)p。392https://franklinpapers.org/framedVolumes.jsp?vol=4&page=392a

G.Peters和J.H.Wilkinson,“关于带旋转的Gauss-Jordan消元法的稳定性”,ACM通讯18, 1975,第20-24页。https://dl.acm.org/doi/10.1145/360569.360653,亦于http://www.stat.uchicago.edu/~lekheng/courses/302/gauss-jordan2.pdf

后记

本·富兰克林也可以为今天的帖子提供图片。

A=富兰克林(8);[~,~,B]=cab(A);地块(B)




发布与MATLAB®R2021a

|

评论

如需留言,请点击在这里登录到您的MathWorks帐户或创建新帐户。