1 2 #include <petsc/private/tsimpl.h> /*I "petscts.h" I*/ 3 4 typedef struct { 5 PetscViewer viewer; 6 } TSTrajectory_Singlefile; 7 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; 12 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(0); 24 } 25 26 static PetscErrorCode TSTrajectoryDestroy_Singlefile(TSTrajectory tj) 27 { 28 TSTrajectory_Singlefile *sf = (TSTrajectory_Singlefile *)tj->data; 29 30 PetscFunctionBegin; 31 PetscCall(PetscViewerDestroy(&sf->viewer)); 32 PetscCall(PetscFree(sf)); 33 PetscFunctionReturn(0); 34 } 35 36 /*MC 37 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 38 39 Level: intermediate 40 41 .seealso: `TSTrajectoryCreate()`, `TS`, `TSTrajectorySetType()` 42 43 M*/ 44 PETSC_EXTERN PetscErrorCode TSTrajectoryCreate_Singlefile(TSTrajectory tj, TS ts) 45 { 46 TSTrajectory_Singlefile *sf; 47 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(0); 56 } 57