Actual source code: singlefile.c
2: #include <petsc/private/tsimpl.h>
4: typedef struct {
5: PetscViewer viewer;
6: } TSTrajectory_Singlefile;
8: static PetscErrorCode TSTrajectorySet_Singlefile(TSTrajectory tj, TS ts, PetscInt stepnum, PetscReal time, Vec X)
9: {
10: TSTrajectory_Singlefile *sf = (TSTrajectory_Singlefile *)tj->data;
11: const char *filename;
13: PetscFunctionBegin;
14: if (stepnum == 0) {
15: PetscCall(PetscViewerCreate(PetscObjectComm((PetscObject)X), &sf->viewer));
16: PetscCall(PetscViewerSetType(sf->viewer, PETSCVIEWERBINARY));
17: PetscCall(PetscViewerFileSetMode(sf->viewer, FILE_MODE_WRITE));
18: PetscCall(PetscObjectGetName((PetscObject)tj, &filename));
19: PetscCall(PetscViewerFileSetName(sf->viewer, filename));
20: }
21: PetscCall(VecView(X, sf->viewer));
22: PetscCall(PetscViewerBinaryWrite(sf->viewer, &time, 1, PETSC_REAL));
23: PetscFunctionReturn(PETSC_SUCCESS);
24: }
26: static PetscErrorCode TSTrajectoryDestroy_Singlefile(TSTrajectory tj)
27: {
28: TSTrajectory_Singlefile *sf = (TSTrajectory_Singlefile *)tj->data;
30: PetscFunctionBegin;
31: PetscCall(PetscViewerDestroy(&sf->viewer));
32: PetscCall(PetscFree(sf));
33: PetscFunctionReturn(PETSC_SUCCESS);
34: }
36: /*MC
37: TSTRAJECTORYSINGLEFILE - Stores all solutions of the ODE/ADE into a single file followed by each timestep.
38: Does not save the intermediate stages in a multistage method
40: Level: intermediate
42: .seealso: [](ch_ts), `TSTrajectoryCreate()`, `TS`, `TSTrajectorySetType()`, `TSTrajectoryType`, `TSTrajectory`
43: M*/
44: PETSC_EXTERN PetscErrorCode TSTrajectoryCreate_Singlefile(TSTrajectory tj, TS ts)
45: {
46: TSTrajectory_Singlefile *sf;
48: PetscFunctionBegin;
49: PetscCall(PetscNew(&sf));
50: tj->data = sf;
51: tj->ops->set = TSTrajectorySet_Singlefile;
52: tj->ops->get = NULL;
53: tj->ops->destroy = TSTrajectoryDestroy_Singlefile;
54: ts->setupcalled = PETSC_TRUE;
55: PetscFunctionReturn(PETSC_SUCCESS);
56: }