Actual source code: ex225.c


  2: static char help[] = "Test Hypre matrix APIs\n";

  4: #include <petscmathypre.h>

  6: int main(int argc, char **args)
  7: {
  8:   Mat         A, B, C;
  9:   PetscReal   err;
 10:   PetscInt    i, j, M = 20;
 11:   PetscMPIInt NP;
 12:   MPI_Comm    comm;
 13:   PetscInt   *rows;

 16:   PetscInitialize(&argc, &args, (char *)0, help);
 17:   comm = PETSC_COMM_WORLD;
 18:   MPI_Comm_size(comm, &NP);
 19:   PetscOptionsGetInt(NULL, NULL, "-M", &M, NULL);
 21:   /* Hypre matrix */
 22:   MatCreate(comm, &B);
 23:   MatSetSizes(B, PETSC_DECIDE, PETSC_DECIDE, M, M);
 24:   MatSetType(B, MATHYPRE);
 25:   MatHYPRESetPreallocation(B, 9, NULL, 9, NULL);

 27:   /* PETSc AIJ matrix */
 28:   MatCreate(comm, &A);
 29:   MatSetSizes(A, PETSC_DECIDE, PETSC_DECIDE, M, M);
 30:   MatSetType(A, MATAIJ);
 31:   MatSeqAIJSetPreallocation(A, 9, NULL);
 32:   MatMPIAIJSetPreallocation(A, 9, NULL, 9, NULL);

 34:   /*Set Values */
 35:   for (i = 0; i < M; i++) {
 36:     PetscInt    cols[]  = {0, 1, 2, 3, 4, 5};
 37:     PetscScalar vals[6] = {0};
 38:     PetscScalar value[] = {100};
 39:     for (j = 0; j < 6; j++) vals[j] = ((PetscReal)j) / NP;

 41:     MatSetValues(B, 1, &i, 6, cols, vals, ADD_VALUES);
 42:     MatSetValues(B, 1, &i, 1, &i, value, ADD_VALUES);
 43:     MatSetValues(A, 1, &i, 6, cols, vals, ADD_VALUES);
 44:     MatSetValues(A, 1, &i, 1, &i, value, ADD_VALUES);
 45:   }

 47:   /* MAT_FLUSH_ASSEMBLY currently not supported */
 48:   MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY);
 49:   MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY);
 50:   MatAssemblyBegin(B, MAT_FINAL_ASSEMBLY);
 51:   MatAssemblyEnd(B, MAT_FINAL_ASSEMBLY);

 53:   /* Compare A and B */
 54:   MatConvert(B, MATAIJ, MAT_INITIAL_MATRIX, &C);
 55:   MatAXPY(C, -1., A, SAME_NONZERO_PATTERN);
 56:   MatNorm(C, NORM_INFINITY, &err);
 58:   MatDestroy(&C);

 60:   /* MatZeroRows */
 61:   PetscMalloc1(M, &rows);
 62:   for (i = 0; i < M; i++) rows[i] = i;
 63:   MatZeroRows(B, M, rows, 10.0, NULL, NULL);
 64:   MatSetOption(A, MAT_KEEP_NONZERO_PATTERN, PETSC_TRUE);
 65:   MatZeroRows(A, M, rows, 10.0, NULL, NULL);
 66:   MatConvert(B, MATAIJ, MAT_INITIAL_MATRIX, &C);
 67:   MatAXPY(C, -1., A, SAME_NONZERO_PATTERN);
 68:   MatNorm(C, NORM_INFINITY, &err);
 70:   MatDestroy(&C);
 71:   PetscFree(rows);

 73:   /* Test MatZeroEntries */
 74:   MatZeroEntries(B);
 75:   MatConvert(B, MATAIJ, MAT_INITIAL_MATRIX, &C);
 76:   MatNorm(C, NORM_INFINITY, &err);
 78:   MatDestroy(&C);

 80:   /* Insert Values */
 81:   for (i = 0; i < M; i++) {
 82:     PetscInt    cols[]  = {0, 1, 2, 3, 4, 5};
 83:     PetscScalar vals[6] = {0};
 84:     PetscScalar value[] = {100};

 86:     for (j = 0; j < 6; j++) vals[j] = ((PetscReal)j) / NP;

 88:     MatSetValues(B, 1, &i, 6, cols, vals, INSERT_VALUES);
 89:     MatSetValues(B, 1, &i, 1, &i, value, INSERT_VALUES);
 90:     MatSetValues(A, 1, &i, 6, cols, vals, INSERT_VALUES);
 91:     MatSetValues(A, 1, &i, 1, &i, value, INSERT_VALUES);
 92:   }

 94:   /* MAT_FLUSH_ASSEMBLY currently not supported */
 95:   MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY);
 96:   MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY);
 97:   MatAssemblyBegin(B, MAT_FINAL_ASSEMBLY);
 98:   MatAssemblyEnd(B, MAT_FINAL_ASSEMBLY);

100:   /* Rows are not sorted with HYPRE so we need an intermediate sort
101:      They use a temporary buffer, so we can sort inplace the const memory */
102:   {
103:     const PetscInt    *idxA, *idxB;
104:     const PetscScalar *vA, *vB;
105:     PetscInt           rstart, rend, nzA, nzB;
106:     PetscInt           cols[] = {0, 1, 2, 3, 4, -5};
107:     PetscInt          *rows;
108:     PetscScalar       *valuesA, *valuesB;
109:     PetscBool          flg;

111:     MatGetOwnershipRange(A, &rstart, &rend);
112:     for (i = rstart; i < rend; i++) {
113:       MatGetRow(A, i, &nzA, &idxA, &vA);
114:       MatGetRow(B, i, &nzB, &idxB, &vB);
116:       PetscSortIntWithScalarArray(nzB, (PetscInt *)idxB, (PetscScalar *)vB);
117:       PetscArraycmp(idxA, idxB, nzA, &flg);
119:       PetscArraycmp(vA, vB, nzA, &flg);
121:       MatRestoreRow(A, i, &nzA, &idxA, &vA);
122:       MatRestoreRow(B, i, &nzB, &idxB, &vB);
123:     }

125:     MatGetOwnershipRange(A, &rstart, &rend);
126:     PetscCalloc3((rend - rstart) * 6, &valuesA, (rend - rstart) * 6, &valuesB, rend - rstart, &rows);
127:     for (i = rstart; i < rend; i++) rows[i - rstart] = i;

129:     MatGetValues(A, rend - rstart, rows, 6, cols, valuesA);
130:     MatGetValues(B, rend - rstart, rows, 6, cols, valuesB);

132:     for (i = 0; i < (rend - rstart); i++) {
133:       PetscArraycmp(valuesA + 6 * i, valuesB + 6 * i, 6, &flg);
135:     }
136:     PetscFree3(valuesA, valuesB, rows);
137:   }

139:   /* Compare A and B */
140:   MatConvert(B, MATAIJ, MAT_INITIAL_MATRIX, &C);
141:   MatAXPY(C, -1., A, SAME_NONZERO_PATTERN);
142:   MatNorm(C, NORM_INFINITY, &err);

145:   MatDestroy(&A);
146:   MatDestroy(&B);
147:   MatDestroy(&C);

149:   PetscFinalize();
150:   return 0;
151: }

153: /*TEST

155:    build:
156:       requires: hypre

158:    test:
159:       suffix: 1
160:       requires: !defined(PETSC_HAVE_HYPRE_DEVICE)

162:    test:
163:       suffix: 2
164:       requires: !defined(PETSC_HAVE_HYPRE_DEVICE)
165:       output_file: output/ex225_1.out
166:       nsize: 2

168: TEST*/