16                                 {
   17 
   19 
   20  try {
   21 
   24    int surface_side_set = 200;
   25    PetscBool flg_vol_block_set;
   26    int vol_block_set = 1;
   29    PetscBool flg_shift;
   30    double shift[] = {0, 0, 0};
   31    int nmax = 3;
   32    int fraction_level = 2;
   33    PetscBool squash_bits = PETSC_TRUE;
   34    PetscBool set_coords = PETSC_TRUE;
   36    int create_surface_side_set = 201;
   37    PetscBool flg_create_surface_side_set;
   38    int nb_ref_cut = 0;
   39    int nb_ref_trim = 0;
   40    PetscBool flg_tol;
   41    double tol[] = {1e-2, 2e-1, 2e-1};
 
   42    int nmax_tol = 3;
   43    PetscBool 
debug = PETSC_FALSE;
 
   44 
   45    PetscOptionsBegin(PETSC_COMM_WORLD, "", "Mesh cut options", "none");
   46    CHKERR PetscOptionsString(
"-my_file", 
"mesh file name", 
"", 
"mesh.h5m",
 
   48    CHKERR PetscOptionsInt(
"-surface_side_set", 
"surface side set", 
"",
 
   49                           surface_side_set, &surface_side_set, PETSC_NULLPTR);
   50    CHKERR PetscOptionsInt(
"-vol_block_set", 
"volume side set", 
"",
 
   51                           vol_block_set, &vol_block_set, &flg_vol_block_set);
   52    CHKERR PetscOptionsInt(
"-edges_block_set", 
"edges block set", 
"",
 
   54    CHKERR PetscOptionsInt(
"-vertex_block_set", 
"vertex block set", 
"",
 
   56    CHKERR PetscOptionsRealArray(
"-shift", 
"shift surface by vector", 
"", shift,
 
   57                                 &nmax, &flg_shift);
   58    CHKERR PetscOptionsInt(
"-fraction_level", 
"fraction of merges merged", 
"",
 
   59                           fraction_level, &fraction_level, PETSC_NULLPTR);
   60    CHKERR PetscOptionsBool(
"-squash_bits", 
"true to squash bits at the end",
 
   61                            "", squash_bits, &squash_bits, PETSC_NULLPTR);
   62    CHKERR PetscOptionsBool(
"-set_coords", 
"true to set coords at the end", 
"",
 
   63                            set_coords, &set_coords, PETSC_NULLPTR);
   64    CHKERR PetscOptionsBool(
"-output_vtk", 
"if true outout vtk file", 
"",
 
   66    CHKERR PetscOptionsInt(
"-create_side_set", 
"crete side set", 
"",
 
   67                           create_surface_side_set, &create_surface_side_set,
   68                           &flg_create_surface_side_set);
   69    CHKERR PetscOptionsInt(
"-nb_ref_cut", 
"nb refinements befor cut", 
"",
 
   70                           nb_ref_cut, &nb_ref_cut, PETSC_NULLPTR);
   71    CHKERR PetscOptionsInt(
"-nb_ref_trim", 
"nb refinements befor trim", 
"",
 
   72                           nb_ref_trim, &nb_ref_trim, PETSC_NULLPTR);
   73    CHKERR PetscOptionsRealArray(
 
   74        "-tol", "tolerances tolCut, tolCutClose, tolTrim, tolTrimClose", "",
   75        tol, &nmax_tol, &flg_tol);
 
   76    CHKERR PetscOptionsBool(
"-debug", 
"debug (produces many VTK files)", 
"",
 
   78    PetscOptionsEnd();
   79 
   82              "*** ERROR -my_file (MESH FILE NEEDED)");
   83    }
   84    if (flg_shift && nmax != 3) {
   86    }
   87 
   88    if (flg_tol && nmax_tol != 3)
   90 
   91    moab::Core mb_instance;
   92    moab::Interface &moab = mb_instance;
   93    ParallelComm *pcomm = ParallelComm::get_pcomm(&moab, 
MYPCOMM_INDEX);
 
   94    if (pcomm == NULL)
   95      pcomm = new ParallelComm(&moab, PETSC_COMM_WORLD);
   96 
   97    const char *option;
   98    option = "";
  100 
  104 
  105    
  108    
  111    
  114 
  117      cout << *it << endl;
  118    }
  119 
  122 
  123    
  128 
  131 
  132    
  133    
  134    
  137    else
  139               "Side set not found %d", surface_side_set);
  140 
  141    
  142    Range fixed_edges, corner_nodes;
 
  145                                                     1, fixed_edges, true);
  148                                                     1, fixed_edges, true);
  152    
  153 
  155    if (flg_vol_block_set) {
  158                                                       3, tets, true);
  159      else
  161                 "Block set %d not found", vol_block_set);
  162    } else
  163      CHKERR moab.get_entities_by_dimension(0, 3, tets, 
false);
 
  164 
  167 
  168    
  169    
  171 
  172    
  175    auto shift_after_ref = [&]() {
  178      mask.set(0);
  179      for (int ll = 1; ll != nb_ref_cut + nb_ref_trim + 1; ++ll)
  180        mask.set(ll);
  182          nb_ref_cut + nb_ref_trim, mask, 
VERBOSE);
 
  184    };
  186 
  187    
  188    double def_position[] = {0, 0, 0};
  190    CHKERR moab.tag_get_handle(
"POSITION", 3, MB_TYPE_DOUBLE, 
th,
 
  191                               MB_TAG_CREAT | MB_TAG_SPARSE, def_position);
  192    
  194 
  195    
  196    int first_bit = 1;
  198                                     tol[1], 
tol[2], fixed_edges, corner_nodes,
 
  199                                     true, false);
  200 
  201    
  202    if (set_coords)
  204 
  205    
  206    if (flg_vol_block_set) {
  211    }
  212 
  213    
  214    if (squash_bits) {
  216      for (int ll = 0; ll != first_bit; ++ll)
  217        shift_mask.set(ll);
  219          first_bit - 1, shift_mask, 
VERBOSE);
 
  220    }
  221 
  225    CHKERR moab.get_adjacencies(cut_mesh->
getSurface(), 1, 
false, surface_edges,
 
  226                                moab::Interface::UNION);
  228    CHKERR moab.delete_entities(surface_edges);
 
  229    CHKERR moab.delete_entities(surface_verts);
 
  230 
  231    if (flg_create_surface_side_set) {
  232      
  234              create_surface_side_set, 
SIDESET))
 
  236      else
  238                  "Warring >>> sideset %d is on the mesh",
  239                  create_surface_side_set);
  240 
  243    }
  244 
  245    CHKERR moab.write_file(
"out.h5m");
 
  246 
  249      CHKERR moab.create_meshset(MESHSET_SET, meshset);
 
  250      if (flg_vol_block_set) {
  253                                                ents, true);
  254        CHKERR moab.add_entities(meshset, ents);
 
  255      } else {
  259      }
  261      CHKERR moab.write_file(
"out.vtk", 
"VTK", 
"", &meshset, 1);
 
  262      CHKERR moab.delete_entities(&meshset, 1);
 
  263    }
  264  }
  266 
  268 
  269  return 0;
  270}
#define MOFEM_LOG_C(channel, severity, format,...)
#define CATCH_ERRORS
Catch errors.
#define MYPCOMM_INDEX
default communicator number PCOMM
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
@ MOFEM_ATOM_TEST_INVALID
@ MOFEM_DATA_INCONSISTENCY
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
#define CHKERR
Inline error check.
MoFEMErrorCode setBitRefLevelByType(const EntityHandle meshset, const EntityType type, const BitRefLevel bit, int verb=QUIET) const
Set the Bit Ref Level By Type object.
MoFEMErrorCode getEntitiesByTypeAndRefLevel(const BitRefLevel bit, const BitRefLevel mask, const EntityType type, const EntityHandle meshset, int verb=0) const
add all ents from ref level given by bit to meshset
MoFEMErrorCode getMeshset(const int ms_id, const unsigned int cubit_bc_type, EntityHandle &meshset) const
get meshset from CUBIT Id and CUBIT type
MoFEMErrorCode addEntitiesToMeshset(const CubitBCType cubit_bc_type, const int ms_id, const Range &ents)
Add entities to CUBIT meshset.
bool checkMeshset(const int ms_id, const CubitBCType cubit_bc_type) const
Check for CUBIT meshset by ID and type.
MoFEMErrorCode addMeshset(const CubitBCType cubit_bc_type, const int ms_id, const std::string name="")
Add CUBIT meshset to manager.
MoFEMErrorCode getEntitiesByDimension(const int ms_id, const unsigned int cubit_bc_type, const int dimension, Range &entities, const bool recursive=true) const
get entities from CUBIT/meshset of a particular entity dimension
#define _IT_CUBITMESHSETS_BY_SET_TYPE_FOR_LOOP_(MESHSET_MANAGER, CUBITBCTYPE, IT)
Iterator that loops over a specific Cubit MeshSet having a particular BC meshset in a moFEM field.
std::bitset< BITREFLEVEL_SIZE > BitRefLevel
Bit structure attached to each entity identifying to what mesh entity is attached.
static MoFEMErrorCode Initialize(int *argc, char ***args, const char file[], const char help[])
Initializes the MoFEM database PETSc, MOAB and MPI.
static MoFEMErrorCode Finalize()
Checks for options to be called at the conclusion of the program.
MoFEMErrorCode copySurface(const Range surface, Tag th=NULL, double *shift=NULL, double *origin=NULL, double *transform=NULL, const std::string save_mesh="")
copy surface entities
MoFEMErrorCode setCoords(Tag th, const BitRefLevel bit=BitRefLevel(), const BitRefLevel mask=BitRefLevel().set())
set coords from tag
const Range & getMergedSurfaces() const
MoFEMErrorCode buildTree()
build tree
const Range & getSurface() const
MoFEMErrorCode setVolume(const Range volume)
set volume entities
MoFEMErrorCode setTagData(Tag th, const BitRefLevel bit=BitRefLevel())
set coords to tag
MoFEMErrorCode cutTrimAndMerge(int &first_bit, const int fraction_level, Tag th, const double tol_geometry, const double tol_cut_close, const double tol_trim_close, Range &fixed_edges, Range &corner_nodes, const bool update_meshsets=false, const bool debug=false)
Cut, trim and merge.
MoFEMErrorCode makeFront(const bool debug=false)
Create front from the surface.
MoFEMErrorCode refineMesh(const int init_bit_level, const int surf_levels, const int front_levels, Range *fixed_edges=nullptr, int verb=QUIET, const bool debug=false)
Refine and set level sets.
Interface for managing meshsets containing materials and boundary conditions.
MoFEMErrorCode getInterface(IFACE *&iface) const
Get interface reference to pointer of interface.