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