斯卡拉矩陣反演

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知識 - 這是很多領域的難以置信的有用的數學。