1 static char help[] = "Simple MOAB example\n\n"; 2 3 #include <petscdmmoab.h> 4 #include "moab/ScdInterface.hpp" 5 6 typedef struct { 7 DM dm; /* DM implementation using the MOAB interface */ 8 PetscLogEvent createMeshEvent; 9 /* Domain and mesh definition */ 10 PetscInt dim; 11 char filename[PETSC_MAX_PATH_LEN]; 12 char tagname[PETSC_MAX_PATH_LEN]; 13 } AppCtx; 14 15 PetscErrorCode ProcessOptions(MPI_Comm comm, AppCtx *options) 16 { 17 PetscErrorCode ierr; 18 PetscBool flg; 19 20 PetscFunctionBegin; 21 ierr = PetscStrcpy(options->filename, "");CHKERRQ(ierr); 22 ierr = PetscStrcpy(options->tagname, "petsc_tag");CHKERRQ(ierr); 23 options->dim = -1; 24 25 ierr = PetscOptionsBegin(comm, "", "MOAB example options", "DMMOAB");CHKERRQ(ierr); 26 ierr = PetscOptionsRangeInt("-dim", "The topological mesh dimension", "ex1.cxx", options->dim, &options->dim, NULL,PETSC_DECIDE,3);CHKERRQ(ierr); 27 ierr = PetscOptionsString("-filename", "The file containing the mesh", "ex1.cxx", options->filename, options->filename, sizeof(options->filename), NULL);CHKERRQ(ierr); 28 ierr = PetscOptionsString("-tagname", "The tag name from which to create a vector", "ex1.cxx", options->tagname, options->tagname, sizeof(options->tagname), &flg);CHKERRQ(ierr); 29 ierr = PetscOptionsEnd();CHKERRQ(ierr); 30 31 ierr = PetscLogEventRegister("CreateMesh", DM_CLASSID, &options->createMeshEvent);CHKERRQ(ierr); 32 PetscFunctionReturn(0); 33 } 34 35 PetscErrorCode CreateMesh(MPI_Comm comm, AppCtx *user, DM *dm) 36 { 37 moab::Interface *iface=NULL; 38 moab::Tag tag=NULL; 39 moab::Tag ltog_tag=NULL; 40 moab::Range range; 41 PetscInt tagsize; 42 moab::ErrorCode merr; 43 PetscErrorCode ierr; 44 45 PetscFunctionBegin; 46 ierr = PetscLogEventBegin(user->createMeshEvent,0,0,0,0);CHKERRQ(ierr); 47 ierr = DMMoabCreateMoab(comm, iface, <og_tag, &range, dm);CHKERRQ(ierr); 48 std::cout << "Created DMMoab using DMMoabCreateMoab." << std::endl; 49 ierr = DMMoabGetInterface(*dm, &iface);CHKERRQ(ierr); 50 51 // load file and get entities of requested or max dimension 52 if (strlen(user->filename) > 0) { 53 merr = iface->load_file(user->filename);MBERRNM(merr); 54 std::cout << "Read mesh from file " << user->filename << std::endl; 55 } 56 else { 57 // make a simple structured mesh 58 moab::ScdInterface *scdi; 59 merr = iface->query_interface(scdi); 60 61 moab::ScdBox *box; 62 merr = scdi->construct_box (moab::HomCoord(0,0,0), moab::HomCoord(5,5,5), NULL, 0, box);MBERRNM(merr); 63 user->dim = 3; 64 merr = iface->set_dimension(user->dim);MBERRNM(merr); 65 std::cout << "Created structured 5x5x5 mesh." << std::endl; 66 } 67 if (-1 == user->dim) { 68 moab::Range tmp_range; 69 merr = iface->get_entities_by_handle(0, tmp_range);MBERRNM(merr); 70 if (tmp_range.empty()) { 71 MBERRNM(moab::MB_FAILURE); 72 } 73 user->dim = iface->dimension_from_handle(*tmp_range.rbegin()); 74 } 75 merr = iface->get_entities_by_dimension(0, user->dim, range);MBERRNM(merr); 76 ierr = DMMoabSetLocalVertices(*dm, &range);CHKERRQ(ierr); 77 78 // get the requested tag and create if necessary 79 std::cout << "Creating tag with name: " << user->tagname << ";\n"; 80 merr = iface->tag_get_handle(user->tagname, 1, moab::MB_TYPE_DOUBLE, tag, moab::MB_TAG_CREAT | moab::MB_TAG_DENSE);MBERRNM(merr); 81 { 82 // initialize new tag with gids 83 std::vector<double> tag_vals(range.size()); 84 merr = iface->tag_get_data(ltog_tag, range, tag_vals.data());MBERRNM(merr); // read them into ints 85 double *dval = tag_vals.data(); int *ival = reinterpret_cast<int*>(dval); // "stretch" them into doubles, from the end 86 for (int i = tag_vals.size()-1; i >= 0; i--) dval[i] = ival[i]; 87 merr = iface->tag_set_data(tag, range, tag_vals.data());MBERRNM(merr); // write them into doubles 88 } 89 merr = iface->tag_get_length(tag, tagsize);MBERRNM(merr); 90 91 ierr = DMSetUp(*dm);CHKERRQ(ierr); 92 93 // create the dmmoab and initialize its data 94 ierr = PetscObjectSetName((PetscObject) *dm, "MOAB mesh");CHKERRQ(ierr); 95 ierr = PetscLogEventEnd(user->createMeshEvent,0,0,0,0);CHKERRQ(ierr); 96 user->dm = *dm; 97 PetscFunctionReturn(0); 98 } 99 100 int main(int argc, char **argv) 101 { 102 AppCtx user; /* user-defined work context */ 103 moab::ErrorCode merr; 104 PetscErrorCode ierr; 105 Vec vec; 106 moab::Interface* mbImpl=NULL; 107 moab::Tag datatag=NULL; 108 109 ierr = PetscInitialize(&argc, &argv, NULL,help);if (ierr) return ierr; 110 ierr = ProcessOptions(PETSC_COMM_WORLD, &user);CHKERRQ(ierr); 111 112 ierr = CreateMesh(PETSC_COMM_WORLD, &user, &user.dm);CHKERRQ(ierr); /* create the MOAB dm and the mesh */ 113 114 ierr = DMMoabGetInterface(user.dm, &mbImpl);CHKERRQ(ierr); 115 merr = mbImpl->tag_get_handle(user.tagname, datatag);MBERRNM(merr); 116 ierr = DMMoabCreateVector(user.dm, datatag, NULL, PETSC_TRUE, PETSC_FALSE,&vec);CHKERRQ(ierr); /* create a vec from user-input tag */ 117 118 std::cout << "Created VecMoab from existing tag." << std::endl; 119 ierr = VecDestroy(&vec);CHKERRQ(ierr); 120 std::cout << "Destroyed VecMoab." << std::endl; 121 ierr = DMDestroy(&user.dm);CHKERRQ(ierr); 122 std::cout << "Destroyed DMMoab." << std::endl; 123 ierr = PetscFinalize(); 124 return ierr; 125 } 126 127 /*TEST 128 129 build: 130 requires: moab 131 132 test: 133 134 TEST*/ 135