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