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 PetscFunctionBegin; 7 PetscCall(PetscViewerCreate(comm, viewer)); 8 PetscCall(PetscViewerSetType(*viewer, PETSCVIEWERBINARY)); 9 PetscCall(PetscViewerFileSetMode(*viewer, FILE_MODE_WRITE)); 10 PetscCall(PetscViewerFileSetName(*viewer, filename)); 11 PetscFunctionReturn(PETSC_SUCCESS); 12 } 13 14 static PetscErrorCode TSTrajectorySet_Visualization(TSTrajectory tj, TS ts, PetscInt stepnum, PetscReal time, Vec X) 15 { 16 PetscViewer viewer; 17 char filename[PETSC_MAX_PATH_LEN]; 18 PetscReal tprev; 19 MPI_Comm comm; 20 21 PetscFunctionBegin; 22 PetscCall(PetscObjectGetComm((PetscObject)ts, &comm)); 23 if (stepnum == 0) { 24 PetscMPIInt rank; 25 PetscCallMPI(MPI_Comm_rank(comm, &rank)); 26 if (rank == 0) { 27 PetscCall(PetscRMTree("Visualization-data")); 28 PetscCall(PetscMkdir("Visualization-data")); 29 } 30 if (tj->names) { 31 PetscViewer bnames; 32 PetscCall(PetscViewerBinaryOpen(comm, "Visualization-data/variablenames", FILE_MODE_WRITE, &bnames)); 33 PetscCall(PetscViewerBinaryWriteStringArray(bnames, (const char *const *)tj->names)); 34 PetscCall(PetscViewerDestroy(&bnames)); 35 } 36 PetscCall(PetscSNPrintf(filename, sizeof(filename), "Visualization-data/SA-%06" PetscInt_FMT ".bin", stepnum)); 37 PetscCall(OutputBIN(comm, filename, &viewer)); 38 if (!tj->transform) { 39 PetscCall(VecView(X, viewer)); 40 } else { 41 Vec XX; 42 PetscCall((*tj->transform)(tj->transformctx, X, &XX)); 43 PetscCall(VecView(XX, viewer)); 44 PetscCall(VecDestroy(&XX)); 45 } 46 PetscCall(PetscViewerBinaryWrite(viewer, &time, 1, PETSC_REAL)); 47 PetscCall(PetscViewerDestroy(&viewer)); 48 PetscFunctionReturn(PETSC_SUCCESS); 49 } 50 PetscCall(PetscSNPrintf(filename, sizeof(filename), "Visualization-data/SA-%06" PetscInt_FMT ".bin", stepnum)); 51 PetscCall(OutputBIN(comm, filename, &viewer)); 52 if (!tj->transform) { 53 PetscCall(VecView(X, viewer)); 54 } else { 55 Vec XX; 56 PetscCall((*tj->transform)(tj->transformctx, X, &XX)); 57 PetscCall(VecView(XX, viewer)); 58 PetscCall(VecDestroy(&XX)); 59 } 60 PetscCall(PetscViewerBinaryWrite(viewer, &time, 1, PETSC_REAL)); 61 62 PetscCall(TSGetPrevTime(ts, &tprev)); 63 PetscCall(PetscViewerBinaryWrite(viewer, &tprev, 1, PETSC_REAL)); 64 65 PetscCall(PetscViewerDestroy(&viewer)); 66 PetscFunctionReturn(PETSC_SUCCESS); 67 } 68 69 /*MC 70 TSTRAJECTORYVISUALIZATION - Stores each solution of the ODE/DAE in a file 71 72 Saves each timestep into a separate file in Visualization-data/SA-%06d.bin 73 74 This version saves only the solutions at each timestep, it does not save the solution at each stage, 75 see `TSTRAJECTORYBASIC` that saves all stages 76 77 $PETSC_DIR/share/petsc/matlab/PetscReadBinaryTrajectory.m and $PETSC_DIR/lib/petsc/bin/PetscBinaryIOTrajectory.py 78 can read in files created with this format into MATLAB and Python. 79 80 Level: intermediate 81 82 .seealso: [](ch_ts), `TSTrajectoryCreate()`, `TS`, `TSTrajectorySetType()`, `TSTrajectoryType`, `TSTrajectorySetVariableNames()`, 83 `TSTrajectoryType`, `TSTrajectory` 84 M*/ 85 PETSC_EXTERN PetscErrorCode TSTrajectoryCreate_Visualization(TSTrajectory tj, TS ts) 86 { 87 PetscFunctionBegin; 88 tj->ops->set = TSTrajectorySet_Visualization; 89 tj->setupcalled = PETSC_TRUE; 90 PetscFunctionReturn(PETSC_SUCCESS); 91 } 92