1 2 #include <petsc/private/tsimpl.h> /*I "petscts.h" I*/ 3 4 typedef struct { 5 PetscViewer viewer; 6 } TSTrajectory_Basic; 7 8 static PetscErrorCode TSTrajectorySet_Basic(TSTrajectory tj,TS ts,PetscInt stepnum,PetscReal time,Vec X) 9 { 10 TSTrajectory_Basic *tjbasic = (TSTrajectory_Basic*)tj->data; 11 char filename[PETSC_MAX_PATH_LEN]; 12 PetscErrorCode ierr; 13 14 PetscFunctionBegin; 15 ierr = PetscSNPrintf(filename,sizeof(filename),tj->dirfiletemplate,stepnum);CHKERRQ(ierr); 16 ierr = PetscViewerFileSetName(tjbasic->viewer,filename);CHKERRQ(ierr); /* this triggers PetscViewer to be set up again */ 17 ierr = PetscViewerSetUp(tjbasic->viewer);CHKERRQ(ierr); 18 ierr = VecView(X,tjbasic->viewer);CHKERRQ(ierr); 19 ierr = PetscViewerBinaryWrite(tjbasic->viewer,&time,1,PETSC_REAL,PETSC_FALSE);CHKERRQ(ierr); 20 if (stepnum && !tj->solution_only) { 21 Vec *Y; 22 PetscReal tprev; 23 PetscInt ns,i; 24 25 ierr = TSGetStages(ts,&ns,&Y);CHKERRQ(ierr); 26 for (i=0;i<ns;i++) { 27 ierr = VecView(Y[i],tjbasic->viewer);CHKERRQ(ierr); 28 } 29 ierr = TSGetPrevTime(ts,&tprev);CHKERRQ(ierr); 30 ierr = PetscViewerBinaryWrite(tjbasic->viewer,&tprev,1,PETSC_REAL,PETSC_FALSE);CHKERRQ(ierr); 31 } 32 PetscFunctionReturn(0); 33 } 34 35 static PetscErrorCode TSTrajectorySetFromOptions_Basic(PetscOptionItems *PetscOptionsObject,TSTrajectory tj) 36 { 37 PetscErrorCode ierr; 38 39 PetscFunctionBegin; 40 ierr = PetscOptionsHead(PetscOptionsObject,"TS trajectory options for Basic type");CHKERRQ(ierr); 41 ierr = PetscOptionsTail();CHKERRQ(ierr); 42 PetscFunctionReturn(0); 43 } 44 45 static PetscErrorCode TSTrajectoryGet_Basic(TSTrajectory tj,TS ts,PetscInt stepnum,PetscReal *t) 46 { 47 PetscViewer viewer; 48 char filename[PETSC_MAX_PATH_LEN]; 49 PetscErrorCode ierr; 50 Vec Sol; 51 52 PetscFunctionBegin; 53 ierr = PetscSNPrintf(filename,sizeof(filename),tj->dirfiletemplate,stepnum);CHKERRQ(ierr); 54 ierr = PetscViewerBinaryOpen(PETSC_COMM_WORLD,filename,FILE_MODE_READ,&viewer);CHKERRQ(ierr); 55 ierr = TSGetSolution(ts,&Sol);CHKERRQ(ierr); 56 ierr = PetscViewerPushFormat(viewer,PETSC_VIEWER_NATIVE);CHKERRQ(ierr); 57 ierr = VecLoad(Sol,viewer);CHKERRQ(ierr); 58 ierr = PetscViewerBinaryRead(viewer,t,1,NULL,PETSC_REAL);CHKERRQ(ierr); 59 if (stepnum && !tj->solution_only) { 60 Vec *Y; 61 PetscInt Nr,i; 62 PetscReal timepre; 63 64 ierr = TSGetStages(ts,&Nr,&Y);CHKERRQ(ierr); 65 for (i=0; i<Nr; i++) { 66 ierr = VecLoad(Y[i],viewer);CHKERRQ(ierr); 67 } 68 ierr = PetscViewerBinaryRead(viewer,&timepre,1,NULL,PETSC_REAL);CHKERRQ(ierr); 69 if (tj->adjoint_solve_mode) { 70 ierr = TSSetTimeStep(ts,-(*t)+timepre);CHKERRQ(ierr); 71 } 72 } 73 ierr = PetscViewerDestroy(&viewer);CHKERRQ(ierr); 74 PetscFunctionReturn(0); 75 } 76 77 static PetscErrorCode TSTrajectorySetUp_Basic(TSTrajectory tj,TS ts) 78 { 79 MPI_Comm comm; 80 PetscMPIInt rank; 81 PetscErrorCode ierr; 82 83 PetscFunctionBegin; 84 ierr = PetscObjectGetComm((PetscObject)tj,&comm);CHKERRQ(ierr); 85 ierr = MPI_Comm_rank(comm,&rank);CHKERRQ(ierr); 86 if (!rank) { 87 const char* dir = tj->dirname; 88 PetscBool flg; 89 90 /* I don't like running PetscRMTree on a directory */ 91 ierr = PetscTestDirectory(dir,'w',&flg);CHKERRQ(ierr); 92 if (!flg) { 93 ierr = PetscTestFile(dir,'r',&flg);CHKERRQ(ierr); 94 if (flg) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_USER,"Specified path is a file - not a dir: %s",dir); 95 ierr = PetscMkdir(dir);CHKERRQ(ierr); 96 } else { 97 ierr = PetscRMTree(tj->dirname);CHKERRQ(ierr); /* avoid having to delete the folder manually */ 98 } 99 } 100 ierr = PetscBarrier((PetscObject)tj);CHKERRQ(ierr); 101 PetscFunctionReturn(0); 102 } 103 104 static PetscErrorCode TSTrajectoryDestroy_Basic(TSTrajectory tj) 105 { 106 TSTrajectory_Basic *tjbasic = (TSTrajectory_Basic*)tj->data; 107 PetscErrorCode ierr; 108 109 PetscFunctionBegin; 110 ierr = PetscViewerDestroy(&tjbasic->viewer);CHKERRQ(ierr); 111 ierr = PetscFree(tjbasic);CHKERRQ(ierr); 112 PetscFunctionReturn(0); 113 } 114 115 /*MC 116 TSTRAJECTORYBASIC - Stores each solution of the ODE/DAE in a file 117 118 Saves each timestep into a seperate file named SA-data/SA-%06d.bin. The file name can be changed. 119 120 This version saves the solutions at all the stages 121 122 $PETSC_DIR/share/petsc/matlab/PetscReadBinaryTrajectory.m can read in files created with this format 123 124 Level: intermediate 125 126 .seealso: TSTrajectoryCreate(), TS, TSTrajectorySetType(), TSTrajectorySetDirname(), TSTrajectorySetFile() 127 128 M*/ 129 PETSC_EXTERN PetscErrorCode TSTrajectoryCreate_Basic(TSTrajectory tj,TS ts) 130 { 131 TSTrajectory_Basic *tjbasic; 132 PetscErrorCode ierr; 133 134 PetscFunctionBegin; 135 ierr = PetscNew(&tjbasic);CHKERRQ(ierr); 136 137 ierr = PetscViewerCreate(PetscObjectComm((PetscObject)tj),&tjbasic->viewer);CHKERRQ(ierr); 138 ierr = PetscViewerSetType(tjbasic->viewer,PETSCVIEWERBINARY);CHKERRQ(ierr); 139 ierr = PetscViewerPushFormat(tjbasic->viewer,PETSC_VIEWER_NATIVE);CHKERRQ(ierr); 140 ierr = PetscViewerFileSetMode(tjbasic->viewer,FILE_MODE_WRITE);CHKERRQ(ierr); 141 tj->data = tjbasic; 142 143 tj->ops->set = TSTrajectorySet_Basic; 144 tj->ops->get = TSTrajectoryGet_Basic; 145 tj->ops->setup = TSTrajectorySetUp_Basic; 146 tj->ops->destroy = TSTrajectoryDestroy_Basic; 147 tj->ops->setfromoptions = TSTrajectorySetFromOptions_Basic; 148 PetscFunctionReturn(0); 149 } 150