1 2 #include <petsc/private/tsimpl.h> /*I "petscts.h" I*/ 3 4 typedef struct { 5 PetscViewer viewer; 6 } TSTrajectory_Singlefile; 7 8 #undef __FUNCT__ 9 #define __FUNCT__ "TSTrajectorySet_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; 21 22 PetscFunctionBeginUser; 23 if (stepnum == 0) { 24 ierr = PetscViewerCreate(PETSC_COMM_WORLD, &sf->viewer);CHKERRQ(ierr); 25 ierr = PetscViewerSetType(sf->viewer, PETSCVIEWERBINARY);CHKERRQ(ierr); 26 ierr = PetscViewerFileSetMode(sf->viewer,FILE_MODE_WRITE);CHKERRQ(ierr); 27 ierr = PetscObjectGetName((PetscObject)jac,&filename);CHKERRQ(ierr); 28 ierr = PetscViewerFileSetName(sf->viewer, filename);CHKERRQ(ierr); 29 } 30 ierr = TSGetTotalSteps(ts,&stepnum);CHKERRQ(ierr); 31 32 ierr = VecView(X,sf->viewer);CHKERRQ(ierr); 33 34 ierr = PetscViewerBinaryWrite(sf->viewer,&time,1,PETSC_REAL,PETSC_FALSE);CHKERRQ(ierr); 35 36 ierr = TSGetStages(ts,&ns,&Y);CHKERRQ(ierr); 37 for (i=0;i<ns;i++) { 38 ierr = VecView(Y[i],sf->viewer);CHKERRQ(ierr); 39 } 40 41 /* tprev is only needed for the adjoint run */ 42 /* 43 ierr = TSGetPrevTime(ts,&tprev);CHKERRQ(ierr); 44 ierr = PetscViewerBinaryWrite(sf->viewer,&tprev,1,PETSC_REAL,PETSC_FALSE);CHKERRQ(ierr); 45 */ 46 PetscFunctionReturn(0); 47 } 48 49 #undef __FUNCT__ 50 #define __FUNCT__ "TSTrajectoryDestroy_Singlefile" 51 PetscErrorCode TSTrajectoryDestroy_Singlefile(TSTrajectory jac) 52 { 53 PetscErrorCode ierr; 54 TSTrajectory_Singlefile *sf = (TSTrajectory_Singlefile*)jac->data; 55 56 PetscFunctionBegin; 57 ierr = PetscViewerDestroy(&sf->viewer);CHKERRQ(ierr); 58 ierr = PetscFree(sf);CHKERRQ(ierr); 59 PetscFunctionReturn(0); 60 } 61 62 /*MC 63 TSTRAJECTORYSINGLEFILE - Stores all solutions of the ODE/ADE into a single file 64 65 Level: intermediate 66 67 .seealso: TSTrajectoryCreate(), TS, TSTrajectorySetType() 68 69 M*/ 70 #undef __FUNCT__ 71 #define __FUNCT__ "TSTrajectoryCreate_Singlefile" 72 PETSC_EXTERN PetscErrorCode TSTrajectoryCreate_Singlefile(TSTrajectory ts) 73 { 74 PetscErrorCode ierr; 75 TSTrajectory_Singlefile *sf; 76 77 PetscFunctionBegin; 78 ierr = PetscNew(&sf);CHKERRQ(ierr); 79 ts->data = sf; 80 ts->ops->set = TSTrajectorySet_Singlefile; 81 ts->ops->get = NULL; 82 ts->ops->destroy = TSTrajectoryDestroy_Singlefile; 83 PetscFunctionReturn(0); 84 } 85