xref: /petsc/src/mat/tests/ex16.c (revision 117ef88edefbfc12e7c19efe87a19a2e1b0acd4f)
1 static char help[] = "Tests MatDenseGetArray() and MatView()/MatLoad() with binary viewers.\n\n";
2 
3 #include <petscmat.h>
4 #include <petscviewer.h>
5 
6 static PetscErrorCode CheckValues(Mat A,PetscBool one)
7 {
8   const PetscScalar *array;
9   PetscInt          M,N,rstart,rend,lda,i,j;
10   PetscErrorCode    ierr;
11 
12   PetscFunctionBegin;
13   ierr = MatDenseGetArrayRead(A,&array);CHKERRQ(ierr);
14   ierr = MatDenseGetLDA(A,&lda);CHKERRQ(ierr);
15   ierr = MatGetSize(A,&M,&N);CHKERRQ(ierr);
16   ierr = MatGetOwnershipRange(A,&rstart,&rend);CHKERRQ(ierr);
17   for (i=rstart; i<rend; i++) {
18     for (j=0; j<N; j++) {
19       PetscInt ii = i - rstart, jj = j;
20       PetscReal v = (PetscReal)(one ? 1 : (1 + i + j*M));
21       PetscReal w = PetscRealPart(array[ii + jj*lda]);
22       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);
23     }
24   }
25   ierr = MatDenseRestoreArrayRead(A,&array);CHKERRQ(ierr);
26   PetscFunctionReturn(0);
27 }
28 
29 #define CheckValuesIJ(A)  CheckValues(A,PETSC_FALSE)
30 #define CheckValuesOne(A) CheckValues(A,PETSC_TRUE)
31 
32 int main(int argc,char **args)
33 {
34   Mat            A;
35   PetscInt       i,j,M = 4,N = 3,rstart,rend;
36   PetscErrorCode ierr;
37   PetscScalar    *array;
38   char           mattype[256];
39   PetscViewer    view;
40 
41   ierr = PetscInitialize(&argc,&args,NULL,help);if (ierr) return ierr;
42   ierr = PetscStrcpy(mattype,MATMPIDENSE);CHKERRQ(ierr);
43   ierr = PetscOptionsGetString(NULL,NULL,"-mat_type",mattype,sizeof(mattype),NULL);CHKERRQ(ierr);
44   /*
45       Create a parallel dense matrix shared by all processors
46   */
47   ierr = MatCreateDense(PETSC_COMM_WORLD,PETSC_DECIDE,PETSC_DECIDE,M,N,NULL,&A);CHKERRQ(ierr);
48   ierr = MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr);
49   ierr = MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr);
50   ierr = MatConvert(A,mattype,MAT_INPLACE_MATRIX,&A);CHKERRQ(ierr);
51   /*
52      Set values into the matrix
53   */
54   for (i=0; i<M; i++) {
55     for (j=0; j<N; j++) {
56       PetscScalar v = (PetscReal)(1 + i + j*M);
57       ierr = MatSetValues(A,1,&i,1,&j,&v,INSERT_VALUES);CHKERRQ(ierr);
58     }
59   }
60   ierr = MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr);
61   ierr = MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr);
62   ierr = MatScale(A,2.0);CHKERRQ(ierr);
63   ierr = MatScale(A,1.0/2.0);CHKERRQ(ierr);
64 
65   /*
66       Store the binary matrix to a file
67   */
68   ierr = PetscViewerBinaryOpen(PETSC_COMM_WORLD, "matrix.dat", FILE_MODE_WRITE, &view);CHKERRQ(ierr);
69   for (i=0; i<2; i++) {
70     ierr = MatView(A,view);CHKERRQ(ierr);
71     ierr = PetscViewerPushFormat(view,PETSC_VIEWER_NATIVE);CHKERRQ(ierr);
72     ierr = MatView(A,view);CHKERRQ(ierr);
73     ierr = PetscViewerPopFormat(view);CHKERRQ(ierr);
74   }
75   ierr = PetscViewerDestroy(&view);CHKERRQ(ierr);
76   ierr = MatDestroy(&A);CHKERRQ(ierr);
77 
78   /*
79       Now reload the matrix and check its values
80   */
81   ierr = PetscViewerBinaryOpen(PETSC_COMM_WORLD,"matrix.dat",FILE_MODE_READ,&view);CHKERRQ(ierr);
82   ierr = MatCreate(PETSC_COMM_WORLD,&A);CHKERRQ(ierr);
83   ierr = MatSetType(A,mattype);CHKERRQ(ierr);
84   for (i=0; i<4; i++) {
85     if (i > 0) {ierr = MatZeroEntries(A);CHKERRQ(ierr);}
86     ierr = MatLoad(A,view);CHKERRQ(ierr);
87     ierr = CheckValuesIJ(A);CHKERRQ(ierr);
88   }
89   ierr = PetscViewerDestroy(&view);CHKERRQ(ierr);
90 
91   ierr = MatGetOwnershipRange(A,&rstart,&rend);CHKERRQ(ierr);
92   ierr = PetscMalloc1((rend-rstart)*N,&array);CHKERRQ(ierr);
93   for (i=0; i<(rend-rstart)*N; i++) array[i] = (PetscReal)1;
94   ierr = MatDensePlaceArray(A,array);CHKERRQ(ierr);
95   ierr = MatScale(A,2.0);CHKERRQ(ierr);
96   ierr = MatScale(A,1.0/2.0);CHKERRQ(ierr);
97   ierr = CheckValuesOne(A);CHKERRQ(ierr);
98   ierr = PetscViewerBinaryOpen(PETSC_COMM_WORLD,"matrix.dat",FILE_MODE_WRITE,&view);CHKERRQ(ierr);
99   ierr = MatView(A,view);CHKERRQ(ierr);
100   ierr = MatDenseResetArray(A);CHKERRQ(ierr);
101   ierr = PetscFree(array);CHKERRQ(ierr);
102   ierr = CheckValuesIJ(A);CHKERRQ(ierr);
103   ierr = PetscViewerBinarySetSkipHeader(view,PETSC_TRUE);CHKERRQ(ierr);
104   ierr = MatView(A,view);CHKERRQ(ierr);
105   ierr = PetscViewerBinarySetSkipHeader(view,PETSC_FALSE);CHKERRQ(ierr);
106   ierr = PetscViewerDestroy(&view);CHKERRQ(ierr);
107   ierr = MatDestroy(&A);CHKERRQ(ierr);
108 
109   ierr = MatCreate(PETSC_COMM_WORLD,&A);CHKERRQ(ierr);
110   ierr = MatSetType(A,mattype);CHKERRQ(ierr);
111   ierr = PetscViewerBinaryOpen(PETSC_COMM_WORLD,"matrix.dat",FILE_MODE_READ,&view);CHKERRQ(ierr);
112   ierr = MatLoad(A,view);CHKERRQ(ierr);
113   ierr = CheckValuesOne(A);CHKERRQ(ierr);
114   ierr = MatZeroEntries(A);CHKERRQ(ierr);
115   ierr = PetscViewerBinarySetSkipHeader(view,PETSC_TRUE);CHKERRQ(ierr);
116   ierr = MatLoad(A,view);CHKERRQ(ierr);
117   ierr = PetscViewerBinarySetSkipHeader(view,PETSC_FALSE);CHKERRQ(ierr);
118   ierr = CheckValuesIJ(A);CHKERRQ(ierr);
119   ierr = PetscViewerDestroy(&view);CHKERRQ(ierr);
120   ierr = MatDestroy(&A);CHKERRQ(ierr);
121 
122   {
123     PetscInt m = PETSC_DECIDE, n = PETSC_DECIDE;
124     ierr = PetscSplitOwnership(PETSC_COMM_WORLD,&m,&M);CHKERRQ(ierr);
125     ierr = PetscSplitOwnership(PETSC_COMM_WORLD,&n,&N);CHKERRQ(ierr);
126     /* TODO: MatCreateDense requires data!=NULL at all processes! */
127     ierr = PetscMalloc1(m*N+1,&array);CHKERRQ(ierr);
128 
129     ierr = PetscViewerBinaryOpen(PETSC_COMM_WORLD,"matrix.dat",FILE_MODE_READ,&view);CHKERRQ(ierr);
130     ierr = MatCreateDense(PETSC_COMM_WORLD,m,n,M,N,array,&A);CHKERRQ(ierr);
131     ierr = MatLoad(A,view);CHKERRQ(ierr);
132     ierr = CheckValuesOne(A);CHKERRQ(ierr);
133     ierr = PetscViewerBinarySetSkipHeader(view,PETSC_TRUE);CHKERRQ(ierr);
134     ierr = MatLoad(A,view);CHKERRQ(ierr);
135     ierr = PetscViewerBinarySetSkipHeader(view,PETSC_FALSE);CHKERRQ(ierr);
136     ierr = CheckValuesIJ(A);CHKERRQ(ierr);
137     ierr = MatDestroy(&A);CHKERRQ(ierr);
138     ierr = PetscViewerDestroy(&view);CHKERRQ(ierr);
139 
140     ierr = MatCreateDense(PETSC_COMM_WORLD,m,n,M,N,array,&A);CHKERRQ(ierr);
141     ierr = CheckValuesIJ(A);CHKERRQ(ierr);
142     ierr = MatDestroy(&A);CHKERRQ(ierr);
143 
144     ierr = PetscFree(array);CHKERRQ(ierr);
145   }
146 
147   ierr = PetscFinalize();
148   return ierr;
149 }
150 
151 
152 /*TEST
153 
154    testset:
155       args: -viewer_binary_mpiio 0
156       output_file: output/ex16.out
157       test:
158         suffix: stdio_1
159         nsize: 1
160         args: -mat_type seqdense
161       test:
162         suffix: stdio_2
163         nsize: 2
164       test:
165         suffix: stdio_3
166         nsize: 3
167       test:
168         suffix: stdio_4
169         nsize: 4
170       test:
171         suffix: stdio_5
172         nsize: 5
173       test:
174         requires: cuda
175         args: -mat_type seqdensecuda
176         suffix: stdio_cuda_1
177         nsize: 1
178       test:
179         requires: cuda
180         args: -mat_type mpidensecuda
181         suffix: stdio_cuda_2
182         nsize: 2
183       test:
184         requires: cuda
185         args: -mat_type mpidensecuda
186         suffix: stdio_cuda_3
187         nsize: 3
188       test:
189         requires: cuda
190         args: -mat_type mpidensecuda
191         suffix: stdio_cuda_4
192         nsize: 4
193       test:
194         requires: cuda
195         args: -mat_type mpidensecuda
196         suffix: stdio_cuda_5
197         nsize: 5
198 
199    testset:
200       requires: mpiio
201       args: -viewer_binary_mpiio 1
202       output_file: output/ex16.out
203       test:
204         suffix: mpiio_1
205         nsize: 1
206       test:
207         suffix: mpiio_2
208         nsize: 2
209       test:
210         suffix: mpiio_3
211         nsize: 3
212       test:
213         suffix: mpiio_4
214         nsize: 4
215       test:
216         suffix: mpiio_5
217         nsize: 5
218       test:
219         requires: cuda
220         args: -mat_type mpidensecuda
221         suffix: mpiio_cuda_1
222         nsize: 1
223       test:
224         requires: cuda
225         args: -mat_type mpidensecuda
226         suffix: mpiio_cuda_2
227         nsize: 2
228       test:
229         requires: cuda
230         args: -mat_type mpidensecuda
231         suffix: mpiio_cuda_3
232         nsize: 3
233       test:
234         requires: cuda
235         args: -mat_type mpidensecuda
236         suffix: mpiio_cuda_4
237         nsize: 4
238       test:
239         requires: cuda
240         args: -mat_type mpidensecuda
241         suffix: mpiio_cuda_5
242         nsize: 5
243 
244 TEST*/
245