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