Actual source code: trajvisualization.c


  2: #include <petsc/private/tsimpl.h>

  4: static PetscErrorCode OutputBIN(MPI_Comm comm, const char *filename, PetscViewer *viewer)
  5: {
  6:   PetscViewerCreate(comm, viewer);
  7:   PetscViewerSetType(*viewer, PETSCVIEWERBINARY);
  8:   PetscViewerFileSetMode(*viewer, FILE_MODE_WRITE);
  9:   PetscViewerFileSetName(*viewer, filename);
 10:   return 0;
 11: }

 13: static PetscErrorCode TSTrajectorySet_Visualization(TSTrajectory tj, TS ts, PetscInt stepnum, PetscReal time, Vec X)
 14: {
 15:   PetscViewer viewer;
 16:   char        filename[PETSC_MAX_PATH_LEN];
 17:   PetscReal   tprev;
 18:   MPI_Comm    comm;

 20:   PetscObjectGetComm((PetscObject)ts, &comm);
 21:   if (stepnum == 0) {
 22:     PetscMPIInt rank;
 23:     MPI_Comm_rank(comm, &rank);
 24:     if (rank == 0) {
 25:       PetscRMTree("Visualization-data");
 26:       PetscMkdir("Visualization-data");
 27:     }
 28:     if (tj->names) {
 29:       PetscViewer bnames;
 30:       PetscViewerBinaryOpen(comm, "Visualization-data/variablenames", FILE_MODE_WRITE, &bnames);
 31:       PetscViewerBinaryWriteStringArray(bnames, (const char *const *)tj->names);
 32:       PetscViewerDestroy(&bnames);
 33:     }
 34:     PetscSNPrintf(filename, sizeof(filename), "Visualization-data/SA-%06" PetscInt_FMT ".bin", stepnum);
 35:     OutputBIN(comm, filename, &viewer);
 36:     if (!tj->transform) {
 37:       VecView(X, viewer);
 38:     } else {
 39:       Vec XX;
 40:       (*tj->transform)(tj->transformctx, X, &XX);
 41:       VecView(XX, viewer);
 42:       VecDestroy(&XX);
 43:     }
 44:     PetscViewerBinaryWrite(viewer, &time, 1, PETSC_REAL);
 45:     PetscViewerDestroy(&viewer);
 46:     return 0;
 47:   }
 48:   PetscSNPrintf(filename, sizeof(filename), "Visualization-data/SA-%06" PetscInt_FMT ".bin", stepnum);
 49:   OutputBIN(comm, filename, &viewer);
 50:   if (!tj->transform) {
 51:     VecView(X, viewer);
 52:   } else {
 53:     Vec XX;
 54:     (*tj->transform)(tj->transformctx, X, &XX);
 55:     VecView(XX, viewer);
 56:     VecDestroy(&XX);
 57:   }
 58:   PetscViewerBinaryWrite(viewer, &time, 1, PETSC_REAL);

 60:   TSGetPrevTime(ts, &tprev);
 61:   PetscViewerBinaryWrite(viewer, &tprev, 1, PETSC_REAL);

 63:   PetscViewerDestroy(&viewer);
 64:   return 0;
 65: }

 67: /*MC
 68:       TSTRAJECTORYVISUALIZATION - Stores each solution of the ODE/DAE in a file

 70:       Saves each timestep into a separate file in Visualization-data/SA-%06d.bin

 72:       This version saves only the solutions at each timestep, it does not save the solution at each stage,
 73:       see `TSTRAJECTORYBASIC` that saves all stages

 75:       $PETSC_DIR/share/petsc/matlab/PetscReadBinaryTrajectory.m and $PETSC_DIR/lib/petsc/bin/PetscBinaryIOTrajectory.py
 76:       can read in files created with this format into MATLAB and Python.

 78:   Level: intermediate

 80: .seealso: [](chapter_ts), `TSTrajectoryCreate()`, `TS`, `TSTrajectorySetType()`, `TSTrajectoryType`, `TSTrajectorySetVariableNames()`,
 81:           `TSTrajectoryType`, `TSTrajectory`
 82: M*/
 83: PETSC_EXTERN PetscErrorCode TSTrajectoryCreate_Visualization(TSTrajectory tj, TS ts)
 84: {
 85:   tj->ops->set    = TSTrajectorySet_Visualization;
 86:   tj->setupcalled = PETSC_TRUE;
 87:   return 0;
 88: }