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