1 2 #include <petsc/private/tsimpl.h> /*I "petscts.h" I*/ 3 4 #undef __FUNCT__ 5 #define __FUNCT__ "OutputBIN" 6 static PetscErrorCode OutputBIN(const char *filename,PetscViewer *viewer) 7 { 8 PetscErrorCode ierr; 9 10 PetscFunctionBegin; 11 ierr = PetscViewerCreate(PETSC_COMM_WORLD,viewer);CHKERRQ(ierr); 12 ierr = PetscViewerSetType(*viewer,PETSCVIEWERBINARY);CHKERRQ(ierr); 13 ierr = PetscViewerFileSetMode(*viewer,FILE_MODE_WRITE);CHKERRQ(ierr); 14 ierr = PetscViewerFileSetName(*viewer,filename);CHKERRQ(ierr); 15 PetscFunctionReturn(0); 16 } 17 18 #undef __FUNCT__ 19 #define __FUNCT__ "TSTrajectorySet_Basic" 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; 27 PetscErrorCode ierr; 28 29 PetscFunctionBegin; 30 ierr = TSGetTotalSteps(ts,&stepnum);CHKERRQ(ierr); 31 if (stepnum == 0) { 32 PetscMPIInt rank; 33 ierr = MPI_Comm_rank(PetscObjectComm((PetscObject)ts),&rank);CHKERRQ(ierr); 34 if (!rank) { 35 ierr = PetscRMTree("SA-data");CHKERRQ(ierr); 36 ierr = PetscMkdir("SA-data");CHKERRQ(ierr); 37 } 38 ierr = PetscSNPrintf(filename,sizeof(filename),"SA-data/SA-%06d.bin",stepnum);CHKERRQ(ierr); 39 ierr = OutputBIN(filename,&viewer);CHKERRQ(ierr); 40 ierr = VecView(X,viewer);CHKERRQ(ierr); 41 ierr = PetscViewerBinaryWrite(viewer,&time,1,PETSC_REAL,PETSC_FALSE);CHKERRQ(ierr); 42 ierr = PetscViewerDestroy(&viewer);CHKERRQ(ierr); 43 PetscFunctionReturn(0); 44 } 45 ierr = PetscSNPrintf(filename,sizeof(filename),"SA-data/SA-%06d.bin",stepnum);CHKERRQ(ierr); 46 ierr = OutputBIN(filename,&viewer);CHKERRQ(ierr); 47 ierr = VecView(X,viewer);CHKERRQ(ierr); 48 ierr = PetscViewerBinaryWrite(viewer,&time,1,PETSC_REAL,PETSC_FALSE);CHKERRQ(ierr); 49 50 ierr = TSGetStages(ts,&ns,&Y);CHKERRQ(ierr); 51 for (i=0;i<ns;i++) { 52 ierr = VecView(Y[i],viewer);CHKERRQ(ierr); 53 } 54 55 ierr = TSGetPrevTime(ts,&tprev);CHKERRQ(ierr); 56 ierr = PetscViewerBinaryWrite(viewer,&tprev,1,PETSC_REAL,PETSC_FALSE);CHKERRQ(ierr); 57 58 ierr = PetscViewerDestroy(&viewer);CHKERRQ(ierr); 59 PetscFunctionReturn(0); 60 } 61 62 #undef __FUNCT__ 63 #define __FUNCT__ "TSTrajectoryGet_Basic" 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]; 71 PetscErrorCode ierr; 72 73 PetscFunctionBegin; 74 ierr = TSGetTotalSteps(ts,&stepnum);CHKERRQ(ierr); 75 ierr = PetscSNPrintf(filename,sizeof filename,"SA-data/SA-%06d.bin",stepnum);CHKERRQ(ierr); 76 ierr = PetscViewerBinaryOpen(PETSC_COMM_WORLD,filename,FILE_MODE_READ,&viewer);CHKERRQ(ierr); 77 78 ierr = TSGetSolution(ts,&Sol);CHKERRQ(ierr); 79 ierr = VecLoad(Sol,viewer);CHKERRQ(ierr); 80 81 ierr = PetscViewerBinaryRead(viewer,t,1,NULL,PETSC_REAL);CHKERRQ(ierr); 82 83 if (stepnum != 0) { 84 ierr = TSGetStages(ts,&Nr,&Y);CHKERRQ(ierr); 85 for (i=0;i<Nr ;i++) { 86 ierr = VecLoad(Y[i],viewer);CHKERRQ(ierr); 87 } 88 ierr = PetscViewerBinaryRead(viewer,&timepre,1,NULL,PETSC_REAL);CHKERRQ(ierr); 89 ierr = TSSetTimeStep(ts,-(*t)+timepre);CHKERRQ(ierr); 90 } 91 92 ierr = PetscViewerDestroy(&viewer);CHKERRQ(ierr); 93 PetscFunctionReturn(0); 94 } 95 96 /*MC 97 TSTRAJECTORYBASIC - Stores each solution of the ODE/ADE in a file 98 99 Level: intermediate 100 101 .seealso: TSTrajectoryCreate(), TS, TSTrajectorySetType() 102 103 M*/ 104 #undef __FUNCT__ 105 #define __FUNCT__ "TSTrajectoryCreate_Basic" 106 PETSC_EXTERN PetscErrorCode TSTrajectoryCreate_Basic(TSTrajectory tj,TS ts) 107 { 108 PetscFunctionBegin; 109 tj->ops->set = TSTrajectorySet_Basic; 110 tj->ops->get = TSTrajectoryGet_Basic; 111 PetscFunctionReturn(0); 112 } 113