Work on UMFPACK solver for multiple rhs.

AMDiS/src/OEMSolver.h

@@ -66,6 +66,7 @@ namespace AMDiS {
 	print_cycle(100),
 	iterations(-1),
 	error(-1),
+	multipleRhs(false),
 	leftPrecon(NULL),
 	rightPrecon(NULL)
     {}
@@ -95,8 +96,9 @@ namespace AMDiS {
      */
     void setLeftPrecon(ITL_BasePreconditioner* p) 
     {
-	if (leftPrecon) delete leftPrecon;
-	leftPrecon= p;
+	if (leftPrecon) 
+	  delete leftPrecon;
+	leftPrecon = p;
     }
 
     /** Set right Preconditioner
@@ -106,8 +108,9 @@ namespace AMDiS {
      */
     void setRightPrecon(ITL_BasePreconditioner* p) 
     {
-	if (rightPrecon) delete rightPrecon;
-	rightPrecon= p;
+	if (rightPrecon) 
+	  delete rightPrecon;
+	rightPrecon = p;
     }
 
     /// Linear System to be solved in the derived class
@@ -248,6 +251,11 @@ namespace AMDiS {
       info = i;
     }
 
+    void setMultipleRhs(bool b)
+    {
+      multipleRhs = b;
+    }
+
     /** \} */
 
   protected:
@@ -278,6 +286,13 @@ namespace AMDiS {
     /// Error code  in last solver (not set by UmfPack)
     int error;
 
+    /** \brief
+     * If true, the solver has to solve multiple right hand sides with the same
+     * matrix. Some solvers, e.g. direct once, may take advanteges from this knowledge,
+     * as they can do the factorization of the matrix only once.
+     */
+    bool multipleRhs;
+
     /// Left preconditioner (ignored by some iterative solvers and by UmfPack)
     ITL_BasePreconditioner *leftPrecon;
 

AMDiS/src/UmfPackSolver.h

@@ -66,8 +66,6 @@ namespace AMDiS {
     {
       GET_PARAMETER(0, name + "->store symbolic", "%d", &store_symbolic);
       GET_PARAMETER(0, name + "->symmetric strategy", "%d", &symmetric_strategy);
-
-      // GET_PARAMETER(0, name + "->multiple rhs", "%d", &multipleRhs);
     }
     
     /// Destructor
@@ -88,15 +86,19 @@ namespace AMDiS {
 	else
 	  solver = new mtl::matrix::umfpack::solver<matrix_type>(A, UMFPACK_STRATEGY_SYMMETRIC);
       } else {
-	if (store_symbolic)
-	  solver->update_numeric();
-	else
-	  solver->update();
+	//	if (!multipleRhs)
+ 	  if (store_symbolic)
+ 	    solver->update_numeric();
+ 	  else
+	    solver->update();
       }
-      
+     
       int code= (*solver)(x, b);
-      mtl::dense_vector<value_type> r(b); r-= A * x; residual= two_norm(r);
+      mtl::dense_vector<value_type> r(b); 
+      r -= A * x; 
+      residual = two_norm(r);
       std::cout << "UmfPackSolver: ||b-Ax|| = " << residual << "\n";
+
       return code;
     }