xref: /petsc/src/dm/impls/moab/tests/ex2.cxx (revision 6a98f8dc3f2c9149905a87dc2e9d0fedaf64e09a)
1 static char help[] = "Create a box mesh with DMMoab and test defining a tag on the mesh\n\n";
2 
3 #include <petscdmmoab.h>
4 
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   PetscBool     simplex;                        /* Use simplex elements */
13   PetscBool     interlace;
14   char          input_file[PETSC_MAX_PATH_LEN];   /* Import mesh from file */
15   char          output_file[PETSC_MAX_PATH_LEN];   /* Output mesh file name */
16   PetscBool     write_output;                        /* Write output mesh and data to file */
17   PetscInt      nfields;         /* Number of fields */
18   char          *fieldnames[PETSC_MAX_PATH_LEN]; /* Name of a defined field on the mesh */
19 } AppCtx;
20 
21 PetscErrorCode ProcessOptions(MPI_Comm comm, AppCtx *options)
22 {
23   PetscErrorCode ierr;
24   PetscBool      flg;
25 
26   PetscFunctionBegin;
27   options->debug             = PETSC_FALSE;
28   options->dim               = 2;
29   options->nele              = 5;
30   options->nfields           = 256;
31   options->simplex           = PETSC_FALSE;
32   options->write_output      = PETSC_FALSE;
33   options->interlace         = PETSC_FALSE;
34   options->input_file[0]     = '\0';
35   ierr = PetscStrcpy(options->output_file,"ex2.h5m");CHKERRQ(ierr);
36 
37   ierr = PetscOptionsBegin(comm, "", "Meshing Problem Options", "DMMOAB");CHKERRQ(ierr);
38   ierr = PetscOptionsBool("-debug", "Enable debug messages", "ex2.cxx", options->debug, &options->debug, NULL);CHKERRQ(ierr);
39   ierr = PetscOptionsBool("-interlace", "Use interlaced arrangement for the field data", "ex2.cxx", options->interlace, &options->interlace, NULL);CHKERRQ(ierr);
40   ierr = PetscOptionsBool("-simplex", "Create simplices instead of tensor product elements", "ex2.cxx", options->simplex, &options->simplex, NULL);CHKERRQ(ierr);
41   ierr = PetscOptionsRangeInt("-dim", "The topological mesh dimension", "ex2.cxx", options->dim, &options->dim, NULL,1,3);CHKERRQ(ierr);
42   ierr = PetscOptionsBoundedInt("-n", "The number of elements in each dimension", "ex2.cxx", options->nele, &options->nele, NULL,1);CHKERRQ(ierr);
43   ierr = PetscOptionsString("-meshfile", "The input mesh file", "ex2.cxx", options->input_file, options->input_file, PETSC_MAX_PATH_LEN, NULL);CHKERRQ(ierr);
44   ierr = PetscOptionsString("-io", "Write out the mesh and solution that is defined on it (Default H5M format)", "ex2.cxx", options->output_file, options->output_file, PETSC_MAX_PATH_LEN, &options->write_output);CHKERRQ(ierr);
45   ierr = PetscOptionsStringArray("-fields", "The list of names of the field variables", "ex2.cxx", options->fieldnames,&options->nfields, &flg);CHKERRQ(ierr);
46   ierr = PetscOptionsEnd();
47 
48   if (options->debug) PetscPrintf(comm, "Total number of fields: %D.\n",options->nfields);
49   if (!flg) { /* if no field names were given by user, assign a default */
50     options->nfields = 1;
51     ierr = PetscStrallocpy("TestEX2Var",&options->fieldnames[0]);CHKERRQ(ierr);
52   }
53 
54   ierr = PetscLogEventRegister("CreateMesh",          DM_CLASSID,   &options->createMeshEvent);CHKERRQ(ierr);
55   PetscFunctionReturn(0);
56 }
57 
58 PetscErrorCode CreateMesh(MPI_Comm comm, AppCtx *user)
59 {
60   PetscInt       i;
61   size_t         len;
62   PetscMPIInt    rank;
63   PetscErrorCode ierr;
64 
65   PetscFunctionBegin;
66   ierr = PetscLogEventBegin(user->createMeshEvent,0,0,0,0);CHKERRQ(ierr);
67   ierr = MPI_Comm_rank(comm, &rank);CHKERRQ(ierr);
68   ierr = PetscStrlen(user->input_file, &len);CHKERRQ(ierr);
69   if (len) {
70     if (user->debug) PetscPrintf(comm, "Loading mesh from file: %s and creating a DM object.\n",user->input_file);
71     ierr = DMMoabLoadFromFile(comm, user->dim, 1, user->input_file, "", &user->dm);CHKERRQ(ierr);
72   }
73   else {
74     if (user->debug) {
75       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);
76     }
77     ierr = DMMoabCreateBoxMesh(comm, user->dim, user->simplex, NULL, user->nele, 1, &user->dm);CHKERRQ(ierr);
78   }
79 
80   if (user->debug) {
81     PetscPrintf(comm, "Setting field names to DM: \n");
82     for (i=0; i<user->nfields; i++)
83       PetscPrintf(comm, "\t Field{%D} = %s.\n",i,user->fieldnames[i]);
84   }
85   ierr     = DMMoabSetFieldNames(user->dm, user->nfields, (const char**)user->fieldnames);CHKERRQ(ierr);
86   ierr     = PetscObjectSetName((PetscObject)user->dm, "Structured Mesh");CHKERRQ(ierr);
87   ierr     = PetscLogEventEnd(user->createMeshEvent,0,0,0,0);CHKERRQ(ierr);
88   PetscFunctionReturn(0);
89 }
90 
91 int main(int argc, char **argv)
92 {
93   AppCtx         user;                 /* user-defined work context */
94   PetscRandom    rctx;
95   Vec            solution;
96   Mat            system;
97   MPI_Comm       comm;
98   PetscInt       i;
99   PetscErrorCode ierr;
100 
101   ierr = PetscInitialize(&argc, &argv, NULL, help);if (ierr) return ierr;
102   comm = PETSC_COMM_WORLD;
103   ierr = ProcessOptions(comm, &user);CHKERRQ(ierr);
104   ierr = CreateMesh(comm, &user);CHKERRQ(ierr);
105 
106   /* set block size */
107   ierr = DMMoabSetBlockSize(user.dm, (user.interlace?user.nfields:1));CHKERRQ(ierr);
108   ierr = DMSetMatType(user.dm,MATAIJ);CHKERRQ(ierr);
109 
110   ierr = DMSetFromOptions(user.dm);CHKERRQ(ierr);
111 
112   /* SetUp the data structures for DMMOAB */
113   ierr = DMSetUp(user.dm);CHKERRQ(ierr);
114 
115   if (user.debug) PetscPrintf(comm, "Creating a global vector defined on DM and setting random data.\n");
116   ierr = DMCreateGlobalVector(user.dm,&solution);CHKERRQ(ierr);
117   ierr = PetscRandomCreate(comm,&rctx);CHKERRQ(ierr);
118   ierr = VecSetRandom(solution,rctx);CHKERRQ(ierr);
119 
120   /* test if matrix allocation for the prescribed matrix type is done correctly */
121   if (user.debug) PetscPrintf(comm, "Creating a global matrix defined on DM with the right block structure.\n");
122   ierr = DMCreateMatrix(user.dm,&system);CHKERRQ(ierr);
123 
124   if (user.write_output) {
125     ierr = DMMoabSetGlobalFieldVector(user.dm, solution);CHKERRQ(ierr);
126     if (user.debug) PetscPrintf(comm, "Output mesh and associated field data to file: %s.\n",user.output_file);
127     ierr = DMMoabOutput(user.dm,(const char*)user.output_file,"");CHKERRQ(ierr);
128   }
129 
130   if (user.nfields) {
131     for(i=0; i<user.nfields; i++) {
132       ierr = PetscFree(user.fieldnames[i]);CHKERRQ(ierr);
133     }
134   }
135   ierr = PetscRandomDestroy(&rctx);CHKERRQ(ierr);
136   ierr = VecDestroy(&solution);CHKERRQ(ierr);
137   ierr = MatDestroy(&system);CHKERRQ(ierr);
138   ierr = DMDestroy(&user.dm);CHKERRQ(ierr);
139   ierr = PetscFinalize();
140   return ierr;
141 }
142 
143 /*TEST
144 
145    build:
146      requires: moab
147 
148    test:
149      args: -debug -fields v1,v2,v3
150 
151 TEST*/
152