4#define OUR_CHK_ERR(a) { { int epetra_err = a; \
5 if (epetra_err != 0) { std::cerr << "Amesos ERROR " << epetra_err << ", " \
6 << __FILE__ << ", line " << __LINE__ << std::endl; \
7relresidual=1e15; return(1);} }\
11#include "Epetra_Comm.h"
12#include "Teuchos_ParameterList.hpp"
14#include "Epetra_CrsMatrix.h"
15#include "Epetra_Map.h"
16#include "Epetra_Vector.h"
17#include "Epetra_LinearProblem.h"
27 const Epetra_Comm &Comm,
30 Teuchos::ParameterList ParamList,
31 Epetra_CrsMatrix *& Amat,
38 int iam = Comm.MyPID() ;
41 const Epetra_Map *RangeMap =
42 transpose?&Amat->OperatorDomainMap():&Amat->OperatorRangeMap() ;
43 const Epetra_Map *DomainMap =
44 transpose?&Amat->OperatorRangeMap():&Amat->OperatorDomainMap() ;
46 Epetra_Vector xexact(*DomainMap);
47 Epetra_Vector x(*DomainMap);
49 Epetra_Vector b(*RangeMap);
50 Epetra_Vector bcheck(*RangeMap);
52 Epetra_Vector difference(*DomainMap);
54 Epetra_LinearProblem Problem;
58 Abase = Afactory.
Create( AmesosClass, Problem ) ;
60 std::string AC = AmesosClass ;
61 if ( AC ==
"Amesos_Mumps" ) {
62 ParamList.set(
"NoDestroy",
true );
66 double relresidual = 0 ;
72 Problem.SetOperator( Amat );
86 xexact.PutScalar(1.0);
91 Amat->Multiply( transpose, xexact, b ) ;
94 std::cout << __FILE__ <<
"::" << __LINE__ <<
"b = " << std::endl ;
95 b.Print( std::cout ) ;
96 std::cout << __FILE__ <<
"::" << __LINE__ <<
"xexact = " << std::endl ;
97 xexact.Print( std::cout ) ;
98 std::cout << __FILE__ <<
"::" << __LINE__ <<
"x = " << std::endl ;
99 x.Print( std::cout ) ;
107 Problem.SetLHS( &x );
108 Problem.SetRHS( &b );
122 Amat->Multiply( transpose, x, bcheck ) ;
127 difference.Update( 1.0, x, -1.0, xexact, 0.0 ) ;
128 difference.Norm2( &norm_diff ) ;
129 x.Norm2( &norm_one ) ;
131 relresidual = norm_diff / norm_one ;
134 if ( relresidual * Rcond > 1e-16 ) {
135 if (
verbose) std::cout << __FILE__ <<
"::"<< __LINE__
136 <<
" norm( x - xexact ) / norm(x) = "
137 << norm_diff /norm_one << std::endl ;
157 const Epetra_Comm &Comm,
160 Teuchos::ParameterList ParamList,
161 Epetra_CrsMatrix *& Amat,
165 std::string AC = AmesosClass ;
167 if (epetra_err != 0) {
168 std::cerr <<
"Amesos ERROR " << epetra_err <<
", " << __FILE__ <<
", line " << __LINE__ << std::endl;
int PartialFactorizationOneStep(const char *AmesosClass, const Epetra_Comm &Comm, bool transpose, bool verbose, Teuchos::ParameterList ParamList, Epetra_CrsMatrix *&Amat, double Rcond, int Steps)
int PartialFactorization(const char *AmesosClass, const Epetra_Comm &Comm, bool transpose, bool verbose, Teuchos::ParameterList ParamList, Epetra_CrsMatrix *&Amat, double Rcond)
Amesos_BaseSolver: A pure virtual class for direct solution of real-valued double-precision operators...
virtual int Solve()=0
Solves A X = B (or AT x = B)
virtual int NumericFactorization()=0
Performs NumericFactorization on the matrix A.
virtual int SymbolicFactorization()=0
Performs SymbolicFactorization on the matrix A.
virtual int SetUseTranspose(bool UseTranspose)=0
If set true, X will be set to the solution of AT X = B (not A X = B)
virtual int SetParameters(Teuchos::ParameterList &ParameterList)=0
Updates internal variables.
Factory for binding a third party direct solver to an Epetra_LinearProblem.
Amesos_BaseSolver * Create(const char *ClassType, const Epetra_LinearProblem &LinearProblem)
Amesos Create method.