Actual source code: singlefile.c
petsc-3.14.6 2021-03-30
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: PetscErrorCode ierr;
12: const char *filename;
15: if (stepnum == 0) {
16: PetscViewerCreate(PetscObjectComm((PetscObject)X),&sf->viewer);
17: PetscViewerSetType(sf->viewer,PETSCVIEWERBINARY);
18: PetscViewerFileSetMode(sf->viewer,FILE_MODE_WRITE);
19: PetscObjectGetName((PetscObject)tj,&filename);
20: PetscViewerFileSetName(sf->viewer,filename);
21: }
22: VecView(X,sf->viewer);
23: PetscViewerBinaryWrite(sf->viewer,&time,1,PETSC_REAL);
24: return(0);
25: }
27: static PetscErrorCode TSTrajectoryDestroy_Singlefile(TSTrajectory tj)
28: {
29: PetscErrorCode ierr;
30: TSTrajectory_Singlefile *sf = (TSTrajectory_Singlefile*)tj->data;
33: PetscViewerDestroy(&sf->viewer);
34: PetscFree(sf);
35: return(0);
36: }
38: /*MC
39: TSTRAJECTORYSINGLEFILE - Stores all solutions of the ODE/ADE into a single file followed by each timestep. Does not save the intermediate stages in a multistage method
41: Level: intermediate
43: .seealso: TSTrajectoryCreate(), TS, TSTrajectorySetType()
45: M*/
46: PETSC_EXTERN PetscErrorCode TSTrajectoryCreate_Singlefile(TSTrajectory tj,TS ts)
47: {
48: PetscErrorCode ierr;
49: TSTrajectory_Singlefile *sf;
52: PetscNew(&sf);
53: tj->data = sf;
54: tj->ops->set = TSTrajectorySet_Singlefile;
55: tj->ops->get = NULL;
56: tj->ops->destroy = TSTrajectoryDestroy_Singlefile;
57: ts->setupcalled = PETSC_TRUE;
58: return(0);
59: }