v0.9.0
Functions | Variables
uniform_mesh_refinement.cpp File Reference

Uniformly refine mesh. More...

#include <MoFEM.hpp>

Go to the source code of this file.

Functions

int main (int argc, char *argv[])
 

Variables

static char help []
 

Detailed Description

Uniformly refine mesh.

Definition in file uniform_mesh_refinement.cpp.

Function Documentation

◆ main()

int main ( int  argc,
char *  argv[] 
)

Definition at line 33 of file uniform_mesh_refinement.cpp.

33  {
34 
35  MoFEM::Core::Initialize(&argc, &argv, (char *)0, help);
36 
37  try {
38 
39  // global variables
40  char mesh_file_name[255];
41  char mesh_out_file[255] = "out.h5m";
42  PetscBool flg_file = PETSC_FALSE;
43 
44  ierr = PetscOptionsBegin(PETSC_COMM_WORLD, "", "none", "none");
45  CHKERRQ(ierr);
46 
47  CHKERR PetscOptionsString("-my_file", "mesh file name", "", "mesh.h5m",
48  mesh_file_name, 255, &flg_file);
49  CHKERR PetscOptionsString("-output_file", "output mesh file name", "",
50  "mesh.h5m", mesh_out_file, 255, PETSC_NULL);
51 
52  ierr = PetscOptionsEnd();
53  CHKERRQ(ierr);
54 
55  moab::Core mb_instance;
56  moab::Interface &moab = mb_instance;
57  ParallelComm *pcomm = ParallelComm::get_pcomm(&moab, MYPCOMM_INDEX);
58  if (pcomm == NULL)
59  pcomm = new ParallelComm(&moab, PETSC_COMM_WORLD);
60 
61  const char *option;
62  option = "";
63  CHKERR moab.load_file(mesh_file_name, 0, option);
64 
65  // Create MoFEM database
66  MoFEM::Core core(moab);
67  MoFEM::Interface &m_field = core;
68 
69  if (flg_file != PETSC_TRUE) {
70  SETERRQ(PETSC_COMM_SELF, 1, "*** ERROR -my_file (MESH FILE NEEDED)");
71  }
72 
73  BitRefManager *bit_ref_manager;
74  CHKERR m_field.getInterface(bit_ref_manager);
75 
76  CHKERR bit_ref_manager->setBitRefLevelByDim(0, 3, BitRefLevel().set(0));
77 
78  Range ents3d;
79  rval = moab.get_entities_by_dimension(0, 3, ents3d, false);
80  Range edges;
81  CHKERR moab.get_adjacencies(ents3d, 1, true, edges, moab::Interface::UNION);
82 
83  EntityHandle meshset_ref_edges;
84  CHKERR moab.create_meshset(MESHSET_SET, meshset_ref_edges);
85  CHKERR moab.add_entities(meshset_ref_edges, edges);
86 
87  MeshRefinement *refine;
88  CHKERR m_field.getInterface(refine);
89 
90  CHKERR refine->add_vertices_in_the_middle_of_edges(meshset_ref_edges,
91  BitRefLevel().set(1));
92  CHKERR refine->refine_TET(0, BitRefLevel().set(1));
93 
94  // update cubit meshsets
95  for (_IT_CUBITMESHSETS_FOR_LOOP_(m_field, ciit)) {
96  EntityHandle cubit_meshset = ciit->meshset;
97  CHKERR bit_ref_manager->updateMeshsetByEntitiesChildren(
98  cubit_meshset, BitRefLevel().set(1), cubit_meshset, MBVERTEX, true);
99  CHKERR bit_ref_manager->updateMeshsetByEntitiesChildren(
100  cubit_meshset, BitRefLevel().set(1), cubit_meshset, MBEDGE, true);
101  CHKERR bit_ref_manager->updateMeshsetByEntitiesChildren(
102  cubit_meshset, BitRefLevel().set(1), cubit_meshset, MBTRI, true);
103  CHKERR bit_ref_manager->updateMeshsetByEntitiesChildren(
104  cubit_meshset, BitRefLevel().set(1), cubit_meshset, MBTET, true);
105  }
106 
107  CHKERR core.getInterface<BitRefManager>()->shiftRightBitRef(
108  1, BitRefLevel().set(), VERBOSE);
109 
110  CHKERR moab.delete_entities(&meshset_ref_edges, 1);
111  CHKERR moab.write_file(mesh_out_file);
112  }
113  CATCH_ERRORS;
114 
116  CHKERRQ(ierr);
117 
118  return 0;
119 }
virtual MoFEMErrorCode refine_TET(const EntityHandle meshset, const BitRefLevel &bit, const bool respect_interface=false, int verb=QUIET, Range *ref_edges=NULL, const bool debug=false)
refine TET in the meshset
virtual MoFEMErrorCode add_vertices_in_the_middle_of_edges(const EntityHandle meshset, const BitRefLevel &bit, const bool recursive=false, int verb=QUIET, EntityHandle start_v=0)
make vertices in the middle of edges in meshset and add them to refinement levels defined by bit
static char help[]
Core (interface) class.
Definition: Core.hpp:50
char mesh_file_name[255]
static MoFEMErrorCode Initialize(int *argc, char ***args, const char file[], const char help[])
Initializes the MoFEM database PETSc, MOAB and MPI.
Definition: Core.cpp:51
static MoFEMErrorCodeGeneric< moab::ErrorCode > rval
Definition: Exceptions.hpp:84
MoFEMErrorCode getInterface(const MOFEMuuid &uuid, IFACE *&iface) const
Get interface by uuid and return reference to pointer of interface.
CHKERRQ(ierr)
MoFEMErrorCode updateMeshsetByEntitiesChildren(const EntityHandle parent, const BitRefLevel &parent_bit, const BitRefLevel &parent_mask, const BitRefLevel &child_bit, const BitRefLevel &child_mask, const EntityHandle child, EntityType child_type, const bool recursive=false, int verb=0)
Get child entities form meshset containing parent entitiesSearch for refined entities of given type w...
Managing BitRefLevels.
Mesh refinement interface.
#define CHKERR
Inline error check.
Definition: definitions.h:596
#define MYPCOMM_INDEX
default communicator number PCOMM
Definition: definitions.h:285
std::bitset< BITREFLEVEL_SIZE > BitRefLevel
Bit structure attached to each entity identifying to what mesh entity is attached.
Definition: Types.hpp:51
#define _IT_CUBITMESHSETS_FOR_LOOP_(MESHSET_MANAGER, IT)
Iterator that loops over all the Cubit MeshSets in a moFEM field.
#define CATCH_ERRORS
Catch errors.
Definition: definitions.h:433
MoFEMErrorCode setBitRefLevelByDim(const EntityHandle meshset, const int dim, const BitRefLevel bit, int verb=QUIET) const
Set the Bit Ref Level By Dim object.
static MoFEMErrorCode Finalize()
Checks for options to be called at the conclusion of the program.
Definition: Core.cpp:61

Variable Documentation

◆ help

char help[]
static
Initial value:
=
"Uniform mesh refinement\n\n"
"Usage example:\n"
"$ ./uniform_mesh_refinement -my_file mesh.h5m -output_file refined_mesh.h5m"
"\n\n"

Definition at line 24 of file uniform_mesh_refinement.cpp.