https://mooseframework.inl.gov
Classes | Public Member Functions | Public Attributes | Protected Types | Protected Member Functions | Protected Attributes | List of all members
PenetrationThread Class Reference

#include <PenetrationThread.h>

Classes

struct  RidgeData
 
struct  RidgeSetData
 

Public Member Functions

 PenetrationThread (SubProblem &subproblem, const MooseMesh &mesh, BoundaryID primary_boundary, BoundaryID secondary_boundary, std::map< dof_id_type, PenetrationInfo *> &penetration_info, bool check_whether_reasonable, bool update_location, Real tangential_tolerance, bool do_normal_smoothing, Real normal_smoothing_distance, PenetrationLocator::NORMAL_SMOOTHING_METHOD normal_smoothing_method, bool use_point_locator, std::vector< std::vector< libMesh::FEBase *>> &fes, libMesh::FEType &fe_type, NearestNodeLocator &nearest_node, const std::map< dof_id_type, std::vector< dof_id_type >> &node_to_elem_map)
 
 PenetrationThread (PenetrationThread &x, Threads::split split)
 
void operator() (const NodeIdRange &range)
 
void join (const PenetrationThread &other)
 

Public Attributes

std::vector< dof_id_type_recheck_secondary_nodes
 List of secondary nodes for which penetration was not detected in the current patch and for which patch has to be updated. More...
 

Protected Types

enum  CompeteInteractionResult { FIRST_WINS, SECOND_WINS, NEITHER_WINS }
 
enum  CommonEdgeResult { NO_COMMON, COMMON_EDGE, COMMON_NODE, EDGE_AND_COMMON_NODE }
 

Protected Member Functions

CompeteInteractionResult competeInteractions (PenetrationInfo *pi1, PenetrationInfo *pi2)
 When interactions are identified between a node and two faces, compete between the faces to determine whether first (pi1) or second (pi2) interaction is stronger. More...
 
CompeteInteractionResult competeInteractionsBothOnFace (PenetrationInfo *pi1, PenetrationInfo *pi2)
 Determine whether first (pi1) or second (pi2) interaction is stronger when it is known that the node projects to both of the two competing faces. More...
 
CommonEdgeResult interactionsOffCommonEdge (PenetrationInfo *pi1, PenetrationInfo *pi2)
 
bool findRidgeContactPoint (libMesh::Point &contact_point, Real &tangential_distance, const Node *&closest_node, unsigned int &index, libMesh::Point &contact_point_ref, std::vector< PenetrationInfo *> &p_info, const unsigned int index1, const unsigned int index2)
 
void getSideCornerNodes (const Elem *side, std::vector< const Node *> &corner_nodes)
 
bool restrictPointToSpecifiedEdgeOfFace (libMesh::Point &p, const Node *&closest_node, const Elem *side, const std::vector< const Node *> &edge_nodes)
 
bool restrictPointToFace (libMesh::Point &p, const Node *&closest_node, const Elem *side)
 
bool isFaceReasonableCandidate (const Elem *primary_elem, const Elem *side, libMesh::FEBase *fe, const libMesh::Point *secondary_point, const Real tangential_tolerance)
 
void smoothNormal (PenetrationInfo *info, std::vector< PenetrationInfo *> &p_info, const Node &node)
 
void getSmoothingFacesAndWeights (PenetrationInfo *info, std::vector< PenetrationInfo *> &edge_face_info, std::vector< Real > &edge_face_weights, std::vector< PenetrationInfo *> &p_info, const Node &secondary_node)
 
void getSmoothingEdgeNodesAndWeights (const libMesh::Point &p, const Elem *side, std::vector< std::vector< const Node *>> &edge_nodes, std::vector< Real > &edge_face_weights)
 
void getInfoForFacesWithCommonNodes (const Node *secondary_node, const std::set< dof_id_type > &elems_to_exclude, const std::vector< const Node *> edge_nodes, std::vector< PenetrationInfo *> &face_info_comm_edge, std::vector< PenetrationInfo *> &p_info)
 
void getInfoForElem (std::vector< PenetrationInfo *> &thisElemInfo, std::vector< PenetrationInfo *> &p_info, const Elem *elem)
 
void createInfoForElem (std::vector< PenetrationInfo *> &thisElemInfo, std::vector< PenetrationInfo *> &p_info, const Node *secondary_node, const Elem *elem, const std::vector< const Node *> &nodes_that_must_be_on_side, const bool check_whether_reasonable=false)
 
void getSidesOnPrimaryBoundary (std::vector< unsigned int > &sides, const Elem *const elem)
 
void computeSlip (libMesh::FEBase &fe, PenetrationInfo &info)
 
void switchInfo (PenetrationInfo *&info, PenetrationInfo *&infoNew)
 

Protected Attributes

SubProblem_subproblem
 
const MooseMesh_mesh
 
BoundaryID _primary_boundary
 
BoundaryID _secondary_boundary
 
std::map< dof_id_type, PenetrationInfo * > & _penetration_info
 
bool _check_whether_reasonable
 
bool _update_location
 
Real _tangential_tolerance
 
bool _do_normal_smoothing
 
Real _normal_smoothing_distance
 
PenetrationLocator::NORMAL_SMOOTHING_METHOD _normal_smoothing_method
 
bool _use_point_locator
 
MooseVariable_nodal_normal_x
 
MooseVariable_nodal_normal_y
 
MooseVariable_nodal_normal_z
 
std::vector< std::vector< libMesh::FEBase * > > & _fes
 
libMesh::FEType_fe_type
 
NearestNodeLocator_nearest_node
 
const std::map< dof_id_type, std::vector< dof_id_type > > & _node_to_elem_map
 
THREAD_ID _tid
 
libMesh::ElemSideBuilder _elem_side_builder
 Helper for building element sides without extraneous allocation. More...
 

Detailed Description

Definition at line 24 of file PenetrationThread.h.

Member Enumeration Documentation

◆ CommonEdgeResult

Enumerator
NO_COMMON 
COMMON_EDGE 
COMMON_NODE 
EDGE_AND_COMMON_NODE 

Definition at line 95 of file PenetrationThread.h.

◆ CompeteInteractionResult

Enumerator
FIRST_WINS 
SECOND_WINS 
NEITHER_WINS 

Definition at line 88 of file PenetrationThread.h.

Constructor & Destructor Documentation

◆ PenetrationThread() [1/2]

PenetrationThread::PenetrationThread ( SubProblem subproblem,
const MooseMesh mesh,
BoundaryID  primary_boundary,
BoundaryID  secondary_boundary,
std::map< dof_id_type, PenetrationInfo *> &  penetration_info,
bool  check_whether_reasonable,
bool  update_location,
Real  tangential_tolerance,
bool  do_normal_smoothing,
Real  normal_smoothing_distance,
PenetrationLocator::NORMAL_SMOOTHING_METHOD  normal_smoothing_method,
bool  use_point_locator,
std::vector< std::vector< libMesh::FEBase *>> &  fes,
libMesh::FEType fe_type,
NearestNodeLocator nearest_node,
const std::map< dof_id_type, std::vector< dof_id_type >> &  node_to_elem_map 
)

Definition at line 98 of file PenetrationThread.C.

115  : _subproblem(subproblem),
116  _mesh(mesh),
117  _primary_boundary(primary_boundary),
118  _secondary_boundary(secondary_boundary),
119  _penetration_info(penetration_info),
120  _check_whether_reasonable(check_whether_reasonable),
121  _update_location(update_location),
122  _tangential_tolerance(tangential_tolerance),
123  _do_normal_smoothing(do_normal_smoothing),
124  _normal_smoothing_distance(normal_smoothing_distance),
125  _normal_smoothing_method(normal_smoothing_method),
126  _use_point_locator(use_point_locator),
127  _nodal_normal_x(NULL),
128  _nodal_normal_y(NULL),
129  _nodal_normal_z(NULL),
130  _fes(fes),
131  _fe_type(fe_type),
132  _nearest_node(nearest_node),
133  _node_to_elem_map(node_to_elem_map)
134 {
135 }
MooseVariable * _nodal_normal_z
NearestNodeLocator & _nearest_node
libMesh::FEType & _fe_type
BoundaryID _primary_boundary
SubProblem & _subproblem
const MooseMesh & _mesh
BoundaryID _secondary_boundary
std::map< dof_id_type, PenetrationInfo * > & _penetration_info
MooseVariable * _nodal_normal_y
const std::map< dof_id_type, std::vector< dof_id_type > > & _node_to_elem_map
MooseVariable * _nodal_normal_x
std::vector< std::vector< libMesh::FEBase * > > & _fes
PenetrationLocator::NORMAL_SMOOTHING_METHOD _normal_smoothing_method

◆ PenetrationThread() [2/2]

PenetrationThread::PenetrationThread ( PenetrationThread x,
Threads::split  split 
)

Definition at line 138 of file PenetrationThread.C.

140  _mesh(x._mesh),
151  _fes(x._fes),
152  _fe_type(x._fe_type),
155 {
156 }
NearestNodeLocator & _nearest_node
libMesh::FEType & _fe_type
BoundaryID _primary_boundary
SubProblem & _subproblem
const MooseMesh & _mesh
BoundaryID _secondary_boundary
std::map< dof_id_type, PenetrationInfo * > & _penetration_info
const std::map< dof_id_type, std::vector< dof_id_type > > & _node_to_elem_map
std::vector< std::vector< libMesh::FEBase * > > & _fes
PenetrationLocator::NORMAL_SMOOTHING_METHOD _normal_smoothing_method

Member Function Documentation

◆ competeInteractions()

PenetrationThread::CompeteInteractionResult PenetrationThread::competeInteractions ( PenetrationInfo pi1,
PenetrationInfo pi2 
)
protected

When interactions are identified between a node and two faces, compete between the faces to determine whether first (pi1) or second (pi2) interaction is stronger.

Parameters
pi1Pointer to the PenetrationInfo object for the first face
pi2Pointer to the PenetrationInfo object for the second face
Returns
Appropriate ComputeInterationResult enum entry identifying which face is a better match

Definition at line 680 of file PenetrationThread.C.

Referenced by operator()().

681 {
682 
684 
686  pi2->_tangential_distance > _tangential_tolerance) // out of tol on both faces
687  result = NEITHER_WINS;
688 
689  else if (pi1->_tangential_distance == 0.0 &&
690  pi2->_tangential_distance > 0.0) // on face 1, off face 2
691  result = FIRST_WINS;
692 
693  else if (pi2->_tangential_distance == 0.0 &&
694  pi1->_tangential_distance > 0.0) // on face 2, off face 1
695  result = SECOND_WINS;
696 
697  else if (pi1->_tangential_distance <= _tangential_tolerance &&
698  pi2->_tangential_distance > _tangential_tolerance) // in face 1 tol, out of face 2 tol
699  result = FIRST_WINS;
700 
701  else if (pi2->_tangential_distance <= _tangential_tolerance &&
702  pi1->_tangential_distance > _tangential_tolerance) // in face 2 tol, out of face 1 tol
703  result = SECOND_WINS;
704 
705  else if (pi1->_tangential_distance == 0.0 && pi2->_tangential_distance == 0.0) // on both faces
706  result = competeInteractionsBothOnFace(pi1, pi2);
707 
708  else if (pi1->_tangential_distance <= _tangential_tolerance &&
709  pi2->_tangential_distance <= _tangential_tolerance) // off but within tol of both faces
710  {
712  if (cer == COMMON_EDGE || cer == COMMON_NODE) // ridge case.
713  {
714  // We already checked for ridges, and it got rejected, so neither face must be valid
715  result = NEITHER_WINS;
716  // mooseError("Erroneously encountered ridge case");
717  }
718  else if (cer == EDGE_AND_COMMON_NODE) // off side of face, off corner of another face. Favor
719  // the off-side face
720  {
721  if (pi1->_off_edge_nodes.size() == pi2->_off_edge_nodes.size())
722  mooseError("Invalid off_edge_nodes counts");
723 
724  else if (pi1->_off_edge_nodes.size() == 2)
725  result = FIRST_WINS;
726 
727  else if (pi2->_off_edge_nodes.size() == 2)
728  result = SECOND_WINS;
729 
730  else
731  mooseError("Invalid off_edge_nodes counts");
732  }
733  else // The node projects to both faces within tangential tolerance.
734  result = competeInteractionsBothOnFace(pi1, pi2);
735  }
736 
737  return result;
738 }
void mooseError(Args &&... args)
Emit an error message with the given stringified, concatenated args and terminate the application...
Definition: MooseError.h:323
std::vector< const Node * > _off_edge_nodes
CompeteInteractionResult competeInteractionsBothOnFace(PenetrationInfo *pi1, PenetrationInfo *pi2)
Determine whether first (pi1) or second (pi2) interaction is stronger when it is known that the node ...
CommonEdgeResult interactionsOffCommonEdge(PenetrationInfo *pi1, PenetrationInfo *pi2)

◆ competeInteractionsBothOnFace()

PenetrationThread::CompeteInteractionResult PenetrationThread::competeInteractionsBothOnFace ( PenetrationInfo pi1,
PenetrationInfo pi2 
)
protected

Determine whether first (pi1) or second (pi2) interaction is stronger when it is known that the node projects to both of the two competing faces.

Parameters
pi1Pointer to the PenetrationInfo object for the first face
pi2Pointer to the PenetrationInfo object for the second face
Returns
Appropriate ComputeInterationResult enum entry identifying which face is a better match

Definition at line 741 of file PenetrationThread.C.

Referenced by competeInteractions().

742 {
744 
745  if (pi1->_distance >= 0.0 && pi2->_distance < 0.0)
746  result = FIRST_WINS; // favor face with positive distance (penetrated) -- first in this case
747 
748  else if (pi2->_distance >= 0.0 && pi1->_distance < 0.0)
749  result = SECOND_WINS; // favor face with positive distance (penetrated) -- second in this case
750 
751  // TODO: This logic below could cause an abrupt jump from one face to the other with small mesh
752  // movement. If there is some way to smooth the transition, we should do it.
754  result = FIRST_WINS; // otherwise, favor the closer face -- first in this case
755 
757  result = SECOND_WINS; // otherwise, favor the closer face -- second in this case
758 
759  else // Equal within tolerance. Favor the one with a smaller element id (for repeatibility)
760  {
761  if (pi1->_elem->id() < pi2->_elem->id())
762  result = FIRST_WINS;
763 
764  else
765  result = SECOND_WINS;
766  }
767 
768  return result;
769 }
MetaPhysicL::DualNumber< V, D, asd > abs(const MetaPhysicL::DualNumber< V, D, asd > &a)
Definition: EigenADReal.h:42
dof_id_type id() const
const Elem * _elem
bool relativeFuzzyLessThan(const T &var1, const T2 &var2, const T3 &tol=libMesh::TOLERANCE *libMesh::TOLERANCE)
Function to check whether a variable is less than another variable within a relative tolerance...
Definition: MooseUtils.h:631

◆ computeSlip()

void PenetrationThread::computeSlip ( libMesh::FEBase fe,
PenetrationInfo info 
)
protected

Definition at line 1389 of file PenetrationThread.C.

Referenced by operator()().

1390 {
1391  // Slip is current projected position of secondary node minus
1392  // original projected position of secondary node
1393  std::vector<Point> points(1);
1394  points[0] = info._starting_closest_point_ref;
1395  const auto & side = _elem_side_builder(*info._starting_elem, info._starting_side_num);
1396  fe.reinit(&side, &points);
1397  const std::vector<Point> & starting_point = fe.get_xyz();
1398  info._incremental_slip = info._closest_point - starting_point[0];
1399  if (info.isCaptured())
1400  {
1401  info._frictional_energy =
1402  info._frictional_energy_old + info._contact_force * info._incremental_slip;
1403  info._accumulated_slip = info._accumulated_slip_old + info._incremental_slip.norm();
1404  }
1405 }
MPI_Info info
virtual void reinit(const Elem *elem, const std::vector< Point > *const pts=nullptr, const std::vector< Real > *const weights=nullptr)=0
virtual_for_inffe const std::vector< Point > & get_xyz() const
libMesh::ElemSideBuilder _elem_side_builder
Helper for building element sides without extraneous allocation.

◆ createInfoForElem()

void PenetrationThread::createInfoForElem ( std::vector< PenetrationInfo *> &  thisElemInfo,
std::vector< PenetrationInfo *> &  p_info,
const Node secondary_node,
const Elem elem,
const std::vector< const Node *> &  nodes_that_must_be_on_side,
const bool  check_whether_reasonable = false 
)
protected

Definition at line 1807 of file PenetrationThread.C.

Referenced by getInfoForFacesWithCommonNodes(), and operator()().

1813 {
1814  const BoundaryInfo & boundary_info = _mesh.getMesh().get_boundary_info();
1815 
1816  for (auto s : elem->side_index_range())
1817  {
1818  if (!boundary_info.has_boundary_id(elem, s, _primary_boundary))
1819  continue;
1820 
1821  // Don't create info for this side if one already exists
1822  bool already_have_info_this_side = false;
1823  for (const auto & pi : thisElemInfo)
1824  if (pi->_side_num == s)
1825  {
1826  already_have_info_this_side = true;
1827  break;
1828  }
1829 
1830  if (already_have_info_this_side)
1831  break;
1832 
1833  const Elem * side = elem->build_side_ptr(s).release();
1834 
1835  // Only continue with creating info for this side if the side contains
1836  // all of the nodes in nodes_that_must_be_on_side
1837  std::vector<const Node *> nodevec;
1838  for (unsigned int ni = 0; ni < side->n_nodes(); ++ni)
1839  nodevec.push_back(side->node_ptr(ni));
1840 
1841  std::sort(nodevec.begin(), nodevec.end());
1842  std::vector<const Node *> common_nodes;
1843  std::set_intersection(nodes_that_must_be_on_side.begin(),
1844  nodes_that_must_be_on_side.end(),
1845  nodevec.begin(),
1846  nodevec.end(),
1847  std::inserter(common_nodes, common_nodes.end()));
1848  if (common_nodes.size() != nodes_that_must_be_on_side.size())
1849  {
1850  delete side;
1851  break;
1852  }
1853 
1854  FEBase * fe_elem = _fes[_tid][elem->dim()];
1855  FEBase * fe_side = _fes[_tid][side->dim()];
1856 
1857  // Optionally check to see whether face is reasonable candidate based on an
1858  // estimate of how closely it is likely to project to the face
1859  if (check_whether_reasonable)
1860  if (!isFaceReasonableCandidate(elem, side, fe_side, secondary_node, _tangential_tolerance))
1861  {
1862  delete side;
1863  break;
1864  }
1865 
1866  Point contact_phys;
1867  Point contact_ref;
1868  Point contact_on_face_ref;
1869  Real distance = 0.;
1870  Real tangential_distance = 0.;
1871  RealGradient normal;
1872  bool contact_point_on_side;
1873  std::vector<const Node *> off_edge_nodes;
1874  std::vector<std::vector<Real>> side_phi;
1875  std::vector<std::vector<RealGradient>> side_grad_phi;
1876  std::vector<RealGradient> dxyzdxi;
1877  std::vector<RealGradient> dxyzdeta;
1878  std::vector<RealGradient> d2xyzdxideta;
1879 
1880  std::unique_ptr<PenetrationInfo> pen_info =
1881  std::make_unique<PenetrationInfo>(elem,
1882  side,
1883  s,
1884  normal,
1885  distance,
1886  tangential_distance,
1887  contact_phys,
1888  contact_ref,
1889  contact_on_face_ref,
1890  off_edge_nodes,
1891  side_phi,
1892  side_grad_phi,
1893  dxyzdxi,
1894  dxyzdeta,
1895  d2xyzdxideta);
1896 
1897  bool search_succeeded = false;
1898  Moose::findContactPoint(*pen_info,
1899  fe_elem,
1900  fe_side,
1901  _fe_type,
1902  *secondary_node,
1903  true,
1905  contact_point_on_side,
1906  search_succeeded);
1907 
1908  // Do not add contact info from failed searches
1909  if (search_succeeded)
1910  {
1911  thisElemInfo.push_back(pen_info.get());
1912  p_info.push_back(pen_info.release());
1913  }
1914  }
1915 }
bool has_boundary_id(const Node *const node, const boundary_id_type id) const
IntRange< unsigned short > side_index_range() const
virtual std::unique_ptr< Elem > build_side_ptr(const unsigned int i)=0
const BoundaryInfo & get_boundary_info() const
Real distance(const Point &p)
libMesh::FEType & _fe_type
BoundaryID _primary_boundary
const MooseMesh & _mesh
MeshBase & getMesh()
Accessor for the underlying libMesh Mesh object.
Definition: MooseMesh.C:3469
virtual unsigned int n_nodes() const=0
void findContactPoint(PenetrationInfo &p_info, libMesh::FEBase *fe_elem, libMesh::FEBase *fe_side, libMesh::FEType &fe_side_type, const libMesh::Point &secondary_point, bool start_with_centroid, const Real tangential_tolerance, bool &contact_point_on_side, bool &search_succeeded)
Finds the closest point (called the contact point) on the primary_elem on side "side" to the secondar...
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
virtual unsigned short dim() const=0
const Node * node_ptr(const unsigned int i) const
std::vector< std::vector< libMesh::FEBase * > > & _fes
bool isFaceReasonableCandidate(const Elem *primary_elem, const Elem *side, libMesh::FEBase *fe, const libMesh::Point *secondary_point, const Real tangential_tolerance)
const Real pi

◆ findRidgeContactPoint()

bool PenetrationThread::findRidgeContactPoint ( libMesh::Point contact_point,
Real tangential_distance,
const Node *&  closest_node,
unsigned int index,
libMesh::Point contact_point_ref,
std::vector< PenetrationInfo *> &  p_info,
const unsigned int  index1,
const unsigned int  index2 
)
protected

Definition at line 826 of file PenetrationThread.C.

Referenced by operator()().

834 {
835  tangential_distance = 0.0;
836  closest_node = NULL;
837  PenetrationInfo * pi1 = p_info[index1];
838  PenetrationInfo * pi2 = p_info[index2];
839  const unsigned sidedim(pi1->_side->dim());
840  mooseAssert(sidedim == pi2->_side->dim(), "Incompatible dimensionalities");
841 
842  // Nodes on faces for the two interactions
843  std::vector<const Node *> side1_nodes;
844  getSideCornerNodes(pi1->_side, side1_nodes);
845  std::vector<const Node *> side2_nodes;
846  getSideCornerNodes(pi2->_side, side2_nodes);
847 
848  std::sort(side1_nodes.begin(), side1_nodes.end());
849  std::sort(side2_nodes.begin(), side2_nodes.end());
850 
851  // Find nodes shared by the two faces
852  std::vector<const Node *> common_nodes;
853  std::set_intersection(side1_nodes.begin(),
854  side1_nodes.end(),
855  side2_nodes.begin(),
856  side2_nodes.end(),
857  std::inserter(common_nodes, common_nodes.end()));
858 
859  if (common_nodes.size() != sidedim)
860  return false;
861 
862  bool found_point1, found_point2;
863  Point closest_coor_ref1(pi1->_closest_point_ref);
864  const Node * closest_node1;
865  found_point1 = restrictPointToSpecifiedEdgeOfFace(
866  closest_coor_ref1, closest_node1, pi1->_side, common_nodes);
867 
868  Point closest_coor_ref2(pi2->_closest_point_ref);
869  const Node * closest_node2;
870  found_point2 = restrictPointToSpecifiedEdgeOfFace(
871  closest_coor_ref2, closest_node2, pi2->_side, common_nodes);
872 
873  if (!found_point1 || !found_point2)
874  return false;
875 
876  // if (sidedim == 2)
877  // {
878  // TODO:
879  // We have the parametric coordinates of the closest intersection point for both faces.
880  // We need to find a point somewhere in the middle of them so there's not an abrupt jump.
881  // Find that point by taking dot products of vector from contact point to secondary node point
882  // with face normal vectors to see which face we're closer to.
883  // }
884 
885  FEBase * fe = NULL;
886  std::vector<Point> points(1);
887 
888  // We have to pick one of the two faces to own the contact point. It doesn't really matter
889  // which one we pick. For repeatibility, pick the face with the lowest index.
890  if (index1 < index2)
891  {
892  fe = _fes[_tid][pi1->_side->dim()];
893  contact_point_ref = closest_coor_ref1;
894  points[0] = closest_coor_ref1;
895  fe->reinit(pi1->_side, &points);
896  index = index1;
897  }
898  else
899  {
900  fe = _fes[_tid][pi2->_side->dim()];
901  contact_point_ref = closest_coor_ref2;
902  points[0] = closest_coor_ref2;
903  fe->reinit(pi2->_side, &points);
904  index = index2;
905  }
906 
907  contact_point = fe->get_xyz()[0];
908 
909  if (sidedim == 2)
910  {
911  if (closest_node1) // point is off the ridge between the two elements
912  {
913  mooseAssert((closest_node1 == closest_node2 || closest_node2 == NULL),
914  "If off edge of ridge, closest node must be the same on both elements");
915  closest_node = closest_node1;
916 
917  RealGradient off_face = *closest_node1 - contact_point;
918  tangential_distance = off_face.norm();
919  }
920  }
921 
922  return true;
923 }
auto norm() const -> decltype(std::norm(Real()))
Data structure used to hold penetration information.
void getSideCornerNodes(const Elem *side, std::vector< const Node *> &corner_nodes)
virtual void reinit(const Elem *elem, const std::vector< Point > *const pts=nullptr, const std::vector< Real > *const weights=nullptr)=0
virtual_for_inffe const std::vector< Point > & get_xyz() const
bool restrictPointToSpecifiedEdgeOfFace(libMesh::Point &p, const Node *&closest_node, const Elem *side, const std::vector< const Node *> &edge_nodes)
const Elem * _side
virtual unsigned short dim() const=0
std::vector< std::vector< libMesh::FEBase * > > & _fes

◆ getInfoForElem()

void PenetrationThread::getInfoForElem ( std::vector< PenetrationInfo *> &  thisElemInfo,
std::vector< PenetrationInfo *> &  p_info,
const Elem elem 
)
protected

Definition at line 1792 of file PenetrationThread.C.

Referenced by getInfoForFacesWithCommonNodes().

1795 {
1796  for (const auto & pi : p_info)
1797  {
1798  if (!pi)
1799  continue;
1800 
1801  if (pi->_elem == elem)
1802  thisElemInfo.push_back(pi);
1803  }
1804 }
const Real pi

◆ getInfoForFacesWithCommonNodes()

void PenetrationThread::getInfoForFacesWithCommonNodes ( const Node secondary_node,
const std::set< dof_id_type > &  elems_to_exclude,
const std::vector< const Node *>  edge_nodes,
std::vector< PenetrationInfo *> &  face_info_comm_edge,
std::vector< PenetrationInfo *> &  p_info 
)
protected

Definition at line 1702 of file PenetrationThread.C.

Referenced by getSmoothingFacesAndWeights().

1708 {
1709  // elems connected to a node on this edge, find one that has the same corners as this, and is not
1710  // the current elem
1711  auto node_to_elem_pair = _node_to_elem_map.find(edge_nodes[0]->id()); // just need one of the
1712  // nodes
1713  mooseAssert(node_to_elem_pair != _node_to_elem_map.end(), "Missing entry in node to elem map");
1714  const std::vector<dof_id_type> & elems_connected_to_node = node_to_elem_pair->second;
1715 
1716  std::vector<const Elem *> elems_connected_to_edge;
1717 
1718  for (unsigned int ecni = 0; ecni < elems_connected_to_node.size(); ecni++)
1719  {
1720  if (elems_to_exclude.find(elems_connected_to_node[ecni]) != elems_to_exclude.end())
1721  continue;
1722  const Elem * elem = _mesh.elemPtr(elems_connected_to_node[ecni]);
1723 
1724  std::vector<const Node *> nodevec;
1725  for (unsigned int ni = 0; ni < elem->n_nodes(); ++ni)
1726  if (elem->is_vertex(ni))
1727  nodevec.push_back(elem->node_ptr(ni));
1728 
1729  std::vector<const Node *> common_nodes;
1730  std::sort(nodevec.begin(), nodevec.end());
1731  std::set_intersection(edge_nodes.begin(),
1732  edge_nodes.end(),
1733  nodevec.begin(),
1734  nodevec.end(),
1735  std::inserter(common_nodes, common_nodes.end()));
1736 
1737  if (common_nodes.size() == edge_nodes.size())
1738  elems_connected_to_edge.push_back(elem);
1739  }
1740 
1741  if (elems_connected_to_edge.size() > 0)
1742  {
1743 
1744  // There are potentially multiple elements that share a common edge
1745  // 2D:
1746  // There can only be one element on the same surface
1747  // 3D:
1748  // If there are two edge nodes, there can only be one element on the same surface
1749  // If there is only one edge node (a corner), there could be multiple elements on the same
1750  // surface
1751  bool allowMultipleNeighbors = false;
1752 
1753  if (elems_connected_to_edge[0]->dim() == 3)
1754  {
1755  if (edge_nodes.size() == 1)
1756  {
1757  allowMultipleNeighbors = true;
1758  }
1759  }
1760 
1761  for (unsigned int i = 0; i < elems_connected_to_edge.size(); ++i)
1762  {
1763  std::vector<PenetrationInfo *> thisElemInfo;
1764  getInfoForElem(thisElemInfo, p_info, elems_connected_to_edge[i]);
1765  if (thisElemInfo.size() > 0 && !allowMultipleNeighbors)
1766  {
1767  if (thisElemInfo.size() > 1)
1768  mooseError(
1769  "Found multiple neighbors to current edge/face on surface when only one is allowed");
1770  face_info_comm_edge.push_back(thisElemInfo[0]);
1771  break;
1772  }
1773 
1775  thisElemInfo, p_info, secondary_node, elems_connected_to_edge[i], edge_nodes);
1776  if (thisElemInfo.size() > 0 && !allowMultipleNeighbors)
1777  {
1778  if (thisElemInfo.size() > 1)
1779  mooseError(
1780  "Found multiple neighbors to current edge/face on surface when only one is allowed");
1781  face_info_comm_edge.push_back(thisElemInfo[0]);
1782  break;
1783  }
1784 
1785  for (unsigned int j = 0; j < thisElemInfo.size(); ++j)
1786  face_info_comm_edge.push_back(thisElemInfo[j]);
1787  }
1788  }
1789 }
virtual Elem * elemPtr(const dof_id_type i)
Definition: MooseMesh.C:3134
void mooseError(Args &&... args)
Emit an error message with the given stringified, concatenated args and terminate the application...
Definition: MooseError.h:323
static constexpr std::size_t dim
This is the dimension of all vector and tensor datastructures used in MOOSE.
Definition: Moose.h:162
void createInfoForElem(std::vector< PenetrationInfo *> &thisElemInfo, std::vector< PenetrationInfo *> &p_info, const Node *secondary_node, const Elem *elem, const std::vector< const Node *> &nodes_that_must_be_on_side, const bool check_whether_reasonable=false)
const MooseMesh & _mesh
virtual unsigned int n_nodes() const=0
const Node * node_ptr(const unsigned int i) const
virtual bool is_vertex(const unsigned int i) const=0
void getInfoForElem(std::vector< PenetrationInfo *> &thisElemInfo, std::vector< PenetrationInfo *> &p_info, const Elem *elem)
const std::map< dof_id_type, std::vector< dof_id_type > > & _node_to_elem_map

◆ getSideCornerNodes()

void PenetrationThread::getSideCornerNodes ( const Elem side,
std::vector< const Node *> &  corner_nodes 
)
protected

Definition at line 926 of file PenetrationThread.C.

Referenced by findRidgeContactPoint().

927 {
928  const ElemType t(side->type());
929  corner_nodes.clear();
930 
931  corner_nodes.push_back(side->node_ptr(0));
932  corner_nodes.push_back(side->node_ptr(1));
933  switch (t)
934  {
935  case EDGE2:
936  case EDGE3:
937  case EDGE4:
938  {
939  break;
940  }
941 
942  case TRI3:
943  case TRI6:
944  case TRI7:
945  {
946  corner_nodes.push_back(side->node_ptr(2));
947  break;
948  }
949 
950  case QUAD4:
951  case QUAD8:
952  case QUAD9:
953  {
954  corner_nodes.push_back(side->node_ptr(2));
955  corner_nodes.push_back(side->node_ptr(3));
956  break;
957  }
958 
959  default:
960  {
961  mooseError("Unsupported face type: ", t);
962  break;
963  }
964  }
965 }
ElemType
QUAD8
void mooseError(Args &&... args)
Emit an error message with the given stringified, concatenated args and terminate the application...
Definition: MooseError.h:323
EDGE4
TRI3
QUAD4
TRI6
EDGE2
const Node * node_ptr(const unsigned int i) const
TRI7
QUAD9
EDGE3
virtual ElemType type() const=0

◆ getSidesOnPrimaryBoundary()

void PenetrationThread::getSidesOnPrimaryBoundary ( std::vector< unsigned int > &  sides,
const Elem *const  elem 
)
protected

◆ getSmoothingEdgeNodesAndWeights()

void PenetrationThread::getSmoothingEdgeNodesAndWeights ( const libMesh::Point p,
const Elem side,
std::vector< std::vector< const Node *>> &  edge_nodes,
std::vector< Real > &  edge_face_weights 
)
protected

Definition at line 1561 of file PenetrationThread.C.

Referenced by getSmoothingFacesAndWeights().

1566 {
1567  const ElemType t(side->type());
1568  const Real & xi = p(0);
1569  const Real & eta = p(1);
1570 
1571  Real smooth_limit = 1.0 - _normal_smoothing_distance;
1572 
1573  switch (t)
1574  {
1575  case EDGE2:
1576  case EDGE3:
1577  case EDGE4:
1578  {
1579  if (xi < -smooth_limit)
1580  {
1581  std::vector<const Node *> en;
1582  en.push_back(side->node_ptr(0));
1583  edge_nodes.push_back(en);
1584  Real fw = 0.5 - (1.0 + xi) / (2.0 * _normal_smoothing_distance);
1585  if (fw > 0.5)
1586  fw = 0.5;
1587  edge_face_weights.push_back(fw);
1588  }
1589  else if (xi > smooth_limit)
1590  {
1591  std::vector<const Node *> en;
1592  en.push_back(side->node_ptr(1));
1593  edge_nodes.push_back(en);
1594  Real fw = 0.5 - (1.0 - xi) / (2.0 * _normal_smoothing_distance);
1595  if (fw > 0.5)
1596  fw = 0.5;
1597  edge_face_weights.push_back(fw);
1598  }
1599  break;
1600  }
1601 
1602  case TRI3:
1603  case TRI6:
1604  case TRI7:
1605  {
1606  if (eta < -smooth_limit)
1607  {
1608  std::vector<const Node *> en;
1609  en.push_back(side->node_ptr(0));
1610  en.push_back(side->node_ptr(1));
1611  edge_nodes.push_back(en);
1612  Real fw = 0.5 - (1.0 + eta) / (2.0 * _normal_smoothing_distance);
1613  if (fw > 0.5)
1614  fw = 0.5;
1615  edge_face_weights.push_back(fw);
1616  }
1617  if ((xi + eta) > smooth_limit)
1618  {
1619  std::vector<const Node *> en;
1620  en.push_back(side->node_ptr(1));
1621  en.push_back(side->node_ptr(2));
1622  edge_nodes.push_back(en);
1623  Real fw = 0.5 - (1.0 - xi - eta) / (2.0 * _normal_smoothing_distance);
1624  if (fw > 0.5)
1625  fw = 0.5;
1626  edge_face_weights.push_back(fw);
1627  }
1628  if (xi < -smooth_limit)
1629  {
1630  std::vector<const Node *> en;
1631  en.push_back(side->node_ptr(2));
1632  en.push_back(side->node_ptr(0));
1633  edge_nodes.push_back(en);
1634  Real fw = 0.5 - (1.0 + xi) / (2.0 * _normal_smoothing_distance);
1635  if (fw > 0.5)
1636  fw = 0.5;
1637  edge_face_weights.push_back(fw);
1638  }
1639  break;
1640  }
1641 
1642  case QUAD4:
1643  case QUAD8:
1644  case QUAD9:
1645  {
1646  if (eta < -smooth_limit)
1647  {
1648  std::vector<const Node *> en;
1649  en.push_back(side->node_ptr(0));
1650  en.push_back(side->node_ptr(1));
1651  edge_nodes.push_back(en);
1652  Real fw = 0.5 - (1.0 + eta) / (2.0 * _normal_smoothing_distance);
1653  if (fw > 0.5)
1654  fw = 0.5;
1655  edge_face_weights.push_back(fw);
1656  }
1657  if (xi > smooth_limit)
1658  {
1659  std::vector<const Node *> en;
1660  en.push_back(side->node_ptr(1));
1661  en.push_back(side->node_ptr(2));
1662  edge_nodes.push_back(en);
1663  Real fw = 0.5 - (1.0 - xi) / (2.0 * _normal_smoothing_distance);
1664  if (fw > 0.5)
1665  fw = 0.5;
1666  edge_face_weights.push_back(fw);
1667  }
1668  if (eta > smooth_limit)
1669  {
1670  std::vector<const Node *> en;
1671  en.push_back(side->node_ptr(2));
1672  en.push_back(side->node_ptr(3));
1673  edge_nodes.push_back(en);
1674  Real fw = 0.5 - (1.0 - eta) / (2.0 * _normal_smoothing_distance);
1675  if (fw > 0.5)
1676  fw = 0.5;
1677  edge_face_weights.push_back(fw);
1678  }
1679  if (xi < -smooth_limit)
1680  {
1681  std::vector<const Node *> en;
1682  en.push_back(side->node_ptr(3));
1683  en.push_back(side->node_ptr(0));
1684  edge_nodes.push_back(en);
1685  Real fw = 0.5 - (1.0 + xi) / (2.0 * _normal_smoothing_distance);
1686  if (fw > 0.5)
1687  fw = 0.5;
1688  edge_face_weights.push_back(fw);
1689  }
1690  break;
1691  }
1692 
1693  default:
1694  {
1695  mooseError("Unsupported face type: ", t);
1696  break;
1697  }
1698  }
1699 }
ElemType
QUAD8
void mooseError(Args &&... args)
Emit an error message with the given stringified, concatenated args and terminate the application...
Definition: MooseError.h:323
EDGE4
TRI3
QUAD4
TRI6
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
EDGE2
const Node * node_ptr(const unsigned int i) const
TRI7
QUAD9
EDGE3
virtual ElemType type() const=0

◆ getSmoothingFacesAndWeights()

void PenetrationThread::getSmoothingFacesAndWeights ( PenetrationInfo info,
std::vector< PenetrationInfo *> &  edge_face_info,
std::vector< Real > &  edge_face_weights,
std::vector< PenetrationInfo *> &  p_info,
const Node secondary_node 
)
protected

Definition at line 1466 of file PenetrationThread.C.

Referenced by smoothNormal().

1471 {
1472  const Elem * side = info->_side;
1473  const Point & p = info->_closest_point_ref;
1474  std::set<dof_id_type> elems_to_exclude;
1475  elems_to_exclude.insert(info->_elem->id());
1476 
1477  std::vector<std::vector<const Node *>> edge_nodes;
1478 
1479  // Get the pairs of nodes along every edge that we are close enough to smooth with
1480  getSmoothingEdgeNodesAndWeights(p, side, edge_nodes, edge_face_weights);
1481  std::vector<Elem *> edge_neighbor_elems;
1482  edge_face_info.resize(edge_nodes.size(), NULL);
1483 
1484  std::vector<unsigned int> edges_without_neighbors;
1485 
1486  for (unsigned int i = 0; i < edge_nodes.size(); ++i)
1487  {
1488  // Sort all sets of edge nodes (needed for comparing edges)
1489  std::sort(edge_nodes[i].begin(), edge_nodes[i].end());
1490 
1491  std::vector<PenetrationInfo *> face_info_comm_edge;
1493  &secondary_node, elems_to_exclude, edge_nodes[i], face_info_comm_edge, p_info);
1494 
1495  if (face_info_comm_edge.size() == 0)
1496  edges_without_neighbors.push_back(i);
1497  else if (face_info_comm_edge.size() > 1)
1498  mooseError("Only one neighbor allowed per edge");
1499  else
1500  edge_face_info[i] = face_info_comm_edge[0];
1501  }
1502 
1503  // Remove edges without neighbors from the vector, starting from end
1504  std::vector<unsigned int>::reverse_iterator rit;
1505  for (rit = edges_without_neighbors.rbegin(); rit != edges_without_neighbors.rend(); ++rit)
1506  {
1507  unsigned int index = *rit;
1508  edge_nodes.erase(edge_nodes.begin() + index);
1509  edge_face_weights.erase(edge_face_weights.begin() + index);
1510  edge_face_info.erase(edge_face_info.begin() + index);
1511  }
1512 
1513  // Handle corner case
1514  if (edge_nodes.size() > 1)
1515  {
1516  if (edge_nodes.size() != 2)
1517  mooseError("Invalid number of smoothing edges");
1518 
1519  // find common node
1520  std::vector<const Node *> common_nodes;
1521  std::set_intersection(edge_nodes[0].begin(),
1522  edge_nodes[0].end(),
1523  edge_nodes[1].begin(),
1524  edge_nodes[1].end(),
1525  std::inserter(common_nodes, common_nodes.end()));
1526 
1527  if (common_nodes.size() != 1)
1528  mooseError("Invalid number of common nodes");
1529 
1530  for (const auto & pinfo : edge_face_info)
1531  elems_to_exclude.insert(pinfo->_elem->id());
1532 
1533  std::vector<PenetrationInfo *> face_info_comm_edge;
1535  &secondary_node, elems_to_exclude, common_nodes, face_info_comm_edge, p_info);
1536 
1537  unsigned int num_corner_neighbors = face_info_comm_edge.size();
1538 
1539  if (num_corner_neighbors > 0)
1540  {
1541  Real fw0 = edge_face_weights[0];
1542  Real fw1 = edge_face_weights[1];
1543 
1544  // Corner weight is product of edge weights. Spread out over multiple neighbors.
1545  Real fw_corner = (fw0 * fw1) / static_cast<Real>(num_corner_neighbors);
1546 
1547  // Adjust original edge weights
1548  edge_face_weights[0] *= (1.0 - fw1);
1549  edge_face_weights[1] *= (1.0 - fw0);
1550 
1551  for (unsigned int i = 0; i < num_corner_neighbors; ++i)
1552  {
1553  edge_face_weights.push_back(fw_corner);
1554  edge_face_info.push_back(face_info_comm_edge[i]);
1555  }
1556  }
1557  }
1558 }
void getSmoothingEdgeNodesAndWeights(const libMesh::Point &p, const Elem *side, std::vector< std::vector< const Node *>> &edge_nodes, std::vector< Real > &edge_face_weights)
MPI_Info info
void mooseError(Args &&... args)
Emit an error message with the given stringified, concatenated args and terminate the application...
Definition: MooseError.h:323
void getInfoForFacesWithCommonNodes(const Node *secondary_node, const std::set< dof_id_type > &elems_to_exclude, const std::vector< const Node *> edge_nodes, std::vector< PenetrationInfo *> &face_info_comm_edge, std::vector< PenetrationInfo *> &p_info)
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real

◆ interactionsOffCommonEdge()

PenetrationThread::CommonEdgeResult PenetrationThread::interactionsOffCommonEdge ( PenetrationInfo pi1,
PenetrationInfo pi2 
)
protected

Definition at line 772 of file PenetrationThread.C.

Referenced by competeInteractions().

773 {
774  CommonEdgeResult common_edge(NO_COMMON);
775  const std::vector<const Node *> & off_edge_nodes1 = pi1->_off_edge_nodes;
776  const std::vector<const Node *> & off_edge_nodes2 = pi2->_off_edge_nodes;
777  const unsigned dim1 = pi1->_side->dim();
778 
779  if (dim1 == 1)
780  {
781  mooseAssert(pi2->_side->dim() == 1, "Incompatible dimensions.");
782  mooseAssert(off_edge_nodes1.size() < 2 && off_edge_nodes2.size() < 2,
783  "off_edge_nodes size should be <2 for 2D contact");
784  if (off_edge_nodes1.size() == 1 && off_edge_nodes2.size() == 1 &&
785  off_edge_nodes1[0] == off_edge_nodes2[0])
786  common_edge = COMMON_EDGE;
787  }
788  else
789  {
790  mooseAssert(dim1 == 2 && pi2->_side->dim() == 2, "Incompatible dimensions.");
791  mooseAssert(off_edge_nodes1.size() < 3 && off_edge_nodes2.size() < 3,
792  "off_edge_nodes size should be <3 for 3D contact");
793  if (off_edge_nodes1.size() == 1)
794  {
795  if (off_edge_nodes2.size() == 1)
796  {
797  if (off_edge_nodes1[0] == off_edge_nodes2[0])
798  common_edge = COMMON_NODE;
799  }
800  else if (off_edge_nodes2.size() == 2)
801  {
802  if (off_edge_nodes1[0] == off_edge_nodes2[0] || off_edge_nodes1[0] == off_edge_nodes2[1])
803  common_edge = EDGE_AND_COMMON_NODE;
804  }
805  }
806  else if (off_edge_nodes1.size() == 2)
807  {
808  if (off_edge_nodes2.size() == 1)
809  {
810  if (off_edge_nodes1[0] == off_edge_nodes2[0] || off_edge_nodes1[1] == off_edge_nodes2[0])
811  common_edge = EDGE_AND_COMMON_NODE;
812  }
813  else if (off_edge_nodes2.size() == 2)
814  {
815  if ((off_edge_nodes1[0] == off_edge_nodes2[0] &&
816  off_edge_nodes1[1] == off_edge_nodes2[1]) ||
817  (off_edge_nodes1[1] == off_edge_nodes2[0] && off_edge_nodes1[0] == off_edge_nodes2[1]))
818  common_edge = COMMON_EDGE;
819  }
820  }
821  }
822  return common_edge;
823 }
const Elem * _side
std::vector< const Node * > _off_edge_nodes
virtual unsigned short dim() const=0

◆ isFaceReasonableCandidate()

bool PenetrationThread::isFaceReasonableCandidate ( const Elem primary_elem,
const Elem side,
libMesh::FEBase fe,
const libMesh::Point secondary_point,
const Real  tangential_tolerance 
)
protected

Definition at line 1323 of file PenetrationThread.C.

Referenced by createInfoForElem().

1328 {
1329  unsigned int dim = primary_elem->dim();
1330 
1331  const std::vector<Point> & phys_point = fe->get_xyz();
1332 
1333  const std::vector<RealGradient> & dxyz_dxi = fe->get_dxyzdxi();
1334  const std::vector<RealGradient> & dxyz_deta = fe->get_dxyzdeta();
1335 
1336  Point ref_point;
1337 
1338  std::vector<Point> points(1); // Default constructor gives us a point at 0,0,0
1339 
1340  fe->reinit(side, &points);
1341 
1342  RealGradient d = *secondary_point - phys_point[0];
1343 
1344  const Real twosqrt2 = 2.8284; // way more precision than we actually need here
1345  Real max_face_length = side->hmax() + twosqrt2 * tangential_tolerance;
1346 
1347  RealVectorValue normal;
1348  if (dim - 1 == 2)
1349  {
1350  normal = dxyz_dxi[0].cross(dxyz_deta[0]);
1351  }
1352  else if (dim - 1 == 1)
1353  {
1354  const Node * const * elem_nodes = primary_elem->get_nodes();
1355  const Point in_plane_vector1 = *elem_nodes[1] - *elem_nodes[0];
1356  const Point in_plane_vector2 = *elem_nodes[2] - *elem_nodes[0];
1357 
1358  Point out_of_plane_normal = in_plane_vector1.cross(in_plane_vector2);
1359  out_of_plane_normal /= out_of_plane_normal.norm();
1360 
1361  normal = dxyz_dxi[0].cross(out_of_plane_normal);
1362  }
1363  else
1364  {
1365  return true;
1366  }
1367  normal /= normal.norm();
1368 
1369  const Real dot(d * normal);
1370 
1371  const RealGradient normcomp = dot * normal;
1372  const RealGradient tangcomp = d - normcomp;
1373 
1374  const Real tangdist = tangcomp.norm();
1375 
1376  // Increase the size of the zone that we consider if the vector from the face
1377  // to the node has a larger normal component
1378  const Real faceExpansionFactor = 2.0 * (1.0 + normcomp.norm() / d.norm());
1379 
1380  bool isReasonableCandidate = true;
1381  if (tangdist > faceExpansionFactor * max_face_length)
1382  {
1383  isReasonableCandidate = false;
1384  }
1385  return isReasonableCandidate;
1386 }
auto norm() const -> decltype(std::norm(Real()))
virtual_for_inffe const std::vector< RealGradient > & get_dxyzdeta() const
static constexpr std::size_t dim
This is the dimension of all vector and tensor datastructures used in MOOSE.
Definition: Moose.h:162
virtual Real hmax() const
const Node *const * get_nodes() const
virtual void reinit(const Elem *elem, const std::vector< Point > *const pts=nullptr, const std::vector< Real > *const weights=nullptr)=0
virtual_for_inffe const std::vector< Point > & get_xyz() const
TypeVector< typename CompareTypes< Real, T2 >::supertype > cross(const TypeVector< T2 > &v) const
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
virtual unsigned short dim() const=0
virtual_for_inffe const std::vector< RealGradient > & get_dxyzdxi() const

◆ join()

void PenetrationThread::join ( const PenetrationThread other)

Definition at line 638 of file PenetrationThread.C.

639 {
641  other._recheck_secondary_nodes.begin(),
642  other._recheck_secondary_nodes.end());
643 }
std::vector< dof_id_type > _recheck_secondary_nodes
List of secondary nodes for which penetration was not detected in the current patch and for which pat...

◆ operator()()

void PenetrationThread::operator() ( const NodeIdRange range)

Definition at line 159 of file PenetrationThread.C.

160 {
161  ParallelUniqueId puid;
162  _tid = puid.id;
163 
164  // Must get the variables every time this is run because _tid can change
165  if (_do_normal_smoothing &&
167  {
168  _nodal_normal_x = &_subproblem.getStandardVariable(_tid, "nodal_normal_x");
169  _nodal_normal_y = &_subproblem.getStandardVariable(_tid, "nodal_normal_y");
170  _nodal_normal_z = &_subproblem.getStandardVariable(_tid, "nodal_normal_z");
171  }
172 
173  const BoundaryInfo & boundary_info = _mesh.getMesh().get_boundary_info();
174  std::unique_ptr<PointLocatorBase> point_locator;
175  if (_use_point_locator)
176  point_locator = _mesh.getPointLocator();
177 
178  for (const auto & node_id : range)
179  {
180  const Node & node = _mesh.nodeRef(node_id);
181 
182  // We're going to get a reference to the pointer for the pinfo for this node
183  // This will allow us to manipulate this pointer without having to go through
184  // the _penetration_info map... meaning this is the only mutex we'll have to do!
185  pinfo_mutex.lock();
188 
189  std::vector<PenetrationInfo *> p_info;
190  bool info_set(false);
191 
192  // See if we already have info about this node
193  if (info)
194  {
195  FEBase * fe_elem = _fes[_tid][info->_elem->dim()];
196  FEBase * fe_side = _fes[_tid][info->_side->dim()];
197 
198  if (!_update_location && (info->_distance >= 0 || info->isCaptured()))
199  {
200  const Point contact_ref = info->_closest_point_ref;
201  bool contact_point_on_side(false);
202 
203  // Secondary position must be the previous contact point
204  // Use the previous reference coordinates
205  std::vector<Point> points(1);
206  points[0] = contact_ref;
207  const std::vector<Point> & secondary_pos = fe_side->get_xyz();
208  bool search_succeeded = false;
209 
211  fe_elem,
212  fe_side,
213  _fe_type,
214  secondary_pos[0],
215  false,
217  contact_point_on_side,
218  search_succeeded);
219 
220  // Restore the original reference coordinates
221  info->_closest_point_ref = contact_ref;
222  // Just calculated as the distance of the contact point off the surface (0). Set to 0 to
223  // avoid round-off.
224  info->_distance = 0.0;
225  info_set = true;
226  }
227  else
228  {
229  Real old_tangential_distance(info->_tangential_distance);
230  bool contact_point_on_side(false);
231  bool search_succeeded = false;
232 
234  fe_elem,
235  fe_side,
236  _fe_type,
237  node,
238  false,
240  contact_point_on_side,
241  search_succeeded);
242 
243  if (contact_point_on_side)
244  {
245  if (info->_tangential_distance <= 0.0) // on the face
246  {
247  info_set = true;
248  }
249  else if (info->_tangential_distance > 0.0 && old_tangential_distance > 0.0)
250  { // off the face but within tolerance, was that way on the last step too
251  if (info->_side->dim() == 2 && info->_off_edge_nodes.size() < 2)
252  { // Closest point on face is on a node rather than an edge. Another
253  // face might be a better candidate.
254  }
255  else
256  {
257  info_set = true;
258  }
259  }
260  }
261  }
262  }
263 
264  if (!info_set)
265  {
266  const Node * closest_node = _nearest_node.nearestNode(node.id());
267 
268  std::vector<dof_id_type> located_elem_ids;
269  const std::vector<dof_id_type> * closest_elems;
270 
271  if (_use_point_locator)
272  {
273  std::set<const Elem *> candidate_elements;
274  (*point_locator)(*closest_node, candidate_elements);
275 
276  if (candidate_elements.empty())
277  mooseError("No proximate elements found at node ",
278  closest_node->id(),
279  " at ",
280  static_cast<const Point &>(*closest_node),
281  " on boundary ",
283  ". This should never happen.");
284 
285  for (const Elem * elem : candidate_elements)
286  {
287  for (auto s : elem->side_index_range())
288  if (boundary_info.has_boundary_id(elem, s, _primary_boundary))
289  {
290  located_elem_ids.push_back(elem->id());
291  break;
292  }
293  }
294 
295  if (located_elem_ids.empty())
296  mooseError("No proximate elements found at node ",
297  closest_node->id(),
298  " at ",
299  static_cast<const Point &>(*closest_node),
300  " on boundary ",
302  " share that boundary. This may happen if the mesh uses the same boundary id "
303  "for a nodeset and an unrelated sideset.");
304 
305  closest_elems = &located_elem_ids;
306  }
307  else
308  {
309  auto node_to_elem_pair = _node_to_elem_map.find(closest_node->id());
310  mooseAssert(node_to_elem_pair != _node_to_elem_map.end(),
311  "Missing entry in node to elem map");
312  closest_elems = &(node_to_elem_pair->second);
313  }
314 
315  for (const auto & elem_id : *closest_elems)
316  {
317  const Elem * elem = _mesh.elemPtr(elem_id);
318 
319  std::vector<PenetrationInfo *> thisElemInfo;
320 
321  std::vector<const Node *> nodesThatMustBeOnSide;
322  // If we have a disconnected mesh, we might not have *any*
323  // nodes that must be on a side we check; we'll rely on
324  // boundary info to find valid sides, then rely on comparing
325  // closest points from each to find the best.
326  //
327  // If we don't have a disconnected mesh, then for maximum
328  // backwards compatibility we're still using the older ridge
329  // and peak detection code, which depends on us ruling out
330  // sides that don't touch closest_node.
331  if (!_use_point_locator)
332  nodesThatMustBeOnSide.push_back(closest_node);
334  thisElemInfo, p_info, &node, elem, nodesThatMustBeOnSide, _check_whether_reasonable);
335  }
336 
337  if (_use_point_locator)
338  {
339  Real min_distance_sq = std::numeric_limits<Real>::max();
340  Point best_point;
341  unsigned int best_i = invalid_uint;
342 
343  // Find closest point in all p_info to the node of interest
344  for (unsigned int i = 0; i < p_info.size(); ++i)
345  {
346  const Point closest_point = closest_point_to_side(node, *p_info[i]->_side);
347  const Real distance_sq = (closest_point - node).norm_sq();
348  if (distance_sq < min_distance_sq)
349  {
350  min_distance_sq = distance_sq;
351  best_point = closest_point;
352  best_i = i;
353  }
354  }
355 
356  p_info[best_i]->_closest_point = best_point;
357  p_info[best_i]->_distance =
358  (p_info[best_i]->_distance >= 0.0 ? 1.0 : -1.0) * std::sqrt(min_distance_sq);
360  mooseError("Normal smoothing not implemented with point locator code");
361  Point normal = (best_point - node).unit();
362  const Real dot = normal * p_info[best_i]->_normal;
363  if (dot < 0)
364  normal *= -1;
365  p_info[best_i]->_normal = normal;
366 
367  switchInfo(info, p_info[best_i]);
368  info_set = true;
369  }
370  else
371  {
372  if (p_info.size() == 1)
373  {
374  if (p_info[0]->_tangential_distance <= _tangential_tolerance)
375  {
376  switchInfo(info, p_info[0]);
377  info_set = true;
378  }
379  }
380  else if (p_info.size() > 1)
381  {
382  // Loop through all pairs of faces, and check for contact on ridge between each face pair
383  std::vector<RidgeData> ridgeDataVec;
384  for (unsigned int i = 0; i + 1 < p_info.size(); ++i)
385  for (unsigned int j = i + 1; j < p_info.size(); ++j)
386  {
387  Point closest_coor;
388  Real tangential_distance(0.0);
389  const Node * closest_node_on_ridge = NULL;
390  unsigned int index = 0;
391  Point closest_coor_ref;
392  bool found_ridge_contact_point = findRidgeContactPoint(closest_coor,
393  tangential_distance,
394  closest_node_on_ridge,
395  index,
396  closest_coor_ref,
397  p_info,
398  i,
399  j);
400  if (found_ridge_contact_point)
401  {
402  RidgeData rpd;
403  rpd._closest_coor = closest_coor;
404  rpd._tangential_distance = tangential_distance;
405  rpd._closest_node = closest_node_on_ridge;
406  rpd._index = index;
407  rpd._closest_coor_ref = closest_coor_ref;
408  ridgeDataVec.push_back(rpd);
409  }
410  }
411 
412  if (ridgeDataVec.size() > 0) // Either find the ridge pair that is the best or find a peak
413  {
414  // Group together ridges for which we are off the edge of a common node.
415  // Those are peaks.
416  std::vector<RidgeSetData> ridgeSetDataVec;
417  for (unsigned int i = 0; i < ridgeDataVec.size(); ++i)
418  {
419  bool foundSetWithMatchingNode = false;
420  for (unsigned int j = 0; j < ridgeSetDataVec.size(); ++j)
421  {
422  if (ridgeDataVec[i]._closest_node != NULL &&
423  ridgeDataVec[i]._closest_node == ridgeSetDataVec[j]._closest_node)
424  {
425  foundSetWithMatchingNode = true;
426  ridgeSetDataVec[j]._ridge_data_vec.push_back(ridgeDataVec[i]);
427  break;
428  }
429  }
430  if (!foundSetWithMatchingNode)
431  {
432  RidgeSetData rsd;
433  rsd._distance = std::numeric_limits<Real>::max();
434  rsd._ridge_data_vec.push_back(ridgeDataVec[i]);
435  rsd._closest_node = ridgeDataVec[i]._closest_node;
436  ridgeSetDataVec.push_back(rsd);
437  }
438  }
439  // Compute distance to each set of ridges
440  for (unsigned int i = 0; i < ridgeSetDataVec.size(); ++i)
441  {
442  if (ridgeSetDataVec[i]._closest_node !=
443  NULL) // Either a peak or off the edge of single ridge
444  {
445  if (ridgeSetDataVec[i]._ridge_data_vec.size() == 1) // off edge of single ridge
446  {
447  if (ridgeSetDataVec[i]._ridge_data_vec[0]._tangential_distance <=
448  _tangential_tolerance) // off within tolerance
449  {
450  ridgeSetDataVec[i]._closest_coor =
451  ridgeSetDataVec[i]._ridge_data_vec[0]._closest_coor;
452  Point contact_point_vec = node - ridgeSetDataVec[i]._closest_coor;
453  ridgeSetDataVec[i]._distance = contact_point_vec.norm();
454  }
455  }
456  else // several ridges join at common node to make peak. The common node is the
457  // contact point
458  {
459  ridgeSetDataVec[i]._closest_coor = *ridgeSetDataVec[i]._closest_node;
460  Point contact_point_vec = node - ridgeSetDataVec[i]._closest_coor;
461  ridgeSetDataVec[i]._distance = contact_point_vec.norm();
462  }
463  }
464  else // on a single ridge
465  {
466  ridgeSetDataVec[i]._closest_coor =
467  ridgeSetDataVec[i]._ridge_data_vec[0]._closest_coor;
468  Point contact_point_vec = node - ridgeSetDataVec[i]._closest_coor;
469  ridgeSetDataVec[i]._distance = contact_point_vec.norm();
470  }
471  }
472  // Find the set of ridges closest to us.
473  unsigned int closest_ridge_set_index(0);
474  Real closest_distance(ridgeSetDataVec[0]._distance);
475  Point closest_point(ridgeSetDataVec[0]._closest_coor);
476  for (unsigned int i = 1; i < ridgeSetDataVec.size(); ++i)
477  {
478  if (ridgeSetDataVec[i]._distance < closest_distance)
479  {
480  closest_ridge_set_index = i;
481  closest_distance = ridgeSetDataVec[i]._distance;
482  closest_point = ridgeSetDataVec[i]._closest_coor;
483  }
484  }
485 
486  if (closest_distance <
487  std::numeric_limits<Real>::max()) // contact point is on the closest ridge set
488  {
489  // find the face in the ridge set with the smallest index, assign that one to the
490  // interaction
491  unsigned int face_index(std::numeric_limits<unsigned int>::max());
492  for (unsigned int i = 0;
493  i < ridgeSetDataVec[closest_ridge_set_index]._ridge_data_vec.size();
494  ++i)
495  {
496  if (ridgeSetDataVec[closest_ridge_set_index]._ridge_data_vec[i]._index < face_index)
497  face_index = ridgeSetDataVec[closest_ridge_set_index]._ridge_data_vec[i]._index;
498  }
499 
500  mooseAssert(face_index < std::numeric_limits<unsigned int>::max(),
501  "face_index invalid");
502 
503  p_info[face_index]->_closest_point = closest_point;
504  p_info[face_index]->_distance =
505  (p_info[face_index]->_distance >= 0.0 ? 1.0 : -1.0) * closest_distance;
506  // Calculate the normal as the vector from the ridge to the point only if we're not
507  // doing normal
508  // smoothing. Normal smoothing will average out the normals on its own.
510  {
511  Point normal(closest_point - node);
512  const Real len(normal.norm());
513  if (len > 0)
514  {
515  normal /= len;
516  }
517  const Real dot(normal * p_info[face_index]->_normal);
518  if (dot < 0)
519  {
520  normal *= -1;
521  }
522  p_info[face_index]->_normal = normal;
523  }
524  p_info[face_index]->_tangential_distance = 0.0;
525 
526  Point closest_point_ref;
527  if (ridgeSetDataVec[closest_ridge_set_index]._ridge_data_vec.size() ==
528  1) // contact with a single ridge rather than a peak
529  {
530  p_info[face_index]->_tangential_distance = ridgeSetDataVec[closest_ridge_set_index]
531  ._ridge_data_vec[0]
532  ._tangential_distance;
533  p_info[face_index]->_closest_point_ref =
534  ridgeSetDataVec[closest_ridge_set_index]._ridge_data_vec[0]._closest_coor_ref;
535  }
536  else
537  { // peak
538  const Node * closest_node_on_face;
539  bool restricted = restrictPointToFace(p_info[face_index]->_closest_point_ref,
540  closest_node_on_face,
541  p_info[face_index]->_side);
542  if (restricted)
543  {
544  if (closest_node_on_face !=
545  ridgeSetDataVec[closest_ridge_set_index]._closest_node)
546  {
547  mooseError("Closest node when restricting point to face != closest node from "
548  "RidgeSetData");
549  }
550  }
551  }
552 
553  FEBase * fe = _fes[_tid][p_info[face_index]->_side->dim()];
554  std::vector<Point> points(1);
555  points[0] = p_info[face_index]->_closest_point_ref;
556  fe->reinit(p_info[face_index]->_side, &points);
557  p_info[face_index]->_side_phi = fe->get_phi();
558  p_info[face_index]->_side_grad_phi = fe->get_dphi();
559  p_info[face_index]->_dxyzdxi = fe->get_dxyzdxi();
560  p_info[face_index]->_dxyzdeta = fe->get_dxyzdeta();
561  p_info[face_index]->_d2xyzdxideta = fe->get_d2xyzdxideta();
562 
563  switchInfo(info, p_info[face_index]);
564  info_set = true;
565  }
566  else
567  { // todo:remove invalid ridge cases so they don't mess up individual face
568  // competition????
569  }
570  }
571 
572  if (!info_set) // contact wasn't on a ridge -- compete individual interactions
573  {
574  unsigned int best(0), i(1);
575  do
576  {
577  CompeteInteractionResult CIResult = competeInteractions(p_info[best], p_info[i]);
578  if (CIResult == FIRST_WINS)
579  {
580  i++;
581  }
582  else if (CIResult == SECOND_WINS)
583  {
584  best = i;
585  i++;
586  }
587  else if (CIResult == NEITHER_WINS)
588  {
589  best = i + 1;
590  i += 2;
591  }
592  } while (i < p_info.size() && best < p_info.size());
593  if (best < p_info.size())
594  {
595  // Ensure final info is within the tangential tolerance
596  if (p_info[best]->_tangential_distance <= _tangential_tolerance)
597  {
598  switchInfo(info, p_info[best]);
599  info_set = true;
600  }
601  }
602  }
603  }
604  }
605  }
606 
607  if (!info_set)
608  {
609  // If penetration is not detected within the saved patch, it is possible
610  // that the secondary node has moved outside the saved patch. So, the patch
611  // for the secondary nodes saved in _recheck_secondary_nodes has to be updated
612  // and the penetration detection has to be re-run on the updated patch.
613 
614  _recheck_secondary_nodes.push_back(node_id);
615 
616  delete info;
617  info = NULL;
618  }
619  else
620  {
621  smoothNormal(info, p_info, node);
622  FEBase * fe = _fes[_tid][info->_side->dim()];
623  computeSlip(*fe, *info);
624  }
625 
626  for (unsigned int j = 0; j < p_info.size(); ++j)
627  {
628  if (p_info[j])
629  {
630  delete p_info[j];
631  p_info[j] = NULL;
632  }
633  }
634  }
635 }
MooseVariable * _nodal_normal_z
auto norm() const -> decltype(std::norm(Real()))
bool has_boundary_id(const Node *const node, const boundary_id_type id) const
const unsigned int invalid_uint
virtual_for_inffe const std::vector< RealGradient > & get_dxyzdeta() const
virtual Elem * elemPtr(const dof_id_type i)
Definition: MooseMesh.C:3134
MPI_Info info
bool findRidgeContactPoint(libMesh::Point &contact_point, Real &tangential_distance, const Node *&closest_node, unsigned int &index, libMesh::Point &contact_point_ref, std::vector< PenetrationInfo *> &p_info, const unsigned int index1, const unsigned int index2)
void mooseError(Args &&... args)
Emit an error message with the given stringified, concatenated args and terminate the application...
Definition: MooseError.h:323
Data structure used to hold penetration information.
NearestNodeLocator & _nearest_node
virtual_for_inffe const std::vector< RealGradient > & get_d2xyzdxideta() const
void createInfoForElem(std::vector< PenetrationInfo *> &thisElemInfo, std::vector< PenetrationInfo *> &p_info, const Node *secondary_node, const Elem *elem, const std::vector< const Node *> &nodes_that_must_be_on_side, const bool check_whether_reasonable=false)
Threads::spin_mutex pinfo_mutex
const BoundaryInfo & get_boundary_info() const
virtual const Node & nodeRef(const dof_id_type i) const
Definition: MooseMesh.C:849
void smoothNormal(PenetrationInfo *info, std::vector< PenetrationInfo *> &p_info, const Node &node)
libMesh::FEType & _fe_type
auto max(const L &left, const R &right)
BoundaryID _primary_boundary
SubProblem & _subproblem
const std::vector< std::vector< OutputGradient > > & get_dphi() const
const MooseMesh & _mesh
dof_id_type id() const
MeshBase & getMesh()
Accessor for the underlying libMesh Mesh object.
Definition: MooseMesh.C:3469
CompeteInteractionResult competeInteractions(PenetrationInfo *pi1, PenetrationInfo *pi2)
When interactions are identified between a node and two faces, compete between the faces to determine...
void computeSlip(libMesh::FEBase &fe, PenetrationInfo &info)
virtual void reinit(const Elem *elem, const std::vector< Point > *const pts=nullptr, const std::vector< Real > *const weights=nullptr)=0
void findContactPoint(PenetrationInfo &p_info, libMesh::FEBase *fe_elem, libMesh::FEBase *fe_side, libMesh::FEType &fe_side_type, const libMesh::Point &secondary_point, bool start_with_centroid, const Real tangential_tolerance, bool &contact_point_on_side, bool &search_succeeded)
Finds the closest point (called the contact point) on the primary_elem on side "side" to the secondar...
std::vector< dof_id_type > _recheck_secondary_nodes
List of secondary nodes for which penetration was not detected in the current patch and for which pat...
virtual MooseVariable & getStandardVariable(const THREAD_ID tid, const std::string &var_name)=0
Returns the variable reference for requested MooseVariable which may be in any system.
virtual_for_inffe const std::vector< Point > & get_xyz() const
bool restrictPointToFace(libMesh::Point &p, const Node *&closest_node, const Elem *side)
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
std::map< dof_id_type, PenetrationInfo * > & _penetration_info
const Node * nearestNode(dof_id_type node_id)
Valid to call this after findNodes() has been called to get a pointer to the nearest node...
CTSub CT_OPERATOR_BINARY CTMul CTCompareLess CTCompareGreater CTCompareEqual _arg template * sqrt(_arg)) *_arg.template D< dtag >()) CT_SIMPLE_UNARY_FUNCTION(tanh
MooseVariable * _nodal_normal_y
auto norm_sq(const T &a) -> decltype(std::norm(a))
virtual std::unique_ptr< libMesh::PointLocatorBase > getPointLocator() const
Proxy function to get a (sub)PointLocator from either the underlying libMesh mesh (default)...
Definition: MooseMesh.C:3754
const std::map< dof_id_type, std::vector< dof_id_type > > & _node_to_elem_map
void switchInfo(PenetrationInfo *&info, PenetrationInfo *&infoNew)
virtual_for_inffe const std::vector< RealGradient > & get_dxyzdxi() const
MooseVariable * _nodal_normal_x
std::vector< std::vector< libMesh::FEBase * > > & _fes
PenetrationLocator::NORMAL_SMOOTHING_METHOD _normal_smoothing_method
const std::vector< std::vector< OutputShape > > & get_phi() const

◆ restrictPointToFace()

bool PenetrationThread::restrictPointToFace ( libMesh::Point p,
const Node *&  closest_node,
const Elem side 
)
protected

Definition at line 1142 of file PenetrationThread.C.

Referenced by operator()().

1143 {
1144  const ElemType t(side->type());
1145  Real & xi = p(0);
1146  Real & eta = p(1);
1147  closest_node = NULL;
1148 
1149  bool off_of_this_face(false);
1150 
1151  switch (t)
1152  {
1153  case EDGE2:
1154  case EDGE3:
1155  case EDGE4:
1156  {
1157  if (xi < -1.0)
1158  {
1159  xi = -1.0;
1160  off_of_this_face = true;
1161  closest_node = side->node_ptr(0);
1162  }
1163  else if (xi > 1.0)
1164  {
1165  xi = 1.0;
1166  off_of_this_face = true;
1167  closest_node = side->node_ptr(1);
1168  }
1169  break;
1170  }
1171 
1172  case TRI3:
1173  case TRI6:
1174  case TRI7:
1175  {
1176  if (eta < 0.0)
1177  {
1178  eta = 0.0;
1179  off_of_this_face = true;
1180  if (xi < 0.5)
1181  {
1182  closest_node = side->node_ptr(0);
1183  if (xi < 0.0)
1184  xi = 0.0;
1185  }
1186  else
1187  {
1188  closest_node = side->node_ptr(1);
1189  if (xi > 1.0)
1190  xi = 1.0;
1191  }
1192  }
1193  else if ((xi + eta) > 1.0)
1194  {
1195  Real delta = (xi + eta - 1.0) / 2.0;
1196  xi -= delta;
1197  eta -= delta;
1198  off_of_this_face = true;
1199  if (xi > 0.5)
1200  {
1201  closest_node = side->node_ptr(1);
1202  if (xi > 1.0)
1203  {
1204  xi = 1.0;
1205  eta = 0.0;
1206  }
1207  }
1208  else
1209  {
1210  closest_node = side->node_ptr(2);
1211  if (xi < 0.0)
1212  {
1213  xi = 0.0;
1214  eta = 1.0;
1215  }
1216  }
1217  }
1218  else if (xi < 0.0)
1219  {
1220  xi = 0.0;
1221  off_of_this_face = true;
1222  if (eta > 0.5)
1223  {
1224  closest_node = side->node_ptr(2);
1225  if (eta > 1.0)
1226  eta = 1.0;
1227  }
1228  else
1229  {
1230  closest_node = side->node_ptr(0);
1231  if (eta < 0.0)
1232  eta = 0.0;
1233  }
1234  }
1235  break;
1236  }
1237 
1238  case QUAD4:
1239  case QUAD8:
1240  case QUAD9:
1241  {
1242  if (eta < -1.0)
1243  {
1244  eta = -1.0;
1245  off_of_this_face = true;
1246  if (xi < 0.0)
1247  {
1248  closest_node = side->node_ptr(0);
1249  if (xi < -1.0)
1250  xi = -1.0;
1251  }
1252  else
1253  {
1254  closest_node = side->node_ptr(1);
1255  if (xi > 1.0)
1256  xi = 1.0;
1257  }
1258  }
1259  else if (xi > 1.0)
1260  {
1261  xi = 1.0;
1262  off_of_this_face = true;
1263  if (eta < 0.0)
1264  {
1265  closest_node = side->node_ptr(1);
1266  if (eta < -1.0)
1267  eta = -1.0;
1268  }
1269  else
1270  {
1271  closest_node = side->node_ptr(2);
1272  if (eta > 1.0)
1273  eta = 1.0;
1274  }
1275  }
1276  else if (eta > 1.0)
1277  {
1278  eta = 1.0;
1279  off_of_this_face = true;
1280  if (xi < 0.0)
1281  {
1282  closest_node = side->node_ptr(3);
1283  if (xi < -1.0)
1284  xi = -1.0;
1285  }
1286  else
1287  {
1288  closest_node = side->node_ptr(2);
1289  if (xi > 1.0)
1290  xi = 1.0;
1291  }
1292  }
1293  else if (xi < -1.0)
1294  {
1295  xi = -1.0;
1296  off_of_this_face = true;
1297  if (eta < 0.0)
1298  {
1299  closest_node = side->node_ptr(0);
1300  if (eta < -1.0)
1301  eta = -1.0;
1302  }
1303  else
1304  {
1305  closest_node = side->node_ptr(3);
1306  if (eta > 1.0)
1307  eta = 1.0;
1308  }
1309  }
1310  break;
1311  }
1312 
1313  default:
1314  {
1315  mooseError("Unsupported face type: ", t);
1316  break;
1317  }
1318  }
1319  return off_of_this_face;
1320 }
ElemType
QUAD8
void mooseError(Args &&... args)
Emit an error message with the given stringified, concatenated args and terminate the application...
Definition: MooseError.h:323
EDGE4
TRI3
QUAD4
TRI6
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
EDGE2
const Node * node_ptr(const unsigned int i) const
TRI7
QUAD9
EDGE3
virtual ElemType type() const=0

◆ restrictPointToSpecifiedEdgeOfFace()

bool PenetrationThread::restrictPointToSpecifiedEdgeOfFace ( libMesh::Point p,
const Node *&  closest_node,
const Elem side,
const std::vector< const Node *> &  edge_nodes 
)
protected

Definition at line 968 of file PenetrationThread.C.

Referenced by findRidgeContactPoint().

972 {
973  const ElemType t = side->type();
974  Real & xi = p(0);
975  Real & eta = p(1);
976  closest_node = NULL;
977 
978  std::vector<unsigned int> local_node_indices;
979  for (const auto & edge_node : edge_nodes)
980  {
981  unsigned int local_index = side->get_node_index(edge_node);
982  if (local_index == libMesh::invalid_uint)
983  mooseError("Side does not contain node");
984  local_node_indices.push_back(local_index);
985  }
986  mooseAssert(local_node_indices.size() == side->dim(),
987  "Number of edge nodes must match side dimensionality");
988  std::sort(local_node_indices.begin(), local_node_indices.end());
989 
990  bool off_of_this_edge = false;
991 
992  switch (t)
993  {
994  case EDGE2:
995  case EDGE3:
996  case EDGE4:
997  {
998  if (local_node_indices[0] == 0)
999  {
1000  if (xi <= -1.0)
1001  {
1002  xi = -1.0;
1003  off_of_this_edge = true;
1004  closest_node = side->node_ptr(0);
1005  }
1006  }
1007  else if (local_node_indices[0] == 1)
1008  {
1009  if (xi >= 1.0)
1010  {
1011  xi = 1.0;
1012  off_of_this_edge = true;
1013  closest_node = side->node_ptr(1);
1014  }
1015  }
1016  else
1017  {
1018  mooseError("Invalid local node indices");
1019  }
1020  break;
1021  }
1022 
1023  case TRI3:
1024  case TRI6:
1025  case TRI7:
1026  {
1027  if ((local_node_indices[0] == 0) && (local_node_indices[1] == 1))
1028  {
1029  if (eta <= 0.0)
1030  {
1031  eta = 0.0;
1032  off_of_this_edge = true;
1033  if (xi < 0.0)
1034  closest_node = side->node_ptr(0);
1035  else if (xi > 1.0)
1036  closest_node = side->node_ptr(1);
1037  }
1038  }
1039  else if ((local_node_indices[0] == 1) && (local_node_indices[1] == 2))
1040  {
1041  if ((xi + eta) > 1.0)
1042  {
1043  Real delta = (xi + eta - 1.0) / 2.0;
1044  xi -= delta;
1045  eta -= delta;
1046  off_of_this_edge = true;
1047  if (xi > 1.0)
1048  closest_node = side->node_ptr(1);
1049  else if (xi < 0.0)
1050  closest_node = side->node_ptr(2);
1051  }
1052  }
1053  else if ((local_node_indices[0] == 0) && (local_node_indices[1] == 2))
1054  {
1055  if (xi <= 0.0)
1056  {
1057  xi = 0.0;
1058  off_of_this_edge = true;
1059  if (eta > 1.0)
1060  closest_node = side->node_ptr(2);
1061  else if (eta < 0.0)
1062  closest_node = side->node_ptr(0);
1063  }
1064  }
1065  else
1066  {
1067  mooseError("Invalid local node indices");
1068  }
1069 
1070  break;
1071  }
1072 
1073  case QUAD4:
1074  case QUAD8:
1075  case QUAD9:
1076  {
1077  if ((local_node_indices[0] == 0) && (local_node_indices[1] == 1))
1078  {
1079  if (eta <= -1.0)
1080  {
1081  eta = -1.0;
1082  off_of_this_edge = true;
1083  if (xi < -1.0)
1084  closest_node = side->node_ptr(0);
1085  else if (xi > 1.0)
1086  closest_node = side->node_ptr(1);
1087  }
1088  }
1089  else if ((local_node_indices[0] == 1) && (local_node_indices[1] == 2))
1090  {
1091  if (xi >= 1.0)
1092  {
1093  xi = 1.0;
1094  off_of_this_edge = true;
1095  if (eta < -1.0)
1096  closest_node = side->node_ptr(1);
1097  else if (eta > 1.0)
1098  closest_node = side->node_ptr(2);
1099  }
1100  }
1101  else if ((local_node_indices[0] == 2) && (local_node_indices[1] == 3))
1102  {
1103  if (eta >= 1.0)
1104  {
1105  eta = 1.0;
1106  off_of_this_edge = true;
1107  if (xi < -1.0)
1108  closest_node = side->node_ptr(3);
1109  else if (xi > 1.0)
1110  closest_node = side->node_ptr(2);
1111  }
1112  }
1113  else if ((local_node_indices[0] == 0) && (local_node_indices[1] == 3))
1114  {
1115  if (xi <= -1.0)
1116  {
1117  xi = -1.0;
1118  off_of_this_edge = true;
1119  if (eta < -1.0)
1120  closest_node = side->node_ptr(0);
1121  else if (eta > 1.0)
1122  closest_node = side->node_ptr(3);
1123  }
1124  }
1125  else
1126  {
1127  mooseError("Invalid local node indices");
1128  }
1129  break;
1130  }
1131 
1132  default:
1133  {
1134  mooseError("Unsupported face type: ", t);
1135  break;
1136  }
1137  }
1138  return off_of_this_edge;
1139 }
ElemType
QUAD8
const unsigned int invalid_uint
unsigned int get_node_index(const Node *node_ptr) const
void mooseError(Args &&... args)
Emit an error message with the given stringified, concatenated args and terminate the application...
Definition: MooseError.h:323
EDGE4
TRI3
QUAD4
TRI6
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
EDGE2
virtual unsigned short dim() const=0
const Node * node_ptr(const unsigned int i) const
TRI7
QUAD9
EDGE3
virtual ElemType type() const=0

◆ smoothNormal()

void PenetrationThread::smoothNormal ( PenetrationInfo info,
std::vector< PenetrationInfo *> &  p_info,
const Node node 
)
protected

Definition at line 1408 of file PenetrationThread.C.

Referenced by operator()().

1411 {
1413  {
1415  {
1416  // If we are within the smoothing distance of any edges or corners, find the
1417  // corner nodes for those edges/corners, and weights from distance to edge/corner
1418  std::vector<Real> edge_face_weights;
1419  std::vector<PenetrationInfo *> edge_face_info;
1420 
1421  getSmoothingFacesAndWeights(info, edge_face_info, edge_face_weights, p_info, node);
1422 
1423  mooseAssert(edge_face_info.size() == edge_face_weights.size(),
1424  "edge_face_info.size() != edge_face_weights.size()");
1425 
1426  if (edge_face_info.size() > 0)
1427  {
1428  // Smooth the normal using the weighting functions for all participating faces.
1429  RealVectorValue new_normal;
1430  Real this_face_weight = 1.0;
1431 
1432  for (unsigned int efwi = 0; efwi < edge_face_weights.size(); ++efwi)
1433  {
1434  PenetrationInfo * npi = edge_face_info[efwi];
1435  if (npi)
1436  new_normal += npi->_normal * edge_face_weights[efwi];
1437 
1438  this_face_weight -= edge_face_weights[efwi];
1439  }
1440  mooseAssert(this_face_weight >= (0.25 - 1e-8),
1441  "Sum of weights of other faces shouldn't exceed 0.75");
1442  new_normal += info->_normal * this_face_weight;
1443 
1444  const Real len = new_normal.norm();
1445  if (len > 0)
1446  new_normal /= len;
1447 
1448  info->_normal = new_normal;
1449  }
1450  }
1452  {
1453  // params.addParam<VariableName>("var_name","description");
1454  // getParam<VariableName>("var_name")
1455  info->_normal(0) = _nodal_normal_x->getValue(info->_side, info->_side_phi);
1456  info->_normal(1) = _nodal_normal_y->getValue(info->_side, info->_side_phi);
1457  info->_normal(2) = _nodal_normal_z->getValue(info->_side, info->_side_phi);
1458  const Real len(info->_normal.norm());
1459  if (len > 0)
1460  info->_normal /= len;
1461  }
1462  }
1463 }
MooseVariable * _nodal_normal_z
auto norm() const -> decltype(std::norm(Real()))
RealVectorValue _normal
MPI_Info info
OutputType getValue(const Elem *elem, const std::vector< std::vector< OutputShape >> &phi) const
Compute the variable value at a point on an element.
Data structure used to hold penetration information.
void getSmoothingFacesAndWeights(PenetrationInfo *info, std::vector< PenetrationInfo *> &edge_face_info, std::vector< Real > &edge_face_weights, std::vector< PenetrationInfo *> &p_info, const Node &secondary_node)
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
MooseVariable * _nodal_normal_y
MooseVariable * _nodal_normal_x
PenetrationLocator::NORMAL_SMOOTHING_METHOD _normal_smoothing_method

◆ switchInfo()

void PenetrationThread::switchInfo ( PenetrationInfo *&  info,
PenetrationInfo *&  infoNew 
)
protected

Definition at line 646 of file PenetrationThread.C.

Referenced by operator()().

647 {
648  mooseAssert(infoNew != NULL, "infoNew object is null");
649  if (info)
650  {
651  infoNew->_starting_elem = info->_starting_elem;
652  infoNew->_starting_side_num = info->_starting_side_num;
653  infoNew->_starting_closest_point_ref = info->_starting_closest_point_ref;
654  infoNew->_incremental_slip = info->_incremental_slip;
655  infoNew->_accumulated_slip = info->_accumulated_slip;
656  infoNew->_accumulated_slip_old = info->_accumulated_slip_old;
657  infoNew->_frictional_energy = info->_frictional_energy;
658  infoNew->_frictional_energy_old = info->_frictional_energy_old;
659  infoNew->_contact_force = info->_contact_force;
660  infoNew->_contact_force_old = info->_contact_force_old;
661  infoNew->_lagrange_multiplier = info->_lagrange_multiplier;
662  infoNew->_lagrange_multiplier_slip = info->_lagrange_multiplier_slip;
663  infoNew->_locked_this_step = info->_locked_this_step;
664  infoNew->_stick_locked_this_step = info->_stick_locked_this_step;
665  infoNew->_mech_status = info->_mech_status;
666  infoNew->_mech_status_old = info->_mech_status_old;
667  }
668  else
669  {
670  infoNew->_starting_elem = infoNew->_elem;
671  infoNew->_starting_side_num = infoNew->_side_num;
673  }
674  delete info;
675  info = infoNew;
676  infoNew = NULL; // Set this to NULL so that we don't delete it (now owned by _penetration_info).
677 }
MPI_Info info
const Elem * _starting_elem
unsigned int _stick_locked_this_step
RealVectorValue _contact_force_old
unsigned int _locked_this_step
const Elem * _elem
unsigned int _starting_side_num
unsigned int _side_num
RealVectorValue _contact_force
MECH_STATUS_ENUM _mech_status
Point _starting_closest_point_ref
MECH_STATUS_ENUM _mech_status_old
RealVectorValue _lagrange_multiplier_slip

Member Data Documentation

◆ _check_whether_reasonable

bool PenetrationThread::_check_whether_reasonable
protected

Definition at line 64 of file PenetrationThread.h.

Referenced by operator()().

◆ _do_normal_smoothing

bool PenetrationThread::_do_normal_smoothing
protected

Definition at line 67 of file PenetrationThread.h.

Referenced by operator()(), and smoothNormal().

◆ _elem_side_builder

libMesh::ElemSideBuilder PenetrationThread::_elem_side_builder
protected

Helper for building element sides without extraneous allocation.

Definition at line 86 of file PenetrationThread.h.

Referenced by computeSlip().

◆ _fe_type

libMesh::FEType& PenetrationThread::_fe_type
protected

Definition at line 77 of file PenetrationThread.h.

Referenced by createInfoForElem(), and operator()().

◆ _fes

std::vector<std::vector<libMesh::FEBase *> >& PenetrationThread::_fes
protected

Definition at line 75 of file PenetrationThread.h.

Referenced by createInfoForElem(), findRidgeContactPoint(), and operator()().

◆ _mesh

const MooseMesh& PenetrationThread::_mesh
protected

◆ _nearest_node

NearestNodeLocator& PenetrationThread::_nearest_node
protected

Definition at line 79 of file PenetrationThread.h.

Referenced by operator()().

◆ _nodal_normal_x

MooseVariable* PenetrationThread::_nodal_normal_x
protected

Definition at line 71 of file PenetrationThread.h.

Referenced by operator()(), and smoothNormal().

◆ _nodal_normal_y

MooseVariable* PenetrationThread::_nodal_normal_y
protected

Definition at line 72 of file PenetrationThread.h.

Referenced by operator()(), and smoothNormal().

◆ _nodal_normal_z

MooseVariable* PenetrationThread::_nodal_normal_z
protected

Definition at line 73 of file PenetrationThread.h.

Referenced by operator()(), and smoothNormal().

◆ _node_to_elem_map

const std::map<dof_id_type, std::vector<dof_id_type> >& PenetrationThread::_node_to_elem_map
protected

Definition at line 81 of file PenetrationThread.h.

Referenced by getInfoForFacesWithCommonNodes(), and operator()().

◆ _normal_smoothing_distance

Real PenetrationThread::_normal_smoothing_distance
protected

Definition at line 68 of file PenetrationThread.h.

Referenced by getSmoothingEdgeNodesAndWeights().

◆ _normal_smoothing_method

PenetrationLocator::NORMAL_SMOOTHING_METHOD PenetrationThread::_normal_smoothing_method
protected

Definition at line 69 of file PenetrationThread.h.

Referenced by operator()(), and smoothNormal().

◆ _penetration_info

std::map<dof_id_type, PenetrationInfo *>& PenetrationThread::_penetration_info
protected

Definition at line 62 of file PenetrationThread.h.

Referenced by operator()().

◆ _primary_boundary

BoundaryID PenetrationThread::_primary_boundary
protected

Definition at line 58 of file PenetrationThread.h.

Referenced by createInfoForElem(), and operator()().

◆ _recheck_secondary_nodes

std::vector<dof_id_type> PenetrationThread::_recheck_secondary_nodes

List of secondary nodes for which penetration was not detected in the current patch and for which patch has to be updated.

Definition at line 52 of file PenetrationThread.h.

Referenced by join(), and operator()().

◆ _secondary_boundary

BoundaryID PenetrationThread::_secondary_boundary
protected

Definition at line 59 of file PenetrationThread.h.

◆ _subproblem

SubProblem& PenetrationThread::_subproblem
protected

Definition at line 55 of file PenetrationThread.h.

Referenced by operator()().

◆ _tangential_tolerance

Real PenetrationThread::_tangential_tolerance
protected

Definition at line 66 of file PenetrationThread.h.

Referenced by competeInteractions(), createInfoForElem(), and operator()().

◆ _tid

THREAD_ID PenetrationThread::_tid
protected

Definition at line 83 of file PenetrationThread.h.

Referenced by createInfoForElem(), findRidgeContactPoint(), and operator()().

◆ _update_location

bool PenetrationThread::_update_location
protected

Definition at line 65 of file PenetrationThread.h.

Referenced by operator()().

◆ _use_point_locator

bool PenetrationThread::_use_point_locator
protected

Definition at line 70 of file PenetrationThread.h.

Referenced by operator()().


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