From 0e974f910924e67d881733207005b69c0bc924be Mon Sep 17 00:00:00 2001
From: Guillaume Bittoun <guillaume.bittoun@gmail.com>
Date: Tue, 20 Jan 2015 21:49:15 +0000
Subject: [PATCH] matrix: using LU-decomposition for matrix inverse (losing
 precision unfortunately)

---
 src/lol/math/matrix.h | 13 ++++---------
 src/t/math/matrix.cpp | 18 +++++++++---------
 2 files changed, 13 insertions(+), 18 deletions(-)

diff --git a/src/lol/math/matrix.h b/src/lol/math/matrix.h
index 3dd94271..4612b11f 100644
--- a/src/lol/math/matrix.h
+++ b/src/lol/math/matrix.h
@@ -612,16 +612,11 @@ mat_t<T, N, N> u_inverse(mat_t<T, N, N> const & U)
 template<typename T, int N>
 mat_t<T, N, N> inverse(mat_t<T, N, N> const &m)
 {
-    mat_t<T, N, N> ret;
+    mat_t<T, N, N> L, U;
+    lu_decomposition(m, L, U);
+
+    mat_t<T, N, N> ret = u_inverse(U) * l_inverse(L);
 
-    T d = determinant(m);
-    if (d)
-    {
-        /* Transpose result matrix on the fly */
-        for (int i = 0; i < N; ++i)
-            for (int j = 0; j < N; ++j)
-                ret[i][j] = cofactor(m, j, i) / d;
-    }
     return ret;
 }
 
diff --git a/src/t/math/matrix.cpp b/src/t/math/matrix.cpp
index 3aec5dc0..63463c0b 100644
--- a/src/t/math/matrix.cpp
+++ b/src/t/math/matrix.cpp
@@ -222,17 +222,17 @@ lolunit_declare_fixture(MatrixTest)
 
         mat3 m2 = m0 * m1;
 
-        lolunit_assert_doubles_equal(m2[0][0], 1.0f, 1e-5);
-        lolunit_assert_doubles_equal(m2[1][0], 0.0f, 1e-5);
-        lolunit_assert_doubles_equal(m2[2][0], 0.0f, 1e-5);
+        lolunit_assert_doubles_equal(m2[0][0], 1.0f, 1e-4);
+        lolunit_assert_doubles_equal(m2[1][0], 0.0f, 1e-4);
+        lolunit_assert_doubles_equal(m2[2][0], 0.0f, 1e-4);
 
-        lolunit_assert_doubles_equal(m2[0][1], 0.0f, 1e-5);
-        lolunit_assert_doubles_equal(m2[1][1], 1.0f, 1e-5);
-        lolunit_assert_doubles_equal(m2[2][1], 0.0f, 1e-5);
+        lolunit_assert_doubles_equal(m2[0][1], 0.0f, 1e-4);
+        lolunit_assert_doubles_equal(m2[1][1], 1.0f, 1e-4);
+        lolunit_assert_doubles_equal(m2[2][1], 0.0f, 1e-4);
 
-        lolunit_assert_doubles_equal(m2[0][2], 0.0f, 1e-5);
-        lolunit_assert_doubles_equal(m2[1][2], 0.0f, 1e-5);
-        lolunit_assert_doubles_equal(m2[2][2], 1.0f, 1e-5);
+        lolunit_assert_doubles_equal(m2[0][2], 0.0f, 1e-4);
+        lolunit_assert_doubles_equal(m2[1][2], 0.0f, 1e-4);
+        lolunit_assert_doubles_equal(m2[2][2], 1.0f, 1e-4);
     }
 
     lolunit_declare_test(Inverse4x4)