1

这是我为作业中遇到的问题编写的 matlab 代码。在 A 和它的转置相乘之后,根据所有同学的代码(不同的代码)给他们的代码(不同的代码),得到的方阵应该有行列式零。为什么我的代码没有给出 c 和 d 的行列式为无穷大

A = rand(500,1500);

b = rand(500,1);

c = (A.')*A;

detc = det(c);

cinv = inv((A.')*A);

d = A*(A.');

detd = det(d);

dinv = inv(A*(A.'));

x1 = (inv((A.')*A))*((A.')*b);

x2 = A.'*((inv(A*(A.')))*b);
4

3 回答 3

3

此行为在的文档的限制部分中进行了解释,并在查找奇异矩阵的行列式小节中进行了举例说明:det

的行列式A很大,尽管它A是单数的。事实上, 的行列式A应该正好为零!的不准确d是由于 LU 分解的 MATLAB® 实现中舍入误差的聚合,该 LU 分解det用于计算行列式。

也就是说,在这种情况下,您可以通过使用同一页面上给出的m 代码实现来产生您想要的结果,但对角元素U按升序排序。考虑示例脚本:

clc();
clear();

A = rand(500,1500);
b = rand(500,1);
c = (A.')*A;

[L,U] = lu(c);
% Since det(L) is always (+/-)1, it doesn't impact anything
diagU = diag(U);
detU1 = prod(diagU);
detU2 = prod(sort(diagU,'descend'));
detU3 = prod(sort(diagU,'ascend'));

fprintf('Minimum: %+9.5e\n',min(abs(diagU)));
fprintf('Maximum: %+9.5e\n',max(abs(diagU)));
fprintf('Determinant:\n');
fprintf('\tNo Sort:          %g\n'  ,detU1);
fprintf('\tDescending Sort:  %g\n'  ,detU2);
fprintf('\tAscending Sort:   %g\n\n',detU3);

这将产生输出:

Minimum: +1.53111e-13
Maximum: +1.72592e+02
Determinant:
    No Sort:          Inf
    Descending Sort:  Inf
    Ascending Sort:   0

请注意,排序的方向很重要,并且由于对角线上不存在Inf真值,因此无排序给出。0降序排序首先看到最大值相乘,显然,它们超过realmax并且永远不会乘以 true 0,这将生成 a NaN。升序排序将所有接近零的对角线值与极少数大的负值聚集在一起(实际上,更稳健的方法会根据幅度进行排序,但这里没有这样做),它们的乘法生成一个真值0(意味着该值低于 IEEE-754 算术中可用的最小非规格化数)产生“正确”结果。

所有这些,正如其他人所暗示的那样,我将引用原始 Matlab 开发人员和 Mathworks 联合创始人 Cleve Moler 的话:

[行列式] 在理论考虑和手工计算中很有用,但不能为稳健的数值软件提供良好的基础。

于 2016-03-14T01:23:40.747 回答
2

好的。因此,det(A'*A) 不为零这一事实并不能很好地表明 A'*A 的(非)奇异性。行列式取决于缩放,并且显然非奇异的矩阵可以具有非常小的行列式。例如,矩阵 1/2 * I_n (其中 I_n 是 nxn 恒等式)具有 (1/2)^n 的行列式,当 n 趋于无穷时,它会(快速)收敛到 0。但是 1/2 * I_n 根本不是单数。

出于这个原因,检查矩阵奇点的最佳方法是条件数。

在你的情况下,在做了一些测试之后

>> A = rand(500, 1500) ;

>> det(A'*A)

ans =

   Inf

您可以看到(计算的)行列式显然非零。但这实际上并不奇怪,它不应该真的打扰你。行列式很难计算,所以是的,它只是舍入误差。如果您想要更好的近似值,可以执行以下操作

>> s = eig(A'*A) ;
>> prod(s)

ans =

     0

在那里,你看到它更接近于零。

另一方面,条件数是矩阵(非)奇异性的更好估计。这里是

>> cond(A'*A)

ans =

   1.4853e+20

而且,由于它远大于 1e+16,因此矩阵显然是奇异的。1e+16 的原因有点繁琐,但主要是由于进行浮点计算时的计算机精度。

于 2016-03-14T00:04:44.173 回答
0

我认为这几乎只是一个四舍五入的问题, Inf 并不意味着您将 Infinity 作为答案,只是您的行列式确实很大并且超出了realmax。正如阿迪尔所说,A*A。生成一个对称矩阵,并且它的行列式应该有一个数值。例如,设置:

A=rand(5,15)

你应该找到 A*A 的 det。只是一个数值。

那么你的朋友是如何得到零的,那么对于大型矩阵来说很容易得到 0 或 inf det(为什么你首先要这样做,我不知道)。所以我认为他们只是得到相同/相似的舍入问题。

于 2016-03-13T14:48:33.727 回答