xref: /petsc/src/mat/tests/ex50.c (revision c4762a1b19cd2af06abeed90e8f9d34fb975dd94) 
1*c4762a1bSJed Brown static char help[] = "Tests MatView()/MatLoad() with binary viewers for SBAIJ matrices.\n\n";
2*c4762a1bSJed Brown 
3*c4762a1bSJed Brown #include <petscmat.h>
4*c4762a1bSJed Brown #include <petscviewer.h>
5*c4762a1bSJed Brown 
6*c4762a1bSJed Brown #include <petsc/private/hashtable.h>
7*c4762a1bSJed Brown static PetscReal MakeValue(PetscInt i,PetscInt j,PetscInt M)
8*c4762a1bSJed Brown {
9*c4762a1bSJed Brown   PetscHash_t h = PetscHashCombine(PetscHashInt(i),PetscHashInt(j));
10*c4762a1bSJed Brown   return (PetscReal) ((h % 5 == 0) ? (1 + i + j*M) : 0);
11*c4762a1bSJed Brown }
12*c4762a1bSJed Brown 
13*c4762a1bSJed Brown static PetscErrorCode CheckValuesAIJ(Mat A)
14*c4762a1bSJed Brown {
15*c4762a1bSJed Brown   PetscInt        M,N,rstart,rend,i,j;
16*c4762a1bSJed Brown   PetscReal       v,w;
17*c4762a1bSJed Brown   PetscScalar     val;
18*c4762a1bSJed Brown   PetscBool       seqsbaij,mpisbaij,sbaij;
19*c4762a1bSJed Brown   PetscErrorCode  ierr;
20*c4762a1bSJed Brown 
21*c4762a1bSJed Brown   PetscFunctionBegin;
22*c4762a1bSJed Brown   ierr = MatGetSize(A,&M,&N);CHKERRQ(ierr);
23*c4762a1bSJed Brown   ierr = MatGetOwnershipRange(A,&rstart,&rend);CHKERRQ(ierr);
24*c4762a1bSJed Brown   ierr = PetscObjectTypeCompare((PetscObject)A,MATSEQSBAIJ,&seqsbaij);CHKERRQ(ierr);
25*c4762a1bSJed Brown   ierr = PetscObjectTypeCompare((PetscObject)A,MATMPISBAIJ,&mpisbaij);CHKERRQ(ierr);
26*c4762a1bSJed Brown   sbaij = (seqsbaij||mpisbaij) ? PETSC_TRUE : PETSC_FALSE;
27*c4762a1bSJed Brown   for (i=rstart; i<rend; i++) {
28*c4762a1bSJed Brown     for (j=(sbaij?i:0); j<N; j++) {
29*c4762a1bSJed Brown       ierr = MatGetValue(A,i,j,&val);CHKERRQ(ierr);
30*c4762a1bSJed Brown       v = MakeValue(i,j,M); w = PetscRealPart(val);
31*c4762a1bSJed Brown       if (PetscAbsReal(v-w) > 0) SETERRQ4(PETSC_COMM_SELF,PETSC_ERR_PLIB,"Matrix entry (%D,%D) should be %g, got %g",i,j,(double)v,(double)w);
32*c4762a1bSJed Brown     }
33*c4762a1bSJed Brown   }
34*c4762a1bSJed Brown   PetscFunctionReturn(0);
35*c4762a1bSJed Brown }
36*c4762a1bSJed Brown 
37*c4762a1bSJed Brown int main(int argc,char **args)
38*c4762a1bSJed Brown {
39*c4762a1bSJed Brown   Mat            A;
40*c4762a1bSJed Brown   PetscInt       M = 24,N = 24,bs = 3;
41*c4762a1bSJed Brown   PetscInt       rstart,rend,i,j;
42*c4762a1bSJed Brown   PetscErrorCode ierr;
43*c4762a1bSJed Brown   PetscViewer    view;
44*c4762a1bSJed Brown 
45*c4762a1bSJed Brown   ierr = PetscInitialize(&argc,&args,NULL,help);if (ierr) return ierr;
46*c4762a1bSJed Brown   /*
47*c4762a1bSJed Brown       Create a parallel SBAIJ matrix shared by all processors
48*c4762a1bSJed Brown   */
49*c4762a1bSJed Brown   ierr = MatCreateSBAIJ(PETSC_COMM_WORLD,
50*c4762a1bSJed Brown                         bs,
51*c4762a1bSJed Brown                         PETSC_DECIDE,PETSC_DECIDE,
52*c4762a1bSJed Brown                         M,N,
53*c4762a1bSJed Brown                         PETSC_DECIDE,NULL,
54*c4762a1bSJed Brown                         PETSC_DECIDE,NULL,
55*c4762a1bSJed Brown                         &A);CHKERRQ(ierr);
56*c4762a1bSJed Brown 
57*c4762a1bSJed Brown   /*
58*c4762a1bSJed Brown       Set values into the matrix
59*c4762a1bSJed Brown   */
60*c4762a1bSJed Brown   ierr = MatGetSize(A,&M,&N);CHKERRQ(ierr);
61*c4762a1bSJed Brown   ierr = MatGetOwnershipRange(A,&rstart,&rend);CHKERRQ(ierr);
62*c4762a1bSJed Brown   for (i=rstart; i<rend; i++) {
63*c4762a1bSJed Brown     for (j=0; j<N; j++) {
64*c4762a1bSJed Brown       PetscReal v = MakeValue(i,j,M);
65*c4762a1bSJed Brown       if (PetscAbsReal(v) > 0) {
66*c4762a1bSJed Brown         ierr = MatSetValue(A,i,j,v,INSERT_VALUES);CHKERRQ(ierr);
67*c4762a1bSJed Brown       }
68*c4762a1bSJed Brown     }
69*c4762a1bSJed Brown   }
70*c4762a1bSJed Brown   ierr = MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr);
71*c4762a1bSJed Brown   ierr = MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr);
72*c4762a1bSJed Brown   ierr = MatViewFromOptions(A,NULL,"-mat_base_view");CHKERRQ(ierr);
73*c4762a1bSJed Brown 
74*c4762a1bSJed Brown   /*
75*c4762a1bSJed Brown       Store the binary matrix to a file
76*c4762a1bSJed Brown   */
77*c4762a1bSJed Brown   ierr = PetscViewerBinaryOpen(PETSC_COMM_WORLD, "matrix.dat", FILE_MODE_WRITE, &view);CHKERRQ(ierr);
78*c4762a1bSJed Brown   for (i=0; i<3; i++) {
79*c4762a1bSJed Brown     ierr = MatView(A,view);CHKERRQ(ierr);
80*c4762a1bSJed Brown   }
81*c4762a1bSJed Brown   ierr = PetscViewerDestroy(&view);CHKERRQ(ierr);
82*c4762a1bSJed Brown   ierr = MatDestroy(&A);CHKERRQ(ierr);
83*c4762a1bSJed Brown 
84*c4762a1bSJed Brown   /*
85*c4762a1bSJed Brown       Now reload the matrix and check its values
86*c4762a1bSJed Brown   */
87*c4762a1bSJed Brown   ierr = PetscViewerBinaryOpen(PETSC_COMM_WORLD,"matrix.dat",FILE_MODE_READ,&view);CHKERRQ(ierr);
88*c4762a1bSJed Brown   ierr = MatCreate(PETSC_COMM_WORLD,&A);CHKERRQ(ierr);
89*c4762a1bSJed Brown   ierr = MatSetType(A,MATSBAIJ);CHKERRQ(ierr);
90*c4762a1bSJed Brown   for (i=0; i<3; i++) {
91*c4762a1bSJed Brown     if (i > 0) {ierr = MatZeroEntries(A);CHKERRQ(ierr);}
92*c4762a1bSJed Brown     ierr = MatLoad(A,view);CHKERRQ(ierr);
93*c4762a1bSJed Brown     ierr = CheckValuesAIJ(A);CHKERRQ(ierr);
94*c4762a1bSJed Brown   }
95*c4762a1bSJed Brown   ierr = PetscViewerDestroy(&view);CHKERRQ(ierr);
96*c4762a1bSJed Brown   ierr = MatViewFromOptions(A,NULL,"-mat_load_view");CHKERRQ(ierr);
97*c4762a1bSJed Brown   ierr = MatDestroy(&A);CHKERRQ(ierr);
98*c4762a1bSJed Brown 
99*c4762a1bSJed Brown   /*
100*c4762a1bSJed Brown       Reload in SEQSBAIJ matrix and check its values
101*c4762a1bSJed Brown   */
102*c4762a1bSJed Brown   ierr = PetscViewerBinaryOpen(PETSC_COMM_SELF,"matrix.dat",FILE_MODE_READ,&view);CHKERRQ(ierr);
103*c4762a1bSJed Brown   ierr = MatCreate(PETSC_COMM_SELF,&A);CHKERRQ(ierr);
104*c4762a1bSJed Brown   ierr = MatSetType(A,MATSEQSBAIJ);CHKERRQ(ierr);
105*c4762a1bSJed Brown   for (i=0; i<3; i++) {
106*c4762a1bSJed Brown     if (i > 0) {ierr = MatZeroEntries(A);CHKERRQ(ierr);}
107*c4762a1bSJed Brown     ierr = MatLoad(A,view);CHKERRQ(ierr);
108*c4762a1bSJed Brown     ierr = CheckValuesAIJ(A);CHKERRQ(ierr);
109*c4762a1bSJed Brown   }
110*c4762a1bSJed Brown   ierr = PetscViewerDestroy(&view);CHKERRQ(ierr);
111*c4762a1bSJed Brown   ierr = MatDestroy(&A);CHKERRQ(ierr);
112*c4762a1bSJed Brown 
113*c4762a1bSJed Brown   /*
114*c4762a1bSJed Brown      Reload in MPISBAIJ matrix and check its values
115*c4762a1bSJed Brown   */
116*c4762a1bSJed Brown   ierr = PetscViewerBinaryOpen(PETSC_COMM_WORLD,"matrix.dat",FILE_MODE_READ,&view);CHKERRQ(ierr);
117*c4762a1bSJed Brown   ierr = MatCreate(PETSC_COMM_WORLD,&A);CHKERRQ(ierr);
118*c4762a1bSJed Brown   ierr = MatSetType(A,MATMPISBAIJ);CHKERRQ(ierr);
119*c4762a1bSJed Brown   for (i=0; i<3; i++) {
120*c4762a1bSJed Brown     if (i > 0) {ierr = MatZeroEntries(A);CHKERRQ(ierr);}
121*c4762a1bSJed Brown     ierr = MatLoad(A,view);CHKERRQ(ierr);
122*c4762a1bSJed Brown     ierr = CheckValuesAIJ(A);CHKERRQ(ierr);
123*c4762a1bSJed Brown   }
124*c4762a1bSJed Brown   ierr = PetscViewerDestroy(&view);CHKERRQ(ierr);
125*c4762a1bSJed Brown   ierr = MatDestroy(&A);CHKERRQ(ierr);
126*c4762a1bSJed Brown 
127*c4762a1bSJed Brown   ierr = PetscFinalize();
128*c4762a1bSJed Brown   return ierr;
129*c4762a1bSJed Brown }
130*c4762a1bSJed Brown 
131*c4762a1bSJed Brown /*TEST
132*c4762a1bSJed Brown 
133*c4762a1bSJed Brown    testset:
134*c4762a1bSJed Brown       args: -viewer_binary_mpiio 0
135*c4762a1bSJed Brown       output_file: output/ex50.out
136*c4762a1bSJed Brown       test:
137*c4762a1bSJed Brown         suffix: stdio_1
138*c4762a1bSJed Brown         nsize: 1
139*c4762a1bSJed Brown       test:
140*c4762a1bSJed Brown         suffix: stdio_2
141*c4762a1bSJed Brown         nsize: 2
142*c4762a1bSJed Brown       test:
143*c4762a1bSJed Brown         suffix: stdio_3
144*c4762a1bSJed Brown         nsize: 3
145*c4762a1bSJed Brown       test:
146*c4762a1bSJed Brown         suffix: stdio_4
147*c4762a1bSJed Brown         nsize: 4
148*c4762a1bSJed Brown       test:
149*c4762a1bSJed Brown         suffix: stdio_5
150*c4762a1bSJed Brown         nsize: 4
151*c4762a1bSJed Brown 
152*c4762a1bSJed Brown    testset:
153*c4762a1bSJed Brown       requires: mpiio
154*c4762a1bSJed Brown       args: -viewer_binary_mpiio 1
155*c4762a1bSJed Brown       output_file: output/ex50.out
156*c4762a1bSJed Brown       test:
157*c4762a1bSJed Brown         suffix: mpiio_1
158*c4762a1bSJed Brown         nsize: 1
159*c4762a1bSJed Brown       test:
160*c4762a1bSJed Brown         suffix: mpiio_2
161*c4762a1bSJed Brown         nsize: 2
162*c4762a1bSJed Brown       test:
163*c4762a1bSJed Brown         suffix: mpiio_3
164*c4762a1bSJed Brown         nsize: 3
165*c4762a1bSJed Brown       test:
166*c4762a1bSJed Brown         suffix: mpiio_4
167*c4762a1bSJed Brown         nsize: 4
168*c4762a1bSJed Brown       test:
169*c4762a1bSJed Brown         suffix: mpiio_5
170*c4762a1bSJed Brown         nsize: 5
171*c4762a1bSJed Brown 
172*c4762a1bSJed Brown TEST*/
173