斯卡拉矩阵反演

Question

呃，是的，我真的需要一个没有创造者眼睛的人的快速输入。根据我的scalacheck测试，这里出现了一些问题......但我并不十分了解它，以便知道错误的位置。

case class Matrix(_1: (Float, Float, Float, Float), _2: (Float, Float, Float, Float), 
        _3: (Float, Float, Float, Float), _4: (Float, Float, Float, Float)) extends Immutable { 
    def invert = { 
    val _11 = _2._2 * _3._3 * _4._4 - _2._2 * _3._4 * _4._3 - _3._2 * _2._3 * _4._4 
     +_3._2 * _2._4 * _4._3 + _4._2 * _2._3 * _3._4 - _4._2 * _2._4 * _3._3 
    val _21 = -_2._1 * _3._3 * _4._4 + _2._1 * _3._4 * _4._3 + _3._1 * _2._3 * _4._4 
     -_3._1 * _2._4 * _4._3 - _4._1 * _2._3 * _3._4 + _4._1 * _2._4 * _3._3 
    val _31 = _2._1 * _3._2 * _4._4 - _2._1 * _3._4 * _4._2 - _3._1 * _2._2 * _4._4 
     +_3._1 * _2._4 * _4._2 + _4._1 * _2._2 * _3._4 - _4._1 * _2._4 * _3._2 
    val _41 = -_2._1 * _3._2 * _4._3 + _2._1 * _3._3 * _4._2 + _3._1 * _2._2 * _4._3 
     -_3._1 * _2._3 * _4._2 - _4._1 * _2._2 * _3._3 + _4._1 * _2._3 * _3._2 
    val _12 = -_1._2 * _3._3 * _4._4 + _1._2 * _3._4 * _4._3 + _3._2 * _1._3 * _4._4 
     -_3._2 * _1._4 * _4._3 - _4._2 * _1._3 * _3._4 + _4._2 * _1._4 * _3._3 
    val _22 = _1._1 * _3._3 * _4._4 - _1._1 * _3._4 * _4._3 - _3._1 * _1._3 * _4._4 
     +_3._1 * _1._4 * _4._3 + _4._1 * _1._3 * _3._4 - _4._1 * _1._4 * _3._3 
    val _32 = -_1._1 * _3._2 * _4._4 + _1._1 * _3._4 * _4._2 + _3._1 * _1._2 * _4._4 
     -_3._1 * _1._4 * _4._2 - _4._1 * _1._2 * _3._4 + _4._1 * _1._4 * _3._2 
    val _42 = _1._1 * _3._2 * _4._3 - _1._1 * _3._3 * _4._2 - _3._1 * _1._2 * _4._3 
     +_3._1 * _1._3 * _4._2 + _4._1 * _1._2 * _3._3 - _4._1 * _1._3 * _3._2 
    val _13 = _1._2 * _2._3 * _4._4 - _1._2 * _2._4 * _4._3 - _2._2 * _1._3 * _4._4 
     +_2._2 * _1._4 * _4._3 + _4._2 * _1._3 * _2._4 - _4._2 * _1._4 * _2._3 
    val _23 = -_1._1 * _2._3 * _4._4 + _1._1 * _2._4 * _4._3 + _2._1 * _1._3 * _4._4 
     -_2._1 * _1._4 * _4._3 - _4._1 * _1._3 * _2._4 + _4._1 * _1._4 * _2._3 
    val _33 = _1._1 * _2._2 * _4._4 - _1._1 * _2._4 * _4._2 - _2._1 * _1._2 * _4._4 
     +_2._1 * _1._4 * _4._2 + _4._1 * _1._2 * _2._4 - _4._1 * _1._4 * _2._2 
    val _43 = -_1._1 * _2._2 * _4._3 + _1._1 * _2._3 * _4._2 + _2._1 * _1._2 * _4._3 
     -_2._1 * _1._3 * _4._2 - _4._1 * _1._2 * _2._3 + _4._1 * _1._3 * _2._2 
    val _14 = -_1._2 * _2._3 * _3._4 + _1._2 * _2._4 * _3._3 + _2._2 * _1._3 * _3._4 
     -_2._2 * _1._4 * _3._3 - _3._2 * _1._3 * _2._4 + _3._2 * _1._4 * _2._3 
    val _24 = _1._1 * _2._3 * _3._4 - _1._1 * _2._4 * _3._3 - _2._1 * _1._3 * _3._4 
     +_2._1 * _1._4 * _3._3 + _3._1 * _1._3 * _2._4 - _3._1 * _1._4 * _2._3 
    val _34 = -_1._1 * _2._2 * _3._4 + _1._1 * _2._4 * _3._2 + _2._1 * _1._2 * _3._4 
     -_2._1 * _1._4 * _3._2 - _3._1 * _1._2 * _2._4 + _3._1 * _1._4 * _2._2 
    val _44 = _1._1 * _2._2 * _3._3 - _1._1 * _2._3 * _3._2 - _2._1 * _1._2 * _3._3 
     +_2._1 * _1._3 * _3._2 + _3._1 * _1._2 * _2._3 - _3._1 * _1._3 * _2._2 

    val det = _1._1 * _11 + _1._2 * _21 + _1._3 * _31 + _1._4 * _41 
    if (det == 0) this 
    else Matrix(
     (_11, _12, _13, _14), 
     (_21, _22, _23, _24), 
     (_31, _32, _33, _34), 
     (_41, _42, _43, _44) 
    ) * (1/det) 
    } 

    def *(f: Float) = Matrix(
    (_1._1 * f, _1._2 * f, _1._3 * f, _1._4 * f), 
    (_2._1 * f, _2._2 * f, _2._3 * f, _2._4 * f), 
    (_3._1 * f, _3._2 * f, _3._3 * f, _3._4 * f), 
    (_4._1 * f, _4._2 * f, _4._3 * f, _4._4 * f) 
) 
} 

另外，我可以加载这个矩阵到OpenGL中，或者我必须先转置它。我总是对这个数学感到困惑。

Answer 1

观看维基百科上的Invertible matrix: Analytic solution。顶部的整组计算将计算矩阵的adjugate，从中计算行列式，然后求逆的矩阵为1/det倍。

整个计算是明确地写出代码中的4×4矩阵，所以如果有存在一个bug会花费一些精力来检查整个事情。维基百科的文章解释了它应该如何工作。

Answer 2

反转矩阵通常是一个坏主意，因为计算可能会出现病态。

如果想要求解方程组，最好使用LU分解和前向 - 后向替换等方法来实现，特别是如果可以重新使用分解来求解右侧向量。

This link显示了一个用于旋转高斯消除的Java示例。

下面是另一个想法：也许您可以在应用程序中使用Java库，如JAMA的继任者Apache Commons Math？

如果您有特殊情况，我建议您将其输入Wolfram Alpha，以便在开始编码之前查看答案。

Answer 3

我很确定Simplex3D实现了这个计算（很可能它在那里正确地完成）。

Answer 4

如果你想玩数值 - 尽一切办法，自己动手 - 你有一些来自Jesper和duffymo的好建议（反转矩阵在实践中没有用 - 查看LU分解）。

然而，如果你只是想完成各项工作^TM直视Scalala和Scalalab。

无论哪种方式，你将需要Linear Algebra background知识 - 这是很多领域的难以置信的有用的数学。