zlm icon indicating copy to clipboard operation
zlm copied to clipboard

Add matrix inverse function

Open whatupdave opened this issue 4 years ago • 3 comments

As an example here's one for a 4x4 matrix

fn inverseMat4(m: Mat4) Mat4 {
    const A2323 = m.fields[2][2] * m.fields[3][3] - m.fields[2][3] * m.fields[3][2];
    const A1323 = m.fields[2][1] * m.fields[3][3] - m.fields[2][3] * m.fields[3][1];
    const A1223 = m.fields[2][1] * m.fields[3][2] - m.fields[2][2] * m.fields[3][1];
    const A0323 = m.fields[2][0] * m.fields[3][3] - m.fields[2][3] * m.fields[3][0];
    const A0223 = m.fields[2][0] * m.fields[3][2] - m.fields[2][2] * m.fields[3][0];
    const A0123 = m.fields[2][0] * m.fields[3][1] - m.fields[2][1] * m.fields[3][0];
    const A2313 = m.fields[1][2] * m.fields[3][3] - m.fields[1][3] * m.fields[3][2];
    const A1313 = m.fields[1][1] * m.fields[3][3] - m.fields[1][3] * m.fields[3][1];
    const A1213 = m.fields[1][1] * m.fields[3][2] - m.fields[1][2] * m.fields[3][1];
    const A2312 = m.fields[1][2] * m.fields[2][3] - m.fields[1][3] * m.fields[2][2];
    const A1312 = m.fields[1][1] * m.fields[2][3] - m.fields[1][3] * m.fields[2][1];
    const A1212 = m.fields[1][1] * m.fields[2][2] - m.fields[1][2] * m.fields[2][1];
    const A0313 = m.fields[1][0] * m.fields[3][3] - m.fields[1][3] * m.fields[3][0];
    const A0213 = m.fields[1][0] * m.fields[3][2] - m.fields[1][2] * m.fields[3][0];
    const A0312 = m.fields[1][0] * m.fields[2][3] - m.fields[1][3] * m.fields[2][0];
    const A0212 = m.fields[1][0] * m.fields[2][2] - m.fields[1][2] * m.fields[2][0];
    const A0113 = m.fields[1][0] * m.fields[3][1] - m.fields[1][1] * m.fields[3][0];
    const A0112 = m.fields[1][0] * m.fields[2][1] - m.fields[1][1] * m.fields[2][0];

    const det = 1 / (m.fields[0][0] * (m.fields[1][1] * A2323 - m.fields[1][2] * A1323 + m.fields[1][3] * A1223) - m.fields[0][1] * (m.fields[1][0] * A2323 - m.fields[1][2] * A0323 + m.fields[1][3] * A0223) + m.fields[0][2] * (m.fields[1][0] * A1323 - m.fields[1][1] * A0323 + m.fields[1][3] * A0123) - m.fields[0][3] * (m.fields[1][0] * A1223 - m.fields[1][1] * A0223 + m.fields[1][2] * A0123));

    return Mat4{
        .fields = [4][4]f32{
            [4]f32{
                det * (m.fields[1][1] * A2323 - m.fields[1][2] * A1323 + m.fields[1][3] * A1223),
                det * -(m.fields[0][1] * A2323 - m.fields[0][2] * A1323 + m.fields[0][3] * A1223),
                det * (m.fields[0][1] * A2313 - m.fields[0][2] * A1313 + m.fields[0][3] * A1213),
                det * -(m.fields[0][1] * A2312 - m.fields[0][2] * A1312 + m.fields[0][3] * A1212),
            },
            [4]f32{
                det * -(m.fields[1][0] * A2323 - m.fields[1][2] * A0323 + m.fields[1][3] * A0223),
                det * (m.fields[0][0] * A2323 - m.fields[0][2] * A0323 + m.fields[0][3] * A0223),
                det * -(m.fields[0][0] * A2313 - m.fields[0][2] * A0313 + m.fields[0][3] * A0213),
                det * (m.fields[0][0] * A2312 - m.fields[0][2] * A0312 + m.fields[0][3] * A0212),
            },
            [4]f32{
                det * (m.fields[1][0] * A1323 - m.fields[1][1] * A0323 + m.fields[1][3] * A0123),
                det * -(m.fields[0][0] * A1323 - m.fields[0][1] * A0323 + m.fields[0][3] * A0123),
                det * (m.fields[0][0] * A1313 - m.fields[0][1] * A0313 + m.fields[0][3] * A0113),
                det * -(m.fields[0][0] * A1312 - m.fields[0][1] * A0312 + m.fields[0][3] * A0112),
            },
            [4]f32{
                det * -(m.fields[1][0] * A1223 - m.fields[1][1] * A0223 + m.fields[1][2] * A0123),
                det * (m.fields[0][0] * A1223 - m.fields[0][1] * A0223 + m.fields[0][2] * A0123),
                det * -(m.fields[0][0] * A1213 - m.fields[0][1] * A0213 + m.fields[0][2] * A0113),
                det * (m.fields[0][0] * A1212 - m.fields[0][1] * A0212 + m.fields[0][2] * A0112),
            },
        },
    };
}

test "inverse" {
    const input = Mat4{
        .fields = [4][4]f32{
            [4]f32{ 4, 0, 0, 0 },
            [4]f32{ 0, 0, 2, 0 },
            [4]f32{ 0, 1, 2, 0 },
            [4]f32{ 1, 0, 0, 1 },
        },
    };

    const expected = Mat4{
        .fields = [4][4]f32{
            [4]f32{ 0.25, 0, 0, 0 },
            [4]f32{ 0, -1, 1, 0 },
            [4]f32{ 0, 0.5, 0, 0 },
            [4]f32{ -0.25, 0, 0, 1 },
        },
    };

    assert(std.meta.eql(inverseMat4(input), expected));
}

Which is ported from the output of this program: https://github.com/willnode/N-Matrix-Programmer

whatupdave avatar Oct 04 '20 22:10 whatupdave

Nice! Could you just do a PR integrating this?

ikskuh avatar Oct 05 '20 07:10 ikskuh

If I can figure out how to generalise it for Mat2,3,4 i'll send a PR.

whatupdave avatar Oct 05 '20 19:10 whatupdave

Thank you!

ikskuh avatar Oct 06 '20 14:10 ikskuh

@MasterQ32 solved by 8500fcacc3050e09c4f068ecf00525ec4b03327b, should be closed.

Flecart avatar Sep 22 '23 19:09 Flecart