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 { 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; 13 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: 33 PetscCall(MatDiagonalSet(A,x,INSERT_VALUES)); 34 break; 35 case 1: 36 for (j=rstart; j<rend; j++) { 37 PetscCall(MatSetValue(A,j,j,one,INSERT_VALUES)); 38 } 39 PetscCall(MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY)); 40 PetscCall(MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY)); 41 break; 42 case 2: 43 for (j=rstart; j<rend; j++) { 44 PetscCall(MatSetValuesRow(A,j,&one)); 45 } 46 PetscCall(MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY)); 47 PetscCall(MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY)); 48 default: 49 break; 50 } 51 52 /* Compute y <- A*x and verify that the difference between y and x is negligible, as it should be since A is the identity. */ 53 PetscCall(MatMult(A,x,y)); 54 PetscCall(VecAXPY(y,negativeone,x)); 55 PetscCall(VecNorm(y,NORM_2,&norm)); 56 if (norm > PETSC_SQRT_MACHINE_EPSILON) { 57 PetscCall(PetscPrintf(PETSC_COMM_WORLD,"Test %" PetscInt_FMT ": Norm of error is %g, but should be near 0.\n",i,(double)norm)); 58 } 59 60 PetscCall(MatDestroy(&A)); 61 PetscCall(VecDestroy(&x)); 62 PetscCall(VecDestroy(&y)); 63 } 64 65 PetscCall(PetscFinalize()); 66 return 0; 67 } 68 69 /*TEST 70 71 test: 72 suffix: aijviennacl_1 73 nsize: 1 74 args: -mat_type aijviennacl 75 requires: viennacl 76 77 test: 78 suffix: aijviennacl_2 79 nsize: 2 80 args: -mat_type aijviennacl 81 requires: viennacl 82 83 test: 84 suffix: aijcusparse_1 85 nsize: 1 86 args: -mat_type aijcusparse 87 requires: cuda 88 89 test: 90 suffix: aijcusparse_2 91 nsize: 2 92 args: -mat_type aijcusparse 93 requires: cuda 94 TEST*/ 95