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