xref: /petsc/src/vec/vec/tutorials/ex6.c (revision 732aec7a18f2199fb53bb9a2f3aef439a834ce31)
1 static char help[] = "Writes an array to a file, then reads an array from a file, then forms a vector.\n\n";
2 
3 /*
4     This uses the low level PetscBinaryWrite() and PetscBinaryRead() to access a binary file. It will not work in parallel!
5 
6     We HIGHLY recommend using instead VecView() and VecLoad() to read and write Vectors in binary format (which also work in parallel). Then you can use
7     share/petsc/matlab/PetscBinaryRead() and share/petsc/matlab/PetscBinaryWrite() to read (or write) the vector into MATLAB.
8 
9     Note this also works for matrices with MatView() and MatLoad().
10 */
11 #include <petscvec.h>
12 
main(int argc,char ** args)13 int main(int argc, char **args)
14 {
15   PetscMPIInt  size;
16   int          fd;
17   PetscInt     i, m = 10, sz;
18   PetscScalar *avec, *array;
19   Vec          vec;
20   PetscViewer  view_out, view_in;
21 
22   PetscFunctionBeginUser;
23   PetscCall(PetscInitialize(&argc, &args, NULL, help));
24   PetscCallMPI(MPI_Comm_size(PETSC_COMM_WORLD, &size));
25   PetscCheck(size == 1, PETSC_COMM_WORLD, PETSC_ERR_WRONG_MPI_SIZE, "This is a uniprocessor example only!");
26 
27   PetscCall(PetscOptionsGetInt(NULL, NULL, "-m", &m, NULL));
28 
29   /* ---------------------------------------------------------------------- */
30   /*          PART 1: Write some data to a file in binary format            */
31   /* ---------------------------------------------------------------------- */
32 
33   /* Allocate array and set values */
34   PetscCall(PetscMalloc1(m, &array));
35   for (i = 0; i < m; i++) array[i] = i * 10.0;
36 
37   /* Open viewer for binary output */
38   PetscCall(PetscViewerBinaryOpen(PETSC_COMM_SELF, "input.dat", FILE_MODE_WRITE, &view_out));
39   PetscCall(PetscViewerBinaryGetDescriptor(view_out, &fd));
40 
41   /* Write binary output */
42   PetscCall(PetscBinaryWrite(fd, &m, 1, PETSC_INT));
43   PetscCall(PetscBinaryWrite(fd, array, m, PETSC_SCALAR));
44 
45   /* Destroy the output viewer and work array */
46   PetscCall(PetscViewerDestroy(&view_out));
47   PetscCall(PetscFree(array));
48 
49   /* ---------------------------------------------------------------------- */
50   /*          PART 2: Read data from file and form a vector                 */
51   /* ---------------------------------------------------------------------- */
52 
53   /* Open input binary viewer */
54   PetscCall(PetscViewerBinaryOpen(PETSC_COMM_SELF, "input.dat", FILE_MODE_READ, &view_in));
55   PetscCall(PetscViewerBinaryGetDescriptor(view_in, &fd));
56 
57   /* Create vector and get pointer to data space */
58   PetscCall(VecCreate(PETSC_COMM_SELF, &vec));
59   PetscCall(VecSetSizes(vec, PETSC_DECIDE, m));
60   PetscCall(VecSetFromOptions(vec));
61   PetscCall(VecGetArray(vec, &avec));
62 
63   /* Read data into vector */
64   PetscCall(PetscBinaryRead(fd, &sz, 1, NULL, PETSC_INT));
65   PetscCheck(sz > 0, PETSC_COMM_SELF, PETSC_ERR_FILE_UNEXPECTED, "Error: Must have array length > 0");
66 
67   PetscCall(PetscPrintf(PETSC_COMM_SELF, "reading data in binary from input.dat, sz =%" PetscInt_FMT " ...\n", sz));
68   PetscCall(PetscBinaryRead(fd, avec, sz, NULL, PETSC_SCALAR));
69 
70   /* View vector */
71   PetscCall(VecRestoreArray(vec, &avec));
72   PetscCall(VecView(vec, PETSC_VIEWER_STDOUT_SELF));
73 
74   /* Free data structures */
75   PetscCall(VecDestroy(&vec));
76   PetscCall(PetscViewerDestroy(&view_in));
77   PetscCall(PetscFinalize());
78   return 0;
79 }
80 
81 /*TEST
82 
83      test:
84 
85 TEST*/
86