v0.13.1
Public Member Functions | Public Attributes | List of all members
FractureMechanics::GetSmoothingElementsSkin Struct Reference
Collaboration diagram for FractureMechanics::GetSmoothingElementsSkin:
[legend]

Public Member Functions

 GetSmoothingElementsSkin (CrackPropagation &cp)
 
MoFEMErrorCode getOptions ()
 
MoFEMErrorCode operator() (Range &fix_material_nodes, const bool fix_small_g=true, const bool debug=false) const
 

Public Attributes

CrackPropagationcP
 
double gcFixThreshold
 

Detailed Description

Definition at line 242 of file CrackPropagation.cpp.

Constructor & Destructor Documentation

◆ GetSmoothingElementsSkin()

FractureMechanics::GetSmoothingElementsSkin::GetSmoothingElementsSkin ( CrackPropagation cp)

Definition at line 245 of file CrackPropagation.cpp.

245 : cP(cp), gcFixThreshold(0.5) {
246 ierr = getOptions();
247 CHKERRABORT(PETSC_COMM_WORLD, ierr);
248 }
static MoFEMErrorCodeGeneric< PetscErrorCode > ierr
Definition: Exceptions.hpp:76

Member Function Documentation

◆ getOptions()

MoFEMErrorCode FractureMechanics::GetSmoothingElementsSkin::getOptions ( )

Definition at line 250 of file CrackPropagation.cpp.

250 {
252 ierr =
253 PetscOptionsBegin(PETSC_COMM_WORLD, "", "Fixing nodes options", "none");
254 CHKERRG(ierr);
255 CHKERR PetscOptionsScalar("-gc_fix_threshold",
256 "fix nodes below given threshold", "",
257 gcFixThreshold, &gcFixThreshold, PETSC_NULL);
258 CHKERR PetscPrintf(PETSC_COMM_WORLD,
259 "### Input parameter: -gc_fix_threshold %6.4e \n",
261
262 ierr = PetscOptionsEnd();
263 CHKERRG(ierr);
265 }
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:346
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:483
#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

◆ operator()()

MoFEMErrorCode FractureMechanics::GetSmoothingElementsSkin::operator() ( Range fix_material_nodes,
const bool  fix_small_g = true,
const bool  debug = false 
) const

Definition at line 267 of file CrackPropagation.cpp.

269 {
270 MoFEM::Interface &m_field = cP.mField;
272
273 Skinner skin(&m_field.get_moab());
274 Range tets_level, prisms_level;
275 CHKERR m_field.getInterface<BitRefManager>()->getEntitiesByTypeAndRefLevel(
276 cP.mapBitLevel["material_domain"], BitRefLevel().set(), MBTET,
277 tets_level);
278 CHKERR m_field.getInterface<BitRefManager>()->getEntitiesByTypeAndRefLevel(
279 cP.mapBitLevel["material_domain"], BitRefLevel().set(), MBPRISM,
280 prisms_level);
281
282 prisms_level = intersect(prisms_level, cP.contactElements);
283 tets_level.merge(prisms_level);
284
285 Range mesh_skin;
286 CHKERR skin.find_skin(0, tets_level, false, mesh_skin);
287
288 mesh_skin = mesh_skin.subset_by_type(MBTRI);
289
290 Range faces_to_remove, contact_faces;
291 faces_to_remove = unite(cP.crackFaces, cP.bodySkin);
292 contact_faces = unite(cP.contactMasterFaces, cP.contactSlaveFaces);
293 faces_to_remove = subtract(faces_to_remove, contact_faces);
294
295 mesh_skin = subtract(mesh_skin, faces_to_remove);
296
297 CHKERR m_field.get_moab().get_connectivity(mesh_skin, fix_material_nodes);
298 // fix fixed edges nodes
299 Range fixed_edges = cP.getInterface<CPMeshCut>()->getFixedEdges();
300 Range fixed_edges_nodes;
301 CHKERR m_field.get_moab().get_connectivity(fixed_edges, fixed_edges_nodes,
302 true);
303
304 Range crack_faces_and_fixed_edges_nodes;
305 CHKERR m_field.get_moab().get_connectivity(
306 cP.crackFaces, crack_faces_and_fixed_edges_nodes, true);
307 crack_faces_and_fixed_edges_nodes =
308 intersect(crack_faces_and_fixed_edges_nodes, fixed_edges_nodes);
309
310 Range crack_faces_edges;
311 CHKERR m_field.get_moab().get_adjacencies(
312 cP.crackFaces, 1, false, crack_faces_edges, moab::Interface::UNION);
313 Range not_fixed_edges = intersect(crack_faces_edges, fixed_edges);
314 Range not_fixed_edges_nodes;
315 CHKERR m_field.get_moab().get_connectivity(not_fixed_edges,
316 not_fixed_edges_nodes, true);
317 crack_faces_and_fixed_edges_nodes =
318 subtract(crack_faces_and_fixed_edges_nodes, not_fixed_edges_nodes);
319 fix_material_nodes.merge(crack_faces_and_fixed_edges_nodes);
320
321 Range not_fixed_edges_skin;
322 CHKERR skin.find_skin(0, not_fixed_edges, false, not_fixed_edges_skin);
323 fix_material_nodes.merge(
324 subtract(not_fixed_edges_skin, cP.crackFrontNodes));
325 if (debug && !m_field.get_comm_rank()) {
326 EntityHandle meshset;
327 CHKERR m_field.get_moab().create_meshset(MESHSET_SET, meshset);
328 CHKERR m_field.get_moab().add_entities(meshset, mesh_skin);
329 CHKERR m_field.get_moab().write_file("fixed_faces.vtk", "VTK", "",
330 &meshset, 1);
331 CHKERR m_field.get_moab().delete_entities(&meshset, 1);
332 }
333
334 // fix corner nodes
335 fix_material_nodes.merge(cP.getInterface<CPMeshCut>()->getCornerNodes());
336
337 // fix nodes on kinematic and static boundary conditions
338 auto bc_fix_nodes = [&m_field,
339 &fix_material_nodes](const EntityHandle meshset) {
341 Range ents;
342 CHKERR m_field.get_moab().get_entities_by_handle(meshset, ents, true);
343 Range nodes;
344 CHKERR m_field.get_moab().get_connectivity(ents, nodes);
345 fix_material_nodes.merge(nodes);
347 };
349 m_field, NODESET | DISPLACEMENTSET, it)) {
350 CHKERR bc_fix_nodes(it->getMeshset());
351 }
353 it)) {
354 CHKERR bc_fix_nodes(it->getMeshset());
355 }
356 if (!cP.isPressureAle) {
357 // fix nodes on surfaces under pressure
359 m_field, SIDESET | PRESSURESET, it)) {
360 CHKERR bc_fix_nodes(it->getMeshset());
361 }
362 }
363
364 // fix nodes on contact surfaces
365 if (cP.fixContactNodes) {
366 Range contact_nodes;
367 CHKERR m_field.get_moab().get_connectivity(cP.contactElements,
368 contact_nodes);
369 fix_material_nodes.merge(contact_nodes);
370 }
371
372 if (!cP.areSpringsAle) {
373 // fix nodes on spring nodes
375 if (bit->getName().compare(0, 9, "SPRING_BC") == 0) {
376 CHKERR bc_fix_nodes(bit->meshset);
377 }
378 }
379 }
380
381 {
382 Range contact_nodes;
383 CHKERR m_field.get_moab().get_connectivity(cP.mortarContactElements,
384 contact_nodes);
385 fix_material_nodes.merge(contact_nodes);
386 }
387
388 // fix nodes on the constrained interface
389 if (!cP.constrainedInterface.empty()) {
390 Range interface_nodes;
391 CHKERR m_field.get_moab().get_connectivity(cP.constrainedInterface,
392 interface_nodes);
393 fix_material_nodes.merge(interface_nodes);
394 }
395
396 auto add_nodes_on_opposite_side_of_prism = [&m_field, &fix_material_nodes,
397 this](bool only_contact_prisms =
398 false) {
400 Range adj_prisms;
401 CHKERR m_field.get_moab().get_adjacencies(
402 fix_material_nodes, 3, false, adj_prisms, moab::Interface::UNION);
403 adj_prisms = adj_prisms.subset_by_type(MBPRISM);
404 if (only_contact_prisms) {
405 adj_prisms = intersect(adj_prisms, cP.contactElements);
406 }
407 for (auto p : adj_prisms) {
408 int num_nodes;
409 const EntityHandle *conn;
410 CHKERR m_field.get_moab().get_connectivity(p, conn, num_nodes, false);
411 Range add_nodes;
412 for (int n = 0; n != 3; ++n) {
413 if (fix_material_nodes.find(conn[n]) != fix_material_nodes.end()) {
414 add_nodes.insert(conn[3 + n]);
415 }
416 if (fix_material_nodes.find(conn[3 + n]) !=
417 fix_material_nodes.end()) {
418 add_nodes.insert(conn[n]);
419 }
420 }
421 fix_material_nodes.merge(add_nodes);
422 }
424 };
425
426 CHKERR add_nodes_on_opposite_side_of_prism();
427
429 if (fix_small_g) {
430
431 if (!cP.mField.get_comm_rank()) {
432 auto copy_map = [&fix_material_nodes](const auto &map) {
433 std::map<EntityHandle, double> map_copy;
434 for (auto &m : map) {
435 if (fix_material_nodes.find(m.first) == fix_material_nodes.end()) {
436 map_copy[m.first] = m.second;
437 }
438 }
439 return map_copy;
440 };
441
442 const std::map<EntityHandle, double> map_j = copy_map(cP.mapJ);
443 const std::map<EntityHandle, double> map_g1 = copy_map(cP.mapG1);
444
445 auto find_max = [](const auto &map) {
446 double max = 0;
447 for (auto &m : map)
448 max = fmax(m.second, max);
449 return max;
450 };
451
452 const double max_j = find_max(map_j);
453 const double max_g1 = find_max(map_g1);
454
455 for (auto &m : map_j) {
456 if ((m.second < gcFixThreshold * max_j) &&
457 (map_g1.at(m.first) < gcFixThreshold * max_g1)) {
458 fix_material_nodes.insert(m.first);
459 } else {
461 }
462 }
463
464 if (map_j.empty() || map_g1.empty() || cP.fractionOfFixedNodes == 0) {
465 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
466 "No points to move on crack front");
467 }
468
469 cP.fractionOfFixedNodes /= map_j.size();
470 }
471
472 Vec vec;
473 CHKERR VecCreateMPI(cP.mField.get_comm(), PETSC_DETERMINE, 1, &vec);
474 if (!cP.mField.get_comm_rank()) {
475 CHKERR VecSetValue(vec, cP.mField.get_comm_rank(),
476 cP.fractionOfFixedNodes, INSERT_VALUES);
477 }
478 CHKERR VecMax(vec, PETSC_NULL, &cP.fractionOfFixedNodes);
479 CHKERR VecDestroy(&vec);
480 }
481
482 CHKERR add_nodes_on_opposite_side_of_prism(true);
483
485 }
static Index< 'p', 3 > p
@ PRESSURESET
Definition: definitions.h:152
@ FORCESET
Definition: definitions.h:151
@ NODESET
Definition: definitions.h:146
@ SIDESET
Definition: definitions.h:147
@ DISPLACEMENTSET
Definition: definitions.h:150
@ BLOCKSET
Definition: definitions.h:148
@ MOFEM_DATA_INCONSISTENCY
Definition: definitions.h:31
static const bool debug
FTensor::Index< 'n', SPACE_DIM > n
FTensor::Index< 'm', SPACE_DIM > m
#define _IT_CUBITMESHSETS_BY_BCDATA_TYPE_FOR_LOOP_(MESHSET_MANAGER, CUBITBCTYPE, IT)
Iterator that loops over a specific Cubit MeshSet in a moFEM field.
#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
const FTensor::Tensor2< T, Dim, Dim > Vec
std::bitset< BITREFLEVEL_SIZE > BitRefLevel
Bit structure attached to each entity identifying to what mesh entity is attached.
Definition: Types.hpp:40
std::map< std::string, BitRefLevel > mapBitLevel
map< EntityHandle, double > mapJ
hashmap of J - release energy at nodes
Range constrainedInterface
Range of faces on the constrained interface.
PetscBool isPressureAle
If true surface pressure is considered in ALE.
map< EntityHandle, double > mapG1
hashmap of g1 - release energy at nodes
PetscBool areSpringsAle
If true surface spring is considered in ALE.
Managing BitRefLevels.
virtual moab::Interface & get_moab()=0
virtual MPI_Comm & get_comm() const =0
virtual int get_comm_rank() const =0
Deprecated interface functions.
MoFEMErrorCode getInterface(IFACE *&iface) const
Get interface refernce to pointer of interface.

Member Data Documentation

◆ cP

CrackPropagation& FractureMechanics::GetSmoothingElementsSkin::cP

Definition at line 243 of file CrackPropagation.cpp.

◆ gcFixThreshold

double FractureMechanics::GetSmoothingElementsSkin::gcFixThreshold

Definition at line 244 of file CrackPropagation.cpp.


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