diff --git a/src/benchmarks/Sofa.Type/Matrix.cpp b/src/benchmarks/Sofa.Type/Matrix.cpp
index 0c8e1e9..f3e57ac 100644
--- a/src/benchmarks/Sofa.Type/Matrix.cpp
+++ b/src/benchmarks/Sofa.Type/Matrix.cpp
@@ -29,6 +29,10 @@ static void BM_Matrix_eigenmat33_assign(benchmark::State& state);
 static void BM_Matrix_typemat3x3f_multTranspose(benchmark::State& state);
 static void BM_Matrix_eigenmat3x3f_multTranspose(benchmark::State& state);
 
+
+static void BM_Matrix_typemat3x3f_noconstexpr_bench(benchmark::State& state);
+static void BM_Matrix_typemat3x3f_constexpr_bench(benchmark::State& state);
+
 constexpr int64_t minSubIterations = 8 << 4;
 constexpr int64_t maxSubIterations = 8 << 6;
 
@@ -58,6 +62,9 @@ BENCHMARK(BM_Matrix_eigenmat3x3f_multTranspose) BMARGS;
 
 #undef BMARGS
 
+BENCHMARK(BM_Matrix_typemat3x3f_noconstexpr_bench)->RangeMultiplier(2)->Ranges({ {minSubIterations, maxSubIterations} })->Unit(benchmark::kMicrosecond);
+BENCHMARK(BM_Matrix_typemat3x3f_constexpr_bench)->RangeMultiplier(2)->Ranges({ {minSubIterations, maxSubIterations} })->Unit(benchmark::kMicrosecond);
+
 void BM_Matrix_typemat3x3f_construct(benchmark::State& state)
 {
     constexpr auto totalsize = maxSubIterations * 9;
@@ -539,3 +546,39 @@ void BM_Matrix_eigenmat3x3f_multTranspose(benchmark::State& state)
         }
     }
 }
+
+void BM_Matrix_typemat3x3f_constexpr_bench(benchmark::State& state)
+{
+    constexpr auto totalsize = 9;
+    constexpr auto values = CompileTimeRandomValuePool<float, totalsize>::get();
+
+    constexpr sofa::type::Mat3x3f mat{
+                sofa::type::Mat3x3f::LineNoInit{values[0], values[1], values[2]},
+                sofa::type::Mat3x3f::LineNoInit{values[3], values[4], values[5]},
+                sofa::type::Mat3x3f::LineNoInit{values[6], values[7], values[8]}
+            };
+
+    for (auto _ : state)
+    {
+        constexpr sofa::type::Mat3x3f res = mat.transposed() * mat.transposed() * mat.transposed() / 3.0f;
+        benchmark::DoNotOptimize(res);
+    }
+}
+
+void BM_Matrix_typemat3x3f_noconstexpr_bench(benchmark::State& state)
+{
+    constexpr auto totalsize = maxSubIterations * 9;
+    const auto values = RandomValuePool<float, totalsize>::get();
+
+    const sofa::type::Mat3x3f mat{
+                sofa::type::Mat3x3f::LineNoInit{values[0], values[1], values[2]},
+                sofa::type::Mat3x3f::LineNoInit{values[3], values[4], values[5]},
+                sofa::type::Mat3x3f::LineNoInit{values[6], values[7], values[8]}
+    };
+
+    for (auto _ : state)
+    {
+        const sofa::type::Mat3x3f res = mat.transposed() * mat.transposed() * mat.transposed() / 3.0f;
+        benchmark::DoNotOptimize(res);
+    }
+}
diff --git a/src/utils/RandomValuePool.h b/src/utils/RandomValuePool.h
index 04f0957..0fb7309 100644
--- a/src/utils/RandomValuePool.h
+++ b/src/utils/RandomValuePool.h
@@ -23,3 +23,89 @@ struct RandomValuePool
     static inline bool bGenerated{ false };
     static inline std::array<Type, Size> v{};
 };
+
+// https://isocpp.org/blog/2017/01/cpp-weekly-episode-44-constexpr-compile-time-randomjason-turner
+template<typename Type>
+struct CompileTimeRandomValue
+{
+    static constexpr auto seed()
+    {
+        std::uint64_t shifted = 0;
+
+        for (const auto c : __TIME__)
+        {
+            shifted <<= 8;
+            shifted |= c;
+        }
+
+        return shifted;
+    }
+
+    struct PCG
+    {
+        struct pcg32_random_t
+        {
+            std::uint64_t state = 0;
+            std::uint64_t inc = seed();
+        };
+
+        pcg32_random_t rng;
+        typedef std::uint32_t result_type;
+
+        constexpr result_type operator()()
+        {
+            return pcg32_random_r();
+        }
+
+        static result_type constexpr min()
+        {
+            return std::numeric_limits<result_type>::min();
+        }
+
+        static result_type constexpr max()
+        {
+            return std::numeric_limits<result_type>::min();
+        }
+
+    private:
+        constexpr std::uint32_t pcg32_random_r()
+        {
+            std::uint64_t oldstate = rng.state;
+            // Advance internal state
+            rng.state = oldstate * 6364136223846793005ULL + (rng.inc | 1);
+            // Calculate output function (XSH RR), uses old state for max ILP
+            std::uint32_t xorshifted = ((oldstate >> 18u) ^ oldstate) >> 27u;
+            std::uint32_t rot = oldstate >> 59u;
+            return (xorshifted >> rot) | (xorshifted << ((-rot) & 31));
+        }
+
+    };
+
+    static constexpr auto get(int count, std::uint64_t s = seed())
+    {
+        PCG pcg;
+        pcg.rng.inc = s;
+
+        while (count > 0) {
+            pcg();
+            --count;
+        }
+
+        return pcg();
+    }
+};
+
+template<typename Type, std::size_t Size>
+struct CompileTimeRandomValuePool
+{
+    static constexpr auto get(Type minValue = static_cast<Type>(0), Type maxValue = static_cast<Type>(100), std::uint64_t seed = CompileTimeRandomValue<Type>::seed()) -> std::array<Type, Size>
+    {
+        std::array<Type, Size> v{};
+        for (auto i = 0; i < Size; i++)
+        {
+            v[i] = CompileTimeRandomValue<Type>::get(i, seed);
+        }
+        return v;
+    }
+};
+