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;

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

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

 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:       MatDiagonalSet(A,x,INSERT_VALUES);
 34:       break;
 35:     case 1:
 36:       for (j=rstart; j<rend; j++) {
 37:         MatSetValue(A,j,j,one,INSERT_VALUES);
 38:       }
 39:       MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY);
 40:       MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY);
 41:       break;
 42:     case 2:
 43:       for (j=rstart; j<rend; j++) {
 44:         MatSetValuesRow(A,j,&one);
 45:       }
 46:       MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY);
 47:       MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY);
 48:     default:
 49:       break;
 50:     }

 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:     MatMult(A,x,y);
 54:     VecAXPY(y,negativeone,x);
 55:     VecNorm(y,NORM_2,&norm);
 56:     if (norm > PETSC_SQRT_MACHINE_EPSILON) {
 57:       PetscPrintf(PETSC_COMM_WORLD,"Test %" PetscInt_FMT ": Norm of error is %g, but should be near 0.\n",i,(double)norm);
 58:     }

 60:     MatDestroy(&A);
 61:     VecDestroy(&x);
 62:     VecDestroy(&y);
 63:   }

 65:   PetscFinalize();
 66:   return 0;
 67: }

 69: /*TEST

 71:    test:
 72:       suffix: aijviennacl_1
 73:       nsize: 1
 74:       args: -mat_type aijviennacl
 75:       requires: viennacl

 77:    test:
 78:       suffix: aijviennacl_2
 79:       nsize: 2
 80:       args: -mat_type aijviennacl
 81:       requires: viennacl

 83:    test:
 84:       suffix: aijcusparse_1
 85:       nsize: 1
 86:       args: -mat_type aijcusparse
 87:       requires: cuda

 89:    test:
 90:       suffix: aijcusparse_2
 91:       nsize: 2
 92:       args: -mat_type aijcusparse
 93:       requires: cuda
 94: TEST*/