Actual source code: ex46.c
1: static char help[] = "Tests PetscViewerBinary VecView()/VecLoad() function correctly when binary header is skipped.\n\n";
3: #include <petscviewer.h>
4: #include <petscvec.h>
6: #define VEC_LEN 10
7: const PetscReal test_values[] = {0.311256, 88.068, 11.077444, 9953.62, 7.345, 64.8943, 3.1458, 6699.95, 0.00084, 0.0647};
9: PetscErrorCode MyVecDump(const char fname[], PetscBool skippheader, PetscBool usempiio, Vec x)
10: {
11: MPI_Comm comm;
12: PetscViewer viewer;
13: PetscBool ismpiio, isskip;
15: PetscFunctionBeginUser;
16: PetscCall(PetscObjectGetComm((PetscObject)x, &comm));
18: PetscCall(PetscViewerCreate(comm, &viewer));
19: PetscCall(PetscViewerSetType(viewer, PETSCVIEWERBINARY));
20: if (skippheader) PetscCall(PetscViewerBinarySetSkipHeader(viewer, PETSC_TRUE));
21: PetscCall(PetscViewerFileSetMode(viewer, FILE_MODE_WRITE));
22: if (usempiio) PetscCall(PetscViewerBinarySetUseMPIIO(viewer, PETSC_TRUE));
23: PetscCall(PetscViewerFileSetName(viewer, fname));
25: PetscCall(VecView(x, viewer));
27: PetscCall(PetscViewerBinaryGetUseMPIIO(viewer, &ismpiio));
28: if (ismpiio) PetscCall(PetscPrintf(comm, "*** PetscViewer[write] using MPI-IO ***\n"));
29: PetscCall(PetscViewerBinaryGetSkipHeader(viewer, &isskip));
30: if (isskip) PetscCall(PetscPrintf(comm, "*** PetscViewer[write] skipping header ***\n"));
32: PetscCall(PetscViewerDestroy(&viewer));
33: PetscFunctionReturn(PETSC_SUCCESS);
34: }
36: PetscErrorCode MyVecLoad(const char fname[], PetscBool skippheader, PetscBool usempiio, Vec x)
37: {
38: MPI_Comm comm;
39: PetscViewer viewer;
40: PetscBool ismpiio, isskip;
42: PetscFunctionBeginUser;
43: PetscCall(PetscObjectGetComm((PetscObject)x, &comm));
45: PetscCall(PetscViewerCreate(comm, &viewer));
46: PetscCall(PetscViewerSetType(viewer, PETSCVIEWERBINARY));
47: if (skippheader) PetscCall(PetscViewerBinarySetSkipHeader(viewer, PETSC_TRUE));
48: PetscCall(PetscViewerFileSetMode(viewer, FILE_MODE_READ));
49: if (usempiio) PetscCall(PetscViewerBinarySetUseMPIIO(viewer, PETSC_TRUE));
50: PetscCall(PetscViewerFileSetName(viewer, fname));
52: PetscCall(VecLoad(x, viewer));
54: PetscCall(PetscViewerBinaryGetSkipHeader(viewer, &isskip));
55: if (isskip) PetscCall(PetscPrintf(comm, "*** PetscViewer[load] skipping header ***\n"));
56: PetscCall(PetscViewerBinaryGetUseMPIIO(viewer, &ismpiio));
57: if (ismpiio) PetscCall(PetscPrintf(comm, "*** PetscViewer[load] using MPI-IO ***\n"));
59: PetscCall(PetscViewerDestroy(&viewer));
60: PetscFunctionReturn(PETSC_SUCCESS);
61: }
63: PetscErrorCode VecFill(Vec x)
64: {
65: PetscInt i, s, e;
67: PetscFunctionBeginUser;
68: PetscCall(VecGetOwnershipRange(x, &s, &e));
69: for (i = s; i < e; i++) PetscCall(VecSetValue(x, i, (PetscScalar)test_values[i], INSERT_VALUES));
70: PetscCall(VecAssemblyBegin(x));
71: PetscCall(VecAssemblyEnd(x));
72: PetscFunctionReturn(PETSC_SUCCESS);
73: }
75: PetscErrorCode VecCompare(Vec a, Vec b)
76: {
77: PetscInt locmin[2], locmax[2];
78: PetscReal min[2], max[2];
79: Vec ref;
81: PetscFunctionBeginUser;
82: PetscCall(VecMin(a, &locmin[0], &min[0]));
83: PetscCall(VecMax(a, &locmax[0], &max[0]));
85: PetscCall(VecMin(b, &locmin[1], &min[1]));
86: PetscCall(VecMax(b, &locmax[1], &max[1]));
88: PetscCall(PetscPrintf(PETSC_COMM_WORLD, "VecCompare\n"));
89: PetscCall(PetscPrintf(PETSC_COMM_WORLD, " min(a) = %+1.2e [loc %" PetscInt_FMT "]\n", (double)min[0], locmin[0]));
90: PetscCall(PetscPrintf(PETSC_COMM_WORLD, " max(a) = %+1.2e [loc %" PetscInt_FMT "]\n", (double)max[0], locmax[0]));
92: PetscCall(PetscPrintf(PETSC_COMM_WORLD, " min(b) = %+1.2e [loc %" PetscInt_FMT "]\n", (double)min[1], locmin[1]));
93: PetscCall(PetscPrintf(PETSC_COMM_WORLD, " max(b) = %+1.2e [loc %" PetscInt_FMT "]\n", (double)max[1], locmax[1]));
95: PetscCall(VecDuplicate(a, &ref));
96: PetscCall(VecCopy(a, ref));
97: PetscCall(VecAXPY(ref, -1.0, b));
98: PetscCall(VecMin(ref, &locmin[0], &min[0]));
99: if (PetscAbsReal(min[0]) > 1.0e-10) {
100: PetscCall(PetscPrintf(PETSC_COMM_WORLD, " ERROR: min(a-b) > 1.0e-10\n"));
101: PetscCall(PetscPrintf(PETSC_COMM_WORLD, " min(a-b) = %+1.10e\n", (double)PetscAbsReal(min[0])));
102: } else {
103: PetscCall(PetscPrintf(PETSC_COMM_WORLD, " min(a-b) < 1.0e-10\n"));
104: }
105: PetscCall(VecDestroy(&ref));
106: PetscFunctionReturn(PETSC_SUCCESS);
107: }
109: PetscErrorCode HeaderlessBinaryRead(const char name[])
110: {
111: int fdes;
112: PetscScalar buffer[VEC_LEN];
113: PetscInt i;
114: PetscMPIInt rank;
115: PetscBool dataverified = PETSC_TRUE;
117: PetscFunctionBeginUser;
118: PetscCallMPI(MPI_Comm_rank(PETSC_COMM_WORLD, &rank));
119: if (rank == 0) {
120: PetscCall(PetscBinaryOpen(name, FILE_MODE_READ, &fdes));
121: PetscCall(PetscBinaryRead(fdes, buffer, VEC_LEN, NULL, PETSC_SCALAR));
122: PetscCall(PetscBinaryClose(fdes));
124: for (i = 0; i < VEC_LEN; i++) {
125: PetscScalar v;
126: v = PetscAbsScalar(test_values[i] - buffer[i]);
127: #if defined(PETSC_USE_COMPLEX)
128: if ((PetscRealPart(v) > 1.0e-10) || (PetscImaginaryPart(v) > 1.0e-10)) {
129: PetscCall(PetscPrintf(PETSC_COMM_SELF, "ERROR: Difference > 1.0e-10 occurred (delta = (%+1.12e,%+1.12e) [loc %" PetscInt_FMT "])\n", (double)PetscRealPart(buffer[i]), (double)PetscImaginaryPart(buffer[i]), i));
130: dataverified = PETSC_FALSE;
131: }
132: #else
133: if (PetscRealPart(v) > 1.0e-10) {
134: PetscCall(PetscPrintf(PETSC_COMM_SELF, "ERROR: Difference > 1.0e-10 occurred (delta = %+1.12e [loc %" PetscInt_FMT "])\n", (double)PetscRealPart(buffer[i]), i));
135: dataverified = PETSC_FALSE;
136: }
137: #endif
138: }
139: if (dataverified) PetscCall(PetscPrintf(PETSC_COMM_SELF, "Headerless read of data verified\n"));
140: }
141: PetscFunctionReturn(PETSC_SUCCESS);
142: }
144: PetscErrorCode TestBinary(void)
145: {
146: Vec x, y;
147: PetscBool skipheader = PETSC_TRUE;
148: PetscBool usempiio = PETSC_FALSE;
150: PetscFunctionBeginUser;
151: PetscCall(VecCreate(PETSC_COMM_WORLD, &x));
152: PetscCall(VecSetSizes(x, PETSC_DECIDE, VEC_LEN));
153: PetscCall(VecSetFromOptions(x));
154: PetscCall(VecFill(x));
155: PetscCall(MyVecDump("xH.pbvec", skipheader, usempiio, x));
157: PetscCall(VecCreate(PETSC_COMM_WORLD, &y));
158: PetscCall(VecSetSizes(y, PETSC_DECIDE, VEC_LEN));
159: PetscCall(VecSetFromOptions(y));
161: PetscCall(MyVecLoad("xH.pbvec", skipheader, usempiio, y));
162: PetscCall(VecCompare(x, y));
164: PetscCall(VecDestroy(&y));
165: PetscCall(VecDestroy(&x));
167: PetscCall(HeaderlessBinaryRead("xH.pbvec"));
168: PetscFunctionReturn(PETSC_SUCCESS);
169: }
171: #if defined(PETSC_HAVE_MPIIO)
172: PetscErrorCode TestBinaryMPIIO(void)
173: {
174: Vec x, y;
175: PetscBool skipheader = PETSC_TRUE;
176: PetscBool usempiio = PETSC_TRUE;
178: PetscFunctionBeginUser;
179: PetscCall(VecCreate(PETSC_COMM_WORLD, &x));
180: PetscCall(VecSetSizes(x, PETSC_DECIDE, VEC_LEN));
181: PetscCall(VecSetFromOptions(x));
182: PetscCall(VecFill(x));
183: PetscCall(MyVecDump("xHmpi.pbvec", skipheader, usempiio, x));
185: PetscCall(VecCreate(PETSC_COMM_WORLD, &y));
186: PetscCall(VecSetSizes(y, PETSC_DECIDE, VEC_LEN));
187: PetscCall(VecSetFromOptions(y));
189: PetscCall(MyVecLoad("xHmpi.pbvec", skipheader, usempiio, y));
190: PetscCall(VecCompare(x, y));
192: PetscCall(VecDestroy(&y));
193: PetscCall(VecDestroy(&x));
195: PetscCall(HeaderlessBinaryRead("xHmpi.pbvec"));
196: PetscFunctionReturn(PETSC_SUCCESS);
197: }
198: #endif
200: int main(int argc, char **args)
201: {
202: PetscBool usempiio = PETSC_FALSE;
204: PetscFunctionBeginUser;
205: PetscCall(PetscInitialize(&argc, &args, (char *)0, help));
206: PetscCall(PetscOptionsGetBool(NULL, NULL, "-usempiio", &usempiio, NULL));
207: if (!usempiio) {
208: PetscCall(TestBinary());
209: } else {
210: #if defined(PETSC_HAVE_MPIIO)
211: PetscCall(TestBinaryMPIIO());
212: #else
213: PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Warning: Executing TestBinaryMPIIO() requires a working MPI-2 implementation\n"));
214: #endif
215: }
216: PetscCall(PetscFinalize());
217: return 0;
218: }
220: /*TEST
222: test:
223: output_file: output/ex46_1_p1.out
225: test:
226: suffix: 2
227: nsize: 6
228: output_file: output/ex46_1_p6.out
230: test:
231: suffix: 3
232: nsize: 12
233: output_file: output/ex46_1_p12.out
235: testset:
236: requires: mpiio
237: args: -usempiio
238: test:
239: suffix: mpiio_1
240: output_file: output/ex46_2_p1.out
241: test:
242: suffix: mpiio_2
243: nsize: 6
244: output_file: output/ex46_2_p6.out
245: test:
246: suffix: mpiio_3
247: nsize: 12
248: output_file: output/ex46_2_p12.out
250: TEST*/