92 const std::string Descriptor,
93 const Epetra_RowMatrix& A,
94 const Epetra_MultiVector& x,
95 const Epetra_MultiVector& b,
96 const Epetra_MultiVector& x_exact)
98 if (A.Comm().MyPID() == 0)
99 std::cout << std::endl <<
"*** " << SolverType <<
" : "
100 << Descriptor <<
" ***" << std::endl;
101 std::vector<double> Norm;
102 int NumVectors = x.NumVectors();
103 Norm.resize(NumVectors);
104 Epetra_MultiVector Ax(x);
105 A.Multiply(
false,x,Ax);
106 Ax.Update(1.0,b,-1.0);
108 bool TestPassed =
false;
109 double TotalNorm = 0.0;
110 for (
int i = 0 ; i < NumVectors ; ++i) {
111 TotalNorm += Norm[i];
113 if (A.Comm().MyPID() == 0)
114 std::cout <<
"||Ax - b|| = " << TotalNorm << std::endl;
115 if (TotalNorm < 1e-5 )
120 Ax.Update (1.0,x,-1.0,x_exact,0.0);
122 for (
int i = 0 ; i < NumVectors ; ++i) {
123 TotalNorm += Norm[i];
125 if (A.Comm().MyPID() == 0)
126 std::cout <<
"||x - x_exact|| = " << TotalNorm << std::endl;
127 if (TotalNorm < 1e-5 )
137 Epetra_RowMatrix& A, Epetra_MultiVector& x_A,
138 Epetra_MultiVector& b_A, Epetra_MultiVector& x_exactA,
139 Epetra_RowMatrix& B, Epetra_MultiVector& x_B,
140 Epetra_MultiVector& b_B, Epetra_MultiVector& x_exactB,
141 Epetra_RowMatrix& C, Epetra_MultiVector& x_C,
142 Epetra_MultiVector& b_C, Epetra_MultiVector& x_exactC)
144 bool TestPassed =
true;
145 Epetra_LinearProblem ProblemA;
146 Epetra_LinearProblem ProblemB;
147 Epetra_LinearProblem ProblemC;
157 ProblemA.SetOperator((Epetra_RowMatrix*)0);
158 ProblemA.SetLHS((Epetra_MultiVector*)0);
159 ProblemA.SetRHS((Epetra_MultiVector*)0);
161 Solver = Factory.Create(SolverType,ProblemA);
163 ProblemA.SetOperator(&A);
164 ProblemA.SetLHS(&x_A);
165 ProblemA.SetRHS(&b_A);
167 std::string ST = SolverType ;
168 if (! ( ST ==
"Amesos_Superludist" ) ) {
171 TestPassed = TestPassed &&
172 CheckError(SolverType,
"Solve() only", A,x_A,b_A,x_exactA) &&
187 ProblemA.SetOperator((Epetra_RowMatrix*)0);
188 ProblemA.SetLHS((Epetra_MultiVector*)0);
189 ProblemA.SetRHS((Epetra_MultiVector*)0);
191 Solver = Factory.Create(SolverType,ProblemA);
193 ProblemA.SetOperator(&A);
195 ProblemA.SetLHS(&x_A);
196 ProblemA.SetRHS(&b_A);
199 TestPassed = TestPassed &&
200 CheckError(SolverType,
"NumFact() + Solve()", A,x_A,b_A,x_exactA) &&
219 ProblemA.SetOperator((Epetra_RowMatrix*)0);
220 ProblemA.SetLHS((Epetra_MultiVector*)0);
221 ProblemA.SetRHS((Epetra_MultiVector*)0);
223 Solver = Factory.Create(SolverType,ProblemA);
225 ProblemA.SetOperator(&A);
233 ProblemA.SetLHS(&x_A);
234 ProblemA.SetRHS(&b_A);
243 TestPassed = TestPassed &&
244 CheckError(SolverType,
"SymFact() + NumFact() + Solve()",
245 A,x_A,b_A,x_exactA) &&
263 ProblemA.SetOperator((Epetra_RowMatrix*)0);
264 ProblemA.SetLHS((Epetra_MultiVector*)0);
265 ProblemA.SetRHS((Epetra_MultiVector*)0);
267 Solver = Factory.Create(SolverType,ProblemA);
269 ProblemA.SetOperator(&A);
271 ProblemA.SetOperator(&B);
273 ProblemA.SetLHS(&x_B);
274 ProblemA.SetRHS(&b_B);
278 TestPassed = TestPassed &&
279 CheckError(SolverType,
"Set A, solve B", B,x_B,b_B,x_exactB) &&
292 ProblemA.SetOperator(&A);
293 ProblemA.SetLHS(&x_A);
294 ProblemA.SetRHS(&b_A);
296 Solver = Factory.Create(SolverType,ProblemA);
298 ProblemA.SetOperator(&C);
300 std::string ST = SolverType ;
301 if (! ( ST ==
"Amesos_Superludist" ) ) {
306 ProblemA.SetLHS(&x_C);
307 ProblemA.SetRHS(&b_C);
310 TestPassed = TestPassed &&
311 CheckError(SolverType,
"Set A, Solve C", C,x_C,b_C,x_exactC) &&
324 ProblemA.SetOperator(&A);
325 ProblemA.SetLHS(&x_A);
326 ProblemA.SetRHS(&b_A);
328 Solver = Factory.Create(SolverType,ProblemA);
330 std::string ST = SolverType ;
331 if (! ( ST ==
"Amesos_Superludist" ) ) {
334 ProblemA.SetOperator(&C);
338 ProblemA.SetLHS(&x_C);
339 ProblemA.SetRHS(&b_C);
342 TestPassed = TestPassed &&
343 CheckError(SolverType,
"Solve A + Solve C", C,x_C,b_C,x_exactC) &&
363 int NumVectors_AB = 7;
364 int NumVectors_C = 13;
367 if ( RUNNING_ON_VALGRIND ) {
375 Teuchos::ParameterList GaleriList;
377 GaleriList.set(
"n", Size_AB * Size_AB);
378 GaleriList.set(
"nx", Size_AB);
379 GaleriList.set(
"ny", Size_AB);
380 Epetra_Map* Map_A = CreateMap(
"Interlaced", Comm, GaleriList);
381 Epetra_CrsMatrix* Matrix_A =
CreateCrsMatrix(
"Recirc2D", Map_A, GaleriList);
382 Epetra_CrsMatrix* Matrix_B =
CreateCrsMatrix(
"Laplace2D", Map_A, GaleriList);
384 GaleriList.set(
"n", Size_C);
385 Epetra_Map* Map_C = CreateMap(
"Interlaced", Comm, GaleriList);
386 Epetra_CrsMatrix* Matrix_C =
CreateCrsMatrix(
"Minij", Map_A, GaleriList);
404 Epetra_MultiVector x_C(C.OperatorDomainMap(),NumVectors_C);
405 Epetra_MultiVector x_exactC(C.OperatorDomainMap(),NumVectors_C);
406 Epetra_MultiVector b_C(C.OperatorRangeMap(),NumVectors_C);
408 C.Multiply(
false,x_exactC,b_C);
410 std::vector<std::string> SolverType;
411 SolverType.push_back(
"Amesos_Klu");
412 SolverType.push_back(
"Amesos_Lapack");
413 SolverType.push_back(
"Amesos_Umfpack");
414 SolverType.push_back(
"Amesos_Superlu");
415 SolverType.push_back(
"Amesos_Superludist");
416 SolverType.push_back(
"Amesos_Mumps");
420 bool TestPassed =
true;
422 for (
unsigned int i = 0 ; i < SolverType.size() ; ++i) {
423 std::string Solver = SolverType[i];
424 if (Factory.Query((
char*)Solver.c_str())) {
425 bool ok =
Test((
char*)Solver.c_str(),
426 A, x_A, b_A, x_exactA,
427 B, x_B, b_B, x_exactB,
428 C, x_C, b_C, x_exactC);
429 TestPassed = TestPassed && ok;
433 if (Comm.MyPID() == 0)
434 std::cout <<
"Solver " << Solver <<
" not available" << std::endl;
445 if (Comm.MyPID() == 0)
446 std::cout << std::endl <<
"TEST PASSED" << std::endl << std::endl;
447 return(EXIT_SUCCESS);
450 if (Comm.MyPID() == 0)
451 std::cout << std::endl <<
"TEST FAILED" << std::endl << std::endl;
int main(int argc, char *argv[])
int SubMain(Epetra_Comm &Comm)
bool TestTiming(const Amesos_BaseSolver *Solver, const Epetra_Comm &Comm)
bool CheckError(const std::string SolverType, const std::string Descriptor, const Epetra_RowMatrix &A, const Epetra_MultiVector &x, const Epetra_MultiVector &b, const Epetra_MultiVector &x_exact)
bool Test(char *SolverType, Epetra_RowMatrix &A, Epetra_MultiVector &x_A, Epetra_MultiVector &b_A, Epetra_MultiVector &x_exactA, Epetra_RowMatrix &B, Epetra_MultiVector &x_B, Epetra_MultiVector &b_B, Epetra_MultiVector &x_exactB, Epetra_RowMatrix &C, Epetra_MultiVector &x_C, Epetra_MultiVector &b_C, Epetra_MultiVector &x_exactC)