Compute the solver only once to get massive speed boost

hspm/src/hspm.cpp

@@ -7,6 +7,8 @@ void HSPM::initHSPM()
     */
     
     this->computeConstantInfluenceCoeffs();
+    solver.compute(A_n);
+    
     this->blVel = Eigen::VectorXd::Zero(N);
     this->U = Eigen::VectorXd::Zero(N);
     this->V = Eigen::VectorXd::Zero(N);

hspm/src/hspm.h

@@ -73,7 +73,7 @@ public:
     Eigen::VectorXd U_wake;
     Eigen::VectorXd V_wake;
 private:
-    
+    Eigen::PartialPivLU<Eigen::MatrixXd> solver;
 };
 
 #endif
\ No newline at end of file

hspm/src/solver.cpp

@@ -28,19 +28,8 @@ void HSPM::solve()
     c += blVel;
 
     // Solve using Eigen's solver
-    Eigen::PartialPivLU<Eigen::MatrixXd> solver1;
-    //solver1.setTolerance(IT_SOLVER_TOLERANCE);
-    s1 = solver1.compute(A_n).solve(b);
-
-    Eigen::PartialPivLU<Eigen::MatrixXd> solver2;
-    //solver2.setTolerance(IT_SOLVER_TOLERANCE);
-    s2 = solver2.compute(A_n).solve(c);
-
-    // Check if the solver converged
-    //if (solver2.info() != Eigen::Success || solver1.info() != Eigen::Success) {
-    //    std::cout << "Iterative solver failed !" << std::endl;
-    //    exit(-1);
-    //}
+    s1 = solver.solve(b);
+    s2 = solver.solve(c);
 
     // q = s1*tau + s2
     // We need Kutta for tau