Actual source code: ex301.c


  2: static char help[] = "Tests for bugs in A->offloadmask consistency for GPU matrices\n\n";

  4: #include <petscmat.h>

  6: int main(int argc, char **args)
  7: {
  8:   Mat         A;
  9:   PetscInt    i, j, rstart, rend, m = 3;
 10:   PetscScalar one = 1.0, zero = 0.0, negativeone = -1.0;
 11:   PetscReal   norm;
 12:   Vec         x, y;

 15:   PetscInitialize(&argc, &args, (char *)0, help);
 16:   PetscOptionsGetInt(NULL, NULL, "-m", &m, NULL);

 18:   for (i = 0; i < 2; i++) {
 19:     /* Create the matrix and set it to contain explicit zero entries on the diagonal. */
 20:     MatCreate(PETSC_COMM_WORLD, &A);
 21:     MatSetSizes(A, PETSC_DECIDE, PETSC_DECIDE, m * m, m * m);
 22:     MatSetFromOptions(A);
 23:     MatSetUp(A);
 24:     MatGetOwnershipRange(A, &rstart, &rend);
 25:     MatCreateVecs(A, &x, &y);
 26:     VecSet(x, one);
 27:     VecSet(y, zero);
 28:     MatDiagonalSet(A, y, INSERT_VALUES);

 30:     /* Now set A to be the identity using various approaches.
 31:      * Note that there may be other approaches that should be added here. */
 32:     switch (i) {
 33:     case 0:
 34:       MatDiagonalSet(A, x, INSERT_VALUES);
 35:       break;
 36:     case 1:
 37:       for (j = rstart; j < rend; j++) MatSetValue(A, j, j, one, INSERT_VALUES);
 38:       MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY);
 39:       MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY);
 40:       break;
 41:     case 2:
 42:       for (j = rstart; j < rend; j++) MatSetValuesRow(A, j, &one);
 43:       MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY);
 44:       MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY);
 45:     default:
 46:       break;
 47:     }

 49:     /* Compute y <- A*x and verify that the difference between y and x is negligible, as it should be since A is the identity. */
 50:     MatMult(A, x, y);
 51:     VecAXPY(y, negativeone, x);
 52:     VecNorm(y, NORM_2, &norm);
 53:     if (norm > PETSC_SQRT_MACHINE_EPSILON) PetscPrintf(PETSC_COMM_WORLD, "Test %" PetscInt_FMT ": Norm of error is %g, but should be near 0.\n", i, (double)norm);

 55:     MatDestroy(&A);
 56:     VecDestroy(&x);
 57:     VecDestroy(&y);
 58:   }

 60:   PetscFinalize();
 61:   return 0;
 62: }

 64: /*TEST

 66:    test:
 67:       suffix: aijviennacl_1
 68:       nsize: 1
 69:       args: -mat_type aijviennacl
 70:       requires: viennacl

 72:    test:
 73:       suffix: aijviennacl_2
 74:       nsize: 2
 75:       args: -mat_type aijviennacl
 76:       requires: viennacl

 78:    test:
 79:       suffix: aijcusparse_1
 80:       nsize: 1
 81:       args: -mat_type aijcusparse
 82:       requires: cuda

 84:    test:
 85:       suffix: aijcusparse_2
 86:       nsize: 2
 87:       args: -mat_type aijcusparse
 88:       requires: cuda
 89: TEST*/