51 Epetra_Map *& readMap,
52 const bool transpose,
const bool distribute,
53 bool& symmetric, Epetra_CrsMatrix *& Matrix ) {
55 Epetra_CrsMatrix * readA = 0;
56 Epetra_Vector * readx = 0;
57 Epetra_Vector * readb = 0;
58 Epetra_Vector * readxexact = 0;
64 FILE *in_file = fopen( in_filename,
"r");
68 filename = &in_filename[1] ;
71 filename = in_filename ;
76 std::string FileName = filename ;
78 int FN_Size = FileName.size() ;
79 std::string LastFiveBytes = FileName.substr( EPETRA_MAX(0,FN_Size-5), FN_Size );
80 std::string LastFourBytes = FileName.substr( EPETRA_MAX(0,FN_Size-4), FN_Size );
82 if ( LastFiveBytes ==
".triU" ) {
84 EPETRA_CHK_ERR( Trilinos_Util_ReadTriples2Epetra( filename,
false, Comm, readMap, readA, readx,
88 if ( LastFiveBytes ==
".triS" ) {
90 EPETRA_CHK_ERR( Trilinos_Util_ReadTriples2Epetra( filename,
true, Comm, readMap, readA, readx,
94 if ( LastFourBytes ==
".mtx" ) {
95 EPETRA_CHK_ERR( Trilinos_Util_ReadMatrixMarket2Epetra( filename, Comm, readMap,
96 readA, readx, readb, readxexact) );
97 FILE* in_file = fopen( filename,
"r");
98 assert (in_file != NULL) ;
99 const int BUFSIZE = 800 ;
100 char buffer[BUFSIZE] ;
101 fgets( buffer, BUFSIZE, in_file ) ;
102 std::string headerline1 = buffer;
104 if ( headerline1.find(
"symmetric") < BUFSIZE ) symmetric =
true;
106 if ( headerline1.find(
"symmetric") != std::string::npos) symmetric =
true;
113 Trilinos_Util_ReadHb2Epetra( filename, Comm, readMap, readA, readx,
115 if ( LastFourBytes ==
".rsa" ) symmetric = true ;
121 if ( readb )
delete readb;
122 if ( readx )
delete readx;
123 if ( readxexact )
delete readxexact;
125 Epetra_CrsMatrix *serialA ;
126 Epetra_CrsMatrix *transposeA;
129 transposeA =
new Epetra_CrsMatrix( Copy, *readMap, 0 );
131 serialA = transposeA ;
138 assert( (
void *) &serialA->Graph() ) ;
139 assert( (
void *) &serialA->RowMap() ) ;
140 assert( serialA->RowMap().SameAs(*readMap) ) ;
144 Epetra_Map DistMap(readMap->NumGlobalElements(), 0, Comm);
147 Epetra_Export exporter( *readMap, DistMap );
149 Epetra_CrsMatrix *Amat =
new Epetra_CrsMatrix( Copy, DistMap, 0 );
150 Amat->Export(*serialA, exporter, Add);
151 assert(Amat->FillComplete()==0);
190int main(
int argc,
char *argv[])
193 MPI_Init(&argc, &argv);
194 Epetra_MpiComm Comm(MPI_COMM_WORLD);
196 Epetra_SerialComm Comm;
199 bool verbose = (Comm.MyPID() == 0);
200 double TotalResidual = 0.0;
207 ParameterList GaleriList;
208 GaleriList.set(
"nx", 4);
209 GaleriList.set(
"ny", 4);
210 GaleriList.set(
"nz", 4 * Comm.NumProc());
211 GaleriList.set(
"mx", 1);
212 GaleriList.set(
"my", 1);
213 GaleriList.set(
"mz", Comm.NumProc());
214 Epetra_Map* Map = CreateMap(
"Cartesian3D", Comm, GaleriList);
222 Epetra_CrsMatrix* Matrix =
CreateCrsMatrix(
"Laplace3D", Map, GaleriList);
224 bool transpose = false ;
225 bool distribute = false ;
227 Epetra_CrsMatrix *Matrix = 0 ;
228 Epetra_Map *Map = 0 ;
229 CreateCrsMatrix(
"ibm.triU", Comm, Map, transpose, distribute, &symmetric, Matrix ) ;
237 Epetra_MultiVector LHS(*Map, 1);
238 Epetra_MultiVector RHS(*Map, 1);
241 Epetra_LinearProblem Problem(Matrix, &LHS, &RHS);
248 List.set(
"PrintTiming",
true);
249 List.set(
"PrintStatus",
true);
250 List.set(
"MaxProcs",Comm.NumProc());
252 std::vector<std::string> SolverType;
253 SolverType.push_back(
"Amesos_Paraklete");
255 Epetra_Time Time(Comm);
263 for (
unsigned int i = 0 ; i < SolverType.size() ; ++i)
266 if (Factory.
Query(SolverType[i]))
272 Matrix->Multiply(
false, LHS, RHS);
279 assert (Solver != 0);
285 Time.ResetStartTime();
297 if (TotalResidual > 1e-9)
304 return(EXIT_SUCCESS);