42 #include "Teuchos_UnitTestHarness.hpp" 43 #include "HYPRE_IJ_mv.h" 46 #include "Epetra_ConfigDefs.h" 47 #include "Epetra_Vector.h" 48 #include "Epetra_RowMatrix.h" 49 #include "Epetra_MultiVector.h" 50 #include "Epetra_CrsMatrix.h" 51 #include "Epetra_Map.h" 54 #include "Epetra_MpiComm.h" 56 #include "Epetra_SerialComm.h" 60 #include "Teuchos_ParameterList.hpp" 61 #include "Teuchos_ParameterEntry.hpp" 62 #include "Teuchos_ParameterListExceptions.hpp" 63 #include "Teuchos_Array.hpp" 72 const double tol = 1E-6;
79 TEST_EQUALITY(Matrix->Filled(),
true);
81 for(
int i = 0; i < Matrix->NumMyRows(); i++){
84 TEST_EQUALITY(entries, 1);
86 Teuchos::Array<double> Values; Values.resize(entries);
87 Teuchos::Array<int> Indices; Indices.resize(entries);
88 ierr += Matrix->
ExtractMyRowCopy(i, entries, numentries, &Values[0], &Indices[0]);
89 TEST_EQUALITY(ierr, 0);
90 TEST_EQUALITY(numentries,1);
91 for(
int j = 0; j < numentries; j++){
92 TEST_FLOATING_EQUALITY(Values[j],1.0,
tol);
93 TEST_EQUALITY(Indices[j],i);
104 Epetra_MultiVector X(Matrix->RowMatrixRowMap(), num_vectors,
true);
106 Epetra_MultiVector Y(Matrix->RowMatrixRowMap(), num_vectors,
true);
108 ierr += Matrix->
Multiply(
false, X, Y);
111 TEST_EQUALITY(ierr, 0);
122 Epetra_MultiVector X(Matrix->RowMatrixRowMap(), num_vectors,
true);
124 Epetra_MultiVector Y1(Matrix->RowMatrixRowMap(), num_vectors,
true);
125 Epetra_MultiVector Y2(Matrix->RowMatrixRowMap(), num_vectors,
true);
127 ierr += Matrix->
Multiply(
false, X, Y1);
128 ierr += TestMat->Multiply(
false, X, Y2);
132 ierr += Matrix->
Multiply(
false, Y1, X);
133 ierr += TestMat->Multiply(
false, Y1, Y2);
137 ierr += Matrix->
Multiply(
false, Y2, X);
138 ierr += TestMat->Multiply(
false, Y2, Y1);
141 TEST_EQUALITY_CONST(ierr, 0);
153 Epetra_MultiVector X(Matrix->RowMatrixRowMap(), num_vectors,
true);
155 Epetra_MultiVector Y1(Matrix->RowMatrixRowMap(), num_vectors,
true);
156 Epetra_MultiVector Y2(Matrix->RowMatrixRowMap(), num_vectors,
true);
158 ierr += Matrix->
Multiply(
true, X, Y1);
159 ierr += TestMat->Multiply(
true, X, Y2);
162 TEST_EQUALITY(ierr, 0);
174 Epetra_Vector X(Matrix->RowMatrixRowMap(),
true);
176 Matrix->NumMyNonzeros();
178 ierr += TestMat->LeftScale(X);
180 TEST_EQUALITY(ierr, 0);
190 Epetra_Vector X(Matrix->RowMatrixRowMap(),
true);
192 Matrix->NumMyNonzeros();
194 ierr += TestMat->RightScale(X);
197 TEST_EQUALITY(ierr, 0);
208 Epetra_Vector X(Matrix->RowMatrixRowMap(),
true);
209 Epetra_Vector Y(TestMat->RowMatrixRowMap(),
true);
211 ierr += Matrix->ExtractDiagonalCopy(X);
212 ierr += TestMat->ExtractDiagonalCopy(Y);
215 TEST_EQUALITY(ierr, 0);
226 Epetra_Vector X(Matrix->RowMatrixRowMap(),
true);
227 Epetra_Vector Y(TestMat->RowMatrixRowMap(),
true);
229 ierr += Matrix->InvRowSums(X);
230 ierr += TestMat->InvRowSums(Y);
233 TEST_EQUALITY(ierr, 0);
244 Epetra_Vector X(Matrix->RowMatrixColMap(),
true);
245 Epetra_Vector Y(TestMat->RowMatrixColMap(),
true);
247 ierr += Matrix->InvColSums(X);
248 ierr += TestMat->InvColSums(Y);
251 TEST_EQUALITY(ierr, 0);
261 double norm1 = Matrix->NormInf();
262 double norm2 = TestMat->NormInf();
264 TEST_FLOATING_EQUALITY(norm1, norm2,
tol);
273 double norm1 = Matrix->NormOne();
274 double norm2 = TestMat->NormOne();
276 TEST_FLOATING_EQUALITY(norm1, norm2,
tol);
286 int nnz1 = Matrix->NumGlobalNonzeros();
287 int nnz2 = TestMat->NumGlobalNonzeros();
289 TEST_EQUALITY(nnz1, nnz2);
299 int rows1 = Matrix->NumGlobalRows();
300 int rows2 = TestMat->NumGlobalRows();
302 TEST_EQUALITY(rows1, rows2);
312 int cols1 = Matrix->NumGlobalCols();
313 int cols2 = TestMat->NumGlobalCols();
315 TEST_EQUALITY(cols1, cols2);
325 int hdiag1 = Matrix->NumGlobalDiagonals();
326 int Ediag2 = TestMat->NumGlobalDiagonals();
328 TEST_EQUALITY(hdiag1, Ediag2);
338 int nnz1 = Matrix->NumMyNonzeros();
339 int nnz2 = TestMat->NumMyNonzeros();
341 TEST_EQUALITY(nnz1, nnz2);
351 int rows1 = Matrix->NumMyRows();
352 int rows2 = TestMat->NumMyRows();
354 TEST_EQUALITY(rows1, rows2);
364 int cols1 = Matrix->NumMyCols();
365 int cols2 = TestMat->NumMyCols();
367 TEST_EQUALITY(cols1, cols2);
377 int diag1 = Matrix->NumMyDiagonals();
378 int diag2 = TestMat->NumMyDiagonals();
380 TEST_EQUALITY(diag1, diag2);
390 int ent1 = Matrix->MaxNumEntries();
391 int ent2 = TestMat->MaxNumEntries();
393 TEST_EQUALITY(ent1, ent2);
404 Epetra_MultiVector X(Matrix->RowMatrixRowMap(), num_vectors,
true);
406 Epetra_MultiVector Y1(Matrix->RowMatrixRowMap(), num_vectors,
true);
407 Epetra_MultiVector Y2(Matrix->RowMatrixRowMap(), num_vectors,
true);
409 TestMat->ApplyInverse(X,Y2);
422 Epetra_MultiVector X1(Matrix->RowMatrixRowMap(), num_vectors,
true);
424 Epetra_MultiVector X2(X1);
427 TestMat->Multiply(
false,X2,X2);
438 Epetra_MultiVector TrueX(Matrix->RowMatrixRowMap(), num_vectors,
true);
440 Epetra_MultiVector RHS(Matrix->RowMatrixRowMap(), num_vectors,
true);
442 Epetra_MultiVector X(Matrix->RowMatrixRowMap(), num_vectors,
true);
464 Matrix->
Solve(
false,
false,
false, RHS, X);
int ExtractMyRowCopy(int MyRow, int Length, int &NumEntries, double *Values, int *Indices) const
Returns a copy of the specified local row in user-provided arrays.
int NumMyRowEntries(int MyRow, int &NumEntries) const
Return the current number of values stored for the specified local row.
int LeftScale(const Epetra_Vector &X)
Scales the EpetraExt_HypreIJMatrix on the left with a Epetra_Vector x.
bool EquivalentVectors(Epetra_MultiVector &Y1, Epetra_MultiVector &Y2, const double tol)
bool EquivalentMatrices(Epetra_RowMatrix &HypreMatrix, Epetra_RowMatrix &CrsMatrix, const double tol)
int ApplyInverse(const Epetra_MultiVector &X, Epetra_MultiVector &Y) const
Returns the result of a Epetra_Operator inverse applied to an Epetra_MultiVector X in Y...
EpetraExt::BlockCrsMatrix: A class for constructing a distributed block matrix.
Epetra_CrsMatrix * newCrsMatrix(EpetraExt_HypreIJMatrix &Matrix)
int RightScale(const Epetra_Vector &X)
Scales the EpetraExt_HypreIJMatrix on the right with a Epetra_Vector x.
TEUCHOS_UNIT_TEST(EpetraExt_hypre, Construct)
int Multiply(bool TransA, const Epetra_MultiVector &X, Epetra_MultiVector &Y) const
Returns the result of a EpetraExt_HypreIJMatrix multiplied by a Epetra_MultiVector X in Y...
EpetraExt_HypreIJMatrix * newHypreMatrix(const int N, const int type)
int Solve(bool Upper, bool Trans, bool UnitDiagonal, const Epetra_MultiVector &X, Epetra_MultiVector &Y) const
Returns the result of a EpetraExt_HypreIJMatrix solving a Epetra_MultiVector X in Y...
int SetParameter(Hypre_Chooser chooser, int(*pt2Func)(HYPRE_Solver, int), int parameter)
Set a parameter that takes a single int.