Actual source code: singlefile.c
petsc-3.6.1 2015-08-06
2: #include <petsc/private/tsimpl.h> /*I "petscts.h" I*/
4: typedef struct {
5: PetscViewer viewer;
6: } TSTrajectory_Singlefile;
10: PetscErrorCode TSTrajectorySet_Singlefile(TSTrajectory jac,TS ts,PetscInt stepnum,PetscReal time,Vec X)
11: {
12: TSTrajectory_Singlefile *sf = (TSTrajectory_Singlefile*)jac->data;
13: PetscInt ns,i;
14: Vec *Y;
15: /* tprev is only needed for the adjoint run */
16: /*
17: PetscReal tprev;
18: */
19: PetscErrorCode ierr;
20: const char *filename;
23: if (stepnum == 0) {
24: PetscViewerCreate(PETSC_COMM_WORLD, &sf->viewer);
25: PetscViewerSetType(sf->viewer, PETSCVIEWERBINARY);
26: PetscViewerFileSetMode(sf->viewer,FILE_MODE_WRITE);
27: PetscObjectGetName((PetscObject)jac,&filename);
28: PetscViewerFileSetName(sf->viewer, filename);
29: }
30: TSGetTotalSteps(ts,&stepnum);
32: VecView(X,sf->viewer);
34: PetscViewerBinaryWrite(sf->viewer,&time,1,PETSC_REAL,PETSC_FALSE);
36: TSGetStages(ts,&ns,&Y);
37: for (i=0;i<ns;i++) {
38: VecView(Y[i],sf->viewer);
39: }
41: /* tprev is only needed for the adjoint run */
42: /*
43: TSGetPrevTime(ts,&tprev);
44: PetscViewerBinaryWrite(sf->viewer,&tprev,1,PETSC_REAL,PETSC_FALSE);
45: */
46: return(0);
47: }
51: PetscErrorCode TSTrajectoryDestroy_Singlefile(TSTrajectory jac)
52: {
53: PetscErrorCode ierr;
54: TSTrajectory_Singlefile *sf = (TSTrajectory_Singlefile*)jac->data;
57: PetscViewerDestroy(&sf->viewer);
58: PetscFree(sf);
59: return(0);
60: }
62: /*MC
63: TSTRAJECTORYSINGLEFILE - Stores all solutions of the ODE/ADE into a single file
65: Level: intermediate
67: .seealso: TSTrajectoryCreate(), TS, TSTrajectorySetType()
69: M*/
72: PETSC_EXTERN PetscErrorCode TSTrajectoryCreate_Singlefile(TSTrajectory ts)
73: {
74: PetscErrorCode ierr;
75: TSTrajectory_Singlefile *sf;
78: PetscNew(&sf);
79: ts->data = sf;
80: ts->ops->set = TSTrajectorySet_Singlefile;
81: ts->ops->get = NULL;
82: ts->ops->destroy = TSTrajectoryDestroy_Singlefile;
83: return(0);
84: }