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: PetscMPIInt rank;
114: PetscBool dataverified = PETSC_TRUE;
116: PetscFunctionBeginUser;
117: PetscCallMPI(MPI_Comm_rank(PETSC_COMM_WORLD, &rank));
118: if (rank == 0) {
119: PetscCall(PetscBinaryOpen(name, FILE_MODE_READ, &fdes));
120: PetscCall(PetscBinaryRead(fdes, buffer, VEC_LEN, NULL, PETSC_SCALAR));
121: PetscCall(PetscBinaryClose(fdes));
123: for (PetscInt i = 0; i < VEC_LEN; i++) {
124: PetscScalar v;
125: v = PetscAbsScalar(test_values[i] - buffer[i]);
126: if (PetscDefined(USE_COMPLEX) && ((PetscRealPart(v) > 1.0e-10) || (PetscImaginaryPart(v) > 1.0e-10))) {
127: 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));
128: dataverified = PETSC_FALSE;
129: } else if (PetscRealPart(v) > 1.0e-10) {
130: PetscCall(PetscPrintf(PETSC_COMM_SELF, "ERROR: Difference > 1.0e-10 occurred (delta = %+1.12e [loc %" PetscInt_FMT "])\n", (double)PetscRealPart(buffer[i]), i));
131: dataverified = PETSC_FALSE;
132: }
133: }
134: if (dataverified) PetscCall(PetscPrintf(PETSC_COMM_SELF, "Headerless read of data verified\n"));
135: }
136: PetscFunctionReturn(PETSC_SUCCESS);
137: }
139: PetscErrorCode TestBinary(void)
140: {
141: Vec x, y;
142: PetscBool skipheader = PETSC_TRUE;
143: PetscBool usempiio = PETSC_FALSE;
145: PetscFunctionBeginUser;
146: PetscCall(VecCreate(PETSC_COMM_WORLD, &x));
147: PetscCall(VecSetSizes(x, PETSC_DECIDE, VEC_LEN));
148: PetscCall(VecSetFromOptions(x));
149: PetscCall(VecFill(x));
150: PetscCall(MyVecDump("xH.pbvec", skipheader, usempiio, x));
152: PetscCall(VecCreate(PETSC_COMM_WORLD, &y));
153: PetscCall(VecSetSizes(y, PETSC_DECIDE, VEC_LEN));
154: PetscCall(VecSetFromOptions(y));
156: PetscCall(MyVecLoad("xH.pbvec", skipheader, usempiio, y));
157: PetscCall(VecCompare(x, y));
159: PetscCall(VecDestroy(&y));
160: PetscCall(VecDestroy(&x));
162: PetscCall(HeaderlessBinaryRead("xH.pbvec"));
163: PetscFunctionReturn(PETSC_SUCCESS);
164: }
166: #if PetscDefined(HAVE_MPIIO)
167: PetscErrorCode TestBinaryMPIIO(void)
168: {
169: Vec x, y;
170: PetscBool skipheader = PETSC_TRUE;
171: PetscBool usempiio = PETSC_TRUE;
173: PetscFunctionBeginUser;
174: PetscCall(VecCreate(PETSC_COMM_WORLD, &x));
175: PetscCall(VecSetSizes(x, PETSC_DECIDE, VEC_LEN));
176: PetscCall(VecSetFromOptions(x));
177: PetscCall(VecFill(x));
178: PetscCall(MyVecDump("xHmpi.pbvec", skipheader, usempiio, x));
180: PetscCall(VecCreate(PETSC_COMM_WORLD, &y));
181: PetscCall(VecSetSizes(y, PETSC_DECIDE, VEC_LEN));
182: PetscCall(VecSetFromOptions(y));
184: PetscCall(MyVecLoad("xHmpi.pbvec", skipheader, usempiio, y));
185: PetscCall(VecCompare(x, y));
187: PetscCall(VecDestroy(&y));
188: PetscCall(VecDestroy(&x));
190: PetscCall(HeaderlessBinaryRead("xHmpi.pbvec"));
191: PetscFunctionReturn(PETSC_SUCCESS);
192: }
193: #endif
195: int main(int argc, char **args)
196: {
197: PetscBool usempiio = PETSC_FALSE;
199: PetscFunctionBeginUser;
200: PetscCall(PetscInitialize(&argc, &args, NULL, help));
201: PetscCall(PetscOptionsGetBool(NULL, NULL, "-usempiio", &usempiio, NULL));
202: if (!usempiio) {
203: PetscCall(TestBinary());
204: } else {
205: #if PetscDefined(HAVE_MPIIO)
206: PetscCall(TestBinaryMPIIO());
207: #else
208: PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Warning: Executing TestBinaryMPIIO() requires a working MPI-2 implementation\n"));
209: #endif
210: }
211: PetscCall(PetscFinalize());
212: return 0;
213: }
215: /*TEST
217: test:
218: output_file: output/ex46_1_p1.out
220: test:
221: suffix: 2
222: nsize: 6
223: output_file: output/ex46_1_p6.out
225: test:
226: suffix: 3
227: nsize: 12
228: output_file: output/ex46_1_p12.out
230: testset:
231: requires: mpiio
232: args: -usempiio
233: test:
234: suffix: mpiio_1
235: output_file: output/ex46_2_p1.out
236: test:
237: suffix: mpiio_2
238: nsize: 6
239: output_file: output/ex46_2_p6.out
240: test:
241: suffix: mpiio_3
242: nsize: 12
243: output_file: output/ex46_2_p12.out
245: TEST*/