numpy中带有分数的线性系统解决方案答案

【问题标题】：Linear system solution with fractions in numpynumpy中带有分数的线性系统解决方案
【发布时间】：2016-01-30 22:33:45
【问题描述】：

我有矩阵 A 和右侧向量 y 以 fractions.Fraction 对象表示：

import random, fractions, numpy as np

A = np.zeros((3, 3), dtype=fractions.Fraction)
y = np.zeros((3, 1), dtype=fractions.Fraction)
for i in range(3):
    for j in range(3):
        A[i, j] = fractions.Fraction(np.random.randint(0, 4), np.random.randint(1, 6))
    y[i] = fractions.Fraction(np.random.randint(0, 4), np.random.randint(1, 6))

我想使用numpy 中提供的函数解决系统A*x = y 并获得以分数对象表示的结果，但不幸的是基本x = np.linalg.solve(A, y) 以标准浮点值返回结果：

>>> np.linalg.solve(A, y)

array([[-1.5245283 ],
       [ 2.36603774],
       [ 0.56352201]])

有没有办法用分数对象得到准确的结果？

编辑

使用 numpy 的内置功能（从 1.10 版开始 - 请参阅 Mad Physicist 的回答），我想做的只是不可行。一个人可以做的是实现他/她自己的基于高斯消元的线性求解器，它依赖于和、减法、乘法和除法，所有这些都是定义明确的，并且可以用分数对象精确执行（只要分子和分母适合在数据类型中，我认为是任意长的）。

如果您真的对此感兴趣，只需自己实现一个求解器，它会简单快速地完成（按照在线众多教程之一进行操作）。我不是很感兴趣，所以我会坚持浮点结果。

【问题讨论】：

您的矩阵总是 3x3 还是只是偶然的？如果它总是 3x3，你可以在最坏的情况下手动编码逆矩阵。
当我尝试运行np.linalg.solve(A, y) 行时，我实际上得到了一个错误。你是怎么让它工作的？ Numpy 给出以下错误：TypeError: No loop matching the specified signature and casting was found for ufunc solve。我在 scipy 中尝试了类似的代码，它给出了ValueError: object arrays are not supported。
矩阵很小，但它的大小是可变的。说不大于 10x10。我也知道皮草肯定不是单一的。
我正在使用 python 3.4.3，numpy 1.9.2，上面的代码运行良好。
我尝试了 Python 3.5.0、Python 2.7.10，同时使用了 numpy 1.10.1 和 scipy 0.16.1。每次都是同样的错误。也许升级让演员们感到厌烦？

标签： python numpy linear-algebra fractions

【解决方案1】：

恕我直言，没有希望。在许多情况下都有效的解决方案：

y = np.zeros(3, dtype=fractions.Fraction)
....
X= np.linalg.solve(A,y)
s=[fractions.Fraction.from_float(x).limit_denominator(6**9) for x in X]
print(s,y==dot(A,s))

它利用解几乎是分子和分母很少的分数的性质，找到它。

【讨论】：

嗯，实际上从系统的大小和实际分母可以分析找到结果的最大分母并使用它而不是 6**9。但我认为最好的方法仍然是在需要时自己实现 LU 分解和后向替换。但是，是的，在阅读了 Mad Physicist 发布的参考资料后，我同意使用 numpy 内置求解器做所有事情是没有希望的。
它并不优雅，但对我来说效果很好，谢谢

【解决方案2】：

根据 python 邮件列表中的this thread，似乎不可能使用纯 numpy 反转有理矩阵。 The response 建议您可以将 sympy 用于大小为 4x4 的有理矩阵。如果由于某种原因您与 numpy 相关，您可以考虑“手动”采用和使用 3x3 矩阵的逆矩阵。可以在http://www.mathsisfun.com/algebra/matrix-inverse-minors-cofactors-adjugate.html 上找到有关如何执行此操作的分步教程，以及有关矩阵求逆的大量其他教程。

【讨论】：