From 4f73690b506dbed07a470f62acf917fc36f9c687 Mon Sep 17 00:00:00 2001 From: CodeMyst Date: Sun, 21 Jun 2020 18:30:49 +0200 Subject: [PATCH] added inverse matrix for mat4 --- source/cosmomyst/math/matrix.d | 80 ++++++++++++++++++++++++++++++++++ 1 file changed, 80 insertions(+) diff --git a/source/cosmomyst/math/matrix.d b/source/cosmomyst/math/matrix.d index 5ebcbff..1d095d8 100644 --- a/source/cosmomyst/math/matrix.d +++ b/source/cosmomyst/math/matrix.d @@ -253,6 +253,80 @@ struct mat(ulong n) if (n >= 2) { 0f, 0f, -1f, 0f); } +/++ + + returns the inverse of the provided matrix, if no inverse can be found it returns a mat4(inf) + +/ +@nogc mat4 mat_inverse(const mat4 a) pure nothrow { + mat4 t = a; + + float det2_01_01 = t[0, 0] * t[1, 1] - t[0, 1] * t[1, 0]; + float det2_01_02 = t[0, 0] * t[1, 2] - t[0, 2] * t[1, 0]; + float det2_01_03 = t[0, 0] * t[1, 3] - t[0, 3] * t[1, 0]; + float det2_01_12 = t[0, 1] * t[1, 2] - t[0, 2] * t[1, 1]; + float det2_01_13 = t[0, 1] * t[1, 3] - t[0, 3] * t[1, 1]; + float det2_01_23 = t[0, 2] * t[1, 3] - t[0, 3] * t[1, 2]; + + float det3_201_012 = t[2, 0] * det2_01_12 - t[2, 1] * det2_01_02 + t[2, 2] * det2_01_01; + float det3_201_013 = t[2, 0] * det2_01_13 - t[2, 1] * det2_01_03 + t[2, 3] * det2_01_01; + float det3_201_023 = t[2, 0] * det2_01_23 - t[2, 2] * det2_01_03 + t[2, 3] * det2_01_02; + float det3_201_123 = t[2, 1] * det2_01_23 - t[2, 2] * det2_01_13 + t[2, 3] * det2_01_12; + + float det = - det3_201_123 * t[3, 0] + det3_201_023 * t[3, 1] - det3_201_013 * t[3, 2] + det3_201_012 * t[3, 3]; + float invDet = 1 / det; + + float det2_03_01 = t[0, 0] * t[3, 1] - t[0, 1] * t[3, 0]; + float det2_03_02 = t[0, 0] * t[3, 2] - t[0, 2] * t[3, 0]; + float det2_03_03 = t[0, 0] * t[3, 3] - t[0, 3] * t[3, 0]; + float det2_03_12 = t[0, 1] * t[3, 2] - t[0, 2] * t[3, 1]; + float det2_03_13 = t[0, 1] * t[3, 3] - t[0, 3] * t[3, 1]; + float det2_03_23 = t[0, 2] * t[3, 3] - t[0, 3] * t[3, 2]; + float det2_13_01 = t[1, 0] * t[3, 1] - t[1, 1] * t[3, 0]; + float det2_13_02 = t[1, 0] * t[3, 2] - t[1, 2] * t[3, 0]; + float det2_13_03 = t[1, 0] * t[3, 3] - t[1, 3] * t[3, 0]; + float det2_13_12 = t[1, 1] * t[3, 2] - t[1, 2] * t[3, 1]; + float det2_13_13 = t[1, 1] * t[3, 3] - t[1, 3] * t[3, 1]; + float det2_13_23 = t[1, 2] * t[3, 3] - t[1, 3] * t[3, 2]; + + float det3_203_012 = t[2, 0] * det2_03_12 - t[2, 1] * det2_03_02 + t[2, 2] * det2_03_01; + float det3_203_013 = t[2, 0] * det2_03_13 - t[2, 1] * det2_03_03 + t[2, 3] * det2_03_01; + float det3_203_023 = t[2, 0] * det2_03_23 - t[2, 2] * det2_03_03 + t[2, 3] * det2_03_02; + float det3_203_123 = t[2, 1] * det2_03_23 - t[2, 2] * det2_03_13 + t[2, 3] * det2_03_12; + + float det3_213_012 = t[2, 0] * det2_13_12 - t[2, 1] * det2_13_02 + t[2, 2] * det2_13_01; + float det3_213_013 = t[2, 0] * det2_13_13 - t[2, 1] * det2_13_03 + t[2, 3] * det2_13_01; + float det3_213_023 = t[2, 0] * det2_13_23 - t[2, 2] * det2_13_03 + t[2, 3] * det2_13_02; + float det3_213_123 = t[2, 1] * det2_13_23 - t[2, 2] * det2_13_13 + t[2, 3] * det2_13_12; + + float det3_301_012 = t[3, 0] * det2_01_12 - t[3, 1] * det2_01_02 + t[3, 2] * det2_01_01; + float det3_301_013 = t[3, 0] * det2_01_13 - t[3, 1] * det2_01_03 + t[3, 3] * det2_01_01; + float det3_301_023 = t[3, 0] * det2_01_23 - t[3, 2] * det2_01_03 + t[3, 3] * det2_01_02; + float det3_301_123 = t[3, 1] * det2_01_23 - t[3, 2] * det2_01_13 + t[3, 3] * det2_01_12; + + mat4 res; + + res[0, 0] = - det3_213_123 * invDet; + res[1, 0] = + det3_213_023 * invDet; + res[2, 0] = - det3_213_013 * invDet; + res[3, 0] = + det3_213_012 * invDet; + + res[0, 1] = + det3_203_123 * invDet; + res[1, 1] = - det3_203_023 * invDet; + res[2, 1] = + det3_203_013 * invDet; + res[3, 1] = - det3_203_012 * invDet; + + res[0, 2] = + det3_301_123 * invDet; + res[1, 2] = - det3_301_023 * invDet; + res[2, 2] = + det3_301_013 * invDet; + res[3, 2] = - det3_301_012 * invDet; + + res[0, 3] = - det3_201_123 * invDet; + res[1, 3] = + det3_201_023 * invDet; + res[2, 3] = - det3_201_013 * invDet; + res[3, 3] = + det3_201_012 * invDet; + + return res; +} + unittest { auto t1 = mat2(2f); auto t2 = mat2(1f, 2f, 3f, 4f); @@ -324,3 +398,9 @@ unittest { auto m1 = mat3(3f); assert(m1 * trans == mat3(3f, 3f, 27f, 3f, 3f, 27f, 3f, 3f, 27f)); } + +unittest { + + auto m1 = mat4(5f, 6f, 6f, 8f, 2f, 2f, 2f, 8f, 6f, 6f, 2f, 8f, 2f, 3f, 6f, 7f); + assert(mat_inverse(m1) * m1 == mat_ident!4()); +}