Actual source code: trajbasic.c

petsc-3.7.3 2016-08-01
Report Typos and Errors
  2: #include <petsc/private/tsimpl.h>        /*I "petscts.h"  I*/

  6: static PetscErrorCode OutputBIN(const char *filename,PetscViewer *viewer)
  7: {

 11:   PetscViewerCreate(PETSC_COMM_WORLD,viewer);
 12:   PetscViewerSetType(*viewer,PETSCVIEWERBINARY);
 13:   PetscViewerFileSetMode(*viewer,FILE_MODE_WRITE);
 14:   PetscViewerFileSetName(*viewer,filename);
 15:   return(0);
 16: }

 20: static PetscErrorCode TSTrajectorySet_Basic(TSTrajectory tj,TS ts,PetscInt stepnum,PetscReal time,Vec X)
 21: {
 22:   PetscViewer    viewer;
 23:   PetscInt       ns,i;
 24:   Vec            *Y;
 25:   char           filename[PETSC_MAX_PATH_LEN];
 26:   PetscReal      tprev;

 30:   TSGetTotalSteps(ts,&stepnum);
 31:   if (stepnum == 0) {
 32:     PetscMPIInt rank;
 33:     MPI_Comm_rank(PetscObjectComm((PetscObject)ts),&rank);
 34:     if (!rank) {
 35:       PetscRMTree("SA-data");
 36:       PetscMkdir("SA-data");
 37:     }
 38:     PetscSNPrintf(filename,sizeof(filename),"SA-data/SA-%06d.bin",stepnum);
 39:     OutputBIN(filename,&viewer);
 40:     VecView(X,viewer);
 41:     PetscViewerBinaryWrite(viewer,&time,1,PETSC_REAL,PETSC_FALSE);
 42:     PetscViewerDestroy(&viewer);
 43:     return(0);
 44:   }
 45:   PetscSNPrintf(filename,sizeof(filename),"SA-data/SA-%06d.bin",stepnum);
 46:   OutputBIN(filename,&viewer);
 47:   VecView(X,viewer);
 48:   PetscViewerBinaryWrite(viewer,&time,1,PETSC_REAL,PETSC_FALSE);

 50:   TSGetStages(ts,&ns,&Y);
 51:   for (i=0;i<ns;i++) {
 52:     VecView(Y[i],viewer);
 53:   }

 55:   TSGetPrevTime(ts,&tprev);
 56:   PetscViewerBinaryWrite(viewer,&tprev,1,PETSC_REAL,PETSC_FALSE);

 58:   PetscViewerDestroy(&viewer);
 59:   return(0);
 60: }

 64: static PetscErrorCode TSTrajectoryGet_Basic(TSTrajectory tj,TS ts,PetscInt stepnum,PetscReal *t)
 65: {
 66:   Vec            Sol,*Y;
 67:   PetscInt       Nr,i;
 68:   PetscViewer    viewer;
 69:   PetscReal      timepre;
 70:   char           filename[PETSC_MAX_PATH_LEN];

 74:   TSGetTotalSteps(ts,&stepnum);
 75:   PetscSNPrintf(filename,sizeof filename,"SA-data/SA-%06d.bin",stepnum);
 76:   PetscViewerBinaryOpen(PETSC_COMM_WORLD,filename,FILE_MODE_READ,&viewer);

 78:   TSGetSolution(ts,&Sol);
 79:   VecLoad(Sol,viewer);

 81:   PetscViewerBinaryRead(viewer,t,1,NULL,PETSC_REAL);

 83:   if (stepnum != 0) {
 84:     TSGetStages(ts,&Nr,&Y);
 85:     for (i=0;i<Nr ;i++) {
 86:       VecLoad(Y[i],viewer);
 87:     }
 88:     PetscViewerBinaryRead(viewer,&timepre,1,NULL,PETSC_REAL);
 89:     TSSetTimeStep(ts,-(*t)+timepre);
 90:   }

 92:   PetscViewerDestroy(&viewer);
 93:   return(0);
 94: }

 96: /*MC
 97:       TSTRAJECTORYBASIC - Stores each solution of the ODE/ADE in a file

 99:   Level: intermediate

101: .seealso:  TSTrajectoryCreate(), TS, TSTrajectorySetType()

103: M*/
106: PETSC_EXTERN PetscErrorCode TSTrajectoryCreate_Basic(TSTrajectory tj,TS ts)
107: {
109:   tj->ops->set  = TSTrajectorySet_Basic;
110:   tj->ops->get  = TSTrajectoryGet_Basic;
111:   return(0);
112: }