1 2 static char help[] = "Tests for bugs in A->offloadmask consistency for GPU matrices\n\n"; 3 4 #include <petscmat.h> 5 6 int main(int argc, char **args) { 7 Mat A; 8 PetscInt i, j, rstart, rend, m = 3; 9 PetscScalar one = 1.0, zero = 0.0, negativeone = -1.0; 10 PetscReal norm; 11 Vec x, y; 12 13 PetscFunctionBeginUser; 14 PetscCall(PetscInitialize(&argc, &args, (char *)0, help)); 15 PetscCall(PetscOptionsGetInt(NULL, NULL, "-m", &m, NULL)); 16 17 for (i = 0; i < 2; i++) { 18 /* Create the matrix and set it to contain explicit zero entries on the diagonal. */ 19 PetscCall(MatCreate(PETSC_COMM_WORLD, &A)); 20 PetscCall(MatSetSizes(A, PETSC_DECIDE, PETSC_DECIDE, m * m, m * m)); 21 PetscCall(MatSetFromOptions(A)); 22 PetscCall(MatSetUp(A)); 23 PetscCall(MatGetOwnershipRange(A, &rstart, &rend)); 24 PetscCall(MatCreateVecs(A, &x, &y)); 25 PetscCall(VecSet(x, one)); 26 PetscCall(VecSet(y, zero)); 27 PetscCall(MatDiagonalSet(A, y, INSERT_VALUES)); 28 29 /* Now set A to be the identity using various approaches. 30 * Note that there may be other approaches that should be added here. */ 31 switch (i) { 32 case 0: PetscCall(MatDiagonalSet(A, x, INSERT_VALUES)); break; 33 case 1: 34 for (j = rstart; j < rend; j++) { PetscCall(MatSetValue(A, j, j, one, INSERT_VALUES)); } 35 PetscCall(MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY)); 36 PetscCall(MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY)); 37 break; 38 case 2: 39 for (j = rstart; j < rend; j++) { PetscCall(MatSetValuesRow(A, j, &one)); } 40 PetscCall(MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY)); 41 PetscCall(MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY)); 42 default: break; 43 } 44 45 /* Compute y <- A*x and verify that the difference between y and x is negligible, as it should be since A is the identity. */ 46 PetscCall(MatMult(A, x, y)); 47 PetscCall(VecAXPY(y, negativeone, x)); 48 PetscCall(VecNorm(y, NORM_2, &norm)); 49 if (norm > PETSC_SQRT_MACHINE_EPSILON) { PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Test %" PetscInt_FMT ": Norm of error is %g, but should be near 0.\n", i, (double)norm)); } 50 51 PetscCall(MatDestroy(&A)); 52 PetscCall(VecDestroy(&x)); 53 PetscCall(VecDestroy(&y)); 54 } 55 56 PetscCall(PetscFinalize()); 57 return 0; 58 } 59 60 /*TEST 61 62 test: 63 suffix: aijviennacl_1 64 nsize: 1 65 args: -mat_type aijviennacl 66 requires: viennacl 67 68 test: 69 suffix: aijviennacl_2 70 nsize: 2 71 args: -mat_type aijviennacl 72 requires: viennacl 73 74 test: 75 suffix: aijcusparse_1 76 nsize: 1 77 args: -mat_type aijcusparse 78 requires: cuda 79 80 test: 81 suffix: aijcusparse_2 82 nsize: 2 83 args: -mat_type aijcusparse 84 requires: cuda 85 TEST*/ 86