逆行列

逆行列（inverse matrix）とは、行列の逆数のようなもの。

主に3DCGの世界では頂点の座標系を「戻す」のによく使われる。

たとえばバーテックスシェーダーでローカル座標系からワールド座標系に変換してしまった頂点位置をピクセルシェーダーでローカル座標系に戻したい場合などに使われる。

3x3行列の場合

$A^{-1}={\frac {1}{{\text{det}}(A)}}{\begin{bmatrix}(a_{11}a_{22}-a_{12}a_{21})&(a_{02}a_{21}-a_{01}a_{22})&(a_{01}a_{12}-a_{02}a_{11})\\(a_{12}a_{20}-a_{10}a_{22})&(a_{00}a_{22}-a_{02}a_{20})&(a_{02}a_{10}-a_{00}a_{12})\\(a_{10}a_{21}-a_{11}a_{20})&(a_{01}a_{20}-a_{00}a_{21})&(a_{00}a_{11}-a_{01}a_{10})\end{bmatrix}}$

C言語風に実装すると以下のような感じ。

bool inverseMatrix(float a[3][3], float inv[3][3]) {
    float det = a[0][0] * (a[1][1] * a[2][2] - a[2][1] * a[1][2]) 
              - a[0][1] * (a[1][0] * a[2][2] - a[1][2] * a[2][0]) 
              + a[0][2] * (a[1][0] * a[2][1] - a[1][1] * a[2][0]);

    if (det == 0) {
        printf("The matrix is singular and cannot be inverted.\n");
        return false;
    }

    float invdet = 1 / det;

    inv[0][0] =  (a[1][1] * a[2][2] - a[2][1] * a[1][2]) * invdet;
    inv[0][1] = -(a[0][1] * a[2][2] - a[0][2] * a[2][1]) * invdet;
    inv[0][2] =  (a[0][1] * a[1][2] - a[0][2] * a[1][1]) * invdet;

    inv[1][0] = -(a[1][0] * a[2][2] - a[1][2] * a[2][0]) * invdet;
    inv[1][1] =  (a[0][0] * a[2][2] - a[0][2] * a[2][0]) * invdet;
    inv[1][2] = -(a[0][0] * a[1][2] - a[1][0] * a[0][2]) * invdet;

    inv[2][0] =  (a[1][0] * a[2][1] - a[2][0] * a[1][1]) * invdet;
    inv[2][1] = -(a[0][0] * a[2][1] - a[2][0] * a[0][1]) * invdet;
    inv[2][2] =  (a[0][0] * a[1][1] - a[1][0] * a[0][1]) * invdet;

    return true;
}

4x4行列の場合

C言語風に実装すると以下のような感じ。

bool inverseMatrix(float a[4][4], float inv[4][4]) {
    float det = a[0][0] * (a[1][1] * (a[2][2] * a[3][3] - a[3][2] * a[2][3]) - a[2][1] * (a[1][2] * a[3][3] - a[3][2] * a[1][3]) + a[3][1] * (a[1][2] * a[2][3] - a[2][2] * a[1][3]))
              - a[0][1] * (a[1][0] * (a[2][2] * a[3][3] - a[3][2] * a[2][3]) - a[2][0] * (a[1][2] * a[3][3] - a[3][2] * a[1][3]) + a[3][0] * (a[1][2] * a[2][3] - a[2][2] * a[1][3]))
              + a[0][2] * (a[1][0] * (a[2][1] * a[3][3] - a[3][1] * a[2][3]) - a[2][0] * (a[1][1] * a[3][3] - a[3][1] * a[1][3]) + a[3][0] * (a[1][1] * a[2][3] - a[2][1] * a[1][3]))
              - a[0][3] * (a[1][0] * (a[2][1] * a[3][2] - a[3][1] * a[2][2]) - a[2][0] * (a[1][1] * a[3][2] - a[3][1] * a[1][2]) + a[3][0] * (a[1][1] * a[2][2] - a[2][1] * a[1][2]));

    if (det == 0) {
        printf("The matrix is singular and cannot be inverted.\n");
        return false;
    }

    float invdet = 1 / det;

    inv[0][0] =  ((a[1][1] * (a[2][2] * a[3][3] - a[3][2] * a[2][3]) - a[2][1] * (a[1][2] * a[3][3] - a[3][2] * a[1][3]) + a[3][1] * (a[1][2] * a[2][3] - a[2][2] * a[1][3])) * invdet;
    inv[0][1] = -((a[1][0] * (a[2][2] * a[3][3] - a[3][2] * a[2][3]) - a[2][0] * (a[1][2] * a[3][3] - a[3][2] * a[1][3]) + a[3][0] * (a[1][2] * a[2][3] - a[2][2] * a[1][3])) * invdet);
    inv[0][2] =  ((a[1][0] * (a[2][1] * a[3][3] - a[3][1] * a[2][3]) - a[2][0] * (a[1][1] * a[3][3] - a[3][1] * a[1][3]) + a[3][0] * (a[1][1] * a[2][3] - a[2][1] * a[1][3])) * invdet;
    inv[0][3] = -((a[1][0] * (a[2][1] * a[3][2] - a[3][1] * a[2][2]) - a[2][0] * (a[1][1] * a[3][2] - a[3][1] * a[1][2]) + a[3][0] * (a[1][1] * a[2][2] - a[2][1] * a[1][2])) * invdet);

    inv[1][0] = -((a[0][1] * (a[2][2] * a[3][3] - a[3][2] * a[2][3]) - a[2][1] * (a[0][2] * a[3][3] - a[3][2] * a[0][3]) + a[3][1] * (a[0][2] * a[2][3] - a[2][2] * a[0][3])) * invdet);
    inv[1][1] =  ((a[0][0] * (a[2][2] * a[3][3] - a[3][2] * a[2][3]) - a[2][0] * (a[0][2] * a[3][3] - a[3][2] * a[0][3]) + a[3][0] * (a[0][2] * a[2][3] - a[2][2] * a[0][3])) * invdet);
    inv[1][2] = -((a[0][0] * (a[2][1] * a[3][3] - a[3][1] * a[2][3]) - a[2][0] * (a[0][1] * a[3][3] - a[3][1] * a[0][3]) + a[3][0] * (a[0][1] * a[2][3] - a[2][1] * a[0][3])) * invdet);
    inv[1][3] =  ((a[0][0] * (a[2][1] * a[3][2] - a[3][1] * a[2][2]) - a[2][0] * (a[0][1] * a[3][2] - a[3][1] * a[0][2]) + a[3][0] * (a[0][1] * a[2][2] - a[2][1] * a[0][2])) * invdet);

    inv[2][0] =  ((a[0][1] * (a[1][2] * a[3][3] - a[3][2] * a[1][3]) - a[1][1] * (a[0][2] * a[3][3] - a[3][2] * a[0][3]) + a[3][1] * (a[0][2] * a[1][3] - a[1][2] * a[0][3])) * invdet);
    inv[2][1] = -((a[0][0] * (a[1][2] * a[3][3] - a[3][2] * a[1][3]) - a[1][0] * (a[0][2] * a[3][3] - a[3][2] * a[0][3]) + a[3][0] * (a[0][2] * a[1][3] - a[1][2] * a[0][3])) * invdet);
    inv[2][2] =  ((a[0][0] * (a[1][1] * a[3][3] - a[3][1] * a[1][3]) - a[1][0] * (a[0][1] * a[3][3] - a[3][1] * a[0][3]) + a[3][0] * (a[0][1] * a[1][3] - a[1][1] * a[0][3])) * invdet);
    inv[2][3] = -((a[0][0] * (a[1][1] * a[3][2] - a[3][1] * a[1][2]) - a[1][0] * (a[0][1] * a[3][2] - a[3][1] * a[0][2]) + a[3][0] * (a[0][1] * a[1][2] - a[1][1] * a[0][2])) * invdet);

    inv[3][0] = -((a[0][1] * (a[1][2] * a[2][3] - a[2][2] * a[1][3]) - a[1][1] * (a[0][2] * a[2][3] - a[2][2] * a[0][3]) + a[2][1] * (a[0][2] * a[1][3] - a[1][2] * a[0][3])) * invdet);
    inv[3][1] =  ((a[0][0] * (a[1][2] * a[2][3] - a[2][2] * a[1][3]) - a[1][0] * (a[0][2] * a[2][3] - a[2][2] * a[0][3]) + a[2][0] * (a[0][2] * a[1][3] - a[1][2] * a[0][3])) * invdet);
    inv[3][2] = -((a[0][0] * (a[1][1] * a[2][3] - a[2][1] * a[1][3]) - a[1][0] * (a[0][1] * a[2][3] - a[2][1] * a[0][3]) + a[2][0] * (a[0][1] * a[1][3] - a[1][1] * a[0][3])) * invdet);
    inv[3][3] =  ((a[0][0] * (a[1][1] * a[2][2] - a[2][1] * a[1][2]) - a[1][0] * (a[0][1] * a[2][2] - a[2][1] * a[0][2]) + a[2][0] * (a[0][1] * a[1][2] - a[1][1] * a[0][2])) * invdet);

    return true;
}

関連項目

• 行列
• 逆行列
• 転置行列