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 659 of file PenetrationThread.C.

Referenced by operator()().

660 {
661 
663 
665  pi2->_tangential_distance > _tangential_tolerance) // out of tol on both faces
666  result = NEITHER_WINS;
667 
668  else if (pi1->_tangential_distance == 0.0 &&
669  pi2->_tangential_distance > 0.0) // on face 1, off face 2
670  result = FIRST_WINS;
671 
672  else if (pi2->_tangential_distance == 0.0 &&
673  pi1->_tangential_distance > 0.0) // on face 2, off face 1
674  result = SECOND_WINS;
675 
676  else if (pi1->_tangential_distance <= _tangential_tolerance &&
677  pi2->_tangential_distance > _tangential_tolerance) // in face 1 tol, out of face 2 tol
678  result = FIRST_WINS;
679 
680  else if (pi2->_tangential_distance <= _tangential_tolerance &&
681  pi1->_tangential_distance > _tangential_tolerance) // in face 2 tol, out of face 1 tol
682  result = SECOND_WINS;
683 
684  else if (pi1->_tangential_distance == 0.0 && pi2->_tangential_distance == 0.0) // on both faces
685  result = competeInteractionsBothOnFace(pi1, pi2);
686 
687  else if (pi1->_tangential_distance <= _tangential_tolerance &&
688  pi2->_tangential_distance <= _tangential_tolerance) // off but within tol of both faces
689  {
691  if (cer == COMMON_EDGE || cer == COMMON_NODE) // ridge case.
692  {
693  // We already checked for ridges, and it got rejected, so neither face must be valid
694  result = NEITHER_WINS;
695  // mooseError("Erroneously encountered ridge case");
696  }
697  else if (cer == EDGE_AND_COMMON_NODE) // off side of face, off corner of another face. Favor
698  // the off-side face
699  {
700  if (pi1->_off_edge_nodes.size() == pi2->_off_edge_nodes.size())
701  mooseError("Invalid off_edge_nodes counts");
702 
703  else if (pi1->_off_edge_nodes.size() == 2)
704  result = FIRST_WINS;
705 
706  else if (pi2->_off_edge_nodes.size() == 2)
707  result = SECOND_WINS;
708 
709  else
710  mooseError("Invalid off_edge_nodes counts");
711  }
712  else // The node projects to both faces within tangential tolerance.
713  result = competeInteractionsBothOnFace(pi1, pi2);
714  }
715 
716  return result;
717 }
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 720 of file PenetrationThread.C.

Referenced by competeInteractions().

721 {
723 
724  if (pi1->_distance >= 0.0 && pi2->_distance < 0.0)
725  result = FIRST_WINS; // favor face with positive distance (penetrated) -- first in this case
726 
727  else if (pi2->_distance >= 0.0 && pi1->_distance < 0.0)
728  result = SECOND_WINS; // favor face with positive distance (penetrated) -- second in this case
729 
730  // TODO: This logic below could cause an abrupt jump from one face to the other with small mesh
731  // movement. If there is some way to smooth the transition, we should do it.
733  result = FIRST_WINS; // otherwise, favor the closer face -- first in this case
734 
736  result = SECOND_WINS; // otherwise, favor the closer face -- second in this case
737 
738  else // Equal within tolerance. Favor the one with a smaller element id (for repeatibility)
739  {
740  if (pi1->_elem->id() < pi2->_elem->id())
741  result = FIRST_WINS;
742 
743  else
744  result = SECOND_WINS;
745  }
746 
747  return result;
748 }
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:629

◆ computeSlip()

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

Definition at line 1368 of file PenetrationThread.C.

Referenced by operator()().

1369 {
1370  // Slip is current projected position of secondary node minus
1371  // original projected position of secondary node
1372  std::vector<Point> points(1);
1373  points[0] = info._starting_closest_point_ref;
1374  const auto & side = _elem_side_builder(*info._starting_elem, info._starting_side_num);
1375  fe.reinit(&side, &points);
1376  const std::vector<Point> & starting_point = fe.get_xyz();
1377  info._incremental_slip = info._closest_point - starting_point[0];
1378  if (info.isCaptured())
1379  {
1380  info._frictional_energy =
1381  info._frictional_energy_old + info._contact_force * info._incremental_slip;
1382  info._accumulated_slip = info._accumulated_slip_old + info._incremental_slip.norm();
1383  }
1384 }
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 1786 of file PenetrationThread.C.

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

1792 {
1793  const BoundaryInfo & boundary_info = _mesh.getMesh().get_boundary_info();
1794 
1795  for (auto s : elem->side_index_range())
1796  {
1797  if (!boundary_info.has_boundary_id(elem, s, _primary_boundary))
1798  continue;
1799 
1800  // Don't create info for this side if one already exists
1801  bool already_have_info_this_side = false;
1802  for (const auto & pi : thisElemInfo)
1803  if (pi->_side_num == s)
1804  {
1805  already_have_info_this_side = true;
1806  break;
1807  }
1808 
1809  if (already_have_info_this_side)
1810  break;
1811 
1812  const Elem * side = elem->build_side_ptr(s).release();
1813 
1814  // Only continue with creating info for this side if the side contains
1815  // all of the nodes in nodes_that_must_be_on_side
1816  std::vector<const Node *> nodevec;
1817  for (unsigned int ni = 0; ni < side->n_nodes(); ++ni)
1818  nodevec.push_back(side->node_ptr(ni));
1819 
1820  std::sort(nodevec.begin(), nodevec.end());
1821  std::vector<const Node *> common_nodes;
1822  std::set_intersection(nodes_that_must_be_on_side.begin(),
1823  nodes_that_must_be_on_side.end(),
1824  nodevec.begin(),
1825  nodevec.end(),
1826  std::inserter(common_nodes, common_nodes.end()));
1827  if (common_nodes.size() != nodes_that_must_be_on_side.size())
1828  {
1829  delete side;
1830  break;
1831  }
1832 
1833  FEBase * fe_elem = _fes[_tid][elem->dim()];
1834  FEBase * fe_side = _fes[_tid][side->dim()];
1835 
1836  // Optionally check to see whether face is reasonable candidate based on an
1837  // estimate of how closely it is likely to project to the face
1838  if (check_whether_reasonable)
1839  if (!isFaceReasonableCandidate(elem, side, fe_side, secondary_node, _tangential_tolerance))
1840  {
1841  delete side;
1842  break;
1843  }
1844 
1845  Point contact_phys;
1846  Point contact_ref;
1847  Point contact_on_face_ref;
1848  Real distance = 0.;
1849  Real tangential_distance = 0.;
1850  RealGradient normal;
1851  bool contact_point_on_side;
1852  std::vector<const Node *> off_edge_nodes;
1853  std::vector<std::vector<Real>> side_phi;
1854  std::vector<std::vector<RealGradient>> side_grad_phi;
1855  std::vector<RealGradient> dxyzdxi;
1856  std::vector<RealGradient> dxyzdeta;
1857  std::vector<RealGradient> d2xyzdxideta;
1858 
1859  std::unique_ptr<PenetrationInfo> pen_info =
1860  std::make_unique<PenetrationInfo>(elem,
1861  side,
1862  s,
1863  normal,
1864  distance,
1865  tangential_distance,
1866  contact_phys,
1867  contact_ref,
1868  contact_on_face_ref,
1869  off_edge_nodes,
1870  side_phi,
1871  side_grad_phi,
1872  dxyzdxi,
1873  dxyzdeta,
1874  d2xyzdxideta);
1875 
1876  bool search_succeeded = false;
1877  Moose::findContactPoint(*pen_info,
1878  fe_elem,
1879  fe_side,
1880  _fe_type,
1881  *secondary_node,
1882  true,
1884  contact_point_on_side,
1885  search_succeeded);
1886 
1887  // Do not add contact info from failed searches
1888  if (search_succeeded)
1889  {
1890  thisElemInfo.push_back(pen_info.get());
1891  p_info.push_back(pen_info.release());
1892  }
1893  }
1894 }
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:3488
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 805 of file PenetrationThread.C.

Referenced by operator()().

813 {
814  tangential_distance = 0.0;
815  closest_node = NULL;
816  PenetrationInfo * pi1 = p_info[index1];
817  PenetrationInfo * pi2 = p_info[index2];
818  const unsigned sidedim(pi1->_side->dim());
819  mooseAssert(sidedim == pi2->_side->dim(), "Incompatible dimensionalities");
820 
821  // Nodes on faces for the two interactions
822  std::vector<const Node *> side1_nodes;
823  getSideCornerNodes(pi1->_side, side1_nodes);
824  std::vector<const Node *> side2_nodes;
825  getSideCornerNodes(pi2->_side, side2_nodes);
826 
827  std::sort(side1_nodes.begin(), side1_nodes.end());
828  std::sort(side2_nodes.begin(), side2_nodes.end());
829 
830  // Find nodes shared by the two faces
831  std::vector<const Node *> common_nodes;
832  std::set_intersection(side1_nodes.begin(),
833  side1_nodes.end(),
834  side2_nodes.begin(),
835  side2_nodes.end(),
836  std::inserter(common_nodes, common_nodes.end()));
837 
838  if (common_nodes.size() != sidedim)
839  return false;
840 
841  bool found_point1, found_point2;
842  Point closest_coor_ref1(pi1->_closest_point_ref);
843  const Node * closest_node1;
844  found_point1 = restrictPointToSpecifiedEdgeOfFace(
845  closest_coor_ref1, closest_node1, pi1->_side, common_nodes);
846 
847  Point closest_coor_ref2(pi2->_closest_point_ref);
848  const Node * closest_node2;
849  found_point2 = restrictPointToSpecifiedEdgeOfFace(
850  closest_coor_ref2, closest_node2, pi2->_side, common_nodes);
851 
852  if (!found_point1 || !found_point2)
853  return false;
854 
855  // if (sidedim == 2)
856  // {
857  // TODO:
858  // We have the parametric coordinates of the closest intersection point for both faces.
859  // We need to find a point somewhere in the middle of them so there's not an abrupt jump.
860  // Find that point by taking dot products of vector from contact point to secondary node point
861  // with face normal vectors to see which face we're closer to.
862  // }
863 
864  FEBase * fe = NULL;
865  std::vector<Point> points(1);
866 
867  // We have to pick one of the two faces to own the contact point. It doesn't really matter
868  // which one we pick. For repeatibility, pick the face with the lowest index.
869  if (index1 < index2)
870  {
871  fe = _fes[_tid][pi1->_side->dim()];
872  contact_point_ref = closest_coor_ref1;
873  points[0] = closest_coor_ref1;
874  fe->reinit(pi1->_side, &points);
875  index = index1;
876  }
877  else
878  {
879  fe = _fes[_tid][pi2->_side->dim()];
880  contact_point_ref = closest_coor_ref2;
881  points[0] = closest_coor_ref2;
882  fe->reinit(pi2->_side, &points);
883  index = index2;
884  }
885 
886  contact_point = fe->get_xyz()[0];
887 
888  if (sidedim == 2)
889  {
890  if (closest_node1) // point is off the ridge between the two elements
891  {
892  mooseAssert((closest_node1 == closest_node2 || closest_node2 == NULL),
893  "If off edge of ridge, closest node must be the same on both elements");
894  closest_node = closest_node1;
895 
896  RealGradient off_face = *closest_node1 - contact_point;
897  tangential_distance = off_face.norm();
898  }
899  }
900 
901  return true;
902 }
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 1771 of file PenetrationThread.C.

Referenced by getInfoForFacesWithCommonNodes().

1774 {
1775  for (const auto & pi : p_info)
1776  {
1777  if (!pi)
1778  continue;
1779 
1780  if (pi->_elem == elem)
1781  thisElemInfo.push_back(pi);
1782  }
1783 }
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 1681 of file PenetrationThread.C.

Referenced by getSmoothingFacesAndWeights().

1687 {
1688  // elems connected to a node on this edge, find one that has the same corners as this, and is not
1689  // the current elem
1690  auto node_to_elem_pair = _node_to_elem_map.find(edge_nodes[0]->id()); // just need one of the
1691  // nodes
1692  mooseAssert(node_to_elem_pair != _node_to_elem_map.end(), "Missing entry in node to elem map");
1693  const std::vector<dof_id_type> & elems_connected_to_node = node_to_elem_pair->second;
1694 
1695  std::vector<const Elem *> elems_connected_to_edge;
1696 
1697  for (unsigned int ecni = 0; ecni < elems_connected_to_node.size(); ecni++)
1698  {
1699  if (elems_to_exclude.find(elems_connected_to_node[ecni]) != elems_to_exclude.end())
1700  continue;
1701  const Elem * elem = _mesh.elemPtr(elems_connected_to_node[ecni]);
1702 
1703  std::vector<const Node *> nodevec;
1704  for (unsigned int ni = 0; ni < elem->n_nodes(); ++ni)
1705  if (elem->is_vertex(ni))
1706  nodevec.push_back(elem->node_ptr(ni));
1707 
1708  std::vector<const Node *> common_nodes;
1709  std::sort(nodevec.begin(), nodevec.end());
1710  std::set_intersection(edge_nodes.begin(),
1711  edge_nodes.end(),
1712  nodevec.begin(),
1713  nodevec.end(),
1714  std::inserter(common_nodes, common_nodes.end()));
1715 
1716  if (common_nodes.size() == edge_nodes.size())
1717  elems_connected_to_edge.push_back(elem);
1718  }
1719 
1720  if (elems_connected_to_edge.size() > 0)
1721  {
1722 
1723  // There are potentially multiple elements that share a common edge
1724  // 2D:
1725  // There can only be one element on the same surface
1726  // 3D:
1727  // If there are two edge nodes, there can only be one element on the same surface
1728  // If there is only one edge node (a corner), there could be multiple elements on the same
1729  // surface
1730  bool allowMultipleNeighbors = false;
1731 
1732  if (elems_connected_to_edge[0]->dim() == 3)
1733  {
1734  if (edge_nodes.size() == 1)
1735  {
1736  allowMultipleNeighbors = true;
1737  }
1738  }
1739 
1740  for (unsigned int i = 0; i < elems_connected_to_edge.size(); ++i)
1741  {
1742  std::vector<PenetrationInfo *> thisElemInfo;
1743  getInfoForElem(thisElemInfo, p_info, elems_connected_to_edge[i]);
1744  if (thisElemInfo.size() > 0 && !allowMultipleNeighbors)
1745  {
1746  if (thisElemInfo.size() > 1)
1747  mooseError(
1748  "Found multiple neighbors to current edge/face on surface when only one is allowed");
1749  face_info_comm_edge.push_back(thisElemInfo[0]);
1750  break;
1751  }
1752 
1754  thisElemInfo, p_info, secondary_node, elems_connected_to_edge[i], edge_nodes);
1755  if (thisElemInfo.size() > 0 && !allowMultipleNeighbors)
1756  {
1757  if (thisElemInfo.size() > 1)
1758  mooseError(
1759  "Found multiple neighbors to current edge/face on surface when only one is allowed");
1760  face_info_comm_edge.push_back(thisElemInfo[0]);
1761  break;
1762  }
1763 
1764  for (unsigned int j = 0; j < thisElemInfo.size(); ++j)
1765  face_info_comm_edge.push_back(thisElemInfo[j]);
1766  }
1767  }
1768 }
virtual Elem * elemPtr(const dof_id_type i)
Definition: MooseMesh.C:3153
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:159
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 905 of file PenetrationThread.C.

Referenced by findRidgeContactPoint().

906 {
907  const ElemType t(side->type());
908  corner_nodes.clear();
909 
910  corner_nodes.push_back(side->node_ptr(0));
911  corner_nodes.push_back(side->node_ptr(1));
912  switch (t)
913  {
914  case EDGE2:
915  case EDGE3:
916  case EDGE4:
917  {
918  break;
919  }
920 
921  case TRI3:
922  case TRI6:
923  case TRI7:
924  {
925  corner_nodes.push_back(side->node_ptr(2));
926  break;
927  }
928 
929  case QUAD4:
930  case QUAD8:
931  case QUAD9:
932  {
933  corner_nodes.push_back(side->node_ptr(2));
934  corner_nodes.push_back(side->node_ptr(3));
935  break;
936  }
937 
938  default:
939  {
940  mooseError("Unsupported face type: ", t);
941  break;
942  }
943  }
944 }
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 1540 of file PenetrationThread.C.

Referenced by getSmoothingFacesAndWeights().

1545 {
1546  const ElemType t(side->type());
1547  const Real & xi = p(0);
1548  const Real & eta = p(1);
1549 
1550  Real smooth_limit = 1.0 - _normal_smoothing_distance;
1551 
1552  switch (t)
1553  {
1554  case EDGE2:
1555  case EDGE3:
1556  case EDGE4:
1557  {
1558  if (xi < -smooth_limit)
1559  {
1560  std::vector<const Node *> en;
1561  en.push_back(side->node_ptr(0));
1562  edge_nodes.push_back(en);
1563  Real fw = 0.5 - (1.0 + xi) / (2.0 * _normal_smoothing_distance);
1564  if (fw > 0.5)
1565  fw = 0.5;
1566  edge_face_weights.push_back(fw);
1567  }
1568  else if (xi > smooth_limit)
1569  {
1570  std::vector<const Node *> en;
1571  en.push_back(side->node_ptr(1));
1572  edge_nodes.push_back(en);
1573  Real fw = 0.5 - (1.0 - xi) / (2.0 * _normal_smoothing_distance);
1574  if (fw > 0.5)
1575  fw = 0.5;
1576  edge_face_weights.push_back(fw);
1577  }
1578  break;
1579  }
1580 
1581  case TRI3:
1582  case TRI6:
1583  case TRI7:
1584  {
1585  if (eta < -smooth_limit)
1586  {
1587  std::vector<const Node *> en;
1588  en.push_back(side->node_ptr(0));
1589  en.push_back(side->node_ptr(1));
1590  edge_nodes.push_back(en);
1591  Real fw = 0.5 - (1.0 + eta) / (2.0 * _normal_smoothing_distance);
1592  if (fw > 0.5)
1593  fw = 0.5;
1594  edge_face_weights.push_back(fw);
1595  }
1596  if ((xi + eta) > smooth_limit)
1597  {
1598  std::vector<const Node *> en;
1599  en.push_back(side->node_ptr(1));
1600  en.push_back(side->node_ptr(2));
1601  edge_nodes.push_back(en);
1602  Real fw = 0.5 - (1.0 - xi - eta) / (2.0 * _normal_smoothing_distance);
1603  if (fw > 0.5)
1604  fw = 0.5;
1605  edge_face_weights.push_back(fw);
1606  }
1607  if (xi < -smooth_limit)
1608  {
1609  std::vector<const Node *> en;
1610  en.push_back(side->node_ptr(2));
1611  en.push_back(side->node_ptr(0));
1612  edge_nodes.push_back(en);
1613  Real fw = 0.5 - (1.0 + xi) / (2.0 * _normal_smoothing_distance);
1614  if (fw > 0.5)
1615  fw = 0.5;
1616  edge_face_weights.push_back(fw);
1617  }
1618  break;
1619  }
1620 
1621  case QUAD4:
1622  case QUAD8:
1623  case QUAD9:
1624  {
1625  if (eta < -smooth_limit)
1626  {
1627  std::vector<const Node *> en;
1628  en.push_back(side->node_ptr(0));
1629  en.push_back(side->node_ptr(1));
1630  edge_nodes.push_back(en);
1631  Real fw = 0.5 - (1.0 + eta) / (2.0 * _normal_smoothing_distance);
1632  if (fw > 0.5)
1633  fw = 0.5;
1634  edge_face_weights.push_back(fw);
1635  }
1636  if (xi > smooth_limit)
1637  {
1638  std::vector<const Node *> en;
1639  en.push_back(side->node_ptr(1));
1640  en.push_back(side->node_ptr(2));
1641  edge_nodes.push_back(en);
1642  Real fw = 0.5 - (1.0 - xi) / (2.0 * _normal_smoothing_distance);
1643  if (fw > 0.5)
1644  fw = 0.5;
1645  edge_face_weights.push_back(fw);
1646  }
1647  if (eta > smooth_limit)
1648  {
1649  std::vector<const Node *> en;
1650  en.push_back(side->node_ptr(2));
1651  en.push_back(side->node_ptr(3));
1652  edge_nodes.push_back(en);
1653  Real fw = 0.5 - (1.0 - eta) / (2.0 * _normal_smoothing_distance);
1654  if (fw > 0.5)
1655  fw = 0.5;
1656  edge_face_weights.push_back(fw);
1657  }
1658  if (xi < -smooth_limit)
1659  {
1660  std::vector<const Node *> en;
1661  en.push_back(side->node_ptr(3));
1662  en.push_back(side->node_ptr(0));
1663  edge_nodes.push_back(en);
1664  Real fw = 0.5 - (1.0 + xi) / (2.0 * _normal_smoothing_distance);
1665  if (fw > 0.5)
1666  fw = 0.5;
1667  edge_face_weights.push_back(fw);
1668  }
1669  break;
1670  }
1671 
1672  default:
1673  {
1674  mooseError("Unsupported face type: ", t);
1675  break;
1676  }
1677  }
1678 }
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 1445 of file PenetrationThread.C.

Referenced by smoothNormal().

1450 {
1451  const Elem * side = info->_side;
1452  const Point & p = info->_closest_point_ref;
1453  std::set<dof_id_type> elems_to_exclude;
1454  elems_to_exclude.insert(info->_elem->id());
1455 
1456  std::vector<std::vector<const Node *>> edge_nodes;
1457 
1458  // Get the pairs of nodes along every edge that we are close enough to smooth with
1459  getSmoothingEdgeNodesAndWeights(p, side, edge_nodes, edge_face_weights);
1460  std::vector<Elem *> edge_neighbor_elems;
1461  edge_face_info.resize(edge_nodes.size(), NULL);
1462 
1463  std::vector<unsigned int> edges_without_neighbors;
1464 
1465  for (unsigned int i = 0; i < edge_nodes.size(); ++i)
1466  {
1467  // Sort all sets of edge nodes (needed for comparing edges)
1468  std::sort(edge_nodes[i].begin(), edge_nodes[i].end());
1469 
1470  std::vector<PenetrationInfo *> face_info_comm_edge;
1472  &secondary_node, elems_to_exclude, edge_nodes[i], face_info_comm_edge, p_info);
1473 
1474  if (face_info_comm_edge.size() == 0)
1475  edges_without_neighbors.push_back(i);
1476  else if (face_info_comm_edge.size() > 1)
1477  mooseError("Only one neighbor allowed per edge");
1478  else
1479  edge_face_info[i] = face_info_comm_edge[0];
1480  }
1481 
1482  // Remove edges without neighbors from the vector, starting from end
1483  std::vector<unsigned int>::reverse_iterator rit;
1484  for (rit = edges_without_neighbors.rbegin(); rit != edges_without_neighbors.rend(); ++rit)
1485  {
1486  unsigned int index = *rit;
1487  edge_nodes.erase(edge_nodes.begin() + index);
1488  edge_face_weights.erase(edge_face_weights.begin() + index);
1489  edge_face_info.erase(edge_face_info.begin() + index);
1490  }
1491 
1492  // Handle corner case
1493  if (edge_nodes.size() > 1)
1494  {
1495  if (edge_nodes.size() != 2)
1496  mooseError("Invalid number of smoothing edges");
1497 
1498  // find common node
1499  std::vector<const Node *> common_nodes;
1500  std::set_intersection(edge_nodes[0].begin(),
1501  edge_nodes[0].end(),
1502  edge_nodes[1].begin(),
1503  edge_nodes[1].end(),
1504  std::inserter(common_nodes, common_nodes.end()));
1505 
1506  if (common_nodes.size() != 1)
1507  mooseError("Invalid number of common nodes");
1508 
1509  for (const auto & pinfo : edge_face_info)
1510  elems_to_exclude.insert(pinfo->_elem->id());
1511 
1512  std::vector<PenetrationInfo *> face_info_comm_edge;
1514  &secondary_node, elems_to_exclude, common_nodes, face_info_comm_edge, p_info);
1515 
1516  unsigned int num_corner_neighbors = face_info_comm_edge.size();
1517 
1518  if (num_corner_neighbors > 0)
1519  {
1520  Real fw0 = edge_face_weights[0];
1521  Real fw1 = edge_face_weights[1];
1522 
1523  // Corner weight is product of edge weights. Spread out over multiple neighbors.
1524  Real fw_corner = (fw0 * fw1) / static_cast<Real>(num_corner_neighbors);
1525 
1526  // Adjust original edge weights
1527  edge_face_weights[0] *= (1.0 - fw1);
1528  edge_face_weights[1] *= (1.0 - fw0);
1529 
1530  for (unsigned int i = 0; i < num_corner_neighbors; ++i)
1531  {
1532  edge_face_weights.push_back(fw_corner);
1533  edge_face_info.push_back(face_info_comm_edge[i]);
1534  }
1535  }
1536  }
1537 }
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 751 of file PenetrationThread.C.

Referenced by competeInteractions().

752 {
753  CommonEdgeResult common_edge(NO_COMMON);
754  const std::vector<const Node *> & off_edge_nodes1 = pi1->_off_edge_nodes;
755  const std::vector<const Node *> & off_edge_nodes2 = pi2->_off_edge_nodes;
756  const unsigned dim1 = pi1->_side->dim();
757 
758  if (dim1 == 1)
759  {
760  mooseAssert(pi2->_side->dim() == 1, "Incompatible dimensions.");
761  mooseAssert(off_edge_nodes1.size() < 2 && off_edge_nodes2.size() < 2,
762  "off_edge_nodes size should be <2 for 2D contact");
763  if (off_edge_nodes1.size() == 1 && off_edge_nodes2.size() == 1 &&
764  off_edge_nodes1[0] == off_edge_nodes2[0])
765  common_edge = COMMON_EDGE;
766  }
767  else
768  {
769  mooseAssert(dim1 == 2 && pi2->_side->dim() == 2, "Incompatible dimensions.");
770  mooseAssert(off_edge_nodes1.size() < 3 && off_edge_nodes2.size() < 3,
771  "off_edge_nodes size should be <3 for 3D contact");
772  if (off_edge_nodes1.size() == 1)
773  {
774  if (off_edge_nodes2.size() == 1)
775  {
776  if (off_edge_nodes1[0] == off_edge_nodes2[0])
777  common_edge = COMMON_NODE;
778  }
779  else if (off_edge_nodes2.size() == 2)
780  {
781  if (off_edge_nodes1[0] == off_edge_nodes2[0] || off_edge_nodes1[0] == off_edge_nodes2[1])
782  common_edge = EDGE_AND_COMMON_NODE;
783  }
784  }
785  else if (off_edge_nodes1.size() == 2)
786  {
787  if (off_edge_nodes2.size() == 1)
788  {
789  if (off_edge_nodes1[0] == off_edge_nodes2[0] || off_edge_nodes1[1] == off_edge_nodes2[0])
790  common_edge = EDGE_AND_COMMON_NODE;
791  }
792  else if (off_edge_nodes2.size() == 2)
793  {
794  if ((off_edge_nodes1[0] == off_edge_nodes2[0] &&
795  off_edge_nodes1[1] == off_edge_nodes2[1]) ||
796  (off_edge_nodes1[1] == off_edge_nodes2[0] && off_edge_nodes1[0] == off_edge_nodes2[1]))
797  common_edge = COMMON_EDGE;
798  }
799  }
800  }
801  return common_edge;
802 }
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 1302 of file PenetrationThread.C.

Referenced by createInfoForElem().

1307 {
1308  unsigned int dim = primary_elem->dim();
1309 
1310  const std::vector<Point> & phys_point = fe->get_xyz();
1311 
1312  const std::vector<RealGradient> & dxyz_dxi = fe->get_dxyzdxi();
1313  const std::vector<RealGradient> & dxyz_deta = fe->get_dxyzdeta();
1314 
1315  Point ref_point;
1316 
1317  std::vector<Point> points(1); // Default constructor gives us a point at 0,0,0
1318 
1319  fe->reinit(side, &points);
1320 
1321  RealGradient d = *secondary_point - phys_point[0];
1322 
1323  const Real twosqrt2 = 2.8284; // way more precision than we actually need here
1324  Real max_face_length = side->hmax() + twosqrt2 * tangential_tolerance;
1325 
1326  RealVectorValue normal;
1327  if (dim - 1 == 2)
1328  {
1329  normal = dxyz_dxi[0].cross(dxyz_deta[0]);
1330  }
1331  else if (dim - 1 == 1)
1332  {
1333  const Node * const * elem_nodes = primary_elem->get_nodes();
1334  const Point in_plane_vector1 = *elem_nodes[1] - *elem_nodes[0];
1335  const Point in_plane_vector2 = *elem_nodes[2] - *elem_nodes[0];
1336 
1337  Point out_of_plane_normal = in_plane_vector1.cross(in_plane_vector2);
1338  out_of_plane_normal /= out_of_plane_normal.norm();
1339 
1340  normal = dxyz_dxi[0].cross(out_of_plane_normal);
1341  }
1342  else
1343  {
1344  return true;
1345  }
1346  normal /= normal.norm();
1347 
1348  const Real dot(d * normal);
1349 
1350  const RealGradient normcomp = dot * normal;
1351  const RealGradient tangcomp = d - normcomp;
1352 
1353  const Real tangdist = tangcomp.norm();
1354 
1355  // Increase the size of the zone that we consider if the vector from the face
1356  // to the node has a larger normal component
1357  const Real faceExpansionFactor = 2.0 * (1.0 + normcomp.norm() / d.norm());
1358 
1359  bool isReasonableCandidate = true;
1360  if (tangdist > faceExpansionFactor * max_face_length)
1361  {
1362  isReasonableCandidate = false;
1363  }
1364  return isReasonableCandidate;
1365 }
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:159
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 617 of file PenetrationThread.C.

618 {
620  other._recheck_secondary_nodes.begin(),
621  other._recheck_secondary_nodes.end());
622 }
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  for (const Elem * elem : candidate_elements)
276  {
277  for (auto s : elem->side_index_range())
278  if (boundary_info.has_boundary_id(elem, s, _primary_boundary))
279  {
280  located_elem_ids.push_back(elem->id());
281  break;
282  }
283  }
284  closest_elems = &located_elem_ids;
285  }
286  else
287  {
288  auto node_to_elem_pair = _node_to_elem_map.find(closest_node->id());
289  mooseAssert(node_to_elem_pair != _node_to_elem_map.end(),
290  "Missing entry in node to elem map");
291  closest_elems = &(node_to_elem_pair->second);
292  }
293 
294  for (const auto & elem_id : *closest_elems)
295  {
296  const Elem * elem = _mesh.elemPtr(elem_id);
297 
298  std::vector<PenetrationInfo *> thisElemInfo;
299 
300  std::vector<const Node *> nodesThatMustBeOnSide;
301  // If we have a disconnected mesh, we might not have *any*
302  // nodes that must be on a side we check; we'll rely on
303  // boundary info to find valid sides, then rely on comparing
304  // closest points from each to find the best.
305  //
306  // If we don't have a disconnected mesh, then for maximum
307  // backwards compatibility we're still using the older ridge
308  // and peak detection code, which depends on us ruling out
309  // sides that don't touch closest_node.
310  if (!_use_point_locator)
311  nodesThatMustBeOnSide.push_back(closest_node);
313  thisElemInfo, p_info, &node, elem, nodesThatMustBeOnSide, _check_whether_reasonable);
314  }
315 
316  if (_use_point_locator)
317  {
318  Real min_distance_sq = std::numeric_limits<Real>::max();
319  Point best_point;
320  unsigned int best_i = invalid_uint;
321 
322  // Find closest point in all p_info to the node of interest
323  for (unsigned int i = 0; i < p_info.size(); ++i)
324  {
325  const Point closest_point = closest_point_to_side(node, *p_info[i]->_side);
326  const Real distance_sq = (closest_point - node).norm_sq();
327  if (distance_sq < min_distance_sq)
328  {
329  min_distance_sq = distance_sq;
330  best_point = closest_point;
331  best_i = i;
332  }
333  }
334 
335  p_info[best_i]->_closest_point = best_point;
336  p_info[best_i]->_distance =
337  (p_info[best_i]->_distance >= 0.0 ? 1.0 : -1.0) * std::sqrt(min_distance_sq);
339  mooseError("Normal smoothing not implemented with point locator code");
340  Point normal = (best_point - node).unit();
341  const Real dot = normal * p_info[best_i]->_normal;
342  if (dot < 0)
343  normal *= -1;
344  p_info[best_i]->_normal = normal;
345 
346  switchInfo(info, p_info[best_i]);
347  info_set = true;
348  }
349  else
350  {
351  if (p_info.size() == 1)
352  {
353  if (p_info[0]->_tangential_distance <= _tangential_tolerance)
354  {
355  switchInfo(info, p_info[0]);
356  info_set = true;
357  }
358  }
359  else if (p_info.size() > 1)
360  {
361  // Loop through all pairs of faces, and check for contact on ridge between each face pair
362  std::vector<RidgeData> ridgeDataVec;
363  for (unsigned int i = 0; i + 1 < p_info.size(); ++i)
364  for (unsigned int j = i + 1; j < p_info.size(); ++j)
365  {
366  Point closest_coor;
367  Real tangential_distance(0.0);
368  const Node * closest_node_on_ridge = NULL;
369  unsigned int index = 0;
370  Point closest_coor_ref;
371  bool found_ridge_contact_point = findRidgeContactPoint(closest_coor,
372  tangential_distance,
373  closest_node_on_ridge,
374  index,
375  closest_coor_ref,
376  p_info,
377  i,
378  j);
379  if (found_ridge_contact_point)
380  {
381  RidgeData rpd;
382  rpd._closest_coor = closest_coor;
383  rpd._tangential_distance = tangential_distance;
384  rpd._closest_node = closest_node_on_ridge;
385  rpd._index = index;
386  rpd._closest_coor_ref = closest_coor_ref;
387  ridgeDataVec.push_back(rpd);
388  }
389  }
390 
391  if (ridgeDataVec.size() > 0) // Either find the ridge pair that is the best or find a peak
392  {
393  // Group together ridges for which we are off the edge of a common node.
394  // Those are peaks.
395  std::vector<RidgeSetData> ridgeSetDataVec;
396  for (unsigned int i = 0; i < ridgeDataVec.size(); ++i)
397  {
398  bool foundSetWithMatchingNode = false;
399  for (unsigned int j = 0; j < ridgeSetDataVec.size(); ++j)
400  {
401  if (ridgeDataVec[i]._closest_node != NULL &&
402  ridgeDataVec[i]._closest_node == ridgeSetDataVec[j]._closest_node)
403  {
404  foundSetWithMatchingNode = true;
405  ridgeSetDataVec[j]._ridge_data_vec.push_back(ridgeDataVec[i]);
406  break;
407  }
408  }
409  if (!foundSetWithMatchingNode)
410  {
411  RidgeSetData rsd;
412  rsd._distance = std::numeric_limits<Real>::max();
413  rsd._ridge_data_vec.push_back(ridgeDataVec[i]);
414  rsd._closest_node = ridgeDataVec[i]._closest_node;
415  ridgeSetDataVec.push_back(rsd);
416  }
417  }
418  // Compute distance to each set of ridges
419  for (unsigned int i = 0; i < ridgeSetDataVec.size(); ++i)
420  {
421  if (ridgeSetDataVec[i]._closest_node !=
422  NULL) // Either a peak or off the edge of single ridge
423  {
424  if (ridgeSetDataVec[i]._ridge_data_vec.size() == 1) // off edge of single ridge
425  {
426  if (ridgeSetDataVec[i]._ridge_data_vec[0]._tangential_distance <=
427  _tangential_tolerance) // off within tolerance
428  {
429  ridgeSetDataVec[i]._closest_coor =
430  ridgeSetDataVec[i]._ridge_data_vec[0]._closest_coor;
431  Point contact_point_vec = node - ridgeSetDataVec[i]._closest_coor;
432  ridgeSetDataVec[i]._distance = contact_point_vec.norm();
433  }
434  }
435  else // several ridges join at common node to make peak. The common node is the
436  // contact point
437  {
438  ridgeSetDataVec[i]._closest_coor = *ridgeSetDataVec[i]._closest_node;
439  Point contact_point_vec = node - ridgeSetDataVec[i]._closest_coor;
440  ridgeSetDataVec[i]._distance = contact_point_vec.norm();
441  }
442  }
443  else // on a single ridge
444  {
445  ridgeSetDataVec[i]._closest_coor =
446  ridgeSetDataVec[i]._ridge_data_vec[0]._closest_coor;
447  Point contact_point_vec = node - ridgeSetDataVec[i]._closest_coor;
448  ridgeSetDataVec[i]._distance = contact_point_vec.norm();
449  }
450  }
451  // Find the set of ridges closest to us.
452  unsigned int closest_ridge_set_index(0);
453  Real closest_distance(ridgeSetDataVec[0]._distance);
454  Point closest_point(ridgeSetDataVec[0]._closest_coor);
455  for (unsigned int i = 1; i < ridgeSetDataVec.size(); ++i)
456  {
457  if (ridgeSetDataVec[i]._distance < closest_distance)
458  {
459  closest_ridge_set_index = i;
460  closest_distance = ridgeSetDataVec[i]._distance;
461  closest_point = ridgeSetDataVec[i]._closest_coor;
462  }
463  }
464 
465  if (closest_distance <
466  std::numeric_limits<Real>::max()) // contact point is on the closest ridge set
467  {
468  // find the face in the ridge set with the smallest index, assign that one to the
469  // interaction
470  unsigned int face_index(std::numeric_limits<unsigned int>::max());
471  for (unsigned int i = 0;
472  i < ridgeSetDataVec[closest_ridge_set_index]._ridge_data_vec.size();
473  ++i)
474  {
475  if (ridgeSetDataVec[closest_ridge_set_index]._ridge_data_vec[i]._index < face_index)
476  face_index = ridgeSetDataVec[closest_ridge_set_index]._ridge_data_vec[i]._index;
477  }
478 
479  mooseAssert(face_index < std::numeric_limits<unsigned int>::max(),
480  "face_index invalid");
481 
482  p_info[face_index]->_closest_point = closest_point;
483  p_info[face_index]->_distance =
484  (p_info[face_index]->_distance >= 0.0 ? 1.0 : -1.0) * closest_distance;
485  // Calculate the normal as the vector from the ridge to the point only if we're not
486  // doing normal
487  // smoothing. Normal smoothing will average out the normals on its own.
489  {
490  Point normal(closest_point - node);
491  const Real len(normal.norm());
492  if (len > 0)
493  {
494  normal /= len;
495  }
496  const Real dot(normal * p_info[face_index]->_normal);
497  if (dot < 0)
498  {
499  normal *= -1;
500  }
501  p_info[face_index]->_normal = normal;
502  }
503  p_info[face_index]->_tangential_distance = 0.0;
504 
505  Point closest_point_ref;
506  if (ridgeSetDataVec[closest_ridge_set_index]._ridge_data_vec.size() ==
507  1) // contact with a single ridge rather than a peak
508  {
509  p_info[face_index]->_tangential_distance = ridgeSetDataVec[closest_ridge_set_index]
510  ._ridge_data_vec[0]
511  ._tangential_distance;
512  p_info[face_index]->_closest_point_ref =
513  ridgeSetDataVec[closest_ridge_set_index]._ridge_data_vec[0]._closest_coor_ref;
514  }
515  else
516  { // peak
517  const Node * closest_node_on_face;
518  bool restricted = restrictPointToFace(p_info[face_index]->_closest_point_ref,
519  closest_node_on_face,
520  p_info[face_index]->_side);
521  if (restricted)
522  {
523  if (closest_node_on_face !=
524  ridgeSetDataVec[closest_ridge_set_index]._closest_node)
525  {
526  mooseError("Closest node when restricting point to face != closest node from "
527  "RidgeSetData");
528  }
529  }
530  }
531 
532  FEBase * fe = _fes[_tid][p_info[face_index]->_side->dim()];
533  std::vector<Point> points(1);
534  points[0] = p_info[face_index]->_closest_point_ref;
535  fe->reinit(p_info[face_index]->_side, &points);
536  p_info[face_index]->_side_phi = fe->get_phi();
537  p_info[face_index]->_side_grad_phi = fe->get_dphi();
538  p_info[face_index]->_dxyzdxi = fe->get_dxyzdxi();
539  p_info[face_index]->_dxyzdeta = fe->get_dxyzdeta();
540  p_info[face_index]->_d2xyzdxideta = fe->get_d2xyzdxideta();
541 
542  switchInfo(info, p_info[face_index]);
543  info_set = true;
544  }
545  else
546  { // todo:remove invalid ridge cases so they don't mess up individual face
547  // competition????
548  }
549  }
550 
551  if (!info_set) // contact wasn't on a ridge -- compete individual interactions
552  {
553  unsigned int best(0), i(1);
554  do
555  {
556  CompeteInteractionResult CIResult = competeInteractions(p_info[best], p_info[i]);
557  if (CIResult == FIRST_WINS)
558  {
559  i++;
560  }
561  else if (CIResult == SECOND_WINS)
562  {
563  best = i;
564  i++;
565  }
566  else if (CIResult == NEITHER_WINS)
567  {
568  best = i + 1;
569  i += 2;
570  }
571  } while (i < p_info.size() && best < p_info.size());
572  if (best < p_info.size())
573  {
574  // Ensure final info is within the tangential tolerance
575  if (p_info[best]->_tangential_distance <= _tangential_tolerance)
576  {
577  switchInfo(info, p_info[best]);
578  info_set = true;
579  }
580  }
581  }
582  }
583  }
584  }
585 
586  if (!info_set)
587  {
588  // If penetration is not detected within the saved patch, it is possible
589  // that the secondary node has moved outside the saved patch. So, the patch
590  // for the secondary nodes saved in _recheck_secondary_nodes has to be updated
591  // and the penetration detection has to be re-run on the updated patch.
592 
593  _recheck_secondary_nodes.push_back(node_id);
594 
595  delete info;
596  info = NULL;
597  }
598  else
599  {
600  smoothNormal(info, p_info, node);
601  FEBase * fe = _fes[_tid][info->_side->dim()];
602  computeSlip(*fe, *info);
603  }
604 
605  for (unsigned int j = 0; j < p_info.size(); ++j)
606  {
607  if (p_info[j])
608  {
609  delete p_info[j];
610  p_info[j] = NULL;
611  }
612  }
613  }
614 }
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:3153
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:3488
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:3773
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 1121 of file PenetrationThread.C.

Referenced by operator()().

1122 {
1123  const ElemType t(side->type());
1124  Real & xi = p(0);
1125  Real & eta = p(1);
1126  closest_node = NULL;
1127 
1128  bool off_of_this_face(false);
1129 
1130  switch (t)
1131  {
1132  case EDGE2:
1133  case EDGE3:
1134  case EDGE4:
1135  {
1136  if (xi < -1.0)
1137  {
1138  xi = -1.0;
1139  off_of_this_face = true;
1140  closest_node = side->node_ptr(0);
1141  }
1142  else if (xi > 1.0)
1143  {
1144  xi = 1.0;
1145  off_of_this_face = true;
1146  closest_node = side->node_ptr(1);
1147  }
1148  break;
1149  }
1150 
1151  case TRI3:
1152  case TRI6:
1153  case TRI7:
1154  {
1155  if (eta < 0.0)
1156  {
1157  eta = 0.0;
1158  off_of_this_face = true;
1159  if (xi < 0.5)
1160  {
1161  closest_node = side->node_ptr(0);
1162  if (xi < 0.0)
1163  xi = 0.0;
1164  }
1165  else
1166  {
1167  closest_node = side->node_ptr(1);
1168  if (xi > 1.0)
1169  xi = 1.0;
1170  }
1171  }
1172  else if ((xi + eta) > 1.0)
1173  {
1174  Real delta = (xi + eta - 1.0) / 2.0;
1175  xi -= delta;
1176  eta -= delta;
1177  off_of_this_face = true;
1178  if (xi > 0.5)
1179  {
1180  closest_node = side->node_ptr(1);
1181  if (xi > 1.0)
1182  {
1183  xi = 1.0;
1184  eta = 0.0;
1185  }
1186  }
1187  else
1188  {
1189  closest_node = side->node_ptr(2);
1190  if (xi < 0.0)
1191  {
1192  xi = 0.0;
1193  eta = 1.0;
1194  }
1195  }
1196  }
1197  else if (xi < 0.0)
1198  {
1199  xi = 0.0;
1200  off_of_this_face = true;
1201  if (eta > 0.5)
1202  {
1203  closest_node = side->node_ptr(2);
1204  if (eta > 1.0)
1205  eta = 1.0;
1206  }
1207  else
1208  {
1209  closest_node = side->node_ptr(0);
1210  if (eta < 0.0)
1211  eta = 0.0;
1212  }
1213  }
1214  break;
1215  }
1216 
1217  case QUAD4:
1218  case QUAD8:
1219  case QUAD9:
1220  {
1221  if (eta < -1.0)
1222  {
1223  eta = -1.0;
1224  off_of_this_face = true;
1225  if (xi < 0.0)
1226  {
1227  closest_node = side->node_ptr(0);
1228  if (xi < -1.0)
1229  xi = -1.0;
1230  }
1231  else
1232  {
1233  closest_node = side->node_ptr(1);
1234  if (xi > 1.0)
1235  xi = 1.0;
1236  }
1237  }
1238  else if (xi > 1.0)
1239  {
1240  xi = 1.0;
1241  off_of_this_face = true;
1242  if (eta < 0.0)
1243  {
1244  closest_node = side->node_ptr(1);
1245  if (eta < -1.0)
1246  eta = -1.0;
1247  }
1248  else
1249  {
1250  closest_node = side->node_ptr(2);
1251  if (eta > 1.0)
1252  eta = 1.0;
1253  }
1254  }
1255  else if (eta > 1.0)
1256  {
1257  eta = 1.0;
1258  off_of_this_face = true;
1259  if (xi < 0.0)
1260  {
1261  closest_node = side->node_ptr(3);
1262  if (xi < -1.0)
1263  xi = -1.0;
1264  }
1265  else
1266  {
1267  closest_node = side->node_ptr(2);
1268  if (xi > 1.0)
1269  xi = 1.0;
1270  }
1271  }
1272  else if (xi < -1.0)
1273  {
1274  xi = -1.0;
1275  off_of_this_face = true;
1276  if (eta < 0.0)
1277  {
1278  closest_node = side->node_ptr(0);
1279  if (eta < -1.0)
1280  eta = -1.0;
1281  }
1282  else
1283  {
1284  closest_node = side->node_ptr(3);
1285  if (eta > 1.0)
1286  eta = 1.0;
1287  }
1288  }
1289  break;
1290  }
1291 
1292  default:
1293  {
1294  mooseError("Unsupported face type: ", t);
1295  break;
1296  }
1297  }
1298  return off_of_this_face;
1299 }
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 947 of file PenetrationThread.C.

Referenced by findRidgeContactPoint().

951 {
952  const ElemType t = side->type();
953  Real & xi = p(0);
954  Real & eta = p(1);
955  closest_node = NULL;
956 
957  std::vector<unsigned int> local_node_indices;
958  for (const auto & edge_node : edge_nodes)
959  {
960  unsigned int local_index = side->get_node_index(edge_node);
961  if (local_index == libMesh::invalid_uint)
962  mooseError("Side does not contain node");
963  local_node_indices.push_back(local_index);
964  }
965  mooseAssert(local_node_indices.size() == side->dim(),
966  "Number of edge nodes must match side dimensionality");
967  std::sort(local_node_indices.begin(), local_node_indices.end());
968 
969  bool off_of_this_edge = false;
970 
971  switch (t)
972  {
973  case EDGE2:
974  case EDGE3:
975  case EDGE4:
976  {
977  if (local_node_indices[0] == 0)
978  {
979  if (xi <= -1.0)
980  {
981  xi = -1.0;
982  off_of_this_edge = true;
983  closest_node = side->node_ptr(0);
984  }
985  }
986  else if (local_node_indices[0] == 1)
987  {
988  if (xi >= 1.0)
989  {
990  xi = 1.0;
991  off_of_this_edge = true;
992  closest_node = side->node_ptr(1);
993  }
994  }
995  else
996  {
997  mooseError("Invalid local node indices");
998  }
999  break;
1000  }
1001 
1002  case TRI3:
1003  case TRI6:
1004  case TRI7:
1005  {
1006  if ((local_node_indices[0] == 0) && (local_node_indices[1] == 1))
1007  {
1008  if (eta <= 0.0)
1009  {
1010  eta = 0.0;
1011  off_of_this_edge = true;
1012  if (xi < 0.0)
1013  closest_node = side->node_ptr(0);
1014  else if (xi > 1.0)
1015  closest_node = side->node_ptr(1);
1016  }
1017  }
1018  else if ((local_node_indices[0] == 1) && (local_node_indices[1] == 2))
1019  {
1020  if ((xi + eta) > 1.0)
1021  {
1022  Real delta = (xi + eta - 1.0) / 2.0;
1023  xi -= delta;
1024  eta -= delta;
1025  off_of_this_edge = true;
1026  if (xi > 1.0)
1027  closest_node = side->node_ptr(1);
1028  else if (xi < 0.0)
1029  closest_node = side->node_ptr(2);
1030  }
1031  }
1032  else if ((local_node_indices[0] == 0) && (local_node_indices[1] == 2))
1033  {
1034  if (xi <= 0.0)
1035  {
1036  xi = 0.0;
1037  off_of_this_edge = true;
1038  if (eta > 1.0)
1039  closest_node = side->node_ptr(2);
1040  else if (eta < 0.0)
1041  closest_node = side->node_ptr(0);
1042  }
1043  }
1044  else
1045  {
1046  mooseError("Invalid local node indices");
1047  }
1048 
1049  break;
1050  }
1051 
1052  case QUAD4:
1053  case QUAD8:
1054  case QUAD9:
1055  {
1056  if ((local_node_indices[0] == 0) && (local_node_indices[1] == 1))
1057  {
1058  if (eta <= -1.0)
1059  {
1060  eta = -1.0;
1061  off_of_this_edge = true;
1062  if (xi < -1.0)
1063  closest_node = side->node_ptr(0);
1064  else if (xi > 1.0)
1065  closest_node = side->node_ptr(1);
1066  }
1067  }
1068  else if ((local_node_indices[0] == 1) && (local_node_indices[1] == 2))
1069  {
1070  if (xi >= 1.0)
1071  {
1072  xi = 1.0;
1073  off_of_this_edge = true;
1074  if (eta < -1.0)
1075  closest_node = side->node_ptr(1);
1076  else if (eta > 1.0)
1077  closest_node = side->node_ptr(2);
1078  }
1079  }
1080  else if ((local_node_indices[0] == 2) && (local_node_indices[1] == 3))
1081  {
1082  if (eta >= 1.0)
1083  {
1084  eta = 1.0;
1085  off_of_this_edge = true;
1086  if (xi < -1.0)
1087  closest_node = side->node_ptr(3);
1088  else if (xi > 1.0)
1089  closest_node = side->node_ptr(2);
1090  }
1091  }
1092  else if ((local_node_indices[0] == 0) && (local_node_indices[1] == 3))
1093  {
1094  if (xi <= -1.0)
1095  {
1096  xi = -1.0;
1097  off_of_this_edge = true;
1098  if (eta < -1.0)
1099  closest_node = side->node_ptr(0);
1100  else if (eta > 1.0)
1101  closest_node = side->node_ptr(3);
1102  }
1103  }
1104  else
1105  {
1106  mooseError("Invalid local node indices");
1107  }
1108  break;
1109  }
1110 
1111  default:
1112  {
1113  mooseError("Unsupported face type: ", t);
1114  break;
1115  }
1116  }
1117  return off_of_this_edge;
1118 }
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 1387 of file PenetrationThread.C.

Referenced by operator()().

1390 {
1392  {
1394  {
1395  // If we are within the smoothing distance of any edges or corners, find the
1396  // corner nodes for those edges/corners, and weights from distance to edge/corner
1397  std::vector<Real> edge_face_weights;
1398  std::vector<PenetrationInfo *> edge_face_info;
1399 
1400  getSmoothingFacesAndWeights(info, edge_face_info, edge_face_weights, p_info, node);
1401 
1402  mooseAssert(edge_face_info.size() == edge_face_weights.size(),
1403  "edge_face_info.size() != edge_face_weights.size()");
1404 
1405  if (edge_face_info.size() > 0)
1406  {
1407  // Smooth the normal using the weighting functions for all participating faces.
1408  RealVectorValue new_normal;
1409  Real this_face_weight = 1.0;
1410 
1411  for (unsigned int efwi = 0; efwi < edge_face_weights.size(); ++efwi)
1412  {
1413  PenetrationInfo * npi = edge_face_info[efwi];
1414  if (npi)
1415  new_normal += npi->_normal * edge_face_weights[efwi];
1416 
1417  this_face_weight -= edge_face_weights[efwi];
1418  }
1419  mooseAssert(this_face_weight >= (0.25 - 1e-8),
1420  "Sum of weights of other faces shouldn't exceed 0.75");
1421  new_normal += info->_normal * this_face_weight;
1422 
1423  const Real len = new_normal.norm();
1424  if (len > 0)
1425  new_normal /= len;
1426 
1427  info->_normal = new_normal;
1428  }
1429  }
1431  {
1432  // params.addParam<VariableName>("var_name","description");
1433  // getParam<VariableName>("var_name")
1434  info->_normal(0) = _nodal_normal_x->getValue(info->_side, info->_side_phi);
1435  info->_normal(1) = _nodal_normal_y->getValue(info->_side, info->_side_phi);
1436  info->_normal(2) = _nodal_normal_z->getValue(info->_side, info->_side_phi);
1437  const Real len(info->_normal.norm());
1438  if (len > 0)
1439  info->_normal /= len;
1440  }
1441  }
1442 }
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 625 of file PenetrationThread.C.

Referenced by operator()().

626 {
627  mooseAssert(infoNew != NULL, "infoNew object is null");
628  if (info)
629  {
630  infoNew->_starting_elem = info->_starting_elem;
631  infoNew->_starting_side_num = info->_starting_side_num;
632  infoNew->_starting_closest_point_ref = info->_starting_closest_point_ref;
633  infoNew->_incremental_slip = info->_incremental_slip;
634  infoNew->_accumulated_slip = info->_accumulated_slip;
635  infoNew->_accumulated_slip_old = info->_accumulated_slip_old;
636  infoNew->_frictional_energy = info->_frictional_energy;
637  infoNew->_frictional_energy_old = info->_frictional_energy_old;
638  infoNew->_contact_force = info->_contact_force;
639  infoNew->_contact_force_old = info->_contact_force_old;
640  infoNew->_lagrange_multiplier = info->_lagrange_multiplier;
641  infoNew->_lagrange_multiplier_slip = info->_lagrange_multiplier_slip;
642  infoNew->_locked_this_step = info->_locked_this_step;
643  infoNew->_stick_locked_this_step = info->_stick_locked_this_step;
644  infoNew->_mech_status = info->_mech_status;
645  infoNew->_mech_status_old = info->_mech_status_old;
646  }
647  else
648  {
649  infoNew->_starting_elem = infoNew->_elem;
650  infoNew->_starting_side_num = infoNew->_side_num;
652  }
653  delete info;
654  info = infoNew;
655  infoNew = NULL; // Set this to NULL so that we don't delete it (now owned by _penetration_info).
656 }
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: