1 2 #include <petsc/private/tsimpl.h> /*I "petscts.h" I*/ 3 4 static PetscErrorCode OutputBIN(MPI_Comm comm,const char *filename,PetscViewer *viewer) 5 { 6 PetscErrorCode ierr; 7 8 PetscFunctionBegin; 9 ierr = PetscViewerCreate(comm,viewer);CHKERRQ(ierr); 10 ierr = PetscViewerSetType(*viewer,PETSCVIEWERBINARY);CHKERRQ(ierr); 11 ierr = PetscViewerFileSetMode(*viewer,FILE_MODE_WRITE);CHKERRQ(ierr); 12 ierr = PetscViewerFileSetName(*viewer,filename);CHKERRQ(ierr); 13 PetscFunctionReturn(0); 14 } 15 16 static PetscErrorCode TSTrajectorySet_Basic(TSTrajectory tj,TS ts,PetscInt stepnum,PetscReal time,Vec X) 17 { 18 PetscViewer viewer; 19 PetscInt ns,i; 20 Vec *Y; 21 char filename[PETSC_MAX_PATH_LEN],template[PETSC_MAX_PATH_LEN]; 22 PetscReal tprev; 23 PetscErrorCode ierr; 24 MPI_Comm comm; 25 26 PetscFunctionBegin; 27 ierr = PetscSNPrintf(template,sizeof(template),"%s%s",tj->dirname,tj->filetemplate);CHKERRQ(ierr); 28 ierr = PetscObjectGetComm((PetscObject)ts,&comm);CHKERRQ(ierr); 29 ierr = TSGetStepNumber(ts,&stepnum);CHKERRQ(ierr); 30 if (stepnum == 0) { 31 PetscMPIInt rank; 32 ierr = MPI_Comm_rank(comm,&rank);CHKERRQ(ierr); 33 if (!rank) { 34 ierr = PetscRMTree(tj->dirname);CHKERRQ(ierr); 35 ierr = PetscMkdir(tj->dirname);CHKERRQ(ierr); 36 } 37 ierr = PetscSNPrintf(filename,sizeof(filename),template,stepnum);CHKERRQ(ierr); 38 ierr = OutputBIN(comm,filename,&viewer);CHKERRQ(ierr); 39 ierr = VecView(X,viewer);CHKERRQ(ierr); 40 ierr = PetscViewerBinaryWrite(viewer,&time,1,PETSC_REAL,PETSC_FALSE);CHKERRQ(ierr); 41 ierr = PetscViewerDestroy(&viewer);CHKERRQ(ierr); 42 PetscFunctionReturn(0); 43 } 44 ierr = PetscSNPrintf(filename,sizeof(filename),template,stepnum);CHKERRQ(ierr); 45 ierr = OutputBIN(comm,filename,&viewer);CHKERRQ(ierr); 46 ierr = VecView(X,viewer);CHKERRQ(ierr); 47 ierr = PetscViewerBinaryWrite(viewer,&time,1,PETSC_REAL,PETSC_FALSE);CHKERRQ(ierr); 48 49 ierr = TSGetStages(ts,&ns,&Y);CHKERRQ(ierr); 50 for (i=0;i<ns;i++) { 51 ierr = VecView(Y[i],viewer);CHKERRQ(ierr); 52 } 53 54 ierr = TSGetPrevTime(ts,&tprev);CHKERRQ(ierr); 55 ierr = PetscViewerBinaryWrite(viewer,&tprev,1,PETSC_REAL,PETSC_FALSE);CHKERRQ(ierr); 56 57 ierr = PetscViewerDestroy(&viewer);CHKERRQ(ierr); 58 PetscFunctionReturn(0); 59 } 60 61 static PetscErrorCode TSTrajectoryGet_Basic(TSTrajectory tj,TS ts,PetscInt stepnum,PetscReal *t) 62 { 63 Vec Sol,*Y; 64 PetscInt Nr,i; 65 PetscViewer viewer; 66 PetscReal timepre; 67 char filename[PETSC_MAX_PATH_LEN],template[PETSC_MAX_PATH_LEN]; 68 PetscErrorCode ierr; 69 70 PetscFunctionBegin; 71 ierr = PetscSNPrintf(template,sizeof(template),"%s%s",tj->dirname,tj->filetemplate);CHKERRQ(ierr); 72 ierr = PetscSNPrintf(filename,sizeof filename,template,stepnum);CHKERRQ(ierr); 73 ierr = PetscViewerBinaryOpen(PETSC_COMM_WORLD,filename,FILE_MODE_READ,&viewer);CHKERRQ(ierr); 74 75 ierr = TSGetSolution(ts,&Sol);CHKERRQ(ierr); 76 ierr = VecLoad(Sol,viewer);CHKERRQ(ierr); 77 78 ierr = PetscViewerBinaryRead(viewer,t,1,NULL,PETSC_REAL);CHKERRQ(ierr); 79 80 if (stepnum != 0) { 81 ierr = TSGetStages(ts,&Nr,&Y);CHKERRQ(ierr); 82 for (i=0;i<Nr ;i++) { 83 ierr = VecLoad(Y[i],viewer);CHKERRQ(ierr); 84 } 85 ierr = PetscViewerBinaryRead(viewer,&timepre,1,NULL,PETSC_REAL);CHKERRQ(ierr); 86 ierr = TSSetTimeStep(ts,-(*t)+timepre);CHKERRQ(ierr); 87 } 88 89 ierr = PetscViewerDestroy(&viewer);CHKERRQ(ierr); 90 PetscFunctionReturn(0); 91 } 92 93 static PetscErrorCode TSTrajectoryDestroy_Basic(TSTrajectory tj) 94 { 95 PetscErrorCode ierr; 96 PetscMPIInt rank; 97 MPI_Comm comm; 98 99 PetscFunctionBegin; 100 if (!tj->keepfiles) { 101 ierr = PetscObjectGetComm((PetscObject)tj,&comm);CHKERRQ(ierr); 102 ierr = MPI_Comm_rank(comm,&rank);CHKERRQ(ierr); 103 if (!rank) { 104 ierr = PetscRMTree(tj->dirname);CHKERRQ(ierr); 105 } 106 } 107 PetscFunctionReturn(0); 108 } 109 110 /*MC 111 TSTRAJECTORYBASIC - Stores each solution of the ODE/DAE in a file 112 113 Saves each timestep into a seperate file named SA-data/SA-%06d.bin. The file name can be changed. 114 115 This version saves the solutions at all the stages 116 117 $PETSC_DIR/share/petsc/matlab/PetscReadBinaryTrajectory.m can read in files created with this format 118 119 Level: intermediate 120 121 .seealso: TSTrajectoryCreate(), TS, TSTrajectorySetType(), TSTrajectorySetDirname(), TSTrajectorySetFile() 122 123 M*/ 124 PETSC_EXTERN PetscErrorCode TSTrajectoryCreate_Basic(TSTrajectory tj,TS ts) 125 { 126 PetscFunctionBegin; 127 tj->keepfiles = PETSC_FALSE; 128 tj->ops->set = TSTrajectorySet_Basic; 129 tj->ops->get = TSTrajectoryGet_Basic; 130 tj->ops->destroy = TSTrajectoryDestroy_Basic; 131 PetscFunctionReturn(0); 132 } 133