v0.13.1
Public Member Functions | Public Attributes | Private Member Functions | Private Attributes | List of all members
FractureMechanics::CPMeshCut Struct Reference

#include <users_modules/fracture_mechanics/src/CPMeshCut.hpp>

Inheritance diagram for FractureMechanics::CPMeshCut:
[legend]
Collaboration diagram for FractureMechanics::CPMeshCut:
[legend]

Public Member Functions

MoFEMErrorCode query_interface (boost::typeindex::type_index type_index, MoFEM::UnknownInterface **iface) const
 Getting interface of core database. More...
 
 CPMeshCut (CrackPropagation &cp)
 
virtual ~CPMeshCut ()
 
MoFEMErrorCode getOptions ()
 Get options from command line. More...
 
MoFEMErrorCode setVolume (const BitRefLevel &bit)
 Set cutting volume from bit ref level. More...
 
MoFEMErrorCode refineMesh (const Range *vol, const bool make_front, const int verb=QUIET, const bool debug=false)
 Refine mesh close to crack surface and crack front. More...
 
MoFEMErrorCode setCutSurfaceFromFile ()
 Set the cut surface from file. More...
 
MoFEMErrorCode copySurface (const std::string save_mesh="")
 
MoFEMErrorCode rebuildCrackSurface (const double factor, const std::string save_mesh="", const int verb=QUIET, bool debug=false)
 
MoFEMErrorCode cutMesh (int &first_bit, const bool debug=false)
 Ref, cut and trim, merge and do TetGen, if option is on. More...
 
const RangegetCuttingSurface () const
 
MoFEMErrorCode findBodySkin (const BitRefLevel bit, const int verb=1, const bool debug=false, std::string file_name="out_body_skin.vtk")
 find body skin More...
 
MoFEMErrorCode findCrack (const BitRefLevel bit, const int verb=1, const bool debug=false)
 get crack front More...
 
MoFEMErrorCode findCrackFromPrisms (const BitRefLevel bit, const int verb=1, const bool debug=false)
 get crack front More...
 
MoFEMErrorCode findContactFromPrisms (const BitRefLevel bit, const int verb=1, const bool debug=false)
 get contact elements More...
 
MoFEMErrorCode splitFaces (const int verb=1, const bool debug=false)
 split crack faces More...
 
MoFEMErrorCode insertContactInterface (const int verb=1, const bool debug=false)
 insert contact interface More...
 
MoFEMErrorCode addMortarContactPrisms (const int verb=1, const bool debug=false)
 insert contact interface More...
 
MoFEMErrorCode refineCrackTip (const int front_id, const int verb=1, const bool debug=false)
 refine elements at crack tip More...
 
MoFEMErrorCode refineAndSplit (const int verb=1, const bool debug=false)
 
MoFEMErrorCode cutRefineAndSplit (const int verb=VERBOSE, bool debug=false)
 Refine, cut, trim, merge and ref front and split. More...
 
MoFEMErrorCode getFrontEdgesAndElements (const BitRefLevel bit, const int verb=1, const bool debug=false)
 get crack front edges and finite elements More...
 
MoFEMErrorCode setFixEdgesAndCorners (const BitRefLevel bit)
 
std::vector< double > & getOrgCoords ()
 
const std::vector< double > & getOrgCoords () const
 
const RangegetFixedEdges () const
 
const RangegetCornerNodes () const
 
const std::string & getCutSurfMeshName () const
 
const int getEdgesBlockSet () const
 
MoFEMErrorCode getInterfacesPtr ()
 
int getCuttingSurfaceSidesetId () const
 
int getSkinOfTheBodyId () const
 
int getCrackSurfaceId () const
 
int getCrackFrontId () const
 
int getContactPrismsBlockSetId () const
 
MoFEMErrorCode clearData ()
 
- Public Member Functions inherited from MoFEM::UnknownInterface
virtual MoFEMErrorCode query_interface (boost::typeindex::type_index type_index, UnknownInterface **iface) const =0
 
template<class IFACE >
MoFEMErrorCode registerInterface (bool error_if_registration_failed=true)
 Register interface. More...
 
template<class IFACE >
MoFEMErrorCode getInterface (IFACE *&iface) const
 Get interface refernce to pointer of interface. More...
 
template<class IFACE >
MoFEMErrorCode getInterface (IFACE **const iface) const
 Get interface pointer to pointer of interface. More...
 
template<class IFACE , typename boost::enable_if< boost::is_pointer< IFACE >, int >::type = 0>
IFACE getInterface () const
 Get interface pointer to pointer of interface. More...
 
template<class IFACE , typename boost::enable_if< boost::is_reference< IFACE >, int >::type = 0>
IFACE getInterface () const
 Get reference to interface. More...
 
template<class IFACE >
IFACE * getInterface () const
 Function returning pointer to interface. More...
 
virtual ~UnknownInterface ()=default
 

Public Attributes

CrackPropagationcP
 

Private Member Functions

MoFEMErrorCode setMeshOrgCoords ()
 
MoFEMErrorCode getMeshOrgCoords ()
 

Private Attributes

PetscBool edgesBlockSetFlg
 
int edgesBlockSet
 
int vertexBlockSet
 
int vertexNodeSet
 
int crackedBodyBlockSetId
 
int contactPrismsBlockSetId
 
CutMeshInterface * cutMeshPtr
 
BitRefManager * bitRefPtr
 
MeshsetsManager * meshsetMngPtr
 
Range fixedEdges
 
Range cornerNodes
 
Range fRont
 
Range sUrface
 
int fractionLevel
 
double tolCut
 
double tolCutClose
 
double tolTrimClose
 
int nbRefBeforCut
 
int nbRefBeforTrim
 
PetscBool debugCut
 
PetscBool removePathologicalFrontTris
 
double sHift [3]
 
int cuttingSurfaceSidesetId
 
int skinOfTheBodyId
 
int crackSurfaceId
 
int crackFrontId
 
std::vector< doubleoriginalCoords
 
std::string cutSurfMeshName
 

Additional Inherited Members

- Static Public Member Functions inherited from MoFEM::UnknownInterface
static MoFEMErrorCode getLibVersion (Version &version)
 Get library version. More...
 
static MoFEMErrorCode getFileVersion (moab::Interface &moab, Version &version)
 Get database major version. More...
 
static MoFEMErrorCode setFileVersion (moab::Interface &moab, Version version=Version(MoFEM_VERSION_MAJOR, MoFEM_VERSION_MINOR, MoFEM_VERSION_BUILD))
 Get database major version. More...
 
static MoFEMErrorCode getInterfaceVersion (Version &version)
 Get database major version. More...
 

Detailed Description

Definition at line 24 of file CPMeshCut.hpp.

Constructor & Destructor Documentation

◆ CPMeshCut()

FractureMechanics::CPMeshCut::CPMeshCut ( CrackPropagation cp)

Definition at line 51 of file CPMeshCut.cpp.

51 : cP(cp) {
53 CHKERRABORT(PETSC_COMM_SELF, ierr);
54
55 try {
56 MoFEM::LogManager::getLog("CPMeshCutSelf");
57 } catch (...) {
58 auto core_log = logging::core::get();
59 core_log->add_sink(
60 LogManager::createSink(LogManager::getStrmSelf(), "CPMeshCutSelf"));
61 LogManager::setLog("CPMeshCutSelf");
62 core_log->add_sink(
63 LogManager::createSink(LogManager::getStrmWorld(), "CPMeshCutWorld"));
64 LogManager::setLog("CPMeshCutWorld");
65 core_log->add_sink(
66 LogManager::createSink(LogManager::getStrmSync(), "CPMeshCutSync"));
67 LogManager::setLog("CPMeshCutSync");
68 MOFEM_LOG_TAG("CPMeshCutSelf", "CPMeshCut");
69 MOFEM_LOG_TAG("CPMeshCutWorld", "CPMeshCut");
70 MOFEM_LOG_TAG("CPMeshCutSync", "CPMeshCut");
71 }
72
74 MOFEM_LOG("CPMeshCutWorld", Sev::noisy) << "CPMeshCutSelf interface created";
75}
#define MOFEM_LOG(channel, severity)
Log.
Definition: LogManager.hpp:301
#define MOFEM_LOG_TAG(channel, tag)
Tag channel.
Definition: LogManager.hpp:332
static LoggerType & getLog(const std::string channel)
Get logger by channel.
Definition: LogManager.cpp:370
#define MOFEM_LOG_FUNCTION()
Set scope.
Definition: LogManager.hpp:318
static MoFEMErrorCodeGeneric< PetscErrorCode > ierr
Definition: Exceptions.hpp:76
MoFEMErrorCode getInterfacesPtr()
Definition: CPMeshCut.cpp:77
CrackPropagation & cP
Definition: CPMeshCut.hpp:35

◆ ~CPMeshCut()

virtual FractureMechanics::CPMeshCut::~CPMeshCut ( )
virtual

Definition at line 37 of file CPMeshCut.hpp.

37{}

Member Function Documentation

◆ addMortarContactPrisms()

MoFEMErrorCode FractureMechanics::CPMeshCut::addMortarContactPrisms ( const int  verb = 1,
const bool  debug = false 
)

insert contact interface

Definition at line 1749 of file CPMeshCut.cpp.

1750 {
1751 MoFEM::Interface &m_field = cP.mField;
1752 moab::Interface &moab = m_field.get_moab();
1754
1755 Range integration_triangles;
1758 cP.mortarContactElements.clear();
1759
1760 Range all_masters, all_slaves;
1761
1763 m_field, all_masters, all_slaves);
1764
1766 BitRefLevel().set(), all_slaves);
1767
1769 BitRefLevel().set(), all_masters);
1770
1771 if (all_slaves.empty() || all_masters.empty())
1773
1774 ContactSearchKdTree contact_search_kd_tree(m_field);
1775
1777 boost::shared_ptr<ContactSearchKdTree::ContactCommonData_multiIndex>(
1779
1780 // create kd_tree with master_surface only
1781 CHKERR contact_search_kd_tree.buildTree(all_masters);
1782
1783 CHKERR contact_search_kd_tree.contactSearchAlgorithm(
1784 all_masters, all_slaves, cP.contactSearchMultiIndexPtr,
1786
1787 Range tris_from_prism;
1788 // find actual masters and slave used
1789 CHKERR m_field.get_moab().get_adjacencies(cP.mortarContactElements, 2, true,
1790 tris_from_prism,
1791 moab::Interface::UNION);
1792
1793 tris_from_prism = tris_from_prism.subset_by_type(MBTRI);
1794 cP.mortarContactMasterFaces = intersect(tris_from_prism, all_masters);
1795 cP.mortarContactSlaveFaces = intersect(tris_from_prism, all_slaves);
1796
1797 EntityHandle meshset_slave_master_prisms;
1798 CHKERR moab.create_meshset(MESHSET_SET, meshset_slave_master_prisms);
1799
1800 CHKERR
1801 moab.add_entities(meshset_slave_master_prisms, cP.mortarContactElements);
1802
1803 CHKERR m_field.getInterface<BitRefManager>()->setBitRefLevelByDim(
1804 meshset_slave_master_prisms, 3, cP.mapBitLevel["spatial_domain"]);
1805
1807 ContactSearchKdTree::Prism_tag>::type::iterator ItMultIndexPrism;
1808
1809 for (Range::iterator it_tri = cP.mortarContactElements.begin();
1810 it_tri != cP.mortarContactElements.end(); ++it_tri) {
1811
1812 ItMultIndexPrism it_mult_index_prism =
1814 .find(*it_tri);
1815
1816 Range range_poly_tris;
1817 range_poly_tris.clear();
1818 range_poly_tris = it_mult_index_prism->get()->commonIntegratedTriangle;
1819 integration_triangles.insert(range_poly_tris.begin(),
1820 range_poly_tris.end());
1821 }
1822
1823 EntityHandle ent_integration_vertex;
1824 CHKERR moab.create_meshset(MESHSET_SET, ent_integration_vertex);
1825
1826 EntityHandle ent_integration_edge;
1827 CHKERR moab.create_meshset(MESHSET_SET, ent_integration_edge);
1828
1829 Range edges, verts;
1830 CHKERR m_field.get_moab().get_adjacencies(integration_triangles, 1, true,
1831 edges, moab::Interface::UNION);
1832 CHKERR moab.add_entities(ent_integration_edge, edges);
1833
1834 CHKERR m_field.get_moab().get_connectivity(integration_triangles, verts,
1835 true);
1836 CHKERR moab.add_entities(ent_integration_vertex, verts);
1837
1838 EntityHandle ent_integration_tris;
1839 CHKERR moab.create_meshset(MESHSET_SET, ent_integration_tris);
1840 CHKERR moab.add_entities(ent_integration_tris, integration_triangles);
1841
1842 CHKERR m_field.getInterface<BitRefManager>()->setBitRefLevelByDim(
1843 ent_integration_vertex, 0, cP.mapBitLevel["spatial_domain"]);
1844 CHKERR m_field.getInterface<BitRefManager>()->setBitRefLevelByDim(
1845 ent_integration_edge, 1, cP.mapBitLevel["spatial_domain"]);
1846 CHKERR m_field.getInterface<BitRefManager>()->setBitRefLevelByDim(
1847 ent_integration_tris, 2, cP.mapBitLevel["spatial_domain"]);
1848
1849 CHKERR moab.delete_entities(&ent_integration_tris, 1);
1850
1851 const std::string output_name = "slave_master_prisms.vtk";
1852 CHKERR moab.write_mesh(output_name.c_str(), &meshset_slave_master_prisms, 1);
1853
1854 // get integration tris of each prism
1855 CHKERR moab.delete_entities(&meshset_slave_master_prisms, 1);
1856
1858}
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:447
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:346
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:416
#define CHKERR
Inline error check.
Definition: definitions.h:535
std::bitset< BITREFLEVEL_SIZE > BitRefLevel
Bit structure attached to each entity identifying to what mesh entity is attached.
Definition: Types.hpp:40
DeprecatedCoreInterface Interface
Definition: Interface.hpp:1955
multi_index_container< boost::shared_ptr< ContactCommonData >, indexed_by< ordered_unique< tag< Prism_tag >, member< ContactCommonData, EntityHandle, &ContactCommonData::pRism > > > > ContactCommonData_multiIndex
std::map< std::string, BitRefLevel > mapBitLevel
boost::shared_ptr< ContactSearchKdTree::ContactCommonData_multiIndex > contactSearchMultiIndexPtr
Managing BitRefLevels.
MoFEMErrorCode filterEntitiesByRefLevel(const BitRefLevel bit, const BitRefLevel mask, Range &ents, int verb=QUIET) const
filter entities by bit ref level
virtual moab::Interface & get_moab()=0
Deprecated interface functions.
MoFEMErrorCode getInterface(IFACE *&iface) const
Get interface refernce to pointer of interface.
static MoFEMErrorCode masterSlaveTrianglesCreation(MoFEM::Interface &m_field, Range &range_surf_master, Range &range_surf_slave, int step=1)

◆ clearData()

MoFEMErrorCode FractureMechanics::CPMeshCut::clearData ( )

Definition at line 2446 of file CPMeshCut.cpp.

2446 {
2448
2449 fixedEdges.clear();
2450 cornerNodes.clear();
2451 fRont.clear();
2452 sUrface.clear();
2453 originalCoords.clear();
2454
2456};
std::vector< double > originalCoords
Definition: CPMeshCut.hpp:210

◆ copySurface()

PetscErrorCode FractureMechanics::CPMeshCut::copySurface ( const std::string  save_mesh = "")

Definition at line 322 of file CPMeshCut.cpp.

322 {
324
325 sUrface.clear();
326 fRont.clear();
327
329 PetscPrintf(PETSC_COMM_WORLD,
330 "Sideset %d with cutting crack surface not found\n",
333 }
335 2, sUrface, true);
336
337 // Copy surface
338 const_cast<Range &>(cutMeshPtr->getSurface()).clear();
339 CHKERR cutMeshPtr->copySurface(sUrface, NULL, sHift, NULL, NULL, save_mesh);
340
342}
@ SIDESET
Definition: definitions.h:147
bool checkMeshset(const int ms_id, const CubitBCType cubit_bc_type) const
check for CUBIT Id and CUBIT type
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
MeshsetsManager * meshsetMngPtr
Definition: CPMeshCut.hpp:187
CutMeshInterface * cutMeshPtr
Definition: CPMeshCut.hpp:185
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
const Range & getSurface() const

◆ cutMesh()

MoFEMErrorCode FractureMechanics::CPMeshCut::cutMesh ( int &  first_bit,
const bool  debug = false 
)

Ref, cut and trim, merge and do TetGen, if option is on.

This is running cut mesh interface function

Parameters
first_bit
debug
Returns
MoFEMErrorCode

Definition at line 1204 of file CPMeshCut.cpp.

1204 {
1205 MoFEM::Interface &m_field = cP.mField;
1207
1208 auto save_mesh = [&](const std::string file, const Range ents) {
1210 if (m_field.get_comm_rank() == 0) {
1211 if (!ents.empty()) {
1212 EntityHandle meshset;
1213 CHKERR m_field.get_moab().create_meshset(MESHSET_SET, meshset);
1214 CHKERR m_field.get_moab().add_entities(meshset, ents);
1215 CHKERR m_field.get_moab().write_file(file.c_str(), "VTK", "", &meshset,
1216 1);
1217 CHKERR m_field.get_moab().delete_entities(&meshset, 1);
1218 }
1219 }
1221 };
1222
1223 auto get_entities_from_contact_meshset = [&](BitRefLevel &bit) {
1225 cP.contactMeshsetFaces.clear();
1226 Range tris;
1228 if (cit->getName().compare(0, 11, "INT_CONTACT") == 0) {
1229 tris.clear();
1230 CHKERR meshsetMngPtr->getEntitiesByDimension(cit->getMeshsetId(),
1231 BLOCKSET, 2, tris, true);
1232 cP.contactMeshsetFaces.merge(tris);
1233 }
1234 }
1238 };
1239
1240 auto get_entities_from_constrained_interface = [&](BitRefLevel &bit) {
1242 cP.constrainedInterface.clear();
1243 Range tris;
1245 if (cit->getName().compare(0, 15, "INT_CONSTRAINED") == 0) {
1246 tris.clear();
1247 CHKERR meshsetMngPtr->getEntitiesByDimension(cit->getMeshsetId(),
1248 BLOCKSET, 2, tris, true);
1249 cP.constrainedInterface.merge(tris);
1250 }
1251 }
1255 };
1256
1257 Range all_constrained_faces;
1258 auto get_all_constrained_faces = [&](BitRefLevel &bit) {
1260 all_constrained_faces.clear();
1261
1262 if (!cP.ignoreContact) {
1263 CHKERR get_entities_from_contact_meshset(bit);
1264 all_constrained_faces.merge(cP.contactMeshsetFaces);
1265 }
1266
1267 CHKERR get_entities_from_constrained_interface(bit);
1268 all_constrained_faces.merge(cP.constrainedInterface);
1269
1271 };
1272
1273 if (debug)
1274 CHKERR save_mesh("crack_surface_before_merge.vtk",
1276
1277 if (debug) {
1278 CHKERR cP.mField.getInterface<BitRefManager>()->writeBitLevelByType(
1279 cP.mapBitLevel["mesh_cut"], BitRefLevel().set(), MBVERTEX,
1280 "verts_on_mesh_cut_level.vtk", "VTK", "", false);
1281 CHKERR cP.mField.getInterface<BitRefManager>()->writeBitLevelByType(
1282 cP.mapBitLevel["mesh_cut"], BitRefLevel().set(), MBEDGE,
1283 "edges_on_mesh_cut_level.vtk", "VTK", "", false);
1284 CHKERR cP.mField.getInterface<BitRefManager>()->writeBitLevelByType(
1285 cP.mapBitLevel["mesh_cut"], BitRefLevel().set(), MBTRI,
1286 "triangles_on_mesh_cut_level.vtk", "VTK", "", false);
1287 CHKERR cP.mField.getInterface<BitRefManager>()->writeBitLevelByType(
1288 cP.mapBitLevel["mesh_cut"], BitRefLevel().set(), MBTET,
1289 "tets_on_mesh_cut_level.vtk", "VTK", "", false);
1290 }
1291
1292 CHKERR get_all_constrained_faces(cP.mapBitLevel["mesh_cut"]);
1293 if (!all_constrained_faces.empty()) {
1294 CHKERR cutMeshPtr->setConstrainSurface(all_constrained_faces);
1295 if (debug)
1296 CHKERR save_mesh("constrained_surface_before_merge.vtk",
1297 all_constrained_faces);
1298 }
1299
1300 // Cut mesh, trim surface and merge bad edges
1302
1305 cornerNodes, true, debug);
1306
1307 BitRefLevel bit_after_merge = BitRefLevel().set(first_bit);
1310 bit_after_merge, const_cast<Range &>(cutMeshPtr->getMergedSurfaces()));
1311
1312 if (debug)
1313 CHKERR save_mesh("crack_surface_after_merge.vtk",
1315
1316 CHKERR get_all_constrained_faces(bit_after_merge);
1317 if (!all_constrained_faces.empty() && debug) {
1318 CHKERR save_mesh("constrained_surface_after_merge.vtk",
1319 all_constrained_faces);
1320 }
1321
1322 auto add_entities_to_crack_meshset = [&]() {
1326 Range crack_faces = cutMeshPtr->getMergedSurfaces();
1327
1328 if (!all_constrained_faces.empty() && crackedBodyBlockSetId < 0) {
1329 SETERRQ(PETSC_COMM_WORLD, MOFEM_NOT_FOUND,
1330 "Block set id of the cracked body is required for a simulation "
1331 "with contact and/or constrained interface, use "
1332 "'-my_cracked_body_block_set' parameter");
1333 }
1334 if (crackedBodyBlockSetId > 0) {
1335 Range cracked_body_tets, cracked_body_tris;
1337 crackedBodyBlockSetId, BLOCKSET, 3, cracked_body_tets, true);
1338 CHKERR m_field.get_moab().get_adjacencies(cracked_body_tets, 2, false,
1339 cracked_body_tris,
1340 moab::Interface::UNION);
1341 crack_faces = intersect(crack_faces, cracked_body_tris);
1342 if (debug)
1343 CHKERR save_mesh("crack_faces_intersected_with_body_block.vtk",
1344 crack_faces);
1345 }
1346
1348 crack_faces);
1350 };
1351
1352 auto clean_copied_cutting_surface_entities = [&]() {
1354 Range surface_verts;
1355 CHKERR m_field.get_moab().get_connectivity(cutMeshPtr->getSurface(),
1356 surface_verts);
1357 Range surface_edges;
1358 CHKERR m_field.get_moab().get_adjacencies(cutMeshPtr->getSurface(), 1,
1359 false, surface_edges,
1360 moab::Interface::UNION);
1361 // Remove entities from meshsets
1362 Range meshsets;
1363 CHKERR m_field.get_moab().get_entities_by_type(0, MBENTITYSET, meshsets,
1364 false);
1365 for (Range::iterator mit = meshsets.begin(); mit != meshsets.end(); mit++) {
1366 CHKERR m_field.get_moab().remove_entities(*mit, cutMeshPtr->getSurface());
1367 CHKERR m_field.get_moab().remove_entities(*mit, surface_edges);
1368 CHKERR m_field.get_moab().remove_entities(*mit, surface_verts);
1369 }
1370
1371 CHKERR m_field.get_moab().delete_entities(cutMeshPtr->getSurface());
1372 CHKERR m_field.get_moab().delete_entities(surface_edges);
1373 CHKERR m_field.get_moab().delete_entities(surface_verts);
1375 };
1376
1377 CHKERR add_entities_to_crack_meshset();
1378 CHKERR clean_copied_cutting_surface_entities();
1379
1381}
@ BLOCKSET
Definition: definitions.h:148
@ MOFEM_NOT_FOUND
Definition: definitions.h:33
static const bool debug
MoFEMErrorCode addMeshset(const CubitBCType cubit_bc_type, const int ms_id, const std::string name="")
add cubit meshset
#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.
auto bit
set bit
MoFEMErrorCode setFixEdgesAndCorners(const BitRefLevel bit)
Definition: CPMeshCut.cpp:1173
Range constrainedInterface
Range of faces on the constrained interface.
virtual int get_comm_rank() const =0
const Range & getMergedSurfaces() const
MoFEMErrorCode removePathologicalFrontTris(const BitRefLevel split_bit, Range &ents)
Remove pathological elements on surface internal front.
MoFEMErrorCode setConstrainSurface(const Range surf)
Set the constrain surface object.
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 addEntitiesToMeshset(const CubitBCType cubit_bc_type, const int ms_id, const Range &ents)
add entities to cubit meshset

◆ cutRefineAndSplit()

MoFEMErrorCode FractureMechanics::CPMeshCut::cutRefineAndSplit ( const int  verb = VERBOSE,
bool  debug = false 
)

Refine, cut, trim, merge and ref front and split.

Note
This functions squashes bits, and do refine of crack front after mesh is cutted. That refinement improves approximation of the fields but not geometry of crack front.
Parameters
verb
debug
Returns
MoFEMErrorCode

Definition at line 2265 of file CPMeshCut.cpp.

2265 {
2266 MoFEM::Interface &m_field = cP.mField;
2267 moab::Interface &moab = m_field.get_moab();
2268 if (debugCut)
2269 debug = true;
2271
2272 // FIXME: If will be needed in future one can cut, and split and then put
2273 // meshes on the right bit ref level.
2274
2275 // Get entities not in database
2276 Range all_free_ents;
2278
2279 // Cut mesh
2280 BitRefLevel bit = cP.mapBitLevel["mesh_cut"];
2281 std::size_t idx = 0;
2282 while (!bit.test(idx))
2283 ++idx;
2284 int first_bit = idx + 1;
2285
2286 CHKERR cutMesh(first_bit, debug);
2287
2288 // Squash bits
2289 BitRefLevel shift_mask;
2290 for (int ll = 0; ll != first_bit + 1; ++ll)
2291 shift_mask.set(ll);
2292 CHKERR bitRefPtr->shiftRightBitRef(first_bit, shift_mask, VERBOSE,
2293 MF_NOT_THROW);
2294
2295 if (debug)
2297 BitRefLevel().set(), MBTET,
2298 "out_mesh_after_cut.vtk", "VTK", "");
2299
2301
2302 if (debug)
2304 cP.mapBitLevel["spatial_domain"], BitRefLevel().set(), MBTET,
2305 "out_mesh_after_split_and_refine.vtk", "VTK", "");
2306
2307 // Delete left obsolete entities
2308
2309 auto get_entities_to_delete = [&]() {
2310 Range ents_not_in_database;
2311 CHKERR moab.get_entities_by_handle(0, ents_not_in_database, false);
2312 ents_not_in_database = subtract(
2313 ents_not_in_database, ents_not_in_database.subset_by_type(MBENTITYSET));
2314 CHKERR bitRefPtr->filterEntitiesNotInDatabase(ents_not_in_database);
2315 return subtract(ents_not_in_database, all_free_ents);
2316 };
2317 Range ents_to_remove = get_entities_to_delete();
2318
2319 auto remove_entities_from_meshsets = [&]() {
2321 // Remove entities from meshsets
2322 Range meshsets;
2323 CHKERR moab.get_entities_by_type(0, MBENTITYSET, meshsets, false);
2324 for (auto m : meshsets)
2325 CHKERR moab.remove_entities(m, ents_to_remove);
2327 };
2328
2329 CHKERR remove_entities_from_meshsets();
2330 CHKERR moab.delete_entities(ents_to_remove);
2331
2333}
@ VERBOSE
Definition: definitions.h:209
@ MF_NOT_THROW
Definition: definitions.h:101
FTensor::Index< 'm', SPACE_DIM > m
MoFEMErrorCode shiftRightBitRef(const int shift, const BitRefLevel mask=BitRefLevel().set(), int verb=DEFAULT_VERBOSITY, MoFEMTypes mf=MF_ZERO) const
right shift bit ref level
MoFEMErrorCode cutMesh(int &first_bit, const bool debug=false)
Ref, cut and trim, merge and do TetGen, if option is on.
Definition: CPMeshCut.cpp:1204
MoFEMErrorCode refineAndSplit(const int verb=1, const bool debug=false)
Definition: CPMeshCut.cpp:2155
MoFEMErrorCode getAllEntitiesNotInDatabase(Range &ents) const
Get all entities not in database.
MoFEMErrorCode writeBitLevelByType(const BitRefLevel bit, const BitRefLevel mask, const EntityType type, const char *file_name, const char *file_type, const char *options, const bool check_for_empty=true) const
Write bit ref level to file.
MoFEMErrorCode filterEntitiesNotInDatabase(Range &ents) const
Get entities not in database.

◆ findBodySkin()

MoFEMErrorCode FractureMechanics::CPMeshCut::findBodySkin ( const BitRefLevel  bit,
const int  verb = 1,
const bool  debug = false,
std::string  file_name = "out_body_skin.vtk" 
)

find body skin

Definition at line 1383 of file CPMeshCut.cpp.

1385 {
1386 MoFEM::Interface &m_field = cP.mField;
1387 moab::Interface &moab = m_field.get_moab();
1389
1390 auto save_mesh = [&](const std::string file, const Range ents) {
1392 if (m_field.get_comm_rank() == 0) {
1393 EntityHandle meshset;
1394 CHKERR moab.create_meshset(MESHSET_SET, meshset);
1395 CHKERR moab.add_entities(meshset, ents);
1396 CHKERR moab.write_file(file.c_str(), "VTK", "", &meshset, 1);
1397 CHKERR moab.delete_entities(&meshset, 1);
1398 }
1400 };
1401
1402 // Find body skin, take it from meshset if is found, if not create meshest,
1403 // and add skin entities
1404 cP.bodySkin.clear();
1406 Range tets, prisms;
1408 MBTET, tets);
1410 MBPRISM, prisms);
1411 prisms = subtract(prisms, cP.contactElements);
1412 prisms = subtract(prisms, cP.mortarContactElements);
1413 Range prisms_faces;
1414 for (auto p : prisms) {
1415 EntityHandle f;
1416 CHKERR moab.side_element(p, 2, 3, f);
1417 prisms_faces.insert(f);
1418 CHKERR moab.side_element(p, 2, 4, f);
1419 prisms_faces.insert(f);
1420 }
1421 Range skin_faces;
1422 Skinner skin(&m_field.get_moab());
1423 CHKERR skin.find_skin(0, tets, false, skin_faces);
1424 skin_faces = subtract(skin_faces, prisms_faces);
1426 EntityHandle meshset;
1428 CHKERR moab.add_entities(meshset, skin_faces);
1429 }
1430
1431 // take skin entities from meshset
1433 cP.bodySkin, true);
1435 cP.bodySkin);
1436
1437 if (verb >= VERY_VERBOSE)
1438 if (m_field.get_comm_rank() == 0)
1439 cerr << cP.bodySkin << endl;
1440
1441 if (debug)
1442 save_mesh(file_name, cP.bodySkin);
1443
1445}
static Index< 'p', 3 > p
@ VERY_VERBOSE
Definition: definitions.h:210
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
auto f
Definition: HenckyOps.hpp:5

◆ findContactFromPrisms()

MoFEMErrorCode FractureMechanics::CPMeshCut::findContactFromPrisms ( const BitRefLevel  bit,
const int  verb = 1,
const bool  debug = false 
)

get contact elements

Definition at line 1987 of file CPMeshCut.cpp.

1989 {
1990 MoFEM::Interface &m_field = cP.mField;
1993
1994 cP.contactElements.clear();
1995 cP.contactSlaveFaces.clear();
1996 cP.contactMasterFaces.clear();
1997
2002
2003 if (!cP.mortarContactElements.empty())
2005
2006 if (!cP.contactElements.empty()) {
2007 EntityHandle tri;
2008 for (auto p : cP.contactElements) {
2009 CHKERR moab.side_element(p, 2, 3, tri);
2010 cP.contactMasterFaces.insert(tri);
2011 CHKERR moab.side_element(p, 2, 4, tri);
2012 cP.contactSlaveFaces.insert(tri);
2013 }
2014 }
2015
2016 Range slave_nodes, master_nodes;
2017 CHKERR moab.get_connectivity(cP.contactSlaveFaces, slave_nodes, false);
2018 CHKERR moab.get_connectivity(cP.contactMasterFaces, master_nodes, false);
2019
2020 if (slave_nodes.size() < master_nodes.size()) {
2023 } else {
2026 }
2027
2029}

◆ findCrack()

MoFEMErrorCode FractureMechanics::CPMeshCut::findCrack ( const BitRefLevel  bit,
const int  verb = 1,
const bool  debug = false 
)

get crack front

Definition at line 1447 of file CPMeshCut.cpp.

1448 {
1449 MoFEM::Interface &m_field = cP.mField;
1452
1453 auto save_mesh = [&](const std::string file, const Range ents) {
1455 if (m_field.get_comm_rank() == 0) {
1456 EntityHandle meshset;
1457 CHKERR moab.create_meshset(MESHSET_SET, meshset);
1458 CHKERR moab.add_entities(meshset, ents);
1459 CHKERR moab.write_file(file.c_str(), "VTK", "", &meshset, 1);
1460 CHKERR moab.delete_entities(&meshset, 1);
1461 }
1463 };
1464
1465 cP.crackFaces.clear();
1467 PetscPrintf(PETSC_COMM_WORLD, "Sideset %d with crack surface not found\n",
1470 }
1472 cP.crackFaces, true);
1474 cP.crackFaces);
1475
1476 // Find crack front
1477 cP.crackFront.clear();
1478
1480 Range skin_faces_edges;
1481 Range constrained_faces =
1483 CHKERR moab.get_adjacencies(unite(cP.bodySkin, constrained_faces), 1, false,
1484 skin_faces_edges, moab::Interface::UNION);
1485 Range crack_front;
1486 Skinner skin(&m_field.get_moab());
1487 CHKERR skin.find_skin(0, cP.crackFaces, false, crack_front);
1488 crack_front = subtract(crack_front, skin_faces_edges);
1490 EntityHandle meshset;
1492 CHKERR moab.add_entities(meshset, crack_front);
1493 }
1494
1496 cP.crackFront, true);
1497
1498 Range crack_faces_edges;
1499 CHKERR moab.get_adjacencies(cP.crackFaces, 1, false, crack_faces_edges,
1500 moab::Interface::UNION);
1501 cP.crackFront = intersect(cP.crackFront, crack_faces_edges);
1502
1503 if (verb >= VERBOSE)
1504 if (m_field.get_comm_rank() == 0)
1505 cerr << cP.crackFront << endl;
1506
1507 if (debug)
1508 save_mesh("out_crack_and_front.vtk", unite(cP.crackFront, cP.crackFaces));
1509
1511}

◆ findCrackFromPrisms()

MoFEMErrorCode FractureMechanics::CPMeshCut::findCrackFromPrisms ( const BitRefLevel  bit,
const int  verb = 1,
const bool  debug = false 
)

get crack front

Definition at line 1513 of file CPMeshCut.cpp.

1515 {
1516 MoFEM::Interface &m_field = cP.mField;
1519
1520 auto save_mesh = [&](const std::string file, const Range ents) {
1522 if (m_field.get_comm_rank() == 0) {
1523 EntityHandle meshset;
1524 CHKERR moab.create_meshset(MESHSET_SET, meshset);
1525 CHKERR moab.add_entities(meshset, ents);
1526 CHKERR moab.write_file(file.c_str(), "VTK", "", &meshset, 1);
1527 CHKERR moab.delete_entities(&meshset, 1);
1528 }
1530 };
1531
1532 Range prisms_level;
1534 MBPRISM, prisms_level);
1535 prisms_level = subtract(prisms_level, cP.contactElements);
1536 prisms_level = subtract(prisms_level, cP.mortarContactElements);
1537
1538 cP.oneSideCrackFaces.clear();
1539 cP.otherSideCrackFaces.clear();
1540 for (auto p : prisms_level) {
1541 EntityHandle face3, face4;
1542 CHKERR m_field.get_moab().side_element(p, 2, 3, face3);
1543 CHKERR m_field.get_moab().side_element(p, 2, 4, face4);
1544 cP.oneSideCrackFaces.insert(face3);
1545 cP.otherSideCrackFaces.insert(face4);
1546 }
1551 cP.crackFaces);
1552
1553 if (cP.otherSideConstrains == PETSC_TRUE)
1555
1556 if (debug)
1557 if (m_field.get_comm_rank() == 0) {
1558 CHKERR save_mesh("out_one_side.vtk", cP.oneSideCrackFaces);
1559 if (!cP.otherSideCrackFaces.empty())
1560 CHKERR save_mesh("out_other_side.vtk", cP.otherSideCrackFaces);
1561 }
1562
1563 if (!cP.otherSideCrackFaces.empty()) {
1564 if (cP.otherSideCrackFaces.size() != cP.oneSideCrackFaces.size())
1565 SETERRQ2(PETSC_COMM_WORLD, MOFEM_DATA_INCONSISTENCY,
1566 "Wrong number of faces on both sides of crack %d != %d",
1568 }
1569
1570 // Find crack front
1571 cP.crackFront.clear();
1572
1574 auto create_front_meshset = [&]() {
1576 Range skin_faces_edges;
1577 CHKERR moab.get_adjacencies(unite(cP.bodySkin, cP.constrainedInterface), 1,
1578 false, skin_faces_edges,
1579 moab::Interface::UNION);
1580 Range crack_front;
1581 Skinner skin(&m_field.get_moab());
1582 CHKERR skin.find_skin(0, cP.oneSideCrackFaces, false, crack_front);
1583 crack_front = subtract(crack_front, skin_faces_edges);
1586 crack_front);
1588 };
1589 CHKERR create_front_meshset();
1591 cP.crackFront, true);
1592
1593 if (verb >= VERBOSE)
1594 if (m_field.get_comm_rank() == 0)
1595 cerr << cP.crackFront << endl;
1596
1597 if (verb >= VERBOSE) {
1598 CHKERR PetscPrintf(m_field.get_comm(),
1599 "number of Crack Surfaces Faces = %d\n",
1600 cP.crackFaces.size());
1601 CHKERR PetscPrintf(m_field.get_comm(),
1602 "number of Crack Surfaces Faces On One Side = %d\n",
1603 cP.oneSideCrackFaces.size());
1604 CHKERR PetscPrintf(m_field.get_comm(),
1605 "number of Crack Surfaces Faces On Other Side = %d\n",
1606 cP.otherSideCrackFaces.size());
1607 CHKERR PetscPrintf(m_field.get_comm(), "number of Crack Front Edges = %d\n",
1608 cP.crackFront.size());
1609 Range prisms_level;
1611 MBPRISM, prisms_level);
1612 CHKERR PetscPrintf(m_field.get_comm(), "number of Prisms = %d\n",
1613 prisms_level.size());
1614 }
1615
1616 if (debug)
1617 save_mesh("out_crack_and_front.vtk", unite(cP.crackFront, cP.crackFaces));
1618
1620}
@ MF_ZERO
Definition: definitions.h:98
@ MOFEM_DATA_INCONSISTENCY
Definition: definitions.h:31
virtual MPI_Comm & get_comm() const =0
MoFEMErrorCode deleteMeshset(const CubitBCType cubit_bc_type, const int ms_id, const MoFEMTypes bh=MF_EXIST)
delete cubit meshset

◆ getContactPrismsBlockSetId()

int FractureMechanics::CPMeshCut::getContactPrismsBlockSetId ( ) const

Definition at line 170 of file CPMeshCut.hpp.

170{ return contactPrismsBlockSetId; };

◆ getCornerNodes()

const Range & FractureMechanics::CPMeshCut::getCornerNodes ( ) const

Definition at line 156 of file CPMeshCut.hpp.

156{ return cornerNodes; }

◆ getCrackFrontId()

int FractureMechanics::CPMeshCut::getCrackFrontId ( ) const

Definition at line 169 of file CPMeshCut.hpp.

169{ return crackFrontId; };

◆ getCrackSurfaceId()

int FractureMechanics::CPMeshCut::getCrackSurfaceId ( ) const

Definition at line 168 of file CPMeshCut.hpp.

168{ return crackSurfaceId; };

◆ getCutSurfMeshName()

const std::string & FractureMechanics::CPMeshCut::getCutSurfMeshName ( ) const

Definition at line 157 of file CPMeshCut.hpp.

157 {
158 return cutSurfMeshName;
159 }

◆ getCuttingSurface()

const Range & FractureMechanics::CPMeshCut::getCuttingSurface ( ) const

Definition at line 96 of file CPMeshCut.hpp.

96{ return sUrface; }

◆ getCuttingSurfaceSidesetId()

int FractureMechanics::CPMeshCut::getCuttingSurfaceSidesetId ( ) const

Definition at line 166 of file CPMeshCut.hpp.

166{ return cuttingSurfaceSidesetId; };

◆ getEdgesBlockSet()

const int FractureMechanics::CPMeshCut::getEdgesBlockSet ( ) const

Definition at line 160 of file CPMeshCut.hpp.

160 {
161 return edgesBlockSetFlg == PETSC_TRUE ? edgesBlockSet : 0;
162 }

◆ getFixedEdges()

const Range & FractureMechanics::CPMeshCut::getFixedEdges ( ) const

Definition at line 155 of file CPMeshCut.hpp.

155{ return fixedEdges; }

◆ getFrontEdgesAndElements()

MoFEMErrorCode FractureMechanics::CPMeshCut::getFrontEdgesAndElements ( const BitRefLevel  bit,
const int  verb = 1,
const bool  debug = false 
)

get crack front edges and finite elements

Definition at line 2335 of file CPMeshCut.cpp.

2337 {
2338 MoFEM::Interface &m_field = cP.mField;
2339 moab::Interface &moab = m_field.get_moab();
2341
2342 Range level_edges;
2344 MBEDGE, level_edges);
2345
2346 cP.crackFrontNodes.clear();
2347 CHKERR moab.get_connectivity(cP.crackFront, cP.crackFrontNodes, true);
2348 cP.crackFrontNodesEdges.clear();
2349 CHKERR moab.get_adjacencies(cP.crackFrontNodes, 1, false,
2350 cP.crackFrontNodesEdges, moab::Interface::UNION);
2352 cP.crackFrontNodesEdges = intersect(cP.crackFrontNodesEdges, level_edges);
2353
2354 Range level_tets;
2356 MBTET, level_tets);
2357
2358 cP.crackFrontElements.clear();
2359 CHKERR moab.get_adjacencies(cP.crackFrontNodes, 3, false,
2360 cP.crackFrontElements, moab::Interface::UNION);
2361 cP.crackFrontElements = intersect(cP.crackFrontElements, level_tets);
2362
2363 if (debug) {
2364 EntityHandle out_meshset;
2365 CHKERR cP.mField.get_moab().create_meshset(MESHSET_SET, out_meshset);
2366 CHKERR cP.mField.get_moab().add_entities(out_meshset, cP.crackFront);
2367 CHKERR cP.mField.get_moab().add_entities(out_meshset, cP.crackFrontNodes);
2368 CHKERR cP.mField.get_moab().write_file("crack_front_nodes.vtk", "VTK", "",
2369 &out_meshset, 1);
2370 CHKERR cP.mField.get_moab().delete_entities(&out_meshset, 1);
2371 }
2372
2373 if (verb >= VERY_VERBOSE) {
2374 cerr << "getFrontEdgesAndElements ";
2375 cerr << "crackFront\n" << cP.crackFront << endl;
2376 cerr << "crackFrontNodes\n" << cP.crackFrontNodes << endl;
2377 cerr << cP.crackFront.size() << " ";
2378 cerr << cP.crackFrontNodes.size() << " ";
2379 cerr << cP.crackFrontNodesEdges.size() << " ";
2380 cerr << cP.crackFrontElements.size() << endl;
2381 }
2382
2383 if (verb >= NOISY) {
2384 cerr << cP.crackFrontNodes << endl;
2385 cerr << cP.crackFrontNodesEdges << endl;
2386 cerr << cP.crackFrontElements << endl;
2387 }
2388
2390}
@ NOISY
Definition: definitions.h:211

◆ getInterfacesPtr()

MoFEMErrorCode FractureMechanics::CPMeshCut::getInterfacesPtr ( )

◆ getMeshOrgCoords()

MoFEMErrorCode FractureMechanics::CPMeshCut::getMeshOrgCoords ( )
private

Definition at line 2420 of file CPMeshCut.cpp.

2420 {
2422
2423 BitRefLevel mask;
2424 for (int ll = 20; ll != BITREFLEVEL_SIZE; ++ll)
2425 mask.set(ll);
2426
2427 // Get node coordinates on original mesh
2428 Range org_verts;
2430 MBVERTEX, org_verts);
2431 originalCoords.resize(org_verts.size() * 3);
2432 CHKERR cP.mField.get_moab().get_coords(org_verts, &*originalCoords.begin());
2433
2434 // Save positions on tag
2435 double def_position[] = {0, 0, 0};
2436 Tag th;
2437 CHKERR cP.mField.get_moab().tag_get_handle("ORG_POSITION", 3, MB_TYPE_DOUBLE,
2438 th, MB_TAG_CREAT | MB_TAG_SPARSE,
2439 def_position);
2440 CHKERR cP.mField.get_moab().tag_set_data(th, org_verts,
2441 &*originalCoords.begin());
2442
2444}
#define BITREFLEVEL_SIZE
max number of refinements
Definition: definitions.h:219

◆ getOptions()

MoFEMErrorCode FractureMechanics::CPMeshCut::getOptions ( )

Get options from command line.

Returns
MoFEMErrorCode

Definition at line 85 of file CPMeshCut.cpp.

85 {
87
88 PetscBool shift_flg;
89 int nmax_shift = 3;
90
91 CHKERR PetscOptionsBegin(PETSC_COMM_WORLD, "", "Mesh cut options", "none");
92
93 edgesBlockSet = 100;
94 CHKERR PetscOptionsInt("-edges_block_set", "edges side set", "",
96 vertexBlockSet = 101;
97 CHKERR PetscOptionsInt("-vertex_block_set", "vertex block set", "",
98 vertexBlockSet, &vertexBlockSet, PETSC_NULL);
99 //to resolve issues with newer cubit
100 vertexNodeSet = 102;
101 CHKERR PetscOptionsInt("-vertex_node_set", "vertex node set", "",
102 vertexNodeSet, &vertexNodeSet, PETSC_NULL);
103 fractionLevel = 2;
104 CHKERR PetscOptionsInt("-fraction_level", "fraction of merges merged", "",
105 fractionLevel, &fractionLevel, PETSC_NULL);
106 tolCut = 1e-2;
107 CHKERR PetscOptionsScalar("-tol_cut", "get tolerance for cut edges", "",
108 tolCut, &tolCut, PETSC_NULL);
109 tolCutClose = 2e-1;
110 CHKERR PetscOptionsScalar("-tol_cut_close", "get tolerance for cut edges", "",
111 tolCutClose, &tolCutClose, PETSC_NULL);
112 tolTrimClose = 2e-1;
113 CHKERR PetscOptionsScalar("-tol_trim_close", "get tolerance for trim edges",
114 "", tolTrimClose, &tolTrimClose, PETSC_NULL);
115
116 debugCut = PETSC_FALSE;
117 CHKERR PetscOptionsBool("-cut_debug_cut", "debug mesh cutting", "", debugCut,
118 &debugCut, NULL);
119
120 removePathologicalFrontTris = PETSC_FALSE;
121 CHKERR PetscOptionsBool(
122 "-remove_pathological_front_tris",
123 "remove triangles which have three nodes on crack front", "",
125
126 MOFEM_LOG_C("CPMeshCutWorld", Sev::inform,
127 "### Input parameter: -fraction_level %d", fractionLevel);
128 MOFEM_LOG_C("CPMeshCutWorld", Sev::inform,
129 "### Input parameter: -tol_cut %6.4e", tolCut);
130 MOFEM_LOG_C("CPMeshCutWorld", Sev::inform,
131 "### Input parameter: -tol_cut_close %6.4e", tolCutClose);
132 MOFEM_LOG_C("CPMeshCutWorld", Sev::inform,
133 "### Input parameter: -tol_trim_close %6.4e", tolTrimClose);
134
135 PetscBool crack_surf_mesh_flg;
136 char crack_surf_mesh_name[255];
137 CHKERR PetscOptionsString(
138 "-my_cut_surface_file", "file with crack surface for initial cut", "",
139 crack_surf_mesh_name, crack_surf_mesh_name, 255, &crack_surf_mesh_flg);
140 if (crack_surf_mesh_flg)
141 cutSurfMeshName = std::string(crack_surf_mesh_name);
142
143 sHift[0] = sHift[1] = sHift[2] = 0;
144 CHKERR PetscOptionsRealArray("-crack_shift", "shift surface by vector", "",
145 sHift, &nmax_shift, &shift_flg);
146
148 CHKERR PetscOptionsInt(
149 "-cut_surface_side_set", "sideset id with cutting surface", "",
151
152 skinOfTheBodyId = 102;
153 CHKERR PetscOptionsInt("-body_skin", "body skin sideset id", "",
154 skinOfTheBodyId, &skinOfTheBodyId, PETSC_NULL);
155 crackSurfaceId = 200;
156 CHKERR PetscOptionsInt("-crack_side_set", "sideset id with crack surface", "",
157 crackSurfaceId, &crackSurfaceId, PETSC_NULL);
158 crackFrontId = 201;
159 CHKERR PetscOptionsInt("-front_side_set", "sideset id with crack front", "",
160 crackFrontId, &crackFrontId, PETSC_NULL);
161
162 nbRefBeforCut = 0;
163 CHKERR PetscOptionsInt("-ref_before_cut", "number of refinements before cut",
164 "", nbRefBeforCut, &nbRefBeforCut, PETSC_NULL);
165 MOFEM_LOG_C("CPMeshCutWorld", Sev::inform,
166 "### Input parameter: -ref_before_cut %d", nbRefBeforCut);
167
168 nbRefBeforTrim = 0;
169 PetscBool flg = PETSC_FALSE;
170 CHKERR PetscOptionsInt("-ref_before_tim", "number of refinements before trim",
171 "", nbRefBeforTrim, &nbRefBeforTrim, &flg);
172 if (flg == PETSC_TRUE) {
173 MOFEM_LOG("CPMeshCutWorld", Sev::inform)
174 << "### Input parameter: -ref_before_tim is still valid but will be "
175 "made obsolete in the next release, please use -ref_before_trim "
176 "instead.";
177 }
178
179 CHKERR PetscOptionsInt("-ref_before_trim",
180 "number of refinements before trim", "",
181 nbRefBeforTrim, &nbRefBeforTrim, PETSC_NULL);
182 MOFEM_LOG_C("CPMeshCutWorld", Sev::inform,
183 "### Input parameter: -ref_before_trim %d", nbRefBeforTrim);
184
186 CHKERR PetscOptionsInt(
187 "-my_cracked_body_block_set", "blockset id of the cracked body", "",
189
190 MOFEM_LOG_C("CPMeshCutWorld", Sev::inform,
191 "### Input parameter: -my_cracked_body_block_set %d",
193
195 CHKERR PetscOptionsInt(
196 "-my_contact_prisms_block_set", "blockset id of contact prisms", "",
198
199 ierr = PetscOptionsEnd();
200 CHKERRG(ierr);
201
202 if (shift_flg && nmax_shift != 3) {
203 SETERRQ(PETSC_COMM_SELF, MOFEM_ATOM_TEST_INVALID, "three values expected");
204 }
205
207}
#define MOFEM_LOG_C(channel, severity, format,...)
Definition: LogManager.hpp:304
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:483
@ MOFEM_ATOM_TEST_INVALID
Definition: definitions.h:40

◆ getOrgCoords() [1/2]

std::vector< double > & FractureMechanics::CPMeshCut::getOrgCoords ( )

Definition at line 152 of file CPMeshCut.hpp.

152{ return originalCoords; }

◆ getOrgCoords() [2/2]

const std::vector< double > & FractureMechanics::CPMeshCut::getOrgCoords ( ) const

Definition at line 153 of file CPMeshCut.hpp.

153{ return originalCoords; }

◆ getSkinOfTheBodyId()

int FractureMechanics::CPMeshCut::getSkinOfTheBodyId ( ) const

Definition at line 167 of file CPMeshCut.hpp.

167{ return skinOfTheBodyId; }

◆ insertContactInterface()

MoFEMErrorCode FractureMechanics::CPMeshCut::insertContactInterface ( const int  verb = 1,
const bool  debug = false 
)

insert contact interface

Definition at line 1860 of file CPMeshCut.cpp.

1861 {
1862 MoFEM::Interface &m_field = cP.mField;
1863 moab::Interface &moab = m_field.get_moab();
1864 PrismInterface *prism_interface;
1866 CHKERR m_field.getInterface(prism_interface);
1867
1868 std::vector<BitRefLevel> bit_levels;
1869 bit_levels.push_back(cP.mapBitLevel["mesh_cut"]);
1870
1871 auto get_cracked_not_in_body_blockset = [&](auto ref_level_meshset) {
1872 Range ref_level_tet;
1873 CHKERR moab.get_entities_by_type(ref_level_meshset, MBTET, ref_level_tet,
1874 0);
1875 Range cracked_body_tets;
1877 crackedBodyBlockSetId, BLOCKSET, 3, cracked_body_tets, true);
1878 return subtract(ref_level_tet, cracked_body_tets);
1879 };
1880
1881 std::size_t idx = 0;
1882 while (!bit_levels.back().test(idx))
1883 ++idx;
1884 int contact_lvl = idx + 1;
1885
1886 // Remove parents, that should be no use anymore
1889
1891 if (cit->getName().compare(0, 11, "INT_CONTACT") == 0) {
1892 CHKERR PetscPrintf(PETSC_COMM_WORLD, "Insert %s (id: %d)\n",
1893 cit->getName().c_str(), cit->getMeshsetId());
1894 EntityHandle cubit_meshset = cit->getMeshset();
1895
1896 // get tet entities form back bit_level
1897 EntityHandle ref_level_meshset;
1898 CHKERR moab.create_meshset(MESHSET_SET, ref_level_meshset);
1900 bit_levels.back(), BitRefLevel().set(), MBTET, ref_level_meshset);
1902 bit_levels.back(), BitRefLevel().set(), MBPRISM, ref_level_meshset);
1903
1904 // get faces and tets to split
1905 CHKERR prism_interface->getSides(
1906 cubit_meshset, bit_levels.back(),
1907 get_cracked_not_in_body_blockset(ref_level_meshset), true, verb);
1908 // set new bit level
1909 bit_levels.push_back(BitRefLevel().set(contact_lvl));
1910 // split faces and tets
1911 CHKERR prism_interface->splitSides(ref_level_meshset, bit_levels.back(),
1912 cubit_meshset, true, true, verb);
1913
1914 // clean meshsets
1915 CHKERR moab.delete_entities(&ref_level_meshset, 1);
1916
1917 // update cubit meshsets
1918 for (_IT_CUBITMESHSETS_FOR_LOOP_(m_field, cubit_it)) {
1919 EntityHandle updated_meshset = cubit_it->meshset;
1921 updated_meshset, bit_levels.front(), BitRefLevel().set(),
1922 bit_levels.back(), bit_levels.back(), updated_meshset, MBMAXTYPE,
1923 true);
1924 }
1925
1926 if (cP.isPartitioned) {
1927 // FIXME check partitioned case
1928 ParallelComm *pcomm = ParallelComm::get_pcomm(&moab, MYPCOMM_INDEX);
1929 Tag part_tag = pcomm->part_tag();
1930 Range tagged_sets;
1931 CHKERR m_field.get_moab().get_entities_by_type_and_tag(
1932 0, MBENTITYSET, &part_tag, NULL, 1, tagged_sets,
1933 moab::Interface::UNION);
1934 for (Range::iterator mit = tagged_sets.begin();
1935 mit != tagged_sets.end(); mit++) {
1936 int part;
1937 CHKERR moab.tag_get_data(part_tag, &*mit, 1, &part);
1939 *mit, bit_levels.front(), BitRefLevel().set(), bit_levels.back(),
1940 bit_levels.back(), *mit, MBTET, true);
1941 Range ents3d;
1942 CHKERR moab.get_entities_by_dimension(*mit, 3, ents3d, true);
1943 for (Range::iterator eit = ents3d.begin(); eit != ents3d.end();
1944 eit++) {
1945 CHKERR moab.tag_set_data(part_tag, &*eit, 1, &part);
1946 }
1947 for (int dd = 0; dd != 3; dd++) {
1948 Range ents;
1949 CHKERR moab.get_adjacencies(ents3d, dd, false, ents,
1950 moab::Interface::UNION);
1951 for (Range::iterator eit = ents.begin(); eit != ents.end(); eit++) {
1952 CHKERR moab.tag_set_data(part_tag, &*eit, 1, &part);
1953 }
1954 }
1955 }
1956 }
1957
1958 if (contact_lvl > 1) {
1959 CHKERR m_field.delete_ents_by_bit_ref(bit_levels.front(),
1960 bit_levels.front(), true);
1961 }
1962 BitRefLevel shift_mask;
1963 for (int ll = contact_lvl - 1; ll != 19; ++ll) {
1964 shift_mask.set(ll);
1965 }
1966 CHKERR m_field.getInterface<BitRefManager>()->shiftRightBitRef(
1967 1, shift_mask, verb);
1968 bit_levels.pop_back();
1969 }
1970 }
1971
1974 }
1976
1977 Range prisms_level;
1979 bit_levels.back(), BitRefLevel().set(), MBPRISM, prisms_level);
1980
1982 prisms_level);
1983
1985}
#define MYPCOMM_INDEX
default communicator number PCOMM
Definition: definitions.h:215
virtual const RefEntity_multiIndex * get_ref_ents() const =0
Get the ref ents object.
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 entities.
#define _IT_CUBITMESHSETS_FOR_LOOP_(MESHSET_MANAGER, IT)
Iterator that loops over all the Cubit MeshSets in a moFEM field.
const Tensor2_symmetric_Expr< const ddTensor0< T, Dim, i, j >, typename promote< T, double >::V, Dim, i, j > dd(const Tensor0< T * > &a, const Index< i, Dim > index1, const Index< j, Dim > index2, const Tensor1< int, Dim > &d_ijk, const Tensor1< double, Dim > &d_xyz)
Definition: ddTensor0.hpp:33
MoFEMErrorCode reconstructMultiIndex(const MI &mi, MO &&mo=Modify_change_nothing())
Template used to reconstruct multi-index.
Definition: Templates.hpp:1432
virtual MoFEMErrorCode delete_ents_by_bit_ref(const BitRefLevel bit, const BitRefLevel mask, const bool remove_parent=false, int verb=DEFAULT_VERBOSITY, MoFEMTypes mf=MF_ZERO)=0
delete entities form mofem and moab database
Create interface from given surface and insert flat prisms in-between.
MoFEMErrorCode getSides(const int msId, const CubitBCType cubit_bc_type, const BitRefLevel mesh_bit_level, const bool recursive, int verb=QUIET)
Store tetrahedra from each side of the interface separately in two child meshsets of the parent meshs...
MoFEMErrorCode splitSides(const EntityHandle meshset, const BitRefLevel &bit, const int msId, const CubitBCType cubit_bc_type, const bool add_interface_entities, const bool recursive=false, int verb=QUIET)
Split nodes and other entities of tetrahedra on both sides of the interface and insert flat prisms in...

◆ query_interface()

MoFEMErrorCode FractureMechanics::CPMeshCut::query_interface ( boost::typeindex::type_index  type_index,
MoFEM::UnknownInterface **  iface 
) const
virtual

Getting interface of core database.

Parameters
uuidunique ID of interface of CPMeshCut solely
ifacereturned pointer to interface
Returns
error code

Implements MoFEM::UnknownInterface.

Definition at line 45 of file CPMeshCut.cpp.

46 {
47 *iface = const_cast<CPMeshCut *>(this);
48 return 0;
49}
CPMeshCut(CrackPropagation &cp)
Definition: CPMeshCut.cpp:51

◆ rebuildCrackSurface()

MoFEMErrorCode FractureMechanics::CPMeshCut::rebuildCrackSurface ( const double  factor,
const std::string  save_mesh = "",
const int  verb = QUIET,
bool  debug = false 
)

Definition at line 344 of file CPMeshCut.cpp.

346 {
347 Skinner skin(&cP.mField.get_moab());
348 MoFEM::Interface &m_field = cP.mField;
349 if (debugCut)
350 debug = true;
352
353 BitRefLevel bit = cP.mapBitLevel["spatial_domain"];
354
356
357 Range bit_edges;
358
359 Tag th_griffith_force;
360 Tag th_griffith_force_projected;
361 Tag th_w;
362 Tag th_front_positions;
363 Tag th_area_change;
364 Tag th_material_force;
365
366 boost::shared_ptr<OrientedBoxTreeTool> tree_body_skin_ptr;
367 tree_body_skin_ptr = boost::shared_ptr<OrientedBoxTreeTool>(
368 new OrientedBoxTreeTool(&cP.mField.get_moab(), "ROOTSETSURF", true));
369 EntityHandle rootset_tree_body_skin;
370
371 /// Maps nodes on old crack surface and new crack surface
372 std::map<EntityHandle, EntityHandle> verts_map;
373
374 // Seed pseudo-random generator with step. Algorithm is deterministic.
375 srand(cP.startStep);
376
377 struct SetFrontPositionsParameters {
378
379 double tolRelSnapToFixedEdges; ///< Relative tolerance to snap edges to
380 ///< fixed edges
381 double tolAbsSnapToFixedEdges; ///< Relative tolerance to snap edges to
382 ///< fixed edges
383 double cornerFactor; ///< extension of the edge beyond corner when node on
384 ///< corner edge
385 double skinFactor; ///< extension of cutting surface when node on the skin
386 double endNodeEdgeDeltaFactor; ///< End node edge extension
387
388 SetFrontPositionsParameters() {
389 tolRelSnapToFixedEdges = 0.1;
390 tolAbsSnapToFixedEdges = 0.;
391 cornerFactor = 0.2;
392 skinFactor = 0.2;
393 endNodeEdgeDeltaFactor = 0.2;
394 ierr = getOptions();
395 CHKERRABORT(PETSC_COMM_WORLD, ierr);
396 }
397
400
401 CHKERR PetscOptionsBegin(PETSC_COMM_WORLD, "cutting_",
402 "Set front positions", "none");
403
404 CHKERR PetscOptionsScalar(
405 "-surf_corner_factor",
406 "extension of the edge beyond corner when node on corer edge", "",
407 cornerFactor, &cornerFactor, PETSC_NULL);
408 CHKERR PetscOptionsScalar(
409 "-surf_skin_factor",
410 "extension of cutting surface when node on the skin", "", skinFactor,
411 &skinFactor, PETSC_NULL);
412 CHKERR PetscOptionsScalar(
413 "-end_node_edge_delta_factor",
414 "extension of the front edge at the node on the surface", "",
415 endNodeEdgeDeltaFactor, &endNodeEdgeDeltaFactor, PETSC_NULL);
416 CHKERR PetscOptionsScalar("-snap_to_fixed_edge_rtol",
417 "snap to fixed edges (relative tolerances)", "",
418 tolRelSnapToFixedEdges, &tolRelSnapToFixedEdges,
419 PETSC_NULL);
420 CHKERR PetscOptionsScalar("-snap_to_fixed_edge_atol",
421 "snap to fixed edges (absolute tolerances)", "",
422 tolAbsSnapToFixedEdges, &tolAbsSnapToFixedEdges,
423 PETSC_NULL);
424
425 MOFEM_LOG_C("CPMeshCutWorld", Sev::inform,
426 "### Input parameter: -cutting_surf_corner_factor %6.4e",
427 cornerFactor);
428 MOFEM_LOG_C("CPMeshCutWorld", Sev::inform,
429 "### Input parameter: -cutting_surf_skin_factor %6.4e",
430 skinFactor);
432 "CPMeshCutWorld", Sev::inform,
433 "### Input parameter: -cutting_end_node_edge_delta_factor %6.4e",
434 endNodeEdgeDeltaFactor);
435 MOFEM_LOG_C("CPMeshCutWorld", Sev::inform,
436 "### Input parameter: -cutting_snap_to_fixed_edge_rtol %6.4e",
437 tolRelSnapToFixedEdges);
438 MOFEM_LOG_C("CPMeshCutWorld", Sev::inform,
439 "### Input parameter: -cutting_snap_to_fixed_edge_atol %6.4e",
440 tolAbsSnapToFixedEdges);
441
442 ierr = PetscOptionsEnd();
443 CHKERRG(ierr);
444
446 }
447 };
448
449 auto save_entities = [&](const std::string name, Range ents) {
451 if (!ents.empty()) {
452 EntityHandle meshset;
453 CHKERR cP.mField.get_moab().create_meshset(MESHSET_SET, meshset);
454 CHKERR cP.mField.get_moab().add_entities(meshset, ents);
455 CHKERR cP.mField.get_moab().write_file(name.c_str(), "VTK", "", &meshset,
456 1);
457 CHKERR cP.mField.get_moab().delete_entities(&meshset, 1);
458 }
460 };
461
462 auto get_tags = [&]() {
464 CHKERR cP.mField.get_moab().tag_get_handle("GRIFFITH_FORCE",
465 th_griffith_force);
466 CHKERR cP.mField.get_moab().tag_get_handle("GRIFFITH_FORCE_PROJECTED",
467 th_griffith_force_projected);
468 CHKERR cP.mField.get_moab().tag_get_handle("W", th_w);
469 CHKERR cP.mField.get_moab().tag_get_handle("MATERIAL_FORCE",
470 th_material_force);
471
472 rval = cP.mField.get_moab().tag_get_handle("FRONT_POSITIONS",
473 th_front_positions);
474 if (rval == MB_SUCCESS)
475 CHKERR cP.mField.get_moab().tag_delete(th_front_positions);
476
477 double def_val[] = {0, 0, 0};
478 CHKERR cP.mField.get_moab().tag_get_handle(
479 "FRONT_POSITIONS", 3, MB_TYPE_DOUBLE, th_front_positions,
480 MB_TAG_CREAT | MB_TAG_SPARSE, def_val);
481
482 rval = cP.mField.get_moab().tag_get_handle("AREA_CHANGE", th_area_change);
483 if (rval == MB_SUCCESS)
484 CHKERR cP.mField.get_moab().tag_delete(th_area_change);
485
486 CHKERR cP.mField.get_moab().tag_get_handle(
487 "AREA_CHANGE", 3, MB_TYPE_DOUBLE, th_area_change,
488 MB_TAG_CREAT | MB_TAG_SPARSE, def_val);
490 };
491
492 auto get_surface = [this, &bit]() {
493 Range surface;
495 surface, true);
496 Range bit_faces;
498 MBTRI, bit_faces);
499 return intersect(surface, bit_faces);
500 };
501
502 auto get_crack_front = [&]() {
503 Range crack_front;
505 crack_front, true);
506 Range bit_edges;
508 MBEDGE, bit_edges);
509 return intersect(crack_front, bit_edges);
510 };
511
512 // blockset where crack propagation takes place
513 auto get_cracked_body_blockset = [&]() {
514 Range cracked_body_tets;
516 crackedBodyBlockSetId, BLOCKSET, 3, cracked_body_tets, true);
517 return cracked_body_tets;
518 };
519
520 // build tree to search distances to body surface
521 auto build_body_skin_tree =
522 [&](boost::shared_ptr<OrientedBoxTreeTool> &tree_body_skin_ptr,
523 EntityHandle &rootset_tree_body_skin) {
525 Range last_bit_tets;
527 BitRefLevel().set(BITREFLEVEL_SIZE - 1), BitRefLevel().set(), MBTET,
528 last_bit_tets);
529 if (crackedBodyBlockSetId > 0)
530 last_bit_tets = intersect(last_bit_tets, get_cracked_body_blockset());
531 Range last_bit_tets_skin;
532 CHKERR skin.find_skin(0, last_bit_tets, false, last_bit_tets_skin);
533 CHKERR tree_body_skin_ptr->build(last_bit_tets_skin,
534 rootset_tree_body_skin);
536 };
537
538 // Set coords
539 auto set_positions = [&](const Range nodes, Tag *th = nullptr) {
541 VectorDouble coords(3 * nodes.size());
542 if (th)
543 CHKERR cP.mField.get_moab().tag_get_data(*th, nodes, &*coords.begin());
544 else
545 CHKERR cP.mField.get_moab().get_coords(nodes, &*coords.begin());
546 CHKERR cP.mField.get_moab().tag_set_data(th_front_positions, nodes,
547 &*coords.begin());
549 };
550
551 auto get_nodes = [this](Range edges) {
552 Range verts;
553 CHKERR cP.mField.get_moab().get_connectivity(edges, verts, true);
554 return verts;
555 };
556
557 auto get_surface_skin = [this, &skin](auto &surface) {
558 Range surface_skin;
559 CHKERR skin.find_skin(0, surface, false, surface_skin);
560 return surface_skin;
561 };
562
563 auto get_crack_front_ends = [this, &skin](Range crack_front) {
564 Range end_nodes;
565 CHKERR skin.find_skin(0, crack_front, false, end_nodes);
566 return end_nodes;
567 };
568
569 auto get_edge_delta = [&](const EntityHandle edge, const EntityHandle node,
572 int num_nodes;
573 const EntityHandle *conn;
574 CHKERR cP.mField.get_moab().get_connectivity(edge, conn, num_nodes, true);
575 VectorDouble3 n0(3), n1(3);
576 if (conn[1] == node) {
577 CHKERR cP.mField.get_moab().get_coords(&conn[0], 1, &n0[0]);
578 CHKERR cP.mField.get_moab().get_coords(&conn[1], 1, &n1[0]);
579 } else {
580 CHKERR cP.mField.get_moab().get_coords(&conn[1], 1, &n0[0]);
581 CHKERR cP.mField.get_moab().get_coords(&conn[0], 1, &n1[0]);
582 }
583 delta = n1 - n0;
585 };
586
587 auto get_node_vec = [](std::vector<double> &v, const int nn) {
588 return getVectorAdaptor(&v[3 * nn], 3);
589 };
590
591 auto check_if_node_is_in_range = [](const EntityHandle node,
592 const Range &ents) {
593 return (ents.find(node) != ents.end());
594 };
595
596 auto check_for_nan = [&](auto &v, std::string msg = "") {
598 for (auto d : v)
599 if (std::isnan(d) || std::isinf(d))
600 if (msg.size()) {
602 "(%s) Not number %3.4e", msg.c_str(), d);
603 } else {
605 "Not number %3.4e", d);
606 }
607
609 };
610
611 auto set_tag = [&](const EntityHandle node, auto coords) {
613 CHKERR check_for_nan(coords, "set_tag");
614 CHKERR cP.mField.get_moab().tag_set_data(th_front_positions, &node, 1,
615 &coords[0]);
617 };
618
619 auto get_str = [](auto v) { return boost::lexical_cast<std::string>(v); };
620
621 auto get_tag_data = [this](Tag th, const Range &ents) {
622 std::vector<double> v(3 * ents.size());
623 CHKERR cP.mField.get_moab().tag_get_data(th, ents, &*v.begin());
624 return v;
625 };
626
627 auto get_coords_vec = [this](const Range &ents) {
628 std::vector<double> v(3 * ents.size());
629 CHKERR cP.mField.get_moab().get_coords(ents, &*v.begin());
630 return v;
631 };
632
633 auto get_prisms = [&]() {
634 Range prisms;
635 CHKERR cP.mField.get_moab().get_adjacencies(get_surface(), 3, false, prisms,
636 moab::Interface::UNION);
637 prisms = prisms.subset_by_type(MBPRISM);
638 return prisms;
639 };
640
641 auto get_prims_surface = [&](Range prisms) {
642 Range prisms_surface;
643 for (auto prism : prisms) {
644 EntityHandle face;
645 CHKERR cP.mField.get_moab().side_element(prism, 2, 4, face);
646 prisms_surface.insert(face);
647 }
648 return prisms_surface;
649 };
650
651 auto get_random_double = []() {
652 return static_cast<double>(rand()) / static_cast<double>(RAND_MAX);
653 };
654
655 // Check if point is outside or inside of body
656 auto is_point_is_outside = [&](const auto &coords) {
657 VectorDouble3 random_ray(3);
658 random_ray[0] =
659 get_random_double() + std::numeric_limits<double>::epsilon();
660 random_ray[1] = get_random_double();
661 random_ray[2] = get_random_double();
662 random_ray /= norm_2(random_ray);
663 std::vector<double> distances_out;
664 std::vector<EntityHandle> facets_out;
665 CHKERR tree_body_skin_ptr->ray_intersect_triangles(
666 distances_out, facets_out, rootset_tree_body_skin, 1e-12, &coords[0],
667 &random_ray[0]);
668 return ((facets_out.size() % 2) == 0);
669 };
670
671 // Check until results is consistent
672 auto is_point_is_outside_with_check = [&](const auto &coords) {
673 bool is_out = is_point_is_outside(coords);
674 bool test_out;
675 do {
676 test_out = is_out;
677 is_out = is_point_is_outside(coords);
678 } while (test_out != is_out);
679 return is_out;
680 };
681
682 auto get_proj_front_disp = [&](auto front_disp,
683 auto griffith_force_projected) {
684 const double mgn_force_projected = norm_2(griffith_force_projected);
685 VectorDouble3 disp(3);
686
687 if (mgn_force_projected > std::numeric_limits<double>::epsilon()) {
688 griffith_force_projected /= mgn_force_projected;
689 noalias(disp) = factor * griffith_force_projected *
690 fmax(0, inner_prod(front_disp, griffith_force_projected));
691 } else {
692 disp.clear();
693 }
694
695 return disp;
696 };
697
698 // get average area increment
699 auto get_ave_a_change = [&](Range &crack_front_nodes,
700 std::vector<double> &area_change) {
701 double ave_a_change = 0;
702 int nn = 0;
703 for (auto node : crack_front_nodes) {
704 auto a_change = get_node_vec(area_change, nn);
705 ave_a_change += norm_2(a_change);
706 }
707 ave_a_change /= crack_front_nodes.size();
708 return ave_a_change;
709 };
710
711 auto get_ray_pout = [&](const auto coords, const auto ray) {
712 std::vector<double> distances_out;
713 std::vector<EntityHandle> facets_out;
714 double ray_length = norm_2(ray);
715 VectorDouble3 unit_ray = ray / ray_length;
716 // send the ray in the direction of the edge
717 CHKERR tree_body_skin_ptr->ray_intersect_triangles(
718 distances_out, facets_out, rootset_tree_body_skin, 1e-12, &coords[0],
719 &unit_ray[0], &ray_length);
720 if (distances_out.empty())
721 return -1.;
722 else
723 return *std::min_element(distances_out.begin(), distances_out.end());
724 };
725
727 cP.mapBitLevel["spatial_domain"], BitRefLevel().set(), MBEDGE, bit_edges);
728
729 SetFrontPositionsParameters front_params;
730
731 CHKERR get_tags();
732 sUrface = get_prims_surface(get_prisms());
733
734 Range surface_skin = get_surface_skin(sUrface);
735 Range crack_front = intersect(get_crack_front(), surface_skin);
736 Range crack_front_nodes = get_nodes(crack_front);
737 Range end_nodes = get_crack_front_ends(crack_front);
738 Range surface_skin_nodes =
739 subtract(subtract(get_nodes(surface_skin), crack_front_nodes),
740 get_nodes(intersect(fixedEdges, surface_skin)));
741 surface_skin_nodes.merge(intersect(cornerNodes, get_nodes(surface_skin)));
742
743 CHKERR GriffithForceElement::SaveGriffithForceOnTags(
744 cP.mField, th_area_change, sUrface, nullptr)();
745 CHKERR set_positions(get_nodes(get_prisms()), nullptr);
746 CHKERR build_body_skin_tree(tree_body_skin_ptr, rootset_tree_body_skin);
747
748 if (debug)
749 CHKERR save_entities("prisms_surface.vtk", get_prims_surface(get_prisms()));
750 if (debug)
751 CHKERR save_entities("surface_skin.vtk", surface_skin);
752 if (debug)
753 CHKERR save_entities("crack_front.vtk", crack_front);
754 if (debug)
755 CHKERR save_entities("crack_front_nodes.vtk", crack_front_nodes);
756 if (debug)
757 CHKERR save_entities("end_nodes.vtk", end_nodes);
758
759 auto surface_snap = [&]() {
762 surface_skin, fixedEdges, front_params.tolRelSnapToFixedEdges,
763 front_params.tolAbsSnapToFixedEdges, nullptr, debug);
765 };
766
767 auto move_front_nodes = [&]() {
769
770 auto w = get_tag_data(th_w, crack_front_nodes);
771 auto griffith_forces_projected =
772 get_tag_data(th_griffith_force_projected, crack_front_nodes);
773 auto area_change = get_tag_data(th_area_change, crack_front_nodes);
774 auto mesh_node_positions = get_coords_vec(crack_front_nodes);
775
776 const double ave_a_change =
777 get_ave_a_change(crack_front_nodes, area_change);
778
779 int nn = 0;
780 for (auto node : crack_front_nodes) {
781
782 auto front_disp = get_node_vec(w, nn);
783 auto force_projected = get_node_vec(griffith_forces_projected, nn);
784 auto a_change = get_node_vec(area_change, nn);
785 auto position = get_node_vec(mesh_node_positions, nn);
786
787 CHKERR check_for_nan(front_disp, "front_disp");
788 CHKERR check_for_nan(force_projected, "force_projected");
789 CHKERR check_for_nan(a_change, "a_change");
790
791 VectorDouble3 disp = get_proj_front_disp(front_disp, force_projected);
792
793 if (debug) {
794 CHKERR check_for_nan(position);
795 CHKERR check_for_nan(disp);
796 }
797
798 VectorDouble3 coords = position + disp;
799
800 // check if nodes on fixed edges
801 Range adj_crack_front_nodes_edges;
802 CHKERR cP.mField.get_moab().get_adjacencies(&node, 1, 1, false,
803 adj_crack_front_nodes_edges,
804 moab::Interface::UNION);
805 Range adj_crack_front_nodes_on_fixed_edges;
806 adj_crack_front_nodes_on_fixed_edges =
807 intersect(adj_crack_front_nodes_edges, fixedEdges);
808 adj_crack_front_nodes_on_fixed_edges =
809 intersect(adj_crack_front_nodes_on_fixed_edges, bit_edges);
810
811 CHKERR check_for_nan(front_disp, "front_disp");
812 CHKERR check_for_nan(coords, "coords");
813 CHKERR check_for_nan(disp, "disp");
814
815 enum Classifier {
816 FREE = 0,
817 FIXED_EDGE = 1 << 0,
818 END_NODE = 1 << 1,
819 FIXED_EDGE_KIND0 = 1 << 2,
820 FIXED_EDGE_KIND1 = 1 << 3
821 };
822
823 unsigned int kind = FREE;
824 if (check_if_node_is_in_range(node, get_nodes(fixedEdges)))
825 kind |= FIXED_EDGE;
826 if (check_if_node_is_in_range(node, end_nodes))
827 kind |= END_NODE;
828
829 if (kind == FREE) {
830
831 // check pushing other nodes, which are not end nodes
832 VectorDouble3 scaled_a_change =
833 ave_a_change * a_change /
834 (norm_2(a_change) + std::numeric_limits<double>::epsilon());
835 if (is_point_is_outside_with_check(coords)) {
836
837 coords += front_params.skinFactor * scaled_a_change;
838 CHKERR check_for_nan(coords,
839 "Pushed through surface (out by factor)");
840
841 MOFEM_LOG("CPMeshCutWorld", Sev::verbose)
842 << "Pushed through surface (out by factor): ";
843
844 } else {
845
846 VectorDouble3 tmp_coords =
847 coords + front_params.skinFactor * scaled_a_change;
848 if (is_point_is_outside_with_check(tmp_coords)) {
849 VectorDouble3 ray =
850 scaled_a_change / (norm_2(scaled_a_change) +
851 std::numeric_limits<double>::epsilon());
852 const double dist = get_ray_pout(coords, scaled_a_change);
853 if (dist > -0)
854 coords += dist * ray;
855
856 coords += front_params.skinFactor * scaled_a_change;
857 CHKERR check_for_nan(
858 coords, "Pushing through surface (out by skin factor)");
859
860 MOFEM_LOG("CPMeshCutWorld", Sev::verbose)
861 << "Pushing through surface (out by skin factor): ";
862 }
863 }
864 }
865
866 if (kind & END_NODE) {
867
868 double min_dist = front_params.skinFactor * ave_a_change;
869
870 if (kind & FIXED_EDGE) {
871
872 VectorDouble3 pt_out(3);
873 CHKERR cP.mField.getInterface<Tools>()->findMinDistanceFromTheEdges(
874 &coords[0], 1, fixedEdges, &min_dist, &pt_out[0]);
875 if (min_dist < front_params.skinFactor * ave_a_change)
876 coords = pt_out;
877
878 CHKERR check_for_nan(coords, "End node on the fixed edge");
879
880 MOFEM_LOG("CPMeshCutWorld", Sev::verbose)
881 << "End node on the fixed edge: ";
882
883 } else {
884 // Find corner edge and check if crack can go through it
885
886 VectorDouble3 pt_out(3);
887 CHKERR cP.mField.getInterface<Tools>()->findMinDistanceFromTheEdges(
888 &coords[0], 1, fixedEdges, &min_dist, &pt_out[0]);
889 if (min_dist < front_params.skinFactor * ave_a_change) {
890 VectorDouble3 ray = pt_out - coords;
891 if (inner_prod(a_change, ray) > 0) {
892 VectorDouble3 scaled_ray =
893 front_params.skinFactor * ave_a_change * ray / norm_2(ray);
894 coords += ray + scaled_ray;
895 CHKERR check_for_nan(
896 coords, "Pushing through surface dist (close to edge)");
897 MOFEM_LOG("CPMeshCutWorld", Sev::verbose)
898 << "Pushing through surface dist (close to edge): ";
899 }
900 }
901 }
902 }
903
904 if (!adj_crack_front_nodes_on_fixed_edges.empty()) {
905
906 const int nb_edges_on_fix_edges =
907 adj_crack_front_nodes_on_fixed_edges.size() - 1;
908 adj_crack_front_nodes_on_fixed_edges =
909 intersect(adj_crack_front_nodes_on_fixed_edges, surface_skin);
910
911 // adj node edge is on fixed edge and on the skin
912 if (!adj_crack_front_nodes_on_fixed_edges.empty()) {
913
914 const int nb_edges_on_fix_edges_and_surface_skin =
915 adj_crack_front_nodes_on_fixed_edges.size();
916 const int diff_nb_edges =
917 nb_edges_on_fix_edges - nb_edges_on_fix_edges_and_surface_skin;
918
919 // If is >- in that case diff_nb_edges = 1
920 // If is >| in that case diff_nb_edges > 1
921
922 if (diff_nb_edges > 1) {
923
924 VectorDouble3 edge_delta;
925 CHKERR get_edge_delta(adj_crack_front_nodes_on_fixed_edges[0], node,
926 edge_delta);
927 coords += front_params.cornerFactor * edge_delta;
928 MOFEM_LOG("CPMeshCutWorld", Sev::verbose)
929 << "Node on fixed edge (crack meeting fixed edge >|): ";
930 CHKERR check_for_nan(
931 coords, "Node on fixed edge (crack meeting fixed edge >|)");
932
933 kind |= FIXED_EDGE_KIND0;
934 }
935
936 } else if ((kind & END_NODE) && (kind & FIXED_EDGE)) {
937
938 // It s =|>
939
940 // adj edge is on fixed edge but that edge is not on the skin
941 Range adj_edges_on_skin =
942 intersect(adj_crack_front_nodes_edges, surface_skin);
943 if (!adj_edges_on_skin.empty()) {
944
945 VectorDouble3 edge_delta;
946 CHKERR get_edge_delta(adj_edges_on_skin[0], node, edge_delta);
947 coords += front_params.cornerFactor * edge_delta;
948 CHKERR check_for_nan(
949 coords, "Node on fixed edge (crack crossing fixed edge =|>)");
950
951 MOFEM_LOG("CPMeshCutWorld", Sev::verbose)
952 << "Node on fixed edge (crack crossing fixed edge =|>): ";
953
954 kind |= FIXED_EDGE_KIND1;
955 }
956 }
957 }
958
959 if (
960
961 (kind & END_NODE) &&
962
963 !(kind & FIXED_EDGE_KIND0) &&
964
965 !(kind & FIXED_EDGE_KIND1)) {
966
967 // is a end node on the skin, just go outside, in the direction
968 // of crack front edged
969 Range adj_edges = intersect(adj_crack_front_nodes_edges, crack_front);
970 VectorDouble3 edge_delta;
971 CHKERR get_edge_delta(adj_edges[0], node, edge_delta);
972 coords += edge_delta * front_params.endNodeEdgeDeltaFactor;
973 CHKERR check_for_nan(coords, "End node on the skin");
974
975 // front_params.skinFactor;
976 MOFEM_LOG("CPMeshCutWorld", Sev::verbose) << "End node on the skin: ";
977 }
978
979 MOFEM_LOG_C("CPMeshCutWorld", Sev::verbose, "Disp front %s at node %lu",
980 get_str(front_disp).c_str(), node);
981
982 CHKERR set_tag(node, coords);
983
984 ++nn;
985 }
986
988 };
989
990 auto smooth_front = [&]() {
992 Range smoothed_nodes = subtract(crack_front_nodes, end_nodes);
993 std::vector<double> new_positions(3 * smoothed_nodes.size());
994 int nn = 0;
995 for (auto node : smoothed_nodes) {
996 Range edges;
997 CHKERR cP.mField.get_moab().get_adjacencies(&node, 1, 1, false, edges);
998 edges = intersect(edges, crack_front);
999 Range edges_nodes;
1000 CHKERR cP.mField.get_moab().get_connectivity(edges, edges_nodes, true);
1001 edges_nodes.erase(node);
1002 if (edges_nodes.size() != 2)
1003 SETERRQ1(PETSC_COMM_WORLD, MOFEM_DATA_INCONSISTENCY,
1004 "Expected two nodes but have %d", edges_nodes.size());
1005 VectorDouble3 node0(3), node1(3), node2(3);
1006 EntityHandle n0 = edges_nodes[0];
1007 EntityHandle n2 = edges_nodes[1];
1008 CHKERR cP.mField.get_moab().tag_get_data(th_front_positions, &n0, 1,
1009 &node0[0]);
1010 CHKERR cP.mField.get_moab().tag_get_data(th_front_positions, &node, 1,
1011 &node1[0]);
1012 CHKERR cP.mField.get_moab().tag_get_data(th_front_positions, &n2, 1,
1013 &node2[0]);
1014 VectorDouble3 delta0 = node0 - node1;
1015 VectorDouble3 delta2 = node2 - node1;
1016 const double l0 = norm_2(delta0);
1017 const double l2 = norm_2(delta2);
1018
1019 if ((l0 / (l0 + l2)) < 0.5)
1020 node1 += 2 * (0.5 - (l0 / (l0 + l2))) * (delta2 / l2) * l0;
1021 else
1022 node1 += 2 * (0.5 - (l2 / (l0 + l2))) * (delta0 / l0) * l2;
1023
1024 VectorAdaptor new_coords = get_node_vec(new_positions, nn);
1025 new_coords = node1;
1026 ++nn;
1027 }
1028 CHKERR cP.mField.get_moab().tag_set_data(th_front_positions, smoothed_nodes,
1029 &*new_positions.begin());
1030
1032 };
1033
1034 auto delete_mofem_meshsets = [&]() {
1036 CHKERR cP.mField.getInterface<MeshsetsManager>()->deleteMeshset(
1038 CHKERR cP.mField.getInterface<MeshsetsManager>()->deleteMeshset(
1040 CHKERR cP.mField.getInterface<MeshsetsManager>()->deleteMeshset(
1043 };
1044
1045 auto build_surface = [&]() {
1047 // Create cutting surface from existing crack
1048 Tag th_interface_side;
1049 CHKERR cP.mField.get_moab().tag_get_handle("INTERFACE_SIDE",
1050 th_interface_side);
1051 Range new_surface;
1052 // Take nodes from top of the prims
1053 for (auto face : sUrface) {
1054
1055 int side;
1056 CHKERR cP.mField.get_moab().tag_get_data(th_interface_side, &face, 1,
1057 &side);
1058 int num_nodes;
1059 const EntityHandle *conn;
1060 CHKERR cP.mField.get_moab().get_connectivity(face, conn, num_nodes, true);
1061
1062 if (debug)
1063 if (!(conn[0] != conn[1] && conn[0] != conn[2] && conn[1] != conn[2]))
1064 SETERRQ(PETSC_COMM_WORLD, MOFEM_DATA_INCONSISTENCY,
1065 "Imposible trinagle");
1066
1067 MatrixDouble3by3 coords(3, 3);
1068 CHKERR cP.mField.get_moab().tag_get_data(th_front_positions, conn, 3,
1069 &coords(0, 0));
1070 if (debug)
1071 CHKERR check_for_nan(coords.data());
1072
1073 EntityHandle new_verts[3] = {0, 0, 0};
1074 for (auto nn : {0, 1, 2}) {
1075 auto node_map = verts_map.find(conn[nn]);
1076 if (node_map == verts_map.end()) {
1077 EntityHandle new_node;
1078 CHKERR cP.mField.get_moab().create_vertex(&coords(nn, 0), new_node);
1079 verts_map[conn[nn]] = new_node;
1080 }
1081 new_verts[nn] = verts_map[conn[nn]];
1082 }
1083 for (auto nn : {0, 1, 2})
1084 if (new_verts[nn] == 0)
1085 SETERRQ(PETSC_COMM_WORLD, MOFEM_DATA_INCONSISTENCY,
1086 "Imposible trinagle");
1087 if (side == 1) {
1088 EntityHandle n0 = new_verts[0];
1089 new_verts[0] = new_verts[1];
1090 new_verts[1] = n0;
1091 }
1092 EntityHandle ele;
1093 CHKERR cP.mField.get_moab().create_element(MBTRI, new_verts, 3, ele);
1094 new_surface.insert(ele);
1095 }
1096 sUrface.swap(new_surface);
1098 };
1099
1100 auto build_front = [&]() {
1102 Range new_surface_skin;
1103 CHKERR skin.find_skin(0, sUrface, false, new_surface_skin);
1104 Range new_front;
1105 // Creating crack front from existing crack
1106 for (auto f : crack_front) {
1107 int num_nodes;
1108 const EntityHandle *conn;
1109 CHKERR cP.mField.get_moab().get_connectivity(f, conn, num_nodes, true);
1110 EntityHandle new_verts[num_nodes];
1111 for (int nn = 0; nn != 2; nn++) {
1112 auto node_map = verts_map.find(conn[nn]);
1113 if (node_map != verts_map.end())
1114 new_verts[nn] = node_map->second;
1115 else
1116 SETERRQ(PETSC_COMM_WORLD, MOFEM_DATA_INCONSISTENCY,
1117 "Expected that vertex should be already created");
1118 }
1119 Range edges;
1120 CHKERR cP.mField.get_moab().get_adjacencies(new_verts, 2, 1, true, edges,
1121 moab::Interface::INTERSECT);
1122 edges = intersect(edges, new_surface_skin);
1123 new_front.merge(edges);
1124 }
1125 fRont.swap(new_front);
1127 };
1128
1129 auto edges_snap = [&]() {
1131 Range new_surface_skin;
1132 CHKERR skin.find_skin(0, sUrface, false, new_surface_skin);
1134 new_surface_skin, fixedEdges, front_params.tolRelSnapToFixedEdges,
1135 front_params.tolAbsSnapToFixedEdges, nullptr, debug);
1137 };
1138
1139 auto delete_meshsets_with_crack_surfaces = [&]() {
1141 BitRefLevel mask;
1142 for (int ll = 0; ll != 19; ++ll)
1143 mask.set(ll);
1146 };
1147
1148 CHKERR surface_snap();
1149 CHKERR move_front_nodes();
1150 CHKERR smooth_front();
1151 CHKERR delete_mofem_meshsets();
1152 CHKERR build_surface();
1153 if (debug)
1154 CHKERR save_entities("set_new_crack_surface.vtk", sUrface);
1155 CHKERR build_front();
1156 CHKERR edges_snap();
1157 if (debug)
1158 CHKERR save_entities("set_new_front.vtk", fRont);
1159 CHKERR delete_meshsets_with_crack_surfaces();
1160
1161 // Set cutting inteface
1162 const_cast<Range &>(cutMeshPtr->getSurface()).clear();
1163 const_cast<Range &>(cutMeshPtr->getFront()).clear();
1166
1167 if (!save_mesh.empty())
1168 CHKERR save_entities(save_mesh, unite(sUrface, fRont));
1169
1171}
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:440
const double v
phase velocity of light in medium (cm/ns)
const Tensor1_Expr< const dTensor0< T, Dim, i >, typename promote< T, double >::V, Dim, i > d(const Tensor0< T * > &a, const Index< i, Dim > index, const Tensor1< int, Dim > &d_ijk, const Tensor1< double, Dim > &d_xyz)
Definition: dTensor0.hpp:27
static MoFEMErrorCodeGeneric< moab::ErrorCode > rval
Definition: Exceptions.hpp:74
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:56
VectorBoundedArray< double, 3 > VectorDouble3
Definition: Types.hpp:92
VectorShallowArrayAdaptor< double > VectorAdaptor
Definition: Types.hpp:115
UBlasVector< double > VectorDouble
Definition: Types.hpp:68
MatrixBoundedArray< double, 9 > MatrixDouble3by3
Definition: Types.hpp:105
auto getVectorAdaptor(T1 ptr, const size_t n)
Get Vector adaptor.
Definition: Templates.hpp:31
double w(const double g, const double t)
static constexpr double delta
MoFEMErrorCode getOptions()
Get options from command line.
Definition: CPMeshCut.cpp:85
MoFEMErrorCode setFront(const Range surface)
const Range & getFront() const
MoFEMErrorCode setSurface(const Range surface)
set surface entities
MoFEMErrorCode snapSurfaceToEdges(const Range surface_edges, const Range fixed_edges, const double rel_tol, const double abs_tol, Tag th=nullptr, const bool debug=false)
Interface for managing meshsets containing materials and boundary conditions.
Auxiliary tools.
Definition: Tools.hpp:19

◆ refineAndSplit()

MoFEMErrorCode FractureMechanics::CPMeshCut::refineAndSplit ( const int  verb = 1,
const bool  debug = false 
)

Definition at line 2155 of file CPMeshCut.cpp.

2155 {
2156 MoFEM::Interface &m_field = cP.mField;
2158
2159 // clear tags data created by previous use of prism interface
2160 Tag th_interface_side;
2161 rval = m_field.get_moab().tag_get_handle("INTERFACE_SIDE", th_interface_side);
2162 if (rval == MB_SUCCESS)
2163 CHKERR m_field.get_moab().tag_delete(th_interface_side);
2164
2165 auto delete_mofem_meshsets = [&]() {
2167 CHKERR cP.mField.getInterface<MeshsetsManager>()->deleteMeshset(
2169 CHKERR cP.mField.getInterface<MeshsetsManager>()->deleteMeshset(
2172 };
2173
2174 BitRefLevel bit0 = cP.mapBitLevel["mesh_cut"];
2175 BitRefLevel bit1 = cP.mapBitLevel["spatial_domain"];
2176
2177 CHKERR delete_mofem_meshsets();
2178 // find body skin
2179 CHKERR findBodySkin(bit0, verb, debug, "out_body_skin_bit0.vtk");
2180 // find/get crack surface
2181 CHKERR findCrack(bit0, verb, debug);
2182
2183 // refine at crack trip
2184 if (cP.refAtCrackTip > 0) {
2185 if (cP.doCutMesh) {
2186 SETERRQ(
2188 "Refinement at the crack front after the crack insertion ('my_ref') "
2189 "should not be used when the mesh cutting is on, use "
2190 "'ref_before_cut' or 'ref_before_trim' instead");
2191 }
2193 CHKERR delete_mofem_meshsets();
2194 CHKERR findBodySkin(bit0, verb, debug, "out_body_skin_bit0_ref.vtk");
2195 CHKERR findCrack(bit0, verb, debug);
2196 }
2197
2198 if (debug)
2199 CHKERR bitRefPtr->writeBitLevelByType(bit0, BitRefLevel().set(), MBTET,
2200 "out_mesh_after_refine.vtk", "VTK",
2201 "");
2202
2203 if (!cP.ignoreContact) {
2204
2205 if (debug) {
2206 EntityHandle meshset_interface;
2208 meshset_interface);
2209 CHKERR cP.mField.get_moab().write_file(
2210 "crack_surface_before_contact_split.vtk", "VTK", "",
2211 &meshset_interface, 1);
2212 }
2213
2214 // split faces to insert contact elements
2216 if (debug) {
2217
2218 EntityHandle meshset_interface;
2220 meshset_interface);
2221 CHKERR cP.mField.get_moab().write_file(
2222 "crack_surface_after_contact_split.vtk", "VTK", "",
2223 &meshset_interface, 1);
2224
2226 cP.mapBitLevel["mesh_cut"], BitRefLevel().set(), MBTET,
2227 "out_mesh_after_contact_split.vtk", "VTK", "");
2229 cP.mapBitLevel["mesh_cut"], BitRefLevel().set(), MBPRISM,
2230 "out_mesh_contact_prisms.vtk", "VTK", "");
2231 }
2232 }
2233
2234 // split faces to create crack
2235 CHKERR splitFaces(verb, debug);
2236
2237 if (!cP.ignoreContact) {
2238 findContactFromPrisms(bit1, verb, debug);
2239 }
2240
2241 if (!cP.ignoreContact) {
2243 }
2244
2245 if (debug) {
2246 CHKERR bitRefPtr->writeBitLevelByType(bit1, BitRefLevel().set(), MBTET,
2247 "out_mesh_after_split.vtk", "VTK",
2248 "");
2249 CHKERR bitRefPtr->writeBitLevelByType(bit1, BitRefLevel().set(), MBPRISM,
2250 "out_mesh_after_split_prisms.vtk",
2251 "VTK", "");
2252 }
2253
2254 // get information about topology
2255 CHKERR delete_mofem_meshsets();
2256
2258 CHKERR findBodySkin(bit1, verb, debug, "out_body_skin_bit1.vtk");
2259 CHKERR findCrackFromPrisms(bit1, verb, debug);
2261
2263}
@ MOFEM_IMPOSIBLE_CASE
Definition: definitions.h:35
MoFEMErrorCode findBodySkin(const BitRefLevel bit, const int verb=1, const bool debug=false, std::string file_name="out_body_skin.vtk")
find body skin
Definition: CPMeshCut.cpp:1383
MoFEMErrorCode findCrackFromPrisms(const BitRefLevel bit, const int verb=1, const bool debug=false)
get crack front
Definition: CPMeshCut.cpp:1513
MoFEMErrorCode findContactFromPrisms(const BitRefLevel bit, const int verb=1, const bool debug=false)
get contact elements
Definition: CPMeshCut.cpp:1987
MoFEMErrorCode refineCrackTip(const int front_id, const int verb=1, const bool debug=false)
refine elements at crack tip
Definition: CPMeshCut.cpp:2031
MoFEMErrorCode insertContactInterface(const int verb=1, const bool debug=false)
insert contact interface
Definition: CPMeshCut.cpp:1860
MoFEMErrorCode findCrack(const BitRefLevel bit, const int verb=1, const bool debug=false)
get crack front
Definition: CPMeshCut.cpp:1447
MoFEMErrorCode addMortarContactPrisms(const int verb=1, const bool debug=false)
insert contact interface
Definition: CPMeshCut.cpp:1749
MoFEMErrorCode getFrontEdgesAndElements(const BitRefLevel bit, const int verb=1, const bool debug=false)
get crack front edges and finite elements
Definition: CPMeshCut.cpp:2335
MoFEMErrorCode splitFaces(const int verb=1, const bool debug=false)
split crack faces
Definition: CPMeshCut.cpp:1622

◆ refineCrackTip()

MoFEMErrorCode FractureMechanics::CPMeshCut::refineCrackTip ( const int  front_id,
const int  verb = 1,
const bool  debug = false 
)

refine elements at crack tip

Definition at line 2031 of file CPMeshCut.cpp.

2032 {
2033 MoFEM::Interface &m_field = cP.mField;
2034 moab::Interface &moab = m_field.get_moab();
2035 MeshRefinement *mesh_refiner_ptr;
2037 CHKERR m_field.getInterface(mesh_refiner_ptr);
2038
2039 BitRefLevel bit = cP.mapBitLevel["mesh_cut"];
2040
2041 std::size_t idx = 0;
2042 while (!bit.test(idx))
2043 ++idx;
2044 BitRefLevel bit_ref = BitRefLevel().set(idx + 1);
2045
2046 // loop to refine at crack tip
2047 for (int ll = 0; ll != cP.refAtCrackTip; ll++) {
2048
2049 // get tets at bit level
2050 Range tets_level; // test at level
2052 MBTET, tets_level);
2053 Range edges_level; // edges at level
2055 MBEDGE, edges_level);
2056 if (!meshsetMngPtr->checkMeshset(front_id, SIDESET)) {
2057 PetscPrintf(PETSC_COMM_WORLD, "Sideset %d with crack front not found\n",
2060 }
2061 Range crack_front; // crack front edges at level
2063 crack_front, true);
2064 crack_front = intersect(crack_front, edges_level);
2065 Range crack_front_nodes; // nodes at crack front
2066 CHKERR moab.get_connectivity(crack_front, crack_front_nodes, true);
2067
2068 Range crack_front_nodes_tets; // tets adjacent to crack front nodes
2069 CHKERR m_field.get_moab().get_adjacencies(crack_front_nodes, 3, false,
2070 crack_front_nodes_tets,
2071 moab::Interface::UNION);
2072 crack_front_nodes_tets = intersect(crack_front_nodes_tets, tets_level);
2073 Range edges_to_refine; // edges which going to be refined
2074 CHKERR m_field.get_moab().get_adjacencies(crack_front_nodes_tets, 1, false,
2075 edges_to_refine,
2076 moab::Interface::UNION);
2077 edges_to_refine = intersect(edges_to_refine, edges_level);
2078 if (edges_to_refine.empty())
2079 continue;
2080
2081 // create vertices in middle of refined edges
2082 CHKERR mesh_refiner_ptr->addVerticesInTheMiddleOfEdges(
2083 edges_to_refine, bit_ref, VERBOSE);
2084 // refine tets
2085 CHKERR mesh_refiner_ptr->refineTets(tets_level, bit_ref, false);
2086
2087 // update meshsets by new entities
2088 // Meshsets are updated such that: if parent is in the meshset then its
2089 // chiled is in the meshset to.
2090 for (_IT_CUBITMESHSETS_FOR_LOOP_(m_field, cubit_it)) {
2091 EntityHandle cubit_meshset = cubit_it->meshset;
2093 cubit_meshset, bit, BitRefLevel().set(), bit_ref, bit_ref,
2094 cubit_meshset, MBMAXTYPE, true);
2095 }
2096
2097 // update meshsets with partition mesh
2098 if (cP.isPartitioned) {
2099 ParallelComm *pcomm = ParallelComm::get_pcomm(&moab, MYPCOMM_INDEX);
2100 Tag part_tag = pcomm->part_tag();
2101 Range tagged_sets;
2102 CHKERR m_field.get_moab().get_entities_by_type_and_tag(
2103 0, MBENTITYSET, &part_tag, NULL, 1, tagged_sets,
2104 moab::Interface::UNION);
2105 for (Range::iterator mit = tagged_sets.begin(); mit != tagged_sets.end();
2106 mit++) {
2107 int part;
2108 CHKERR moab.tag_get_data(part_tag, &*mit, 1, &part);
2110 *mit, bit, BitRefLevel().set(), bit_ref, bit_ref, *mit, MBTET,
2111 true);
2112 Range ents3d;
2113 CHKERR moab.get_entities_by_dimension(*mit, 3, ents3d, true);
2114 for (Range::iterator eit = ents3d.begin(); eit != ents3d.end(); eit++) {
2115 CHKERR moab.tag_set_data(part_tag, &*eit, 1, &part);
2116 }
2117 for (int dd = 0; dd != 3; dd++) {
2118 Range ents;
2119 CHKERR moab.get_adjacencies(ents3d, dd, false, ents,
2120 moab::Interface::UNION);
2121 for (Range::iterator eit = ents.begin(); eit != ents.end(); eit++) {
2122 CHKERR moab.tag_set_data(part_tag, &*eit, 1, &part);
2123 }
2124 }
2125 }
2126 }
2127
2128 // shift bits, if ents in on bit level 0, and only on that bit it is deleted
2129 BitRefLevel shift_mask;
2130 for (int ll = 0; ll != 19; ++ll)
2131 shift_mask.set(ll);
2132 CHKERR bitRefPtr->shiftRightBitRef(1, shift_mask);
2133
2134 if (debug) {
2135 if (m_field.get_comm_rank() == 0) {
2136 std::string name;
2137 name = "out_crack_surface_level_" +
2138 boost::lexical_cast<std::string>(ll) + ".vtk";
2139 Range crack_faces;
2141 crack_faces, true);
2142 EntityHandle meshset;
2143 CHKERR m_field.get_moab().create_meshset(MESHSET_SET, meshset);
2144 CHKERR m_field.get_moab().add_entities(meshset, crack_faces);
2145 CHKERR m_field.get_moab().write_file(name.c_str(), "VTK", "", &meshset,
2146 1);
2147 CHKERR m_field.get_moab().delete_entities(&meshset, 1);
2148 }
2149 }
2150 }
2151
2153}
Mesh refinement interface.
MoFEMErrorCode refineTets(const EntityHandle meshset, const BitRefLevel &bit, int verb=QUIET, const bool debug=false)
refine TET in the meshset
MoFEMErrorCode addVerticesInTheMiddleOfEdges(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

◆ refineMesh()

MoFEMErrorCode FractureMechanics::CPMeshCut::refineMesh ( const Range vol,
const bool  make_front,
const int  verb = QUIET,
const bool  debug = false 
)

Refine mesh close to crack surface and crack front.

Parameters
verb
volpointer to volumes to be refined (can be set to null if set somwere else)
make_frontmake crack front from skin surface
update_frontupdate crack surface information only in vicinity of crack front
debug
Returns
MoFEMErrorCode

Definition at line 218 of file CPMeshCut.cpp.

219 {
221
222 if (vol) {
224 *vol, BitRefLevel().set(BITREFLEVEL_SIZE - 1), true, verb);
226 } else
227 // Set cutted front
229
231
232 if (make_front)
233 // make crack front from mesh skin
235
236 auto save_entities = [&](const std::string name, Range ents) {
238 if (!ents.empty()) {
239 EntityHandle meshset;
240 CHKERR cP.mField.get_moab().create_meshset(MESHSET_SET, meshset);
241 CHKERR cP.mField.get_moab().add_entities(meshset, ents);
242 CHKERR cP.mField.get_moab().write_file(name.c_str(), "VTK", "", &meshset,
243 1);
244 CHKERR cP.mField.get_moab().delete_entities(&meshset, 1);
245 }
247 };
248
249 // build tree
251 // refine and calculate level sets
252 BitRefLevel bit = cP.mapBitLevel["mesh_cut"];
253 std::size_t idx = 0;
254 while (!bit.test(idx))
255 ++idx;
257 verb, debug);
258
259 auto set_last_bit_ref = [&]() {
262 BitRefLevel bits;
263 for (int ll = 1; ll <= nbRefBeforCut + nbRefBeforTrim + 2; ++ll)
264 bits.set(ll);
265 BitRefLevel mask = bits;
267 BitRefLevel().set(20), true,
269 Range ents;
270 CHKERR bitRefPtr->getEntitiesByRefLevel(bits, mask, ents, verb);
271 ents = subtract(ents, ents.subset_by_type(MBENTITYSET));
272 CHKERR bitRefPtr->addBitRefLevel(ents, BitRefLevel().set(20));
273 }
275 };
276 CHKERR set_last_bit_ref();
277
278 auto shift_after_ref = [&]() {
281 BitRefLevel mask;
282 for (int ll = 0; ll <= nbRefBeforCut + nbRefBeforTrim + 2; ++ll)
283 mask.set(ll);
285 verb, MF_NOT_THROW);
286 }
288 };
289 CHKERR shift_after_ref();
290
291 // save coordinates to be later restored
293
294 if (debug)
296 BitRefLevel().set(), MBTET, "all_after_ref_bits.vtk", "VTK", "");
297
299}
@ QUIET
Definition: definitions.h:208
MoFEMErrorCode setBitRefLevel(const Range &ents, const BitRefLevel bit, const bool only_tets=true, int verb=0, Range *adj_ents_ptr=nullptr) const
add entities to database and set bit ref level
MoFEMErrorCode addBitRefLevel(const Range &ents, const BitRefLevel bit, int verb=QUIET) const
add bit ref level to ref entity
MoFEMErrorCode getEntitiesByRefLevel(const BitRefLevel bit, const BitRefLevel mask, const EntityHandle meshset, const int verb=QUIET) const
add all ents from ref level given by bit to meshset
MoFEMErrorCode setMeshOrgCoords()
Definition: CPMeshCut.cpp:2392
MoFEMErrorCode getMeshOrgCoords()
Definition: CPMeshCut.cpp:2420
MoFEMErrorCode setVolume(const BitRefLevel &bit)
Set cutting volume from bit ref level.
Definition: CPMeshCut.cpp:209
MoFEMErrorCode writeEntitiesAllBitLevelsByType(const BitRefLevel mask, const EntityType type, const char *file_name, const char *file_type, const char *options)
Write all entities by bit levels and type.
MoFEMErrorCode buildTree()
build tree
MoFEMErrorCode setVolume(const Range volume)
set volume entities
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.

◆ setCutSurfaceFromFile()

MoFEMErrorCode FractureMechanics::CPMeshCut::setCutSurfaceFromFile ( )

Set the cut surface from file.

Load cutting mesh from file.

Returns
MoFEMErrorCode

Definition at line 301 of file CPMeshCut.cpp.

301 {
303 if (!cutSurfMeshName.empty()) {
304 EntityHandle meshset_from_file;
306 meshset_from_file);
307 // FIXME:
308 // if its cubit file (.cub) it will show the error :
309 // [0]MOAB ERROR: -m-------------------- Error Message
310 // ------------------------------------ [0]MOAB ERROR: Non-geometric entity
311 // provided! [0]MOAB ERROR: set_sense() line 741 in src/GeomTopoTool.cpp
312 // Failed to restore topology
313 //
314 // but the program will continue
315 CHKERR cP.mField.get_moab().load_file(cutSurfMeshName.c_str(),
316 &meshset_from_file, "");
317 cutSurfMeshName.clear();
318 }
320}

◆ setFixEdgesAndCorners()

MoFEMErrorCode FractureMechanics::CPMeshCut::setFixEdgesAndCorners ( const BitRefLevel  bit)

Definition at line 1173 of file CPMeshCut.cpp.

1173 {
1175 fixedEdges.clear();
1176 cornerNodes.clear();
1177
1180 fixedEdges, true);
1181
1184 SETERRQ(PETSC_COMM_SELF, MOFEM_INVALID_DATA,
1185 "There are both BLOCKSET and NODESET with ID equal to input "
1186 "parameter -vertex_block_set");
1187
1190 cornerNodes, true);
1191
1194 cornerNodes, true);
1195
1197 fixedEdges);
1199 cornerNodes);
1200
1202}
@ NODESET
Definition: definitions.h:146
@ MOFEM_INVALID_DATA
Definition: definitions.h:36

◆ setMeshOrgCoords()

MoFEMErrorCode FractureMechanics::CPMeshCut::setMeshOrgCoords ( )
private

Definition at line 2392 of file CPMeshCut.cpp.

2392 {
2394
2395 BitRefLevel mask;
2396 for (int ll = 20; ll != BITREFLEVEL_SIZE; ++ll)
2397 mask.set(ll);
2398
2399 // Get original mesh
2400 Range org_verts;
2402 MBVERTEX, org_verts);
2403
2404 // Get original positions from tag
2405 Tag th;
2406 rval = cP.mField.get_moab().tag_get_handle("ORG_POSITION", th);
2407 if (rval == MB_SUCCESS) {
2408 originalCoords.resize(org_verts.size() * 3);
2409 CHKERR cP.mField.get_moab().tag_get_data(th, org_verts,
2410 &*originalCoords.begin());
2411 // Set cords
2412 CHKERR cP.mField.get_moab().set_coords(org_verts, &*originalCoords.begin());
2413 } else if (rval != MB_TAG_NOT_FOUND)
2415 "Can not get tag handle");
2416
2418}

◆ setVolume()

MoFEMErrorCode FractureMechanics::CPMeshCut::setVolume ( const BitRefLevel bit)

Set cutting volume from bit ref level.

Parameters
bit
Returns
MoFEMErrorCode

Definition at line 209 of file CPMeshCut.cpp.

209 {
211 Range bit_tets;
213 MBTET, bit_tets);
214 CHKERR cutMeshPtr->setVolume(bit_tets);
216}

◆ splitFaces()

MoFEMErrorCode FractureMechanics::CPMeshCut::splitFaces ( const int  verb = 1,
const bool  debug = false 
)

split crack faces

Definition at line 1622 of file CPMeshCut.cpp.

1622 {
1623 MoFEM::Interface &m_field = cP.mField;
1624 moab::Interface &moab = m_field.get_moab();
1625 PrismInterface *prism_interface;
1627 CHKERR m_field.getInterface(prism_interface);
1628
1629 BitRefLevel bit = cP.mapBitLevel["mesh_cut"];
1630 BitRefLevel bit_split = cP.mapBitLevel["spatial_domain"];
1631
1632 // split sides
1633 EntityHandle meshset_interface;
1635 PetscPrintf(PETSC_COMM_WORLD, "Sideset %d with crack surface not found\n",
1637 Range ents;
1639 CHKERR bitRefPtr->addBitRefLevel(ents, bit_split);
1641 } else {
1642
1643 // clear tags data created by previous use of prism interface
1644 Tag th_interface_side;
1645 rval =
1646 m_field.get_moab().tag_get_handle("INTERFACE_SIDE", th_interface_side);
1647 if (rval == MB_SUCCESS)
1648 rval = m_field.get_moab().tag_delete(th_interface_side);
1649
1650 // create meshset and insert entities from bit level to it
1651 EntityHandle bit_meshset;
1652 CHKERR m_field.get_moab().create_meshset(MESHSET_SET, bit_meshset);
1653
1655 MBTET, bit_meshset);
1657 MBPRISM, bit_meshset);
1658
1659 auto save_mesh = [&](const std::string file, const EntityHandle meshset) {
1661 if (m_field.get_comm_rank() == 0)
1662 CHKERR moab.write_file(file.c_str(), "VTK", "", &meshset, 1);
1664 };
1665
1667 meshset_interface);
1668 if (debug)
1669 save_mesh("crack_surface_before_split.vtk", meshset_interface);
1670
1671 CHKERR prism_interface->getSides(meshset_interface, bit, true, QUIET);
1672 CHKERR prism_interface->splitSides(bit_meshset, bit_split,
1673 meshset_interface, true, true, QUIET);
1674 // add refined ent to cubit meshsets
1675 for (_IT_CUBITMESHSETS_FOR_LOOP_(m_field, cubit_it)) {
1676 EntityHandle update_meshset = cubit_it->meshset;
1678 update_meshset, bit, BitRefLevel().set(), bit_split, bit_split,
1679 update_meshset, MBMAXTYPE, true);
1680 }
1681 CHKERR moab.delete_entities(&bit_meshset, 1);
1682
1683 if (cP.isPartitioned) {
1684 ParallelComm *pcomm = ParallelComm::get_pcomm(&moab, MYPCOMM_INDEX);
1685 Tag part_tag = pcomm->part_tag();
1686 Range tagged_sets;
1687 CHKERR m_field.get_moab().get_entities_by_type_and_tag(
1688 0, MBENTITYSET, &part_tag, NULL, 1, tagged_sets,
1689 moab::Interface::UNION);
1690 for (Range::iterator mit = tagged_sets.begin(); mit != tagged_sets.end();
1691 mit++) {
1692 int part;
1693 CHKERR moab.tag_get_data(part_tag, &*mit, 1, &part);
1695 *mit, bit, BitRefLevel().set(), bit_split, bit_split, *mit, MBTET,
1696 true);
1697 Range ents3d;
1698 CHKERR moab.get_entities_by_dimension(*mit, 3, ents3d, true);
1699 for (Range::iterator eit = ents3d.begin(); eit != ents3d.end(); eit++) {
1700 CHKERR moab.tag_set_data(part_tag, &*eit, 1, &part);
1701 }
1702 for (int dd = 0; dd != 3; dd++) {
1703 Range ents;
1704 CHKERR moab.get_adjacencies(ents3d, dd, false, ents,
1705 moab::Interface::UNION);
1706 for (Range::iterator eit = ents.begin(); eit != ents.end(); eit++) {
1707 CHKERR moab.tag_set_data(part_tag, &*eit, 1, &part);
1708 }
1709 }
1710 }
1711 }
1712
1713 if (debug) {
1714 Range tets_level;
1716 bit_split, BitRefLevel().set(), MBTET, tets_level);
1717 Range prisms_level;
1719 bit_split, BitRefLevel().set(), MBPRISM, prisms_level);
1720 EntityHandle meshset_level;
1721 CHKERR moab.create_meshset(MESHSET_SET, meshset_level);
1722 CHKERR moab.add_entities(meshset_level, tets_level);
1723 CHKERR moab.write_file("out_tets_level.vtk", "VTK", "", &meshset_level,
1724 1);
1725 CHKERR moab.delete_entities(&meshset_level, 1);
1726 EntityHandle meshset_prims_level;
1727 CHKERR moab.create_meshset(MESHSET_SET, meshset_prims_level);
1728 CHKERR moab.add_entities(meshset_prims_level, prisms_level);
1729 CHKERR moab.write_file("out_prisms_level.vtk", "VTK", "",
1730 &meshset_prims_level, 1);
1731 CHKERR moab.delete_entities(&meshset_prims_level, 1);
1732 }
1733 }
1734 if (debug) {
1735 Range tets_level;
1736 CHKERR bitRefPtr->getEntitiesByTypeAndRefLevel(bit_split, bit_split, MBTET,
1737 tets_level);
1738 EntityHandle meshset_level;
1739 CHKERR moab.create_meshset(MESHSET_SET, meshset_level);
1740 CHKERR moab.add_entities(meshset_level, tets_level);
1741 CHKERR moab.write_file("out_tets_level_only.vtk", "VTK", "", &meshset_level,
1742 1);
1743 CHKERR moab.delete_entities(&meshset_level, 1);
1744 }
1745
1747}

Member Data Documentation

◆ bitRefPtr

BitRefManager* FractureMechanics::CPMeshCut::bitRefPtr
private

Definition at line 186 of file CPMeshCut.hpp.

◆ contactPrismsBlockSetId

int FractureMechanics::CPMeshCut::contactPrismsBlockSetId
private

Definition at line 183 of file CPMeshCut.hpp.

◆ cornerNodes

Range FractureMechanics::CPMeshCut::cornerNodes
private

Definition at line 190 of file CPMeshCut.hpp.

◆ cP

CrackPropagation& FractureMechanics::CPMeshCut::cP

Definition at line 35 of file CPMeshCut.hpp.

◆ crackedBodyBlockSetId

int FractureMechanics::CPMeshCut::crackedBodyBlockSetId
private

Definition at line 182 of file CPMeshCut.hpp.

◆ crackFrontId

int FractureMechanics::CPMeshCut::crackFrontId
private

Definition at line 208 of file CPMeshCut.hpp.

◆ crackSurfaceId

int FractureMechanics::CPMeshCut::crackSurfaceId
private

Definition at line 207 of file CPMeshCut.hpp.

◆ cutMeshPtr

CutMeshInterface* FractureMechanics::CPMeshCut::cutMeshPtr
private

Definition at line 185 of file CPMeshCut.hpp.

◆ cutSurfMeshName

std::string FractureMechanics::CPMeshCut::cutSurfMeshName
private

Definition at line 211 of file CPMeshCut.hpp.

◆ cuttingSurfaceSidesetId

int FractureMechanics::CPMeshCut::cuttingSurfaceSidesetId
private

Definition at line 205 of file CPMeshCut.hpp.

◆ debugCut

PetscBool FractureMechanics::CPMeshCut::debugCut
private

Definition at line 201 of file CPMeshCut.hpp.

◆ edgesBlockSet

int FractureMechanics::CPMeshCut::edgesBlockSet
private

Definition at line 179 of file CPMeshCut.hpp.

◆ edgesBlockSetFlg

PetscBool FractureMechanics::CPMeshCut::edgesBlockSetFlg
private

Definition at line 178 of file CPMeshCut.hpp.

◆ fixedEdges

Range FractureMechanics::CPMeshCut::fixedEdges
private

Definition at line 189 of file CPMeshCut.hpp.

◆ fractionLevel

int FractureMechanics::CPMeshCut::fractionLevel
private

Definition at line 194 of file CPMeshCut.hpp.

◆ fRont

Range FractureMechanics::CPMeshCut::fRont
private

Definition at line 191 of file CPMeshCut.hpp.

◆ meshsetMngPtr

MeshsetsManager* FractureMechanics::CPMeshCut::meshsetMngPtr
private

Definition at line 187 of file CPMeshCut.hpp.

◆ nbRefBeforCut

int FractureMechanics::CPMeshCut::nbRefBeforCut
private

Definition at line 198 of file CPMeshCut.hpp.

◆ nbRefBeforTrim

int FractureMechanics::CPMeshCut::nbRefBeforTrim
private

Definition at line 199 of file CPMeshCut.hpp.

◆ originalCoords

std::vector<double> FractureMechanics::CPMeshCut::originalCoords
private

Definition at line 210 of file CPMeshCut.hpp.

◆ removePathologicalFrontTris

PetscBool FractureMechanics::CPMeshCut::removePathologicalFrontTris
private

Definition at line 202 of file CPMeshCut.hpp.

◆ sHift

double FractureMechanics::CPMeshCut::sHift[3]
private

Definition at line 204 of file CPMeshCut.hpp.

◆ skinOfTheBodyId

int FractureMechanics::CPMeshCut::skinOfTheBodyId
private

Definition at line 206 of file CPMeshCut.hpp.

◆ sUrface

Range FractureMechanics::CPMeshCut::sUrface
private

Definition at line 192 of file CPMeshCut.hpp.

◆ tolCut

double FractureMechanics::CPMeshCut::tolCut
private

Definition at line 195 of file CPMeshCut.hpp.

◆ tolCutClose

double FractureMechanics::CPMeshCut::tolCutClose
private

Definition at line 196 of file CPMeshCut.hpp.

◆ tolTrimClose

double FractureMechanics::CPMeshCut::tolTrimClose
private

Definition at line 197 of file CPMeshCut.hpp.

◆ vertexBlockSet

int FractureMechanics::CPMeshCut::vertexBlockSet
private

Definition at line 180 of file CPMeshCut.hpp.

◆ vertexNodeSet

int FractureMechanics::CPMeshCut::vertexNodeSet
private

Definition at line 181 of file CPMeshCut.hpp.


The documentation for this struct was generated from the following files: