Actual source code: ex3.cxx

petsc-3.11.4 2019-09-28
Report Typos and Errors
  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: {

 26:   options->debug             = PETSC_FALSE;
 27:   options->nlevels           = 1;
 28:   options->nghost            = 1;
 29:   options->dim               = 2;
 30:   options->nele              = 5;
 31:   options->degree            = 2;
 32:   options->simplex           = PETSC_FALSE;
 33:   options->write_output      = PETSC_FALSE;
 34:   options->input_file[0]     = '\0';
 35:   PetscStrcpy(options->output_file,"ex3.h5m");

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

 49:   PetscLogEventRegister("CreateMesh",          DM_CLASSID,   &options->createMeshEvent);
 50:   return(0);
 51: };

 53: PetscErrorCode CreateMesh(MPI_Comm comm, AppCtx *user)
 54: {
 55:   size_t         len;
 56:   PetscMPIInt    rank;

 60:   PetscLogEventBegin(user->createMeshEvent,0,0,0,0);
 61:   MPI_Comm_rank(comm, &rank);
 62:   PetscStrlen(user->input_file, &len);
 63:   if (len) {
 64:     if (user->debug) PetscPrintf(comm, "Loading mesh from file: %s and creating the coarse level DM object.\n",user->input_file);
 65:     DMMoabLoadFromFile(comm, user->dim, user->nghost, user->input_file, "", &user->dm);
 66:   }
 67:   else {
 68:     if (user->debug) {
 69:       PetscPrintf(comm, "Creating a %D-dimensional structured %s mesh of %Dx%Dx%D in memory and creating a DM object.\n",user->dim,(user->simplex?"simplex":"regular"),user->nele,user->nele,user->nele);
 70:     }
 71:     DMMoabCreateBoxMesh(comm, user->dim, user->simplex, NULL, user->nele, user->nghost, &user->dm);
 72:   }

 74:   PetscObjectSetName((PetscObject)user->dm, "Coarse Mesh");
 75:   PetscLogEventEnd(user->createMeshEvent,0,0,0,0);
 76:   return(0);
 77: }

 79: int main(int argc, char **argv)
 80: {
 81:   AppCtx         user;                 /* user-defined work context */
 82:   MPI_Comm       comm;
 83:   PetscInt       i;
 84:   Mat            R;
 85:   DM            *dmhierarchy;
 86:   PetscInt      *degrees;
 87:   PetscBool      createR;

 90:   PetscInitialize(&argc, &argv, NULL, help);if (ierr) return ierr;
 91:   comm = PETSC_COMM_WORLD;
 92:   createR = PETSC_FALSE;

 94:   ProcessOptions(comm, &user);
 95:   CreateMesh(comm, &user);
 96:   DMSetFromOptions(user.dm);

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

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

104:   // coarsest grid = 0
105:   // finest grid = nlevels
106:   dmhierarchy[0] = user.dm;
107:   PetscObjectReference((PetscObject)user.dm);

109:   if (user.nlevels) {
110:     PetscMalloc1(user.nlevels, &degrees);
111:     for (i=0; i < user.nlevels; i++) degrees[i] = user.degree;
112:     if (user.debug) PetscPrintf(comm, "Generate the MOAB mesh hierarchy with %D levels.\n", user.nlevels);
113:     DMMoabGenerateHierarchy(user.dm,user.nlevels,degrees);

115:     PetscBool usehierarchy=PETSC_FALSE;
116:     if (usehierarchy) {
117:       DMRefineHierarchy(user.dm,user.nlevels,&dmhierarchy[1]);
118:     }
119:     else {
120:       if (user.debug) {
121:         PetscPrintf(PETSC_COMM_WORLD, "Level %D\n", 0);
122:         DMView(user.dm, 0);
123:       }
124:       for (i=1; i<=user.nlevels; i++) {
125:         if (user.debug) PetscPrintf(PETSC_COMM_WORLD, "Level %D\n", i);
126:         DMRefine(dmhierarchy[i-1],MPI_COMM_NULL,&dmhierarchy[i]);
127:         if (createR) {
128:           DMCreateInterpolation(dmhierarchy[i-1],dmhierarchy[i],&R,NULL);
129:         }
130:         if (user.debug) {
131:           DMView(dmhierarchy[i], 0);
132:           if (createR) {
133:             MatView(R,0);
134:           }
135:         }
136:         /* Solvers could now set operator "R" to the multigrid PC object for level i 
137:             PCMGSetInterpolation(pc,i,R)
138:         */
139:         if (createR) {
140:           MatDestroy(&R);
141:         }
142:       }
143:     }
144:   }

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

151:   for (i=0; i<=user.nlevels; i++) {
152:     DMDestroy(&dmhierarchy[i]);
153:   }
154:   PetscFree(degrees);
155:   PetscFree(dmhierarchy);
156:   DMDestroy(&user.dm);
157:   PetscFinalize();
158:   return 0;
159: }

161: /*TEST

163:      build:
164:        requires: moab

166:      test:
167:        args: -debug -n 2 -dim 2 -levels 2 -simplex
168:        filter:  grep -v "DM_0x*"

170:      test:
171:        args: -debug -n 2 -dim 3 -levels 2
172:        filter:  grep -v "DM_0x*"
173:        suffix: 1_2

175:      test:
176:        args: -debug -n 2 -dim 3 -ghost 1 -levels 2
177:        filter:  grep -v "DM_0x*"
178:        nsize: 2
179:        suffix: 2_1

181: TEST*/