simplex_interpolator: fixing bug on simplex matrix generation

src/lol/math/simplex_interpolator.h

@@ -127,22 +127,22 @@ protected:
         /*  Matrix coordinate transformation to skew simplex referential is done
             by inversing the base matrix M which is written as follows: 
 
-             M = | a b b b … |        M^(-1) = | c d d d … |
-                 | b a b b … |                 | d c d d … |
-                 | b b a b … |                 | d d c d … |
-                 | …         |                 | …         |
+            M = | a b b b … |        M^(-1) = | c d d d … |
+                | b a b b … |                 | d c d d … |
+                | b b a b … |                 | d d c d … |
+                | …         |                 | …         |
 
             where a and b are computed below ↴
         */
 
-        T b = (1 - sqrt((N+1)))/(N * N * N));
-        T a = b + sqrt((N+1)/N);
+        T b = (1 - sqrt((T)(N+1)))/sqrt((T)(N * N * N));
+        T a = b + sqrt((N+1)/(T)N);
 
         // determinant of matrix M
         T det = a * (a + (N-2) * b) - (b * b) * (N-1);
 
         T c = (a + (N-2)*b) / det;
-        T d = b / det;
+        T d = -b / det;
 
         for (int i = 0; i < N; ++i)
         {

src/t/math/simplex_interpolator.cpp

@@ -12,17 +12,65 @@
 
 #include <lol/engine-internal.h>
 
+#include <lol/math/simplex_interpolator.h>
+
 #include <lolunit.h>
 
 namespace lol
 {
 
+template<int N, typename T = float>
+class test_interpolator : public simplex_interpolator<N, T>
+{
+public:
+
+    test_interpolator(arraynd<N, T> const & samples) :
+        simplex_interpolator<N, T>(samples)
+    {
+    }
+
+
+    void DumpMatrix()
+    {
+        std::cout << std::endl;
+        for (int i = 0; i < N; ++i)
+        {
+            for (int j = 0; j < N; ++j)
+            {
+                std::cout << this->m_base[i][j] << ", ";
+            }
+            std::cout << std::endl;
+        }
+
+        std::cout << std::endl;
+
+        for (int i = 0; i < N; ++i)
+        {
+            for (int j = 0; j < N; ++j)
+            {
+                std::cout << this->m_base_inverse[i][j] << ", ";
+            }
+            std::cout << std::endl;
+        }
+    }
+};
+
 lolunit_declare_fixture(SimplexInterpolatorTest)
 {
     void SetUp() {}
 
     void TearDown() {}
 
+    lolunit_declare_test(CoordinateMatrixTest)
+    {
+        test_interpolator<2> s2({{1.f}});
+        s2.DumpMatrix();
+
+        test_interpolator<3> s3(arraynd<3, float>({{{1.f}}}));
+        s3.DumpMatrix();
+    }
+
+#if 0
     lolunit_declare_test(FloatGridPoints2D1x1)
     {
         simplex_interpolator<2> s({{1.f}});
@@ -151,6 +199,7 @@ private:
 
         return (float)i * vi + (float)j * vj;
     }
+#endif
 };
 
 }