不稳定希尔伯特矩阵单位向量答案

【问题标题】：Instability Hilbert Matrix unit vector不稳定希尔伯特矩阵单位向量
【发布时间】：2016-02-09 13:03:30
【问题描述】：

令H(n) 为n 阶希尔伯特矩阵。

让e = (0,...,0,1) - 单位向量。

让e_im := H(n) * e.

我用一些计算机代数系统解决了H(n) * x_1 = e_im。

设r = (r_1,...,r_n) 为小的随机向量。

让(e + r)_im := H(n) * (e + r).

我用一些计算机代数系统解决了H(n) * x_2 = (e + r)_im。

为什么剩余的|| H(n) * x_1 - e_im ||这么小，而|| H(n) * x_2 - (e + r)_im||却这么大？

我使用 numpy 和 scipy.linalg，我的代码：

H = scipy.linalg.hilbert(500)

e = numpy.zeros((500, 1))

e[499] = 1

e_im = H.dot(e)

x_1 = scipy.linalg.solve(H, e_im)

r = 0.0001 * numpy.random.rand(500, 1)

e_plusr_im = e + r

x_2 = scipy.linalg.solve(H, e_plusr_im)

Residials = [scipy.linalg.norm(H.dot(x_1) - b_1, 2), scipy.linalg.norm(H.dot(x_2) - b_2, 2)]

【问题讨论】：

您应该显示一些代码以提供帮助。这可能取决于所使用的算法类型。
@Warren Weckesser，是的

标签： numpy matrix scipy linear-algebra numerical-methods

【解决方案1】：

这是众所周知的数值线性代数现象，（参见大多数线性代数课程资料）。

矩阵的condition numberkappa(A) = ||A|| ||A^-1|| 说明在一般情况下，矩阵求逆放大了多少误差。这里：

>>> import scipy.linalg
>>> import numpy as np
>>> H = scipy.linalg.hilbert(500)
>>> np.linalg.cond(H)
4.6335026663215786e+20

根据经验，（确定性）浮点误差对于 64 位浮点数的相对幅度为 1e-16。条件数是如此之大，以至于产生的误差的相对幅度 > 1，所以你可能会失去所有的精度，除非你很幸运并且解决方案和中间计算具有精确的浮点表示。如果你手动添加噪音，它也会被放大。

【讨论】：