Use Eigen on CPUs.

BlueBrain · Nov 7, 2024 · 1f8f3c1 · 1f8f3c1
1 parent 25da002
commit 1f8f3c1
Showing 1 changed file with 6 additions and 2 deletions.
diff --git a/src/solver/newton/newton.hpp b/src/solver/newton/newton.hpp
@@ -95,6 +95,9 @@ EIGEN_DEVICE_FUNC int newton_solver(Eigen::Matrix<double, N, 1>& X,
         if (is_converged(X, J, F, eps)) {
             return iter;
         }
+        Eigen::Matrix<double, N, 1> X_solve;
+
+#if defined(CORENEURON_ENABLE_GPU)
         // In Eigen the default storage order is ColMajor.
         // Crout's implementation requires matrices stored in RowMajor order (C-style arrays).
         // Therefore, the transposeInPlace is critical such that the data() method to give the rows
@@ -107,9 +110,10 @@ EIGEN_DEVICE_FUNC int newton_solver(Eigen::Matrix<double, N, 1>& X,
         // Check if J is singular
         if (nmodl::crout::Crout<double>(N, J.data(), pivot.data(), rowmax.data()) < 0) {
             return -1;
-        }
-        Eigen::Matrix<double, N, 1> X_solve;
         nmodl::crout::solveCrout<double>(N, J.data(), F.data(), X_solve.data(), pivot.data());
+#else
+        X_solve = J.partialPivLu().solve(F);
+#endif
         X -= X_solve;
     }
     // If we fail to converge after max_iter iterations, return -1