完成枢轴和阿达玛矩阵

几年来,我们认为阿达玛矩阵在完全旋转的高斯消去中显示了最大的元素增长。我们错了。

内容

完成旋转增长因子

我想继续从我的上一篇博客高斯消去过程中元素的生长。假设我们正在解一个线性方程组,它包含一个阶为$n的矩阵$ a。设$a_{i,j}^{(k)}$表示消元步骤$k$之后的矩阵元素。回想一下生长因子对于我们使用的旋转过程,是数量$$ \rho_n = {\max_{i,j,k} |a_{i,j}^{(k)}| \over \max_{i,j} |a_{i,j}|} $$完全旋转的目的是通过在约简的每一步使用行和列交换来控制生长因子,将整个未约简矩阵中最大的元素移到枢轴位置。在1962年对高斯消去的舍入误差分析中,J. H. Wilkinson证明了完全旋转的生长因子满足$$ \rho_n \le g(n) $$生长函数= $$ g(n) = (n \ 2 \ 3^{1/2} 4^{1/3} \cdots n^{1/ 1/(n-1)})^{1/2} \约1.8 \ n^{1/2} n^{1/4 \log{n}} $$威尔金森的分析清楚地表明,实际的增长永远不可能接近这么大。我们看到部分轴转的生长因子是$2^{n-1}$。对于完全旋转,它要小得多。为了对现在的情况有个大致的了解,可以将$g(n)$表示为MATLAB的一行程序
G = @(n)√(n*prod((2:n).^(1:(n-1))))))
g = function_handle值:@ (n) sqrt (n * prod ((2: n)。^ (1. / (1:(n - 1)))))
查看$n = 100$
g (100)
Ans = 3.5703e+03
这里$g(n) \约35n$。

阿达玛矩阵

雅克·阿达马是一位法国数学家,生于1865年,卒于1963年。他在许多领域都有贡献,从数论到偏微分方程。以他的名字命名的矩阵有+1和-1的项,并且行和列相互正交。由阿达玛矩阵的行张成的$n$维平行四边形在所有由以1为界的矩阵张成的平行四边形中具有最大的可能体积。我们对MATLAB世界中的阿达玛矩阵很感兴趣,因为它们是阿达玛变换的基础,而阿达玛变换与傅里叶变换密切相关。它们也适用于纠错码和统计上的方差估计。MATLAB已经有了阿达玛函数。这是8 × 8的阿达玛矩阵。
H = hadamard(8)
H = 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
检查这些行是否互相垂直。
H的* H
Ans = 8 0 0 0 0 0 0 0 0 0 8 0 0 0 0 0 0 0 0 0 0 8 0 0 0 0 0 0 0 0 0 0 0 8 0 0 0 0 0 0 0 0 0 0 8 0 0 0 0 0 0 0 0 0 0 0 8 0 0 0 0 0 0 0 0 0 0 0 0 8
在高斯消去过程中,行的正交性导致元素增长。如果我们要求LU因式分解H在美国,实际上并没有发生转向,但如果完全转向,就会产生对现任总统有利的结果。
[L,U] = lu(H)
L = 1 0 0 0 0 0 0 0 1 1 0 0 0 0 0 0 1 0 1 0 0 0 0 0 1 1 1 1 0 0 0 0 0 0 0 1 0 0 0 1 1 0 0 1 1 0 0 1 0 1 0 1 0 1 0 1 1 1 1 1 1 1 1 U = 1 1 1 1 1 1 1 1 0 2 0 2 0 2 0 2 0 0 2 2 0 0 2 2 0 0 0 4 0 0 0 4 0 0 0 0 2 2 2 2 0 0 0 0 0 4 0 4 0 0 0 0 0 0 4 4 0 0 0 0 0 0 0 8
那么,这是一个矩阵$n = 8$,其中$$ \rho_n = n $$

(n)是否等于n ?

近30年来,所有人都认为阿达玛矩阵是极端情况,实矩阵完全轴转的生长因子$\rho_n$等于$n$。在20世纪60年代和70年代的不同时期,我个人做了一些数值实验,寻找具有大枢轴增长的矩阵。我当时拥有的计算机将我的矩阵的顺序限制在几十个。我从未发现任何$\rho_n > n$的案例。加州雪佛龙研究所的伦纳德·托恩海姆(Leonard Tornheim)就是研究这个问题的人之一。他还做了很多数值实验。他发现了一个复杂的3 × 3矩阵$\rho_3 > 3$。他还证明了对于实矩阵$\rho_3 = 2 \ 1/4$。1968年,Colin Cryer提出了一个官方猜想:对于实矩阵,$\rho_n$等于$n$。简·戴和布莱恩·彼得森在1988年的一篇论文中很好地调查了当时关于完全旋转的元素生长的情况。

尼克·古尔德的惊喜

然后,在1991年,来自牛津的尼克·古尔德(Nick Gould)宣布发现了真正的矩阵,其完整的枢轴$\rho_n$大于$n$。对于一个$n$ × - $n$的矩阵,其中$n$在13到16之间,Nick建立了一个大的、稀疏的非线性优化问题,大约涉及$n^3/3$变量,并用他和他的同事Andy Conn和Philippe Toint开发的LANCELOT包来解决它。他的大部分论文都是关于一个13 × 13的矩阵,其中$$ \rho_{13} = 13.0205 $,这刚刚超过$\rho_n = n$猜想。但他也发现了一些稍大一点的矩阵,可以做得更好。$$ 16 × 16矩阵特别引人注目,因为它不是一个阿达玛矩阵。所以,据我所知,这就是今天的情况。没人对大型实验感兴趣。我不知道$\rho_{100}$或$\rho_{1000}$可能是什么。当然,我们真的不需要知道,因为我们在实践中没有使用完全旋转。

P ^ 8

这让我想起了一个我喜欢讲的故事。我不确定这是不是真的。彼德·布辛格是斯坦福大学数学与计算机科学专业的一名研究生当时我正处于60年代初。他与Gene Golub合作编写了第一个发表的使用Householder反射的QR矩阵分解的Algol程序。彼得离开斯坦福大学,加入了贝尔实验室乔·特劳布的团队,该团队正在开发数学软件的波特库。他们从其他人那里引进程序,在贝尔实验室的机器上运行程序,彻底测试软件,决定哪些代码是最好的,并生成一个统一的库。彼得负责矩阵计算。他检查了所有进口的线性方程求解器,包括我的。他认为没有一个是令人满意的,因为这个旋转增长的问题。所以他自己写了一个程序来进行部分转轴,并监测生长。 If the growth was too large, the program switched to complete pivoting. At the time, Bell Labs was trying to patent software, so Peter named his new subroutine "P^8". This stood for "Peter's Pseudo Perfect Partial Pivot Picker, Patent Pending". In my opinion, that is one of the great all time subroutine names. The experience inspired Peter to go to law school at night and change jobs. He switched to the Bell Labs patent department. I've lost contact with Peter, so I can't ask him if he really did name the routine P^8. If he didn't, he should have.

车绕轴旋转

我在上一篇文章中提到过Les Foster,他发现了部分旋转的指数增长的例子,他提倡他所谓的“车旋转”。沿着列向下搜索最大的元素,就像在部分旋转中一样,但也沿着行搜索最后一个元素。这介于部分旋转和完全旋转之间。Les能够证明一些关于支点生长的结果。但该算法不能很好地推广到块形式。

92号骑士团的阿达玛德。

我见证了阿达玛矩阵历史上的一个里程碑。正交性条件很容易暗示,一个阿达玛矩阵的阶为$n$必须是1、2或4的倍数。但问题仍然存在,对于每个$k$,是否存在阶为$n = 4k$的阿达玛矩阵?这是数论中的一个开放性问题。创建一些大小的阿达玛矩阵是相当容易的。一种生成2次幂的阿达玛矩阵的好方法是递归。
H = 1;N = 2,4,8,...H = [H H H -h];结束
如果一个而且B是阿达玛矩阵,那么也是
克隆亚麻(A, B)
到1961年,有了这些和其他类似的结构,人们知道了如何为$n \le 100$的所有阶(除了$n = 92$)的4的倍构造阿达玛矩阵。1961年,我在JPL(加州理工学院的喷气推进实验室)做暑期工。我的办公室在一个临时拖车里,我的拖车伙伴是伦恩·鲍默特。莱恩自豪地给我看了一张他刚做的彩色图片,他打算把它作为杂志的封面科学美国人.它是一个92 × 92的矩阵,由23 × 23的块组成,分别代表+1和-1的明暗交替单元格。封面上没有这张图片科学美国人但是我的那本保存了很长时间。莱恩在喷气推进实验室的机器上做了一个计算,填充了最后一个小于100的值n美元。他与同事Sol Golomb(南加州大学教授)和Marshall Hall, Jr.(加州理工学院教授)的论文发表在著名的《AMS公报》上。原来我上过一门关于差分集的课,这门课涉及到生成矩阵的数学,是加州理工学院霍尔教授的。这是一个MATLAB图Baumert92.你可以从这个链接
H = Baumert92;pcolor92 (H);
让我们检查一下我们是否得到了期望的支点增长。
[L,U] = lu(H);unn = U(end,end)
Unn = 92.0000

参考文献

关于Trefethen和Bau的书的参考应该包含在我之前关于部分枢轴的文章中。第22讲在实践中解释了部分枢轴的稳定性。Lloyd N. Trefethen和David Bau,数值线性代数< http://bookstore.siam.org/ot50《中国医学杂志》,1997,362pp。Peter Businger,“高斯消去的数值稳定性监测”,< http://link.springer.com/article/10.1007/BF02165006>,数字数学16,4,1971,360-361。简·戴和布莱恩·彼得森,“高斯消去法的增长”,< http://www.jstor.org/stable/2322755>,美国数学月刊,95,6,1988,489-513。尼克·古尔德,"完全旋转高斯消去的增长"< http://www.numerical.rl.ac.uk/people/nimg/pubs/Goul91_simax.pdf>矩阵分析与应用学报12,1991,351-364。莱斯利·福斯特,“基于车枢轴的高斯消去的生长因子和效率”,< http://www.math.sjsu.edu/福斯特/ gerpwithcor.pdf[j] .数学与工程学报,2017,31(3):397 - 397。Leslie Foster,“LURP:基于车旋转的高斯消去”,matlabcentral / fileexchange / 37721 -车-旋转, MATLAB中央文件交换,2012。莱昂纳德·鲍默特、s·w·戈罗姆和小马歇尔·霍尔,《92阶阿达玛矩阵的发现》,< http://dx.doi.org/10.1090%2fs0002 - 9904 - 1962 - 10761 - 7>,美国数学学会公报68(1962),237-238。

已发布与MATLAB®R2018a

|

评论

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