Actual source code: ex3.cxx

  1: static char help[] = "Create a box mesh with DMMoab and test defining a tag on the mesh\n\n";

  3: #include <petscdmmoab.h>

  5: typedef struct {
  6:   DM            dm;    /* DM implementation using the MOAB interface */
  7:   PetscBool     debug; /* The debugging level */
  8:   PetscLogEvent createMeshEvent;
  9:   /* Domain and mesh definition */
 10:   PetscInt  dim;                             /* The topological mesh dimension */
 11:   PetscInt  nele;                            /* Elements in each dimension */
 12:   PetscInt  degree;                          /* Degree of refinement */
 13:   PetscBool simplex;                         /* Use simplex elements */
 14:   PetscInt  nlevels;                         /* Number of levels in mesh hierarchy */
 15:   PetscInt  nghost;                          /* Number of ghost layers in the mesh */
 16:   char      input_file[PETSC_MAX_PATH_LEN];  /* Import mesh from file */
 17:   char      output_file[PETSC_MAX_PATH_LEN]; /* Output mesh file name */
 18:   PetscBool write_output;                    /* Write output mesh and data to file */
 19: } AppCtx;

 21: PetscErrorCode ProcessOptions(MPI_Comm comm, AppCtx *options)
 22: {
 23:   options->debug         = PETSC_FALSE;
 24:   options->nlevels       = 1;
 25:   options->nghost        = 1;
 26:   options->dim           = 2;
 27:   options->nele          = 5;
 28:   options->degree        = 2;
 29:   options->simplex       = PETSC_FALSE;
 30:   options->write_output  = PETSC_FALSE;
 31:   options->input_file[0] = '\0';
 32:   PetscStrcpy(options->output_file, "ex3.h5m");

 34:   PetscOptionsBegin(comm, "", "Uniform Mesh Refinement Options", "DMMOAB");
 35:   PetscOptionsBool("-debug", "Enable debug messages", "ex2.cxx", options->debug, &options->debug, NULL);
 36:   PetscOptionsRangeInt("-dim", "The topological mesh dimension", "ex3.cxx", options->dim, &options->dim, NULL, 0, 3);
 37:   PetscOptionsBoundedInt("-n", "The number of elements in each dimension", "ex3.cxx", options->nele, &options->nele, NULL, 1);
 38:   PetscOptionsBoundedInt("-levels", "Number of levels in the hierarchy", "ex3.cxx", options->nlevels, &options->nlevels, NULL, 0);
 39:   PetscOptionsBoundedInt("-degree", "Number of degrees at each level of refinement", "ex3.cxx", options->degree, &options->degree, NULL, 0);
 40:   PetscOptionsBoundedInt("-ghost", "Number of ghost layers in the mesh", "ex3.cxx", options->nghost, &options->nghost, NULL, 0);
 41:   PetscOptionsBool("-simplex", "Create simplices instead of tensor product elements", "ex3.cxx", options->simplex, &options->simplex, NULL);
 42:   PetscOptionsString("-input", "The input mesh file", "ex3.cxx", options->input_file, options->input_file, sizeof(options->input_file), NULL);
 43:   PetscOptionsString("-io", "Write out the mesh and solution that is defined on it (Default H5M format)", "ex3.cxx", options->output_file, options->output_file, sizeof(options->output_file), &options->write_output);
 44:   PetscOptionsEnd();

 46:   PetscLogEventRegister("CreateMesh", DM_CLASSID, &options->createMeshEvent);
 47:   return 0;
 48: };

 50: PetscErrorCode CreateMesh(MPI_Comm comm, AppCtx *user)
 51: {
 52:   size_t      len;
 53:   PetscMPIInt rank;

 55:   PetscLogEventBegin(user->createMeshEvent, 0, 0, 0, 0);
 56:   MPI_Comm_rank(comm, &rank);
 57:   PetscStrlen(user->input_file, &len);
 58:   if (len) {
 59:     if (user->debug) PetscPrintf(comm, "Loading mesh from file: %s and creating the coarse level DM object.\n", user->input_file);
 60:     DMMoabLoadFromFile(comm, user->dim, user->nghost, user->input_file, "", &user->dm);
 61:   } else {
 62:     if (user->debug)
 63:       PetscCall(PetscPrintf(comm, "Creating a %" PetscInt_FMT "-dimensional structured %s mesh of %" PetscInt_FMT "x%" PetscInt_FMT "x%" PetscInt_FMT " in memory and creating a DM object.\n", user->dim, (user->simplex ? "simplex" : "regular"), user->nele,
 64:                             user->nele, user->nele));
 65:     DMMoabCreateBoxMesh(comm, user->dim, user->simplex, NULL, user->nele, user->nghost, &user->dm);
 66:   }
 67:   PetscObjectSetName((PetscObject)user->dm, "Coarse Mesh");
 68:   PetscLogEventEnd(user->createMeshEvent, 0, 0, 0, 0);
 69:   return 0;
 70: }

 72: int main(int argc, char **argv)
 73: {
 74:   AppCtx    user; /* user-defined work context */
 75:   MPI_Comm  comm;
 76:   PetscInt  i;
 77:   Mat       R;
 78:   DM       *dmhierarchy;
 79:   PetscInt *degrees;
 80:   PetscBool createR;

 83:   PetscInitialize(&argc, &argv, NULL, help);
 84:   comm    = PETSC_COMM_WORLD;
 85:   createR = PETSC_FALSE;

 87:   ProcessOptions(comm, &user);
 88:   CreateMesh(comm, &user);
 89:   DMSetFromOptions(user.dm);

 91:   /* SetUp the data structures for DMMOAB */
 92:   DMSetUp(user.dm);

 94:   PetscMalloc(sizeof(DM) * (user.nlevels + 1), &dmhierarchy);
 95:   for (i = 0; i <= user.nlevels; i++) dmhierarchy[i] = NULL;

 97:   // coarsest grid = 0
 98:   // finest grid = nlevels
 99:   dmhierarchy[0] = user.dm;
100:   PetscObjectReference((PetscObject)user.dm);

102:   if (user.nlevels) {
103:     PetscMalloc1(user.nlevels, &degrees);
104:     for (i = 0; i < user.nlevels; i++) degrees[i] = user.degree;
105:     if (user.debug) PetscPrintf(comm, "Generate the MOAB mesh hierarchy with %" PetscInt_FMT " levels.\n", user.nlevels);
106:     DMMoabGenerateHierarchy(user.dm, user.nlevels, degrees);

108:     PetscBool usehierarchy = PETSC_FALSE;
109:     if (usehierarchy) DMRefineHierarchy(user.dm, user.nlevels, &dmhierarchy[1]);
110:     else {
111:       if (user.debug) {
112:         PetscPrintf(PETSC_COMM_WORLD, "Level %" PetscInt_FMT "\n", 0);
113:         DMView(user.dm, 0);
114:       }
115:       for (i = 1; i <= user.nlevels; i++) {
116:         if (user.debug) PetscPrintf(PETSC_COMM_WORLD, "Level %" PetscInt_FMT "\n", i);
117:         DMRefine(dmhierarchy[i - 1], MPI_COMM_NULL, &dmhierarchy[i]);
118:         if (createR) DMCreateInterpolation(dmhierarchy[i - 1], dmhierarchy[i], &R, NULL);
119:         if (user.debug) {
120:           DMView(dmhierarchy[i], 0);
121:           if (createR) MatView(R, 0);
122:         }
123:         /* Solvers could now set operator "R" to the multigrid PC object for level i
124:             PCMGSetInterpolation(pc,i,R)
125:         */
126:         if (createR) MatDestroy(&R);
127:       }
128:     }
129:   }

131:   if (user.write_output) {
132:     if (user.debug) PetscPrintf(comm, "Output mesh hierarchy to file: %s.\n", user.output_file);
133:     DMMoabOutput(dmhierarchy[user.nlevels], (const char *)user.output_file, "");
134:   }

136:   for (i = 0; i <= user.nlevels; i++) DMDestroy(&dmhierarchy[i]);
137:   PetscFree(degrees);
138:   PetscFree(dmhierarchy);
139:   DMDestroy(&user.dm);
140:   PetscFinalize();
141:   return 0;
142: }

144: /*TEST

146:      build:
147:        requires: moab

149:      test:
150:        args: -debug -n 2 -dim 2 -levels 2 -simplex
151:        filter:  grep -v "DM_0x"

153:      test:
154:        args: -debug -n 2 -dim 3 -levels 2
155:        filter:  grep -v "DM_0x"
156:        suffix: 1_2

158:      test:
159:        args: -debug -n 2 -dim 3 -ghost 1 -levels 2
160:        filter:  grep -v "DM_0x"
161:        nsize: 2
162:        suffix: 2_1

164: TEST*/