计算36x36矩阵的Jordan形式时出现了错误的解释

Question

我一直在尝试计算一个由仅有三个不同的元素1、1/2和0组成的36x36矩阵的Jordan标准形。该矩阵是一个概率转移矩阵，因此，根据这些条目，该矩阵显然是稀疏的。

我遇到的问题是：每当我尝试计算时，我的计算机似乎会崩溃。

[V, J] = jordan(A),

或者

[V, J] = jordan(sym(A)),

我收到以下错误信息：

使用mupadmex时出错 MuPAD命令中的错误：相似矩阵太大。 在sym / mupadmexnout（第1546行）中出错         out = mupadmex（fcn，args {:}）; 在sym / jordan（第32行）中出错         [Vsym，Jsym] = mupadmexnout（'symobj :: jordan'，A，'All'）;

我在MATLAB帮助文件中读到，计算Jordan形式对扰动非常敏感。然而，因为矩阵的所有条目都是整数或整数比率，所以我认为我的计算不会有问题。

我的问题是：

1.如何解释我收到的错误输出？ 2.我收到的错误是否可以解决？ 3.如果错误无法解决，是否有替代方法（Matlab中的函数）可以尝试计算Jordan形式？

- mrY

Answer 1

1) 我如何解释我收到的错误输出？

问题在于Matlab使用符号计算来评估Jordan形式。这就是它要求您提供有理数的原因。当我们考虑数值编程时，一个36x36的矩阵非常小，但（我不确定）对于符号编程来说，这个大小可能很大。

2) 为什么Matlab没有用于数值评估Jordan形式的工具箱？

问题在于这种评估是数值不稳定的。请参见维基百科中的示例。基本上，任何具有多个特征值（共享相同块）的矩阵的扰动都可以导致这些特征值在所需的Jordan形式的分离块中变得不同。

3) 如果错误无法解决，是否有替代方法（Matlab中的函数）可尝试计算Jordan形式？

我认为Matlab没有数值函数来解决这个任务。

我不知道您正在寻找什么类型的应用程序...话虽如此，一个（非常普遍的）选项是评估Schur分解（两种转换都将矩阵转换为上三角分解），这是数值稳定的。它使用一个酉相似变换。Matlab的schur函数实现了这一点。

另请参见Math.StackExchange问题：Jordan和Schur分解有什么区别？

- DanielTheRocketMan

---

好答案。+1。希望您不介意我做了一些清理。 - horchler

我认为问题的提出者已经更新了问题，所以你的第二点可能已经过时了。 - horchler

Jordan形式在数值上非常不稳定，以至于它已被提议作为无法作弊的数字基准（PDF）的基础。 - horchler

@horchler 有趣，我不知道那个！ - DanielTheRocketMan

Answer 2

我主要回答你问题的第三部分：计算Jordan形式的替代方法。

根据您想要评估的最大矩阵以及可能使用的Matlab版本，是的，您可以符号地计算Jordan形式及其相似变换。

我之前有一个相关问题，所以您的问题让我再次研究了一下。在我的问题中，我不能使用jordan函数来计算大于82x82的测试矩阵的广义特征向量，这在Matlab R2013a（我现在使用的R2014b相同）中是不可能的。

更新：R2015a+似乎已经修复了82x82以上的矩阵的问题，至少对于下面300x300的测试矩阵而言是如此。以下是R2015a之前的解决方法。

这个技巧的关键在于使用不同（可能是更新的）MuPAD函数来计算Jordan形式和相关的相似变换：linalg::jordanForm，而不是symobj::jordan。因为你甚至不能计算一个36×36的Jordan形式，我想知道你是否使用了一个更早的版本。所以我不能保证这在你的Matlab版本中能够工作。这是必要的代码：

% My test matrix
A = -eye(40);
A(1,2:end-1) = -2;
A(1,end) = -1/2;
A(2,2) = 1/2;

% [V,J] = jordan(A); doesn't work for A larger than 82-by-82
P = feval(symengine,'linalg::jordanForm',sym(A),'All');
J = double(P(1));
V = double(P(2));

我使用这种方法可以计算出一个300 x 300的测试矩阵的J和V，尽管在我的MacBook Pro上需要大约一分钟。这仍然是一种符号方法，对于许多符号操作来说，超过100 x 100可能被认为是相当大的。旧版本的Symbolic Math工具箱效率甚至更低。

顺便说一下，如果您实际上不需要相似变换V，那么您可能仍然能够使用J = jordan(A);而不会出现错误。

- horchler