Actual source code: dmmbio.cxx
1: #include <petsc/private/dmmbimpl.h>
2: #include <petscdmmoab.h>
4: static PetscErrorCode DMMoab_GetWriteOptions_Private(PetscInt fsetid, PetscInt numproc, PetscInt dim, MoabWriteMode mode, PetscInt dbglevel, const char *dm_opts, const char *extra_opts, const char **write_opts)
5: {
6: char *wopts;
7: char wopts_par[PETSC_MAX_PATH_LEN];
8: char wopts_parid[PETSC_MAX_PATH_LEN];
9: char wopts_dbg[PETSC_MAX_PATH_LEN];
10: PetscFunctionBegin;
12: PetscCall(PetscMalloc1(PETSC_MAX_PATH_LEN, &wopts));
13: PetscCall(PetscMemzero(&wopts_par, PETSC_MAX_PATH_LEN));
14: PetscCall(PetscMemzero(&wopts_parid, PETSC_MAX_PATH_LEN));
15: PetscCall(PetscMemzero(&wopts_dbg, PETSC_MAX_PATH_LEN));
17: // do parallel read unless only one processor
18: if (numproc > 1) {
19: PetscCall(PetscSNPrintf(wopts_par, PETSC_MAX_PATH_LEN, "PARALLEL=%s;", MoabWriteModes[mode]));
20: if (fsetid >= 0) PetscCall(PetscSNPrintf(wopts_parid, PETSC_MAX_PATH_LEN, "PARALLEL_COMM=%d;", fsetid));
21: }
23: if (dbglevel) PetscCall(PetscSNPrintf(wopts_dbg, PETSC_MAX_PATH_LEN, "CPUTIME;DEBUG_IO=%d;", dbglevel));
25: PetscCall(PetscSNPrintf(wopts, PETSC_MAX_PATH_LEN, "%s%s%s%s%s", wopts_par, wopts_parid, wopts_dbg, (extra_opts ? extra_opts : ""), (dm_opts ? dm_opts : "")));
26: *write_opts = wopts;
27: PetscFunctionReturn(PETSC_SUCCESS);
28: }
30: /*@C
31: DMMoabOutput - Output the solution vectors that are stored in the DMMoab object as tags
32: along with the complete mesh data structure in the native H5M or VTK format. The H5M output file
33: can be visualized directly with Paraview (if compiled with appropriate plugin) or converted
34: with MOAB/tools/mbconvert to a VTK or Exodus file.
36: This routine can also be used for check-pointing purposes to store a complete history of
37: the solution along with any other necessary data to restart computations.
39: Collective
41: Input Parameters:
42: + dm - the discretization manager object containing solution in MOAB tags.
43: . filename - the name of the output file: e.g., poisson.h5m
44: - usrwriteopts - the parallel write options needed for serializing a MOAB mesh database. Can be NULL.
45: Reference (Parallel Mesh Initialization: http://ftp.mcs.anl.gov/pub/fathom/moab-docs/contents.html#fivetwo)
47: Level: intermediate
49: .seealso: `DMMoabLoadFromFile()`, `DMMoabSetGlobalFieldVector()`
50: @*/
51: PetscErrorCode DMMoabOutput(DM dm, const char *filename, const char *usrwriteopts)
52: {
53: DM_Moab *dmmoab;
54: const char *writeopts;
55: PetscBool isftype;
56: moab::ErrorCode merr;
58: PetscFunctionBegin;
60: dmmoab = (DM_Moab *)(dm)->data;
62: PetscCall(PetscStrendswith(filename, "h5m", &isftype));
64: /* add mesh loading options specific to the DM */
65: if (isftype) {
66: #ifdef MOAB_HAVE_MPI
67: PetscCall(DMMoab_GetWriteOptions_Private(dmmoab->pcomm->get_id(), dmmoab->pcomm->size(), dmmoab->dim, dmmoab->write_mode, dmmoab->rw_dbglevel, dmmoab->extra_write_options, usrwriteopts, &writeopts));
68: #else
69: PetscCall(DMMoab_GetWriteOptions_Private(0, 1, dmmoab->dim, dmmoab->write_mode, dmmoab->rw_dbglevel, dmmoab->extra_write_options, usrwriteopts, &writeopts));
70: #endif
71: PetscCall(PetscInfo(dm, "Writing file %s with options: %s\n", filename, writeopts));
72: } else {
73: writeopts = NULL;
74: }
76: /* output file, using parallel write */
77: merr = dmmoab->mbiface->write_file(filename, NULL, writeopts, &dmmoab->fileset, 1);
78: MBERRVM(dmmoab->mbiface, "Writing output of DMMoab failed.", merr);
79: PetscCall(PetscFree(writeopts));
80: PetscFunctionReturn(PETSC_SUCCESS);
81: }