https://mooseframework.inl.gov
Classes | Public Types | Public Member Functions | Public Attributes | Static Public Attributes | Private Member Functions | Private Attributes | Friends | List of all members
AutomaticMortarGeneration Class Reference

This class is a container/interface for the objects involved in automatic generation of mortar spaces. More...

#include <AutomaticMortarGeneration.h>

Inheritance diagram for AutomaticMortarGeneration:
[legend]

Classes

struct  MsmSubdomainStats
 Statistics for one primary-secondary subdomain pair. More...
 

Public Types

using MortarFilterIter = std::unordered_map< dof_id_type, std::set< Elem *, CompareDofObjectsByID > >::const_iterator
 

Public Member Functions

 AutomaticMortarGeneration (MooseApp &app, MeshBase &mesh_in, const std::pair< BoundaryID, BoundaryID > &boundary_key, const std::pair< SubdomainID, SubdomainID > &subdomain_key, bool on_displaced, bool periodic, const bool debug, const bool correct_edge_dropping, const Real minimum_projection_angle, const Mortar3DSubpatchPlane mortar_3d_subpatch_plane, const MortarSegmentTriangulationMode triangulation_mode, const bool triangulate_triangles, const Mortar3DQuadraturePointMapping mortar_3d_qp_mapping=Mortar3DQuadraturePointMapping::NORMAL_PROJECTION)
 Must be constructed with a reference to the Mesh we are generating mortar spaces for. More...
 
void buildNodeToElemMaps ()
 Once the secondary_requested_boundary_ids and primary_requested_boundary_ids containers have been filled in, call this function to build node-to-Elem maps for the lower-dimensional elements. More...
 
void computeNodalGeometry ()
 Computes and stores the nodal normal/tangent vectors in a local data structure instead of using the ExplicitSystem/NumericVector approach. More...
 
void projectSecondaryNodes ()
 Project secondary nodes (find xi^(2) values) to the closest points on the primary surface. More...
 
void projectPrimaryNodes ()
 (Inverse) project primary nodes to the points on the secondary surface where they would have come from (find (xi^(1) values)). More...
 
void buildMortarSegmentMesh ()
 Builds the mortar segment mesh once the secondary and primary node projections have been completed. More...
 
void buildMortarSegmentMesh3d ()
 Builds the mortar segment mesh once the secondary and primary node projections have been completed. More...
 
std::vector< MsmSubdomainStatscomputeMsmStatistics ()
 Computes mortar segment mesh statistics and returns one entry per subdomain pair. More...
 
void msmStatistics ()
 Prints mortar segment mesh statistics to console (calls computeMsmStatistics internally) More...
 
void clear ()
 Clears the mortar segment mesh and accompanying data structures. More...
 
void meshChanged ()
 Invalidates the cached MSM node/element ID starting offset so that the next call to buildMortarSegmentMesh3d() recomputes it via allgather. More...
 
bool onDisplaced () const
 returns whether this object is on the displaced mesh More...
 
std::vector< PointgetNodalNormals (const Elem &secondary_elem) const
 
std::array< MooseUtils::SemidynamicVector< Point, 9 >, 2 > getNodalTangents (const Elem &secondary_elem) const
 Compute the two nodal tangents, which are built on-the-fly. More...
 
std::map< unsigned int, unsigned intgetSecondaryIpToLowerElementMap (const Elem &lower_secondary_elem) const
 Compute on-the-fly mapping from secondary interior parent nodes to lower dimensional nodes. More...
 
std::map< unsigned int, unsigned intgetPrimaryIpToLowerElementMap (const Elem &primary_elem, const Elem &primary_elem_ip, const Elem &lower_secondary_elem) const
 Compute on-the-fly mapping from primary interior parent nodes to its corresponding lower dimensional nodes. More...
 
std::vector< PointgetNormals (const Elem &secondary_elem, const std::vector< Point > &xi1_pts) const
 Compute the normals at given reference points on a secondary element. More...
 
std::vector< PointgetNormals (const Elem &secondary_elem, const std::vector< Real > &oned_xi1_pts) const
 Compute the normals at given reference points on a secondary element. More...
 
const ElemgetSecondaryLowerdElemFromSecondaryElem (dof_id_type secondary_elem_id) const
 Return lower dimensional secondary element given its interior parent. More...
 
void computeInactiveLMNodes ()
 Get list of secondary nodes that don't contribute to interaction with any primary element. More...
 
void computeIncorrectEdgeDroppingInactiveLMNodes ()
 Computes inactive secondary nodes when incorrect edge dropping behavior is enabled (any node touching a partially or fully dropped element is dropped) More...
 
void computeInactiveLMElems ()
 Get list of secondary elems without any corresponding primary elements. More...
 
const std::unordered_map< dof_id_type, std::unordered_set< dof_id_type > > & mortarInterfaceCoupling () const
 
const std::pair< BoundaryID, BoundaryID > & primarySecondaryBoundaryIDPair () const
 
const MeshBasemortarSegmentMesh () const
 
const std::unordered_map< const Elem *, MortarSegmentInfo > & mortarSegmentMeshElemToInfo () const
 
int dim () const
 
Mortar3DQuadraturePointMapping mortar3DQpMapping () const
 Return the 3D mortar quadrature-point mapping method. More...
 
const MortarSegmentReferencePointsmortarSegmentReferencePoints (const Elem &mortar_segment_elem) const
 Return the parent-face reference coordinates for a mortar segment. More...
 
const std::unordered_set< const Node * > & getInactiveLMNodes () const
 
const std::unordered_set< const Elem * > & getInactiveLMElems () const
 
bool incorrectEdgeDropping () const
 
std::vector< MortarFilterItersecondariesToMortarSegments (const Node &node) const
 
const std::unordered_map< dof_id_type, std::set< Elem *, CompareDofObjectsByID > > & secondariesToMortarSegments () const
 
const std::set< SubdomainID > & secondaryIPSubIDs () const
 
const std::set< SubdomainID > & primaryIPSubIDs () const
 
const std::unordered_map< dof_id_type, std::vector< const Elem * > > & nodesToSecondaryElem () const
 
void initOutput ()
 initialize mortar-mesh based output More...
 

Public Attributes

const ConsoleStream _console
 An instance of helper class to write streams to the Console objects. More...
 

Static Public Attributes

static const std::string system_name
 The name of the nodal normals system. More...
 

Private Member Functions

void outputMortarMesh ()
 Write the mortar segment mesh to exodus. More...
 
std::string mortarInterfaceName () const
 
void buildCouplingInformation ()
 build the _mortar_interface_coupling data More...
 
void projectSecondaryNodesSinglePair (SubdomainID lower_dimensional_primary_subdomain_id, SubdomainID lower_dimensional_secondary_subdomain_id)
 Helper function responsible for projecting secondary nodes onto primary elements for a single primary/secondary pair. More...
 
void projectPrimaryNodesSinglePair (SubdomainID lower_dimensional_primary_subdomain_id, SubdomainID lower_dimensional_secondary_subdomain_id)
 Helper function used internally by AutomaticMortarGeneration::project_primary_nodes(). More...
 
void householderOrthogolization (const Point &normal, Point &tangent_one, Point &tangent_two) const
 Householder orthogonalization procedure to obtain proper basis for tangent and binormal vectors. More...
 
bool processAlignedNodes (const Node &secondary_node, const Node &primary_node, const std::vector< const Elem *> *secondary_node_neighbors, const std::vector< const Elem *> *primary_node_neighbors, const VectorValue< Real > &nodal_normal, const Elem &candidate_element, std::set< const Elem *> &rejected_element_candidates)
 Process aligned nodes. More...
 

Private Attributes

MooseApp_app
 The Moose app. More...
 
MeshBase_mesh
 Reference to the mesh stored in equation_systems. More...
 
std::set< BoundaryID_secondary_requested_boundary_ids
 The boundary ids corresponding to all the secondary surfaces. More...
 
std::set< BoundaryID_primary_requested_boundary_ids
 The boundary ids corresponding to all the primary surfaces. More...
 
std::vector< std::pair< BoundaryID, BoundaryID > > _primary_secondary_boundary_id_pairs
 A list of primary/secondary boundary id pairs corresponding to each side of the mortar interface. More...
 
std::unordered_map< dof_id_type, std::vector< const Elem * > > _nodes_to_secondary_elem_map
 Map from nodes to connected lower-dimensional elements on the secondary/primary subdomains. More...
 
std::unordered_map< dof_id_type, std::vector< const Elem * > > _nodes_to_primary_elem_map
 
std::unordered_map< std::pair< const Node *, const Elem * >, std::pair< Real, const Elem * > > _secondary_node_and_elem_to_xi2_primary_elem
 Similar to the map above, but associates a (Secondary Node, Secondary Elem) pair to a (xi^(2), primary Elem) pair. More...
 
std::map< std::tuple< dof_id_type, const Node *, const Elem * >, std::pair< Real, const Elem * > > _primary_node_and_elem_to_xi1_secondary_elem
 Same type of container, but for mapping (Primary Node ID, Primary Node, Primary Elem) -> (xi^(1), Secondary Elem) where they are inverse-projected along the nodal normal direction. More...
 
std::unique_ptr< MeshBase_mortar_segment_mesh
 1D Mesh of mortar segment elements which gets built by the call to build_mortar_segment_mesh(). More...
 
std::unordered_map< const Elem *, MortarSegmentInfo_msm_elem_to_info
 Map between Elems in the mortar segment mesh and their info structs. More...
 
std::unordered_map< const Elem *, unsigned int_lower_elem_to_side_id
 Keeps track of the mapping between lower-dimensional elements and the side_id of the interior_parent which they are. More...
 
std::vector< std::pair< SubdomainID, SubdomainID > > _primary_secondary_subdomain_id_pairs
 A list of primary/secondary subdomain id pairs corresponding to each side of the mortar interface. More...
 
std::set< SubdomainID_secondary_boundary_subdomain_ids
 The secondary/primary lower-dimensional boundary subdomain ids are the secondary/primary boundary ids. More...
 
std::set< SubdomainID_primary_boundary_subdomain_ids
 
std::unordered_map< dof_id_type, std::unordered_set< dof_id_type > > _mortar_interface_coupling
 Used by the AugmentSparsityOnInterface functor to determine whether a given Elem is coupled to any others across the gap, and to explicitly set up the dependence between interior_parent() elements on the secondary side and their lower-dimensional sides which are on the interface. More...
 
std::unordered_map< const Node *, Point_secondary_node_to_nodal_normal
 Container for storing the nodal normal vector associated with each secondary node. More...
 
std::unordered_map< const Node *, std::array< Point, 2 > > _secondary_node_to_hh_nodal_tangents
 Container for storing the nodal tangent/binormal vectors associated with each secondary node (Householder approach). More...
 
std::unordered_map< dof_id_type, const Elem * > _secondary_element_to_secondary_lowerd_element
 Map from full dimensional secondary element id to lower dimensional secondary element. More...
 
std::unordered_set< const Node * > _inactive_local_lm_nodes
 
std::unordered_set< const Elem * > _inactive_local_lm_elems
 List of inactive lagrange multiplier nodes (for elemental variables) More...
 
std::unordered_map< dof_id_type, std::set< Elem *, CompareDofObjectsByID > > _secondary_elems_to_mortar_segments
 We maintain a mapping from lower-dimensional secondary elements in the original mesh to (sets of) elements in mortar_segment_mesh. More...
 
std::set< SubdomainID_secondary_ip_sub_ids
 All the secondary interior parent subdomain IDs associated with the mortar mesh. More...
 
std::set< SubdomainID_primary_ip_sub_ids
 All the primary interior parent subdomain IDs associated with the mortar mesh. More...
 
const bool _debug
 Whether to print debug output. More...
 
const bool _on_displaced
 Whether this object is on the displaced mesh. More...
 
const bool _periodic
 Whether this object will be generating a mortar segment mesh for periodic constraints. More...
 
const bool _distributed
 Whether the mortar segment mesh is distributed. More...
 
Real _newton_tolerance = 1e-12
 Newton solve tolerance for node projections. More...
 
Real _xi_tolerance = 1e-6
 Tolerance for checking projection xi values. More...
 
const bool _correct_edge_dropping
 Flag to enable regressed treatment of edge dropping where all LM DoFs on edge dropping element are strongly set to 0. More...
 
const Real _minimum_projection_angle
 Parameter to control which angle (in degrees) is admissible for the creation of mortar segments. More...
 
const Mortar3DSubpatchPlane _mortar_3d_subpatch_plane
 Method used to define the local projection planes for 3D secondary subpatches. More...
 
const MortarSegmentTriangulationMode _triangulation_mode
 Triangulation mode used for clipped 3D mortar polygons. More...
 
const bool _triangulate_triangles
 Whether already-triangular clipped polygons should still be centroid-subdivided. More...
 
const Mortar3DQuadraturePointMapping _mortar_3d_qp_mapping
 Method used to map 3D mortar segment quadrature points to primary and secondary faces. More...
 
std::unordered_map< const Elem *, MortarSegmentReferencePoints_msm_elem_to_reference_points
 Reference-coordinate data used only by the reference-interpolation mapping mode. More...
 
std::unique_ptr< InputParameters_output_params
 Storage for the input parameters used by the mortar nodal geometry output. More...
 
std::optional< dof_id_type_msm_node_id_start
 Cached per-rank starting ID for 3D MSM nodes/elements. More...
 
std::unordered_set< dof_id_type_projected_secondary_nodes
 Debugging container for printing information about fraction of successful projections for secondary nodes. More...
 
std::unordered_set< dof_id_type_failed_secondary_node_projections
 Secondary nodes that failed to project. More...
 

Friends

class MortarNodalGeometryOutput
 
class AugmentSparsityOnInterface
 

Detailed Description

This class is a container/interface for the objects involved in automatic generation of mortar spaces.

Definition at line 68 of file AutomaticMortarGeneration.h.

Member Typedef Documentation

◆ MortarFilterIter

using AutomaticMortarGeneration::MortarFilterIter = std::unordered_map<dof_id_type, std::set<Elem *, CompareDofObjectsByID> >::const_iterator

Definition at line 358 of file AutomaticMortarGeneration.h.

Constructor & Destructor Documentation

◆ AutomaticMortarGeneration()

AutomaticMortarGeneration::AutomaticMortarGeneration ( MooseApp app,
MeshBase mesh_in,
const std::pair< BoundaryID, BoundaryID > &  boundary_key,
const std::pair< SubdomainID, SubdomainID > &  subdomain_key,
bool  on_displaced,
bool  periodic,
const bool  debug,
const bool  correct_edge_dropping,
const Real  minimum_projection_angle,
const Mortar3DSubpatchPlane  mortar_3d_subpatch_plane,
const MortarSegmentTriangulationMode  triangulation_mode,
const bool  triangulate_triangles,
const Mortar3DQuadraturePointMapping  mortar_3d_qp_mapping = Mortar3DQuadraturePointMapping::NORMAL_PROJECTION 
)

Must be constructed with a reference to the Mesh we are generating mortar spaces for.

Definition at line 389 of file AutomaticMortarGeneration.C.

403  : ConsoleStreamInterface(app),
404  _app(app),
405  _mesh(mesh_in),
406  _debug(debug),
407  _on_displaced(on_displaced),
408  _periodic(periodic),
409  // 3D mortar always builds the mortar segment mesh distributedly (each rank adds only its local
410  // secondary elements). For 2D, we ghost the entire mortar interface when displaced, so
411  // displaced meshes are always replicated; otherwise follow the parent mesh.
413  _correct_edge_dropping(correct_edge_dropping),
414  _minimum_projection_angle(minimum_projection_angle),
415  _mortar_3d_subpatch_plane(mortar_3d_subpatch_plane),
416  _triangulation_mode(triangulation_mode),
417  _triangulate_triangles(triangulate_triangles),
418  _mortar_3d_qp_mapping(mortar_3d_qp_mapping)
419 {
420  _primary_secondary_boundary_id_pairs.push_back(boundary_key);
421  _primary_requested_boundary_ids.insert(boundary_key.first);
422  _secondary_requested_boundary_ids.insert(boundary_key.second);
423  _primary_secondary_subdomain_id_pairs.push_back(subdomain_key);
424  _primary_boundary_subdomain_ids.insert(subdomain_key.first);
425  _secondary_boundary_subdomain_ids.insert(subdomain_key.second);
426 
427  if (_distributed)
429  std::make_unique<DistributedMesh>(_mesh.comm(), _mesh.spatial_dimension());
430  else
432  std::make_unique<ReplicatedMesh>(_mesh.comm(), _mesh.spatial_dimension());
433 }
std::set< SubdomainID > _primary_boundary_subdomain_ids
std::vector< std::pair< BoundaryID, BoundaryID > > _primary_secondary_boundary_id_pairs
A list of primary/secondary boundary id pairs corresponding to each side of the mortar interface...
std::set< SubdomainID > _secondary_boundary_subdomain_ids
The secondary/primary lower-dimensional boundary subdomain ids are the secondary/primary boundary ids...
const bool _periodic
Whether this object will be generating a mortar segment mesh for periodic constraints.
const bool _distributed
Whether the mortar segment mesh is distributed.
const Parallel::Communicator & comm() const
const bool _debug
Whether to print debug output.
const bool _triangulate_triangles
Whether already-triangular clipped polygons should still be centroid-subdivided.
const bool _correct_edge_dropping
Flag to enable regressed treatment of edge dropping where all LM DoFs on edge dropping element are st...
ConsoleStreamInterface(MooseApp &app)
A class for providing a helper stream object for writting message to all the Output objects...
std::set< BoundaryID > _primary_requested_boundary_ids
The boundary ids corresponding to all the primary surfaces.
const MortarSegmentTriangulationMode _triangulation_mode
Triangulation mode used for clipped 3D mortar polygons.
const bool _on_displaced
Whether this object is on the displaced mesh.
unsigned int spatial_dimension() const
const Mortar3DSubpatchPlane _mortar_3d_subpatch_plane
Method used to define the local projection planes for 3D secondary subpatches.
virtual bool is_replicated() const
unsigned int mesh_dimension() const
const Real _minimum_projection_angle
Parameter to control which angle (in degrees) is admissible for the creation of mortar segments...
MeshBase & _mesh
Reference to the mesh stored in equation_systems.
std::vector< std::pair< SubdomainID, SubdomainID > > _primary_secondary_subdomain_id_pairs
A list of primary/secondary subdomain id pairs corresponding to each side of the mortar interface...
std::set< BoundaryID > _secondary_requested_boundary_ids
The boundary ids corresponding to all the secondary surfaces.
const Mortar3DQuadraturePointMapping _mortar_3d_qp_mapping
Method used to map 3D mortar segment quadrature points to primary and secondary faces.
MooseApp & _app
The Moose app.
std::unique_ptr< MeshBase > _mortar_segment_mesh
1D Mesh of mortar segment elements which gets built by the call to build_mortar_segment_mesh().

Member Function Documentation

◆ buildCouplingInformation()

void AutomaticMortarGeneration::buildCouplingInformation ( )
private

build the _mortar_interface_coupling data

Definition at line 1713 of file AutomaticMortarGeneration.C.

Referenced by buildMortarSegmentMesh(), and buildMortarSegmentMesh3d().

1714 {
1715  std::unordered_map<processor_id_type, std::vector<std::pair<dof_id_type, dof_id_type>>>
1716  coupling_info;
1717 
1718  // Loop over the msm_elem_to_info object and build a bi-directional
1719  // multimap from secondary elements to the primary Elems which they are
1720  // coupled to and vice-versa. This is used in the
1721  // AugmentSparsityOnInterface functor to determine whether a given
1722  // secondary Elem is coupled across the mortar interface to a primary
1723  // element.
1724  for (const auto & pr : _msm_elem_to_info)
1725  {
1726  const Elem * secondary_elem = pr.second.secondary_elem;
1727  const Elem * primary_elem = pr.second.primary_elem;
1728 
1729  // LowerSecondary
1730  coupling_info[secondary_elem->processor_id()].emplace_back(
1731  secondary_elem->id(), secondary_elem->interior_parent()->id());
1732  if (secondary_elem->processor_id() != _mesh.processor_id())
1733  // We want to keep information for nonlocal lower-dimensional secondary element point
1734  // neighbors for mortar nodal aux kernels
1735  _mortar_interface_coupling[secondary_elem->id()].insert(
1736  secondary_elem->interior_parent()->id());
1737 
1738  // LowerPrimary
1739  coupling_info[secondary_elem->processor_id()].emplace_back(
1740  secondary_elem->id(), primary_elem->interior_parent()->id());
1741  if (secondary_elem->processor_id() != _mesh.processor_id())
1742  // We want to keep information for nonlocal lower-dimensional secondary element point
1743  // neighbors for mortar nodal aux kernels
1744  _mortar_interface_coupling[secondary_elem->id()].insert(
1745  primary_elem->interior_parent()->id());
1746 
1747  // Lower-LowerDimensionalPrimary
1748  coupling_info[secondary_elem->processor_id()].emplace_back(secondary_elem->id(),
1749  primary_elem->id());
1750  if (secondary_elem->processor_id() != _mesh.processor_id())
1751  // We want to keep information for nonlocal lower-dimensional secondary element point
1752  // neighbors for mortar nodal aux kernels
1753  _mortar_interface_coupling[secondary_elem->id()].insert(primary_elem->id());
1754 
1755  // SecondaryLower
1756  coupling_info[secondary_elem->interior_parent()->processor_id()].emplace_back(
1757  secondary_elem->interior_parent()->id(), secondary_elem->id());
1758 
1759  // SecondaryPrimary
1760  coupling_info[secondary_elem->interior_parent()->processor_id()].emplace_back(
1761  secondary_elem->interior_parent()->id(), primary_elem->interior_parent()->id());
1762 
1763  // PrimaryLower
1764  coupling_info[primary_elem->interior_parent()->processor_id()].emplace_back(
1765  primary_elem->interior_parent()->id(), secondary_elem->id());
1766 
1767  // PrimarySecondary
1768  coupling_info[primary_elem->interior_parent()->processor_id()].emplace_back(
1769  primary_elem->interior_parent()->id(), secondary_elem->interior_parent()->id());
1770  }
1771 
1772  // Push the coupling information
1773  auto action_functor =
1774  [this](processor_id_type,
1775  const std::vector<std::pair<dof_id_type, dof_id_type>> & coupling_info)
1776  {
1777  for (auto [i, j] : coupling_info)
1778  _mortar_interface_coupling[i].insert(j);
1779  };
1780  TIMPI::push_parallel_vector_data(_mesh.comm(), coupling_info, action_functor);
1781 }
const Elem * interior_parent() const
std::unordered_map< const Elem *, MortarSegmentInfo > _msm_elem_to_info
Map between Elems in the mortar segment mesh and their info structs.
const Parallel::Communicator & comm() const
void push_parallel_vector_data(const Communicator &comm, MapToVectors &&data, const ActionFunctor &act_on_data)
uint8_t processor_id_type
dof_id_type id() const
MeshBase & _mesh
Reference to the mesh stored in equation_systems.
processor_id_type processor_id() const
std::unordered_map< dof_id_type, std::unordered_set< dof_id_type > > _mortar_interface_coupling
Used by the AugmentSparsityOnInterface functor to determine whether a given Elem is coupled to any ot...
processor_id_type processor_id() const

◆ buildMortarSegmentMesh()

void AutomaticMortarGeneration::buildMortarSegmentMesh ( )

Builds the mortar segment mesh once the secondary and primary node projections have been completed.

Inputs:

  • mesh
  • primary_node_and_elem_to_xi1_secondary_elem
  • secondary_node_and_elem_to_xi2_primary_elem
  • nodes_to_primary_elem_map

Outputs:

  • mortar_segment_mesh
  • msm_elem_to_info

Defined in the file build_mortar_segment_mesh.C.

This was a change to how inactive LM DoFs are handled. Now mortar segment elements are not used in assembly if there is no corresponding primary element and inactive LM DoFs (those with no contribution to an active primary element) are zeroed.

Definition at line 654 of file AutomaticMortarGeneration.C.

Referenced by MortarInterfaceWarehouse::update().

655 {
656  using std::abs;
657 
658  dof_id_type local_id_index = 0;
659  std::size_t node_unique_id_offset = 0;
660 
661  // Create an offset by the maximum number of mortar segment elements that can be created *plus*
662  // the number of lower-dimensional secondary subdomain elements. Recall that the number of mortar
663  // segments created is a function of node projection, *and* that if we split elems we will delete
664  // that elem which has already taken a unique id
665  for (const auto & pr : _primary_secondary_boundary_id_pairs)
666  {
667  const auto primary_bnd_id = pr.first;
668  const auto secondary_bnd_id = pr.second;
669  const auto num_primary_nodes =
670  std::distance(_mesh.bid_nodes_begin(primary_bnd_id), _mesh.bid_nodes_end(primary_bnd_id));
671  const auto num_secondary_nodes = std::distance(_mesh.bid_nodes_begin(secondary_bnd_id),
672  _mesh.bid_nodes_end(secondary_bnd_id));
673  mooseAssert(num_primary_nodes,
674  "There are no primary nodes on boundary ID "
675  << primary_bnd_id << ". Does that bondary ID even exist on the mesh?");
676  mooseAssert(num_secondary_nodes,
677  "There are no secondary nodes on boundary ID "
678  << secondary_bnd_id << ". Does that bondary ID even exist on the mesh?");
679 
680  node_unique_id_offset += num_primary_nodes + 2 * num_secondary_nodes;
681  }
682 
683  // 1.) Add all lower-dimensional secondary side elements as the "initial" mortar segments.
684  for (MeshBase::const_element_iterator el = _mesh.active_elements_begin(),
685  end_el = _mesh.active_elements_end();
686  el != end_el;
687  ++el)
688  {
689  const Elem * secondary_elem = *el;
690 
691  // If this is not one of the lower-dimensional secondary side elements, go on to the next one.
692  if (!this->_secondary_boundary_subdomain_ids.count(secondary_elem->subdomain_id()))
693  continue;
694 
695  std::vector<Node *> new_nodes;
696  for (MooseIndex(secondary_elem->n_nodes()) n = 0; n < secondary_elem->n_nodes(); ++n)
697  {
698  new_nodes.push_back(_mortar_segment_mesh->add_point(
699  secondary_elem->point(n), secondary_elem->node_id(n), secondary_elem->processor_id()));
700  Node * const new_node = new_nodes.back();
701  new_node->set_unique_id(new_node->id() + node_unique_id_offset);
702  }
703 
704  std::unique_ptr<Elem> new_elem;
705  if (secondary_elem->default_order() == SECOND)
706  new_elem = std::make_unique<Edge3>();
707  else
708  new_elem = std::make_unique<Edge2>();
709 
710  new_elem->processor_id() = secondary_elem->processor_id();
711  new_elem->subdomain_id() = secondary_elem->subdomain_id();
712  new_elem->set_id(local_id_index++);
713  new_elem->set_unique_id(new_elem->id());
714 
715  for (MooseIndex(new_elem->n_nodes()) n = 0; n < new_elem->n_nodes(); ++n)
716  new_elem->set_node(n, new_nodes[n]);
717 
718  Elem * new_elem_ptr = _mortar_segment_mesh->add_elem(new_elem.release());
719 
720  // The xi^(1) values for this mortar segment are initially -1 and 1.
721  MortarSegmentInfo msinfo;
722  msinfo.xi1_a = -1;
723  msinfo.xi1_b = +1;
724  msinfo.secondary_elem = secondary_elem;
725 
726  auto new_container_it0 = _secondary_node_and_elem_to_xi2_primary_elem.find(
727  std::make_pair(secondary_elem->node_ptr(0), secondary_elem)),
728  new_container_it1 = _secondary_node_and_elem_to_xi2_primary_elem.find(
729  std::make_pair(secondary_elem->node_ptr(1), secondary_elem));
730 
731  bool new_container_node0_found =
732  (new_container_it0 != _secondary_node_and_elem_to_xi2_primary_elem.end()),
733  new_container_node1_found =
734  (new_container_it1 != _secondary_node_and_elem_to_xi2_primary_elem.end());
735 
736  const Elem * node0_primary_candidate = nullptr;
737  const Elem * node1_primary_candidate = nullptr;
738 
739  if (new_container_node0_found)
740  {
741  const auto & xi2_primary_elem_pair = new_container_it0->second;
742  msinfo.xi2_a = xi2_primary_elem_pair.first;
743  node0_primary_candidate = xi2_primary_elem_pair.second;
744  }
745 
746  if (new_container_node1_found)
747  {
748  const auto & xi2_primary_elem_pair = new_container_it1->second;
749  msinfo.xi2_b = xi2_primary_elem_pair.first;
750  node1_primary_candidate = xi2_primary_elem_pair.second;
751  }
752 
753  // If both node0 and node1 agree on the primary element they are
754  // projected into, then this mortar segment fits entirely within
755  // a single primary element, and we can go ahead and set the
756  // msinfo.primary_elem pointer now.
757  if (node0_primary_candidate == node1_primary_candidate)
758  msinfo.primary_elem = node0_primary_candidate;
759 
760  // Associate this MSM elem with the MortarSegmentInfo.
761  _msm_elem_to_info.emplace(new_elem_ptr, msinfo);
762 
763  // Maintain the mapping between secondary elems and mortar segment elems contained within them.
764  // Initially, only the original secondary_elem is present.
765  _secondary_elems_to_mortar_segments[secondary_elem->id()].insert(new_elem_ptr);
766  }
767 
768  // 2.) Insert new nodes from primary side and split mortar segments as necessary.
769  for (const auto & pr : _primary_node_and_elem_to_xi1_secondary_elem)
770  {
771  auto key = pr.first;
772  auto val = pr.second;
773 
774  const Node * primary_node = std::get<1>(key);
775  Real xi1 = val.first;
776  const Elem * secondary_elem = val.second;
777 
778  // If this is an aligned node, we don't need to do anything.
779  if (abs(abs(xi1) - 1.) < _xi_tolerance)
780  continue;
781 
782  auto && order = secondary_elem->default_order();
783 
784  // Determine physical location of new point to be inserted.
785  Point new_pt(0);
786  for (MooseIndex(secondary_elem->n_nodes()) n = 0; n < secondary_elem->n_nodes(); ++n)
787  new_pt += Moose::fe_lagrange_1D_shape(order, n, xi1) * secondary_elem->point(n);
788 
789  // Find the current mortar segment that will have to be split.
790  auto & mortar_segment_set = _secondary_elems_to_mortar_segments[secondary_elem->id()];
791  Elem * current_mortar_segment = nullptr;
792  MortarSegmentInfo * info = nullptr;
793 
794  for (const auto & mortar_segment_candidate : mortar_segment_set)
795  {
796  try
797  {
798  info = &_msm_elem_to_info.at(mortar_segment_candidate);
799  }
800  catch (std::out_of_range &)
801  {
802  mooseError("MortarSegmentInfo not found for the mortar segment candidate");
803  }
804  if (info->xi1_a <= xi1 && xi1 <= info->xi1_b)
805  {
806  current_mortar_segment = mortar_segment_candidate;
807  break;
808  }
809  }
810 
811  // Make sure we found one.
812  if (current_mortar_segment == nullptr)
813  mooseError("Unable to find appropriate mortar segment during linear search!");
814 
815  // If node lands on endpoint of segment, don't split.
816  // Jacob: This condition was getting missed by the < comparison a few lines above. To fix it I
817  // just made it <= and put this condition in to handle equality different. It probably could be
818  // done with a tolerance but the the toleranced equality is already handled later when we drop
819  // segments with small volume.
820  if (info->xi1_a == xi1 || xi1 == info->xi1_b)
821  continue;
822 
823  const auto new_id = _mortar_segment_mesh->max_node_id();
824  mooseAssert(_mortar_segment_mesh->comm().verify(new_id),
825  "new_id must be the same on all processes");
826  Node * const new_node =
827  _mortar_segment_mesh->add_point(new_pt, new_id, secondary_elem->processor_id());
828  new_node->set_unique_id(new_id + node_unique_id_offset);
829 
830  // Reconstruct the nodal normal at xi1. This will help us
831  // determine the orientation of the primary elems relative to the
832  // new mortar segments.
833  const Point normal = getNormals(*secondary_elem, std::vector<Real>({xi1}))[0];
834 
835  // Get the set of primary_node neighbors.
836  if (this->_nodes_to_primary_elem_map.find(primary_node->id()) ==
837  this->_nodes_to_primary_elem_map.end())
838  mooseError("We should already have built this primary node to elem pair!");
839  const std::vector<const Elem *> & primary_node_neighbors =
840  this->_nodes_to_primary_elem_map[primary_node->id()];
841 
842  // Sanity check
843  if (primary_node_neighbors.size() == 0 || primary_node_neighbors.size() > 2)
844  mooseError("We must have either 1 or 2 primary side nodal neighbors, but we had ",
845  primary_node_neighbors.size());
846 
847  // Primary Elem pointers which we will eventually assign to the
848  // mortar segments being created. We start by assuming
849  // primary_node_neighbor[0] is on the "left" and
850  // primary_node_neighbor[1]/"nothing" is on the "right" and then
851  // swap them if that's not the case.
852  const Elem * left_primary_elem = primary_node_neighbors[0];
853  const Elem * right_primary_elem =
854  (primary_node_neighbors.size() == 2) ? primary_node_neighbors[1] : nullptr;
855 
857 
858  // Storage for z-component of cross products for determining
859  // orientation.
860  std::array<Real, 2> secondary_node_cps;
861  std::vector<Real> primary_node_cps(primary_node_neighbors.size());
862 
863  // Store z-component of left and right secondary node cross products with the nodal normal.
864  for (unsigned int nid = 0; nid < 2; ++nid)
865  secondary_node_cps[nid] = normal.cross(secondary_elem->point(nid) - new_pt)(2);
866 
867  for (MooseIndex(primary_node_neighbors) mnn = 0; mnn < primary_node_neighbors.size(); ++mnn)
868  {
869  const Elem * primary_neigh = primary_node_neighbors[mnn];
870  Point opposite = (primary_neigh->node_ptr(0) == primary_node) ? primary_neigh->point(1)
871  : primary_neigh->point(0);
872  Point cp = normal.cross(opposite - new_pt);
873  primary_node_cps[mnn] = cp(2);
874  }
875 
876  // We will verify that only 1 orientation is actually valid.
877  bool orientation1_valid = false, orientation2_valid = false;
878 
879  if (primary_node_neighbors.size() == 2)
880  {
881  // 2 primary neighbor case
882  orientation1_valid = (secondary_node_cps[0] * primary_node_cps[0] > 0.) &&
883  (secondary_node_cps[1] * primary_node_cps[1] > 0.);
884 
885  orientation2_valid = (secondary_node_cps[0] * primary_node_cps[1] > 0.) &&
886  (secondary_node_cps[1] * primary_node_cps[0] > 0.);
887  }
888  else if (primary_node_neighbors.size() == 1)
889  {
890  // 1 primary neighbor case
891  orientation1_valid = (secondary_node_cps[0] * primary_node_cps[0] > 0.);
892  orientation2_valid = (secondary_node_cps[1] * primary_node_cps[0] > 0.);
893  }
894  else
895  mooseError("Invalid primary node neighbors size ", primary_node_neighbors.size());
896 
897  // Verify that both orientations are not simultaneously valid/invalid. If they are not, then we
898  // are going to throw an exception instead of erroring out since we can easily reach this point
899  // if we have one bad linear solve. It's better in general to catch the error and then try a
900  // smaller time-step
901  if (orientation1_valid && orientation2_valid)
902  throw MooseException(
903  "AutomaticMortarGeneration: Both orientations cannot simultaneously be valid.");
904 
905  // We are going to treat the case where both orientations are invalid as a case in which we
906  // should not be splitting the mortar mesh to incorporate primary mesh elements.
907  // In practice, this case has appeared for very oblique projections, so we assume these cases
908  // will not be considered in mortar thermomechanical contact.
909  if (!orientation1_valid && !orientation2_valid)
910  {
911  mooseDoOnce(mooseWarning(
912  "AutomaticMortarGeneration: Unable to determine valid secondary-primary orientation. "
913  "Consequently we will consider projection of the primary node invalid and not split the "
914  "mortar segment. "
915  "This situation can indicate there are very oblique projections between primary (mortar) "
916  "and secondary (non-mortar) surfaces for a good problem set up. It can also mean your "
917  "time step is too large. This message is only printed once."));
918  continue;
919  }
920 
921  // Make an Elem on the left
922  std::unique_ptr<Elem> new_elem_left;
923  if (order == SECOND)
924  new_elem_left = std::make_unique<Edge3>();
925  else
926  new_elem_left = std::make_unique<Edge2>();
927 
928  new_elem_left->processor_id() = current_mortar_segment->processor_id();
929  new_elem_left->subdomain_id() = current_mortar_segment->subdomain_id();
930  new_elem_left->set_id(local_id_index++);
931  new_elem_left->set_unique_id(new_elem_left->id());
932  new_elem_left->set_node(0, current_mortar_segment->node_ptr(0));
933  new_elem_left->set_node(1, new_node);
934 
935  // Make an Elem on the right
936  std::unique_ptr<Elem> new_elem_right;
937  if (order == SECOND)
938  new_elem_right = std::make_unique<Edge3>();
939  else
940  new_elem_right = std::make_unique<Edge2>();
941 
942  new_elem_right->processor_id() = current_mortar_segment->processor_id();
943  new_elem_right->subdomain_id() = current_mortar_segment->subdomain_id();
944  new_elem_right->set_id(local_id_index++);
945  new_elem_right->set_unique_id(new_elem_right->id());
946  new_elem_right->set_node(0, new_node);
947  new_elem_right->set_node(1, current_mortar_segment->node_ptr(1));
948 
949  if (order == SECOND)
950  {
951  // left
952  Point left_interior_point(0);
953  Real left_interior_xi = (xi1 + info->xi1_a) / 2;
954 
955  // This is eta for the current mortar segment that we're splitting
956  Real current_left_interior_eta =
957  (2. * left_interior_xi - info->xi1_a - info->xi1_b) / (info->xi1_b - info->xi1_a);
958 
959  for (MooseIndex(current_mortar_segment->n_nodes()) n = 0;
960  n < current_mortar_segment->n_nodes();
961  ++n)
962  left_interior_point += Moose::fe_lagrange_1D_shape(order, n, current_left_interior_eta) *
963  current_mortar_segment->point(n);
964 
965  const auto new_interior_left_id = _mortar_segment_mesh->max_node_id();
966  mooseAssert(_mortar_segment_mesh->comm().verify(new_interior_left_id),
967  "new_id must be the same on all processes");
968  Node * const new_interior_node_left = _mortar_segment_mesh->add_point(
969  left_interior_point, new_interior_left_id, new_elem_left->processor_id());
970  new_elem_left->set_node(2, new_interior_node_left);
971  new_interior_node_left->set_unique_id(new_interior_left_id + node_unique_id_offset);
972 
973  // right
974  Point right_interior_point(0);
975  Real right_interior_xi = (xi1 + info->xi1_b) / 2;
976  // This is eta for the current mortar segment that we're splitting
977  Real current_right_interior_eta =
978  (2. * right_interior_xi - info->xi1_a - info->xi1_b) / (info->xi1_b - info->xi1_a);
979 
980  for (MooseIndex(current_mortar_segment->n_nodes()) n = 0;
981  n < current_mortar_segment->n_nodes();
982  ++n)
983  right_interior_point += Moose::fe_lagrange_1D_shape(order, n, current_right_interior_eta) *
984  current_mortar_segment->point(n);
985 
986  const auto new_interior_id_right = _mortar_segment_mesh->max_node_id();
987  mooseAssert(_mortar_segment_mesh->comm().verify(new_interior_id_right),
988  "new_id must be the same on all processes");
989  Node * const new_interior_node_right = _mortar_segment_mesh->add_point(
990  right_interior_point, new_interior_id_right, new_elem_right->processor_id());
991  new_elem_right->set_node(2, new_interior_node_right);
992  new_interior_node_right->set_unique_id(new_interior_id_right + node_unique_id_offset);
993  }
994 
995  // If orientation 2 was valid, swap the left and right primaries.
996  if (orientation2_valid)
997  std::swap(left_primary_elem, right_primary_elem);
998 
999  // Now that we know left_primary_elem and right_primary_elem, we can determine left_xi2 and
1000  // right_xi2.
1001  if (left_primary_elem)
1002  left_xi2 = (primary_node == left_primary_elem->node_ptr(0)) ? -1 : +1;
1003  if (right_primary_elem)
1004  right_xi2 = (primary_node == right_primary_elem->node_ptr(0)) ? -1 : +1;
1005 
1006  // Grab the MortarSegmentInfo object associated with this
1007  // segment. We can use "at()" here since we want this to fail if
1008  // current_mortar_segment is not found... Since we're going to
1009  // erase this entry from the map momentarily, we make an actual
1010  // copy rather than grabbing a reference.
1011  auto msm_it = _msm_elem_to_info.find(current_mortar_segment);
1012  if (msm_it == _msm_elem_to_info.end())
1013  mooseError("MortarSegmentInfo not found for current_mortar_segment.");
1014  MortarSegmentInfo current_msinfo = msm_it->second;
1015 
1016  // add_left
1017  {
1018  Elem * msm_new_elem = _mortar_segment_mesh->add_elem(new_elem_left.release());
1019 
1020  // Create new MortarSegmentInfo objects for new_elem_left
1021  MortarSegmentInfo new_msinfo_left;
1022 
1023  // The new MortarSegmentInfo info objects inherit their "outer"
1024  // information from current_msinfo and the rest is determined by
1025  // the Node being inserted.
1026  new_msinfo_left.xi1_a = current_msinfo.xi1_a;
1027  new_msinfo_left.xi2_a = current_msinfo.xi2_a;
1028  new_msinfo_left.secondary_elem = secondary_elem;
1029  new_msinfo_left.xi1_b = xi1;
1030  new_msinfo_left.xi2_b = left_xi2;
1031  new_msinfo_left.primary_elem = left_primary_elem;
1032 
1033  // Add new msinfo objects to the map.
1034  _msm_elem_to_info.emplace(msm_new_elem, new_msinfo_left);
1035 
1036  // We need to insert new_elem_left in
1037  // the mortar_segment_set for this secondary_elem.
1038  mortar_segment_set.insert(msm_new_elem);
1039  }
1040 
1041  // add_right
1042  {
1043  Elem * msm_new_elem = _mortar_segment_mesh->add_elem(new_elem_right.release());
1044 
1045  // Create new MortarSegmentInfo objects for new_elem_right
1046  MortarSegmentInfo new_msinfo_right;
1047 
1048  new_msinfo_right.xi1_b = current_msinfo.xi1_b;
1049  new_msinfo_right.xi2_b = current_msinfo.xi2_b;
1050  new_msinfo_right.secondary_elem = secondary_elem;
1051  new_msinfo_right.xi1_a = xi1;
1052  new_msinfo_right.xi2_a = right_xi2;
1053  new_msinfo_right.primary_elem = right_primary_elem;
1054 
1055  _msm_elem_to_info.emplace(msm_new_elem, new_msinfo_right);
1056 
1057  mortar_segment_set.insert(msm_new_elem);
1058  }
1059 
1060  // Erase the MortarSegmentInfo object for current_mortar_segment from the map.
1061  _msm_elem_to_info.erase(msm_it);
1062 
1063  // current_mortar_segment must be erased from the
1064  // mortar_segment_set since it has now been split.
1065  mortar_segment_set.erase(current_mortar_segment);
1066 
1067  // The original mortar segment has been split, so erase it from
1068  // the mortar segment mesh.
1069  _mortar_segment_mesh->delete_elem(current_mortar_segment);
1070  }
1071 
1072  // Remove all MSM elements without a primary contribution
1078  for (auto msm_elem : _mortar_segment_mesh->active_element_ptr_range())
1079  {
1080  MortarSegmentInfo & msinfo = libmesh_map_find(_msm_elem_to_info, msm_elem);
1081  Elem * primary_elem = const_cast<Elem *>(msinfo.primary_elem);
1082  if (primary_elem == nullptr || abs(msinfo.xi2_a) > 1.0 + TOLERANCE ||
1083  abs(msinfo.xi2_b) > 1.0 + TOLERANCE)
1084  {
1085  // Erase from secondary to msms map
1086  auto it = _secondary_elems_to_mortar_segments.find(msinfo.secondary_elem->id());
1087  mooseAssert(it != _secondary_elems_to_mortar_segments.end(),
1088  "We should have found the element");
1089  auto & msm_set = it->second;
1090  msm_set.erase(msm_elem);
1091  // We may be creating nodes with only one element neighbor where before this removal there
1092  // were two. But the nodal normal used in computations will reflect the two-neighbor geometry.
1093  // For a lower-d secondary mesh corner, that will imply the corner node will have a tilted
1094  // normal vector (same for tangents) despite the mortar segment mesh not including its
1095  // vertical neighboring element. It is the secondary element neighbors (not mortar segment
1096  // mesh neighbors) that determine the nodal normal field.
1097  if (msm_set.empty())
1099 
1100  // Erase msinfo
1101  _msm_elem_to_info.erase(msm_elem);
1102 
1103  // Remove element from mortar segment mesh
1104  _mortar_segment_mesh->delete_elem(msm_elem);
1105  }
1106  else
1107  {
1110  }
1111  }
1112 
1113  std::unordered_set<Node *> msm_connected_nodes;
1114 
1115  // Deleting elements may produce isolated nodes.
1116  // Loops for identifying and removing such nodes from mortar segment mesh.
1117  for (const auto & element : _mortar_segment_mesh->element_ptr_range())
1118  for (auto & n : element->node_ref_range())
1119  msm_connected_nodes.insert(&n);
1120 
1121  for (const auto & node : _mortar_segment_mesh->node_ptr_range())
1122  if (!msm_connected_nodes.count(node))
1123  _mortar_segment_mesh->delete_node(node);
1124 
1125 #ifdef DEBUG
1126  // Verify that all segments without primary contribution have been deleted
1127  for (auto msm_elem : _mortar_segment_mesh->active_element_ptr_range())
1128  {
1129  const MortarSegmentInfo & msinfo = libmesh_map_find(_msm_elem_to_info, msm_elem);
1130  mooseAssert(msinfo.primary_elem != nullptr,
1131  "All mortar segment elements should have valid "
1132  "primary element.");
1133  }
1134 #endif
1135 
1136  _mortar_segment_mesh->cache_elem_data();
1137 
1138  // (Optionally) Write the mortar segment mesh to file for inspection
1139  if (_debug)
1140  outputMortarMesh();
1141 
1143 }
MetaPhysicL::DualNumber< V, D, asd > abs(const MetaPhysicL::DualNumber< V, D, asd > &a)
Definition: EigenADReal.h:50
std::vector< std::pair< BoundaryID, BoundaryID > > _primary_secondary_boundary_id_pairs
A list of primary/secondary boundary id pairs corresponding to each side of the mortar interface...
const Elem * interior_parent() const
MPI_Info info
std::set< SubdomainID > _secondary_boundary_subdomain_ids
The secondary/primary lower-dimensional boundary subdomain ids are the secondary/primary boundary ids...
void mooseError(Args &&... args)
Emit an error message with the given stringified, concatenated args and terminate the application...
Definition: MooseError.h:311
std::unordered_map< dof_id_type, std::set< Elem *, CompareDofObjectsByID > > _secondary_elems_to_mortar_segments
We maintain a mapping from lower-dimensional secondary elements in the original mesh to (sets of) ele...
void swap(std::vector< T > &data, const std::size_t idx0, const std::size_t idx1, const libMesh::Parallel::Communicator &comm)
Swap function for serial or distributed vector of data.
Definition: Shuffle.h:495
void outputMortarMesh()
Write the mortar segment mesh to exodus.
void mooseWarning(Args &&... args)
Emit a warning message with the given stringified, concatenated args.
Definition: MooseError.h:345
std::unordered_map< const Elem *, MortarSegmentInfo > _msm_elem_to_info
Map between Elems in the mortar segment mesh and their info structs.
void buildCouplingInformation()
build the _mortar_interface_coupling data
static const Real invalid_xi
SECOND
const bool _debug
Whether to print debug output.
const Elem * primary_elem
std::set< SubdomainID > _secondary_ip_sub_ids
All the secondary interior parent subdomain IDs associated with the mortar mesh.
dof_id_type id() const
virtual unsigned int n_nodes() const=0
std::map< std::tuple< dof_id_type, const Node *, const Elem * >, std::pair< Real, const Elem * > > _primary_node_and_elem_to_xi1_secondary_elem
Same type of container, but for mapping (Primary Node ID, Primary Node, Primary Elem) -> (xi^(1)...
std::unordered_map< dof_id_type, std::vector< const Elem * > > _nodes_to_primary_elem_map
TypeVector< typename CompareTypes< Real, T2 >::supertype > cross(const TypeVector< T2 > &v) const
Real _xi_tolerance
Tolerance for checking projection xi values.
const Elem * secondary_elem
Provides a way for users to bail out of the current solve.
void set_unique_id(unique_id_type new_id)
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
std::unordered_map< std::pair< const Node *, const Elem * >, std::pair< Real, const Elem * > > _secondary_node_and_elem_to_xi2_primary_elem
Similar to the map above, but associates a (Secondary Node, Secondary Elem) pair to a (xi^(2)...
Holds xi^(1), xi^(2), and other data for a given mortar segment.
subdomain_id_type subdomain_id() const
const Node * node_ptr(const unsigned int i) const
T fe_lagrange_1D_shape(const Order order, const unsigned int i, const T &xi)
MeshBase & _mesh
Reference to the mesh stored in equation_systems.
virtual Order default_order() const=0
std::set< SubdomainID > _primary_ip_sub_ids
All the primary interior parent subdomain IDs associated with the mortar mesh.
processor_id_type processor_id() const
dof_id_type node_id(const unsigned int i) const
std::vector< Point > getNormals(const Elem &secondary_elem, const std::vector< Point > &xi1_pts) const
Compute the normals at given reference points on a secondary element.
const Point & point(const unsigned int i) const
std::unique_ptr< MeshBase > _mortar_segment_mesh
1D Mesh of mortar segment elements which gets built by the call to build_mortar_segment_mesh().
uint8_t dof_id_type

◆ buildMortarSegmentMesh3d()

void AutomaticMortarGeneration::buildMortarSegmentMesh3d ( )

Builds the mortar segment mesh once the secondary and primary node projections have been completed.

Inputs:

  • mesh

Outputs:

  • mortar_segment_mesh
  • msm_elem_to_info

Step 1: Build mortar segments for all secondary elements

Step 1.1: Linearize secondary face elements

For first order face elements (Tri3 and Quad4) elements are simply linearized around center For second order (Tri6 and Quad9) and third order (Tri7) face elements, elements are sub-divided into four first order elements then each of the sub-elements is linearized around their respective centers For Quad8 elements, they are sub-divided into one quad and four triangle elements and each sub-element is linearized around their respective centers

Step 1.2: Coarse screening using a k-d tree to find nodes on the primary interface that are 'close to' a center point of the secondary element.

Step 1.3: Loop through primary candidate nodes, create mortar segments

Once an element with non-trivial projection onto secondary element identified, switch to breadth-first search (drop all current candidates and add only neighbors of elements with non-trivial overlap)

Step 1.3.2: Sub-divide primary element candidate, then project onto secondary sub-elements, perform polygon clipping, and triangulate to form mortar segments

Step 1.3.3: Create mortar segments and add to mortar segment mesh

Definition at line 1161 of file AutomaticMortarGeneration.C.

Referenced by MortarInterfaceWarehouse::update().

1162 {
1163  const bool use_reference_interpolation =
1164  _mortar_3d_qp_mapping == Mortar3DQuadraturePointMapping::REFERENCE_INTERPOLATION;
1165 
1166  // Add an integer flag to mortar segment mesh to keep track of which subelem
1167  // of second order primal elements mortar segments correspond to
1168  auto secondary_sub_elem = _mortar_segment_mesh->add_elem_integer("secondary_sub_elem");
1169  auto primary_sub_elem = _mortar_segment_mesh->add_elem_integer("primary_sub_elem");
1170 
1171  // Assign globally unique node/element IDs via an exclusive prefix scan: each rank's bound is
1172  // local_secondary_sub_elems * visible_primary_sub_elems * 9, where 9 is the maximum nodes a
1173  // single secondary/primary sub-element pair can produce (8-vertex clipped polygon + center).
1174  // The result is cached and invalidated by meshChanged(), so the allgather only runs on topology
1175  // changes, not on every displaced-mesh residual update.
1176  if (!_msm_node_id_start.has_value())
1177  {
1178  dof_id_type local_secondary_sub_elems = 0, visible_primary_sub_elems = 0;
1179  for (const auto & [primary_sub_id, secondary_sub_id] : _primary_secondary_subdomain_id_pairs)
1180  {
1181  for (const auto * const el :
1183  local_secondary_sub_elems += el->n_sub_elem();
1184  for (const auto * const el : _mesh.active_subdomain_elements_ptr_range(primary_sub_id))
1185  visible_primary_sub_elems += el->n_sub_elem();
1186  }
1187  const dof_id_type per_rank_bound = local_secondary_sub_elems * visible_primary_sub_elems * 9;
1188  std::vector<dof_id_type> per_rank_bounds;
1189  _mesh.comm().allgather(per_rank_bound, per_rank_bounds);
1190  dof_id_type start = 0;
1191  for (const auto r : make_range(_mesh.processor_id()))
1192  start += per_rank_bounds[r];
1193  _msm_node_id_start = start;
1194  }
1195  dof_id_type next_node_id = *_msm_node_id_start;
1196  // Element IDs use the same starting offset: node and element IDs are separately numbered, and
1197  // element count per clip (n triangles) is always <= node count (n+1), so per_rank_bound covers
1198  // both.
1199  dof_id_type next_elem_id = next_node_id;
1200 
1201  // Loop through mortar secondary and primary pairs to create mortar segment mesh between each
1202  for (const auto & pr : _primary_secondary_subdomain_id_pairs)
1203  {
1204  const auto primary_subd_id = pr.first;
1205  const auto secondary_subd_id = pr.second;
1206 
1207  // Build k-d tree for use in Step 1.2 for primary interface coarse screening
1208  NanoflannMeshSubdomainAdaptor<3> mesh_adaptor(_mesh, primary_subd_id);
1209  subdomain_kd_tree_t kd_tree(
1210  3, mesh_adaptor, nanoflann::KDTreeSingleIndexAdaptorParams(/*max leaf=*/10));
1211 
1212  // Construct the KD tree.
1213  kd_tree.buildIndex();
1214 
1215  // Return the unoriented geometric normal of a linearized subpatch. These expressions are the
1216  // TRI3 and QUAD4 mapping tangents evaluated at the reference center, equivalent to evaluating
1217  // the first-order finite-element normal there without constructing a temporary element.
1218  auto get_sub_elem_geometric_normal = [](const std::vector<Point> & nodes)
1219  {
1220  Point dxdxi;
1221  Point dxdeta;
1222  if (nodes.size() == 3)
1223  {
1224  dxdxi = nodes[1] - nodes[0];
1225  dxdeta = nodes[2] - nodes[0];
1226  }
1227  else if (nodes.size() == 4)
1228  {
1229  // Bilinear center tangents define one normal for the full quad instead of selecting one of
1230  // the two diagonal triangle normals.
1231  dxdxi = 0.25 * (nodes[1] + nodes[2] - nodes[0] - nodes[3]);
1232  dxdeta = 0.25 * (nodes[2] + nodes[3] - nodes[0] - nodes[1]);
1233  }
1234  else
1235  mooseError("GEOMETRIC_NORMAL 3D mortar subpatch plane construction only supports "
1236  "triangular and quadrilateral subpatches, but received ",
1237  nodes.size(),
1238  " nodes.");
1239 
1240  Point geometric_normal = dxdxi.cross(dxdeta);
1241  const auto normal_norm = geometric_normal.norm();
1242  // The cross product has units of area, so compare it with the product of tangent lengths.
1243  // Their ratio is the sine of the included angle and is independent of the mesh length scale.
1244  if (normal_norm <= TOLERANCE * dxdxi.norm() * dxdeta.norm())
1245  mooseError("GEOMETRIC_NORMAL 3D mortar subpatch plane construction encountered a "
1246  "degenerate subpatch.");
1247 
1248  geometric_normal /= normal_norm;
1249  return geometric_normal;
1250  };
1251 
1255  for (MeshBase::const_element_iterator el = _mesh.active_local_elements_begin(),
1256  end_el = _mesh.active_local_elements_end();
1257  el != end_el;
1258  ++el)
1259  {
1260  const Elem * secondary_side_elem = *el;
1261 
1262  const Real secondary_volume = secondary_side_elem->volume();
1263 
1264  // If this Elem is not in the current secondary subdomain, go on to the next one.
1265  if (secondary_side_elem->subdomain_id() != secondary_subd_id)
1266  continue;
1267 
1268  auto [secondary_elem_to_msm_map_it, insertion_happened] =
1269  _secondary_elems_to_mortar_segments.emplace(secondary_side_elem->id(),
1270  std::set<Elem *, CompareDofObjectsByID>{});
1271  libmesh_ignore(insertion_happened);
1272  auto & secondary_to_msm_element_set = secondary_elem_to_msm_map_it->second;
1273 
1274  std::vector<std::unique_ptr<MortarSegmentHelper>> mortar_segment_helper(
1275  secondary_side_elem->n_sub_elem());
1276  const auto nodal_normals = getNodalNormals(*secondary_side_elem);
1277 
1288  for (auto sel : make_range(secondary_side_elem->n_sub_elem()))
1289  {
1290  // Get indices of sub-element nodes in element
1291  const auto sub_elem_nodes =
1292  Moose::Mortar::getMortarSubElementNodeIndices(*secondary_side_elem, sel);
1293 
1294  // Secondary sub-element center, normal, and nodes
1295  Point center;
1296  Point normal;
1297  std::vector<Point> nodes(sub_elem_nodes.size());
1298 
1299  // Collect the sub-element points and evaluate its center and averaged nodal normal.
1300  for (auto iv : make_range(sub_elem_nodes.size()))
1301  {
1302  const auto n = sub_elem_nodes[iv];
1303  nodes[iv] = secondary_side_elem->point(n);
1304  center += secondary_side_elem->point(n);
1305  normal += nodal_normals[n];
1306  }
1307  center /= sub_elem_nodes.size();
1308  normal = normal.unit();
1309 
1310  if (_mortar_3d_subpatch_plane == Mortar3DSubpatchPlane::GEOMETRIC_NORMAL)
1311  {
1312  const Point averaged_normal = normal;
1313  normal = get_sub_elem_geometric_normal(nodes);
1314  if (normal * averaged_normal < 0)
1315  normal *= -1;
1316  }
1317 
1318  if (use_reference_interpolation)
1319  {
1320  std::vector<Point> sub_elem_reference_points;
1321  sub_elem_reference_points.reserve(sub_elem_nodes.size());
1322  for (const auto node_index : sub_elem_nodes)
1323  sub_elem_reference_points.push_back(secondary_side_elem->master_point(node_index));
1324 
1325  mortar_segment_helper[sel] =
1326  std::make_unique<MortarSegmentHelper>(std::move(nodes),
1327  std::move(sub_elem_reference_points),
1328  center,
1329  normal,
1332  }
1333  else
1334  mortar_segment_helper[sel] = std::make_unique<MortarSegmentHelper>(
1335  std::move(nodes), center, normal, _triangulation_mode, _triangulate_triangles);
1336  }
1337 
1343  // Search point for performing Nanoflann (k-d tree) searches.
1344  // In each case we use the center point of the original element (not sub-elements for second
1345  // order elements). This is to do search for all sub-elements simultaneously
1346  std::array<Real, 3> query_pt;
1347  Point center_point;
1348  switch (secondary_side_elem->type())
1349  {
1350  case TRI3:
1351  case QUAD4:
1352  center_point = mortar_segment_helper[0]->center();
1353  query_pt = {{center_point(0), center_point(1), center_point(2)}};
1354  break;
1355  case TRI6:
1356  case TRI7:
1357  center_point = mortar_segment_helper[1]->center();
1358  query_pt = {{center_point(0), center_point(1), center_point(2)}};
1359  break;
1360  case QUAD8:
1361  center_point = mortar_segment_helper[4]->center();
1362  query_pt = {{center_point(0), center_point(1), center_point(2)}};
1363  break;
1364  case QUAD9:
1365  center_point = secondary_side_elem->point(8);
1366  query_pt = {{center_point(0), center_point(1), center_point(2)}};
1367  break;
1368  default:
1369  mooseError(
1370  "Face element type: ", secondary_side_elem->type(), "not supported for 3D mortar");
1371  }
1372 
1373  // The number of results we want to get. These results will only be used to find
1374  // a single element with non-trivial overlap, after an element is identified a breadth
1375  // first search is done on neighbors
1376  const std::size_t num_results = 3;
1377 
1378  // Initialize result_set and do the search.
1379  std::vector<size_t> ret_index(num_results);
1380  std::vector<Real> out_dist_sqr(num_results);
1381  nanoflann::KNNResultSet<Real> result_set(num_results);
1382  result_set.init(&ret_index[0], &out_dist_sqr[0]);
1383  kd_tree.findNeighbors(result_set, &query_pt[0], nanoflann::SearchParameters());
1384 
1385  // Initialize list of processed primary elements, we don't want to revisit processed elements
1386  std::set<const Elem *, CompareDofObjectsByID> processed_primary_elems;
1387 
1388  // Initialize candidate set and flag for switching between coarse screening and breadth-first
1389  // search
1390  bool primary_elem_found = false;
1391  std::set<const Elem *, CompareDofObjectsByID> primary_elem_candidates;
1392  const bool use_geometric_subpatch_normals =
1393  _mortar_3d_subpatch_plane == Mortar3DSubpatchPlane::GEOMETRIC_NORMAL;
1394  // In geometric mode the projection-angle cutoff also rejects near-orthogonal subpatch pairs.
1395  // The absolute dot product below keeps opposing primary/secondary orientations admissible.
1396  const Real minimum_subpatch_normal_alignment =
1397  use_geometric_subpatch_normals ? std::sin(_minimum_projection_angle * libMesh::pi / 180.0)
1398  : 0.0;
1399 
1400  // Loop candidate nodes (returned by Nanoflann) and add all adjoining elems to candidate set
1401  for (auto r : make_range(result_set.size()))
1402  {
1403  // Verify that the squared distance we compute is the same as nanoflann's
1404  mooseAssert(abs((_mesh.point(ret_index[r]) - center_point).norm_sq() - out_dist_sqr[r]) <=
1405  TOLERANCE,
1406  "Lower-dimensional element squared distance verification failed.");
1407 
1408  // Get list of elems connected to node
1409  std::vector<const Elem *> & node_elems =
1410  this->_nodes_to_primary_elem_map.at(static_cast<dof_id_type>(ret_index[r]));
1411 
1412  // Uniquely add elems to candidate set
1413  for (auto elem : node_elems)
1414  primary_elem_candidates.insert(elem);
1415  }
1416 
1424  while (!primary_elem_candidates.empty())
1425  {
1426  const Elem * primary_elem_candidate = *primary_elem_candidates.begin();
1427 
1428  // If we've already processed this candidate, we don't need to check it again.
1429  if (processed_primary_elems.count(primary_elem_candidate))
1430  {
1431  primary_elem_candidates.erase(primary_elem_candidate);
1432  continue;
1433  }
1434 
1435  // Initialize set of nodes used to construct mortar segment elements
1436  std::vector<Point> nodal_points;
1437 
1438  // Initialize map from mortar segment elements to nodes
1439  std::vector<std::vector<unsigned int>> elem_to_node_map;
1440 
1441  // Initialize list of secondary and primary sub-elements that formed each mortar segment
1442  std::vector<std::pair<unsigned int, unsigned int>> sub_elem_map;
1443  std::vector<std::array<Point, 3>> elem_to_secondary_reference_points;
1444  std::vector<std::array<Point, 3>> elem_to_primary_reference_points;
1445 
1450  for (auto p_el : make_range(primary_elem_candidate->n_sub_elem()))
1451  {
1452  // Get nodes of primary sub-elements
1453  const auto sub_elem_nodes =
1454  Moose::Mortar::getMortarSubElementNodeIndices(*primary_elem_candidate, p_el);
1455 
1456  // Get list of primary sub-element vertex nodes
1457  std::vector<Point> primary_sub_elem(sub_elem_nodes.size());
1458  for (auto iv : make_range(sub_elem_nodes.size()))
1459  {
1460  const auto n = sub_elem_nodes[iv];
1461  primary_sub_elem[iv] = primary_elem_candidate->point(n);
1462  }
1463  Point primary_sub_elem_normal;
1464  if (use_geometric_subpatch_normals)
1465  primary_sub_elem_normal = get_sub_elem_geometric_normal(primary_sub_elem);
1466 
1467  std::vector<Point> sub_elem_reference_points;
1468  if (use_reference_interpolation)
1469  {
1470  sub_elem_reference_points.reserve(sub_elem_nodes.size());
1471  for (const auto node_index : sub_elem_nodes)
1472  sub_elem_reference_points.push_back(primary_elem_candidate->master_point(node_index));
1473  }
1474 
1475  // Loop through secondary sub-elements
1476  for (auto s_el : make_range(secondary_side_elem->n_sub_elem()))
1477  {
1478  // Nearby primary candidates can include adjacent corner faces. Those faces may clip to
1479  // numerical slivers, which we do not consider valid face-to-face mortar pairs for this
1480  // search.
1481  if (use_geometric_subpatch_normals &&
1482  std::abs(primary_sub_elem_normal * mortar_segment_helper[s_el]->normal()) <
1483  minimum_subpatch_normal_alignment)
1484  continue;
1485 
1486  // Mortar segment helpers were defined for each secondary sub-element, they will:
1487  // 1. Project primary sub-element onto linearized secondary sub-element
1488  // 2. Clip projected primary sub-element against secondary sub-element
1489  // 3. Triangulate clipped polygon to form mortar segments
1490  //
1491  // Mortar segment helpers append a list of mortar segment nodes and connectivities that
1492  // can be directly used to build mortar segments
1493  const auto segments_before_helper = elem_to_node_map.size();
1494  if (use_reference_interpolation)
1495  mortar_segment_helper[s_el]->getMortarSegments(primary_sub_elem,
1496  sub_elem_reference_points,
1497  nodal_points,
1498  elem_to_node_map,
1499  elem_to_secondary_reference_points,
1500  elem_to_primary_reference_points,
1501  TOLERANCE * secondary_volume);
1502  else
1503  mortar_segment_helper[s_el]->getMortarSegments(
1504  primary_sub_elem, nodal_points, elem_to_node_map);
1505 
1506  // Keep track of which secondary and primary sub-elements created segment
1507  for (auto i = segments_before_helper; i < elem_to_node_map.size(); ++i)
1508  sub_elem_map.push_back(std::make_pair(s_el, p_el));
1509  }
1510  }
1511 
1512  // Mark primary element as processed and remove from candidate list
1513  processed_primary_elems.insert(primary_elem_candidate);
1514  primary_elem_candidates.erase(primary_elem_candidate);
1515 
1516  // If overlap of polygons was non-trivial (created mortar segment elements)
1517  if (!elem_to_node_map.empty())
1518  {
1519  if (sub_elem_map.size() != elem_to_node_map.size())
1520  mooseError("The mortar segment subpatch map is not aligned with the mortar segment "
1521  "connectivity map.");
1522  if (use_reference_interpolation &&
1523  (elem_to_secondary_reference_points.size() != elem_to_node_map.size() ||
1524  elem_to_primary_reference_points.size() != elem_to_node_map.size()))
1525  mooseError("The mortar segment reference-point maps are not aligned with the mortar "
1526  "segment connectivity map.");
1527 
1528  // Only overlap polygons large enough to become mortar segments may switch the candidate
1529  // search to breadth first.
1530  bool seed_breadth_first_search = false;
1531  std::vector<bool> retained_mortar_segments(elem_to_node_map.size(), false);
1532  for (const auto el : index_range(elem_to_node_map))
1533  {
1534  const auto & node_map = elem_to_node_map[el];
1535  if (node_map.size() != 3)
1536  mooseError(
1537  "Active mortar segments only supports TRI elements, 3 nodes expected but: ",
1538  node_map.size(),
1539  " provided.");
1540 
1541  const Point e1 = nodal_points[node_map[1]] - nodal_points[node_map[0]];
1542  const Point e2 = nodal_points[node_map[2]] - nodal_points[node_map[0]];
1543  retained_mortar_segments[el] =
1544  0.5 * e1.cross(e2).norm() / secondary_volume >= TOLERANCE;
1545  seed_breadth_first_search = seed_breadth_first_search || retained_mortar_segments[el];
1546  }
1547 
1548  if (seed_breadth_first_search)
1549  {
1550  // If this is the first element with a qualifying overlap, set flag. Candidates will
1551  // now be neighbors of elements that had qualifying overlap.
1552  if (!primary_elem_found)
1553  {
1554  primary_elem_found = true;
1555  primary_elem_candidates.clear();
1556  }
1557 
1558  // Add neighbors to candidate list
1559  for (auto neighbor : primary_elem_candidate->neighbor_ptr_range())
1560  {
1561  // If not valid or not on lower dimensional secondary subdomain, skip
1562  if (neighbor == nullptr || neighbor->subdomain_id() != primary_subd_id)
1563  continue;
1564  // If already processed, skip
1565  if (processed_primary_elems.count(neighbor))
1566  continue;
1567  // Otherwise, add to candidates
1568  primary_elem_candidates.insert(neighbor);
1569  }
1570  }
1571 
1575  std::vector<Node *> new_nodes;
1576  // Clipping can append points for triangles later rejected by the area tolerance. Add only
1577  // points referenced by retained triangles so the mortar mesh has no orphan nodes.
1578  std::vector<bool> retained_nodes(nodal_points.size(), false);
1579  for (const auto el : index_range(elem_to_node_map))
1580  if (retained_mortar_segments[el])
1581  for (const auto node : elem_to_node_map[el])
1582  retained_nodes[node] = true;
1583 
1584  new_nodes.resize(nodal_points.size(), nullptr);
1585  for (const auto node : index_range(nodal_points))
1586  if (retained_nodes[node])
1587  new_nodes[node] = _mortar_segment_mesh->add_point(
1588  nodal_points[node], next_node_id++, secondary_side_elem->processor_id());
1589 
1590  // Loop through triangular elements in map
1591  for (auto el : index_range(elem_to_node_map))
1592  {
1593  if (!retained_mortar_segments[el])
1594  continue;
1595 
1596  std::unique_ptr<Elem> new_elem;
1597  if (elem_to_node_map[el].size() == 3)
1598  new_elem = std::make_unique<Tri3>();
1599  else
1600  mooseError("Active mortar segments only supports TRI elements, 3 nodes expected "
1601  "but: ",
1602  elem_to_node_map[el].size(),
1603  " provided.");
1604 
1605  new_elem->processor_id() = secondary_side_elem->processor_id();
1606  new_elem->subdomain_id() = secondary_side_elem->subdomain_id();
1607  new_elem->set_id(next_elem_id++);
1608 
1609  // Attach newly created nodes
1610  for (auto i : index_range(elem_to_node_map[el]))
1611  new_elem->set_node(i, new_nodes[elem_to_node_map[el][i]]);
1612 
1613  // If element is smaller than tolerance, don't add to msm
1614  if (new_elem->volume() / secondary_volume < TOLERANCE)
1615  continue;
1616 
1617  // Add elements to mortar segment mesh
1618  Elem * msm_new_elem = _mortar_segment_mesh->add_elem(new_elem.release());
1619 
1620  msm_new_elem->set_extra_integer(secondary_sub_elem, sub_elem_map[el].first);
1621  msm_new_elem->set_extra_integer(primary_sub_elem, sub_elem_map[el].second);
1622 
1623  // Fill out mortar segment info
1624  MortarSegmentInfo msinfo;
1625  msinfo.secondary_elem = secondary_side_elem;
1626  msinfo.primary_elem = primary_elem_candidate;
1627 
1628  // Associate this MSM elem with the MortarSegmentInfo.
1629  _msm_elem_to_info.emplace(msm_new_elem, msinfo);
1630 
1631  // Store reference data only for retained segments.
1632  if (use_reference_interpolation)
1633  {
1634  MortarSegmentReferencePoints reference_points{elem_to_secondary_reference_points[el],
1635  elem_to_primary_reference_points[el]};
1636  _msm_elem_to_reference_points.emplace(msm_new_elem, reference_points);
1637  }
1638 
1639  // Add this mortar segment to the secondary elem to mortar segment map
1640  secondary_to_msm_element_set.insert(msm_new_elem);
1641 
1643  // Unlike for 2D, we always have a primary when building the mortar mesh so we don't
1644  // have to check for null
1646  }
1647  }
1648  // End loop through primary element candidates
1649  }
1650 
1651  if (use_geometric_subpatch_normals)
1652  {
1653  // A geometric corner filter may intentionally leave individual subpatches uncovered. Warn
1654  // only when the complete secondary element failed to produce a retained segment.
1655  if (secondary_to_msm_element_set.empty())
1656  mooseDoOnce(
1657  mooseWarning("Some secondary elements on mortar interface were unable to identify"
1658  " a corresponding primary element; this may be expected depending on"
1659  " problem geometry but may indicate a failure of the element search"
1660  " or projection"));
1661  }
1662  else
1663  for (auto sel : make_range(secondary_side_elem->n_sub_elem()))
1664  if (mortar_segment_helper[sel]->remainder() == 1.0)
1665  mooseDoOnce(
1666  mooseWarning("Some secondary elements on mortar interface were unable to identify"
1667  " a corresponding primary element; this may be expected depending on"
1668  " problem geometry but may indicate a failure of the element search"
1669  " or projection"));
1670 
1671  if (secondary_to_msm_element_set.empty())
1672  _secondary_elems_to_mortar_segments.erase(secondary_elem_to_msm_map_it);
1673  } // End loop through secondary elements
1674  } // End loop through mortar constraint pairs
1675 
1676  mooseAssert(!use_reference_interpolation ||
1678  "Mortar segment info and reference-point maps must remain aligned.");
1679 
1680  _mortar_segment_mesh->cache_elem_data();
1681 
1682  // The mesh was built distributedly (each rank owns only its local elements), so mark it
1683  // as such so MeshSerializer correctly gathers it to proc 0 for Exodus output.
1684  _mortar_segment_mesh->set_distributed();
1685 
1686  // Output mortar segment mesh
1687  if (_debug)
1688  {
1689  // If element is not triangular, increment subdomain id
1690  // (ExodusII does not support mixed element types in a single subdomain)
1691  for (const auto msm_el : _mortar_segment_mesh->active_local_element_ptr_range())
1692  if (msm_el->type() != TRI3)
1693  msm_el->subdomain_id()++;
1694 
1695  outputMortarMesh();
1696 
1697  // Undo increment
1698  for (const auto msm_el : _mortar_segment_mesh->active_local_element_ptr_range())
1699  if (msm_el->type() != TRI3)
1700  msm_el->subdomain_id()--;
1701  }
1702 
1704 
1705  // Print mortar segment mesh statistics
1706  if (_debug)
1707  {
1708  msmStatistics();
1709  }
1710 }
MetaPhysicL::DualNumber< V, D, asd > abs(const MetaPhysicL::DualNumber< V, D, asd > &a)
Definition: EigenADReal.h:50
void allgather(const T &send_data, std::vector< T, A > &recv_data) const
CTSub CT_OPERATOR_BINARY CTMul CTCompareLess CTCompareGreater CTCompareEqual _arg template * sin(_arg) *_arg.template D< dtag >()) CT_SIMPLE_UNARY_FUNCTION(tan
auto norm() const
QUAD8
const Elem * interior_parent() const
std::vector< unsigned int > getMortarSubElementNodeIndices(const Elem &parent_elem, unsigned int sub_elem)
Return the node indices for a first-order sub-element of a parent face.
Definition: MortarUtils.C:27
void mooseError(Args &&... args)
Emit an error message with the given stringified, concatenated args and terminate the application...
Definition: MooseError.h:311
std::unordered_map< dof_id_type, std::set< Elem *, CompareDofObjectsByID > > _secondary_elems_to_mortar_segments
We maintain a mapping from lower-dimensional secondary elements in the original mesh to (sets of) ele...
auto norm_sq(const T &a)
void outputMortarMesh()
Write the mortar segment mesh to exodus.
void mooseWarning(Args &&... args)
Emit a warning message with the given stringified, concatenated args.
Definition: MooseError.h:345
std::unordered_map< const Elem *, MortarSegmentInfo > _msm_elem_to_info
Map between Elems in the mortar segment mesh and their info structs.
const Parallel::Communicator & comm() const
void msmStatistics()
Prints mortar segment mesh statistics to console (calls computeMsmStatistics internally) ...
void buildCouplingInformation()
build the _mortar_interface_coupling data
std::vector< Point > getNodalNormals(const Elem &secondary_elem) const
Special adaptor that works with subdomains of the Mesh.
TRI3
QUAD4
const bool _debug
Whether to print debug output.
const Elem * primary_elem
const bool _triangulate_triangles
Whether already-triangular clipped polygons should still be centroid-subdivided.
std::set< SubdomainID > _secondary_ip_sub_ids
All the secondary interior parent subdomain IDs associated with the mortar mesh.
std::optional< dof_id_type > _msm_node_id_start
Cached per-rank starting ID for 3D MSM nodes/elements.
TypeVector< Real > unit() const
void libmesh_ignore(const Args &...)
dof_id_type id() const
TRI6
std::unordered_map< dof_id_type, std::vector< const Elem * > > _nodes_to_primary_elem_map
TypeVector< typename CompareTypes< Real, T2 >::supertype > cross(const TypeVector< T2 > &v) const
nanoflann::KDTreeSingleIndexAdaptor< subdomain_adatper_t, NanoflannMeshSubdomainAdaptor< 3 >, 3 > subdomain_kd_tree_t
std::unordered_map< const Elem *, MortarSegmentReferencePoints > _msm_elem_to_reference_points
Reference-coordinate data used only by the reference-interpolation mapping mode.
const Elem * secondary_elem
const MortarSegmentTriangulationMode _triangulation_mode
Triangulation mode used for clipped 3D mortar polygons.
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
Holds xi^(1), xi^(2), and other data for a given mortar segment.
virtual SimpleRange< element_iterator > active_local_subdomain_elements_ptr_range(subdomain_id_type sid)=0
subdomain_id_type subdomain_id() const
Parent-face reference coordinates associated with the vertices of one triangular mortar segment...
const Mortar3DSubpatchPlane _mortar_3d_subpatch_plane
Method used to define the local projection planes for 3D secondary subpatches.
virtual Point master_point(const unsigned int i) const=0
virtual Real volume() const
IntRange< T > make_range(T beg, T end)
TRI7
unsigned int level ElemType type std::set< subdomain_id_type > ss processor_id_type pid unsigned int level std::set< subdomain_id_type > virtual ss SimpleRange< element_iterator > active_subdomain_elements_ptr_range(subdomain_id_type sid)=0
const Real _minimum_projection_angle
Parameter to control which angle (in degrees) is admissible for the creation of mortar segments...
QUAD9
MeshBase & _mesh
Reference to the mesh stored in equation_systems.
virtual const Point & point(const dof_id_type i) const=0
std::vector< std::pair< SubdomainID, SubdomainID > > _primary_secondary_subdomain_id_pairs
A list of primary/secondary subdomain id pairs corresponding to each side of the mortar interface...
SimpleRange< NeighborPtrIter > neighbor_ptr_range()
virtual unsigned int n_sub_elem() const=0
processor_id_type processor_id() const
std::set< SubdomainID > _primary_ip_sub_ids
All the primary interior parent subdomain IDs associated with the mortar mesh.
SearchParams SearchParameters
processor_id_type processor_id() const
virtual ElemType type() const=0
const Mortar3DQuadraturePointMapping _mortar_3d_qp_mapping
Method used to map 3D mortar segment quadrature points to primary and secondary faces.
const Point & point(const unsigned int i) const
std::unique_ptr< MeshBase > _mortar_segment_mesh
1D Mesh of mortar segment elements which gets built by the call to build_mortar_segment_mesh().
auto index_range(const T &sizable)
void set_extra_integer(const unsigned int index, const dof_id_type value)
uint8_t dof_id_type
const Real pi

◆ buildNodeToElemMaps()

void AutomaticMortarGeneration::buildNodeToElemMaps ( )

Once the secondary_requested_boundary_ids and primary_requested_boundary_ids containers have been filled in, call this function to build node-to-Elem maps for the lower-dimensional elements.

Definition at line 505 of file AutomaticMortarGeneration.C.

Referenced by MortarInterfaceWarehouse::update().

506 {
508  mooseError(
509  "Must specify secondary and primary boundary ids before building node-to-elem maps.");
510 
511  // Construct nodes_to_secondary_elem_map
512  for (const auto & secondary_elem :
513  as_range(_mesh.active_elements_begin(), _mesh.active_elements_end()))
514  {
515  // If this is not one of the lower-dimensional secondary side elements, go on to the next one.
516  if (!this->_secondary_boundary_subdomain_ids.count(secondary_elem->subdomain_id()))
517  continue;
518 
519  for (const auto & nd : secondary_elem->node_ref_range())
520  {
521  std::vector<const Elem *> & vec = _nodes_to_secondary_elem_map[nd.id()];
522  vec.push_back(secondary_elem);
523  }
524  }
525 
526  // Construct nodes_to_primary_elem_map
527  for (const auto & primary_elem :
528  as_range(_mesh.active_elements_begin(), _mesh.active_elements_end()))
529  {
530  // If this is not one of the lower-dimensional primary side elements, go on to the next one.
531  if (!this->_primary_boundary_subdomain_ids.count(primary_elem->subdomain_id()))
532  continue;
533 
534  for (const auto & nd : primary_elem->node_ref_range())
535  {
536  std::vector<const Elem *> & vec = _nodes_to_primary_elem_map[nd.id()];
537  vec.push_back(primary_elem);
538  }
539  }
540 }
std::set< SubdomainID > _primary_boundary_subdomain_ids
std::set< SubdomainID > _secondary_boundary_subdomain_ids
The secondary/primary lower-dimensional boundary subdomain ids are the secondary/primary boundary ids...
void mooseError(Args &&... args)
Emit an error message with the given stringified, concatenated args and terminate the application...
Definition: MooseError.h:311
std::unordered_map< dof_id_type, std::vector< const Elem * > > _nodes_to_primary_elem_map
SimpleRange< IndexType > as_range(const std::pair< IndexType, IndexType > &p)
std::set< BoundaryID > _primary_requested_boundary_ids
The boundary ids corresponding to all the primary surfaces.
MeshBase & _mesh
Reference to the mesh stored in equation_systems.
std::set< BoundaryID > _secondary_requested_boundary_ids
The boundary ids corresponding to all the secondary surfaces.
std::unordered_map< dof_id_type, std::vector< const Elem * > > _nodes_to_secondary_elem_map
Map from nodes to connected lower-dimensional elements on the secondary/primary subdomains.

◆ clear()

void AutomaticMortarGeneration::clear ( )

Clears the mortar segment mesh and accompanying data structures.

Definition at line 466 of file AutomaticMortarGeneration.C.

Referenced by MortarInterfaceWarehouse::update().

467 {
469  _mortar_segment_mesh->clear();
474  _msm_elem_to_info.clear();
475  _lower_elem_to_side_id.clear();
481  _secondary_ip_sub_ids.clear();
482  _primary_ip_sub_ids.clear();
485 }
std::unordered_map< dof_id_type, std::set< Elem *, CompareDofObjectsByID > > _secondary_elems_to_mortar_segments
We maintain a mapping from lower-dimensional secondary elements in the original mesh to (sets of) ele...
std::unordered_map< const Elem *, MortarSegmentInfo > _msm_elem_to_info
Map between Elems in the mortar segment mesh and their info structs.
std::unordered_set< dof_id_type > _projected_secondary_nodes
Debugging container for printing information about fraction of successful projections for secondary n...
std::set< SubdomainID > _secondary_ip_sub_ids
All the secondary interior parent subdomain IDs associated with the mortar mesh.
std::unordered_map< dof_id_type, const Elem * > _secondary_element_to_secondary_lowerd_element
Map from full dimensional secondary element id to lower dimensional secondary element.
std::map< std::tuple< dof_id_type, const Node *, const Elem * >, std::pair< Real, const Elem * > > _primary_node_and_elem_to_xi1_secondary_elem
Same type of container, but for mapping (Primary Node ID, Primary Node, Primary Elem) -> (xi^(1)...
std::unordered_map< dof_id_type, std::vector< const Elem * > > _nodes_to_primary_elem_map
std::unordered_map< const Elem *, MortarSegmentReferencePoints > _msm_elem_to_reference_points
Reference-coordinate data used only by the reference-interpolation mapping mode.
std::unordered_set< dof_id_type > _failed_secondary_node_projections
Secondary nodes that failed to project.
std::unordered_map< std::pair< const Node *, const Elem * >, std::pair< Real, const Elem * > > _secondary_node_and_elem_to_xi2_primary_elem
Similar to the map above, but associates a (Secondary Node, Secondary Elem) pair to a (xi^(2)...
std::unordered_map< const Node *, std::array< Point, 2 > > _secondary_node_to_hh_nodal_tangents
Container for storing the nodal tangent/binormal vectors associated with each secondary node (Househo...
std::unordered_map< const Node *, Point > _secondary_node_to_nodal_normal
Container for storing the nodal normal vector associated with each secondary node.
std::unordered_map< dof_id_type, std::unordered_set< dof_id_type > > _mortar_interface_coupling
Used by the AugmentSparsityOnInterface functor to determine whether a given Elem is coupled to any ot...
std::set< SubdomainID > _primary_ip_sub_ids
All the primary interior parent subdomain IDs associated with the mortar mesh.
std::unique_ptr< MeshBase > _mortar_segment_mesh
1D Mesh of mortar segment elements which gets built by the call to build_mortar_segment_mesh().
std::unordered_map< dof_id_type, std::vector< const Elem * > > _nodes_to_secondary_elem_map
Map from nodes to connected lower-dimensional elements on the secondary/primary subdomains.
std::unordered_map< const Elem *, unsigned int > _lower_elem_to_side_id
Keeps track of the mapping between lower-dimensional elements and the side_id of the interior_parent ...

◆ computeInactiveLMElems()

void AutomaticMortarGeneration::computeInactiveLMElems ( )

Get list of secondary elems without any corresponding primary elements.

Used to enforce zero values on inactive DoFs of elemental variables.

Definition at line 2084 of file AutomaticMortarGeneration.C.

Referenced by MortarInterfaceWarehouse::update().

2085 {
2086  // Mark all active secondary elements
2087  std::unordered_set<const Elem *> active_local_elems;
2088 
2089  //****
2090  // Note that in 3D our trick to check whether an element has edge dropping needs loose tolerances
2091  // since the mortar segments are on the linearized element and comparing the volume of the
2092  // linearized element does not have the same volume as the warped element
2093  const Real tol = (dim() == 3) ? 0.1 : TOLERANCE;
2094 
2095  std::unordered_map<const Elem *, Real> active_volume;
2096 
2097  // Compute fraction of elements with corresponding primary elements
2099  for (const auto msm_elem : _mortar_segment_mesh->active_local_element_ptr_range())
2100  {
2101  const MortarSegmentInfo & msinfo = _msm_elem_to_info.at(msm_elem);
2102  const Elem * secondary_elem = msinfo.secondary_elem;
2103 
2104  active_volume[secondary_elem] += msm_elem->volume();
2105  }
2106  //****
2107 
2108  for (const auto msm_elem : _mortar_segment_mesh->active_local_element_ptr_range())
2109  {
2110  const MortarSegmentInfo & msinfo = _msm_elem_to_info.at(msm_elem);
2111  const Elem * secondary_elem = msinfo.secondary_elem;
2112 
2113  //****
2115  if (abs(active_volume[secondary_elem] / secondary_elem->volume() - 1.0) > tol)
2116  continue;
2117  //****
2118 
2119  active_local_elems.insert(secondary_elem);
2120  }
2121 
2122  // Take complement of active elements in active local subdomain to get inactive local elements
2123  _inactive_local_lm_elems.clear();
2124  for (const auto & pr : _primary_secondary_subdomain_id_pairs)
2125  for (const auto el : _mesh.active_local_subdomain_elements_ptr_range(
2126  /*secondary_subd_id*/ pr.second))
2127  if (active_local_elems.find(el) == active_local_elems.end())
2128  _inactive_local_lm_elems.insert(el);
2129 }
MetaPhysicL::DualNumber< V, D, asd > abs(const MetaPhysicL::DualNumber< V, D, asd > &a)
Definition: EigenADReal.h:50
std::unordered_set< const Elem * > _inactive_local_lm_elems
List of inactive lagrange multiplier nodes (for elemental variables)
std::unordered_map< const Elem *, MortarSegmentInfo > _msm_elem_to_info
Map between Elems in the mortar segment mesh and their info structs.
const bool _correct_edge_dropping
Flag to enable regressed treatment of edge dropping where all LM DoFs on edge dropping element are st...
const Elem * secondary_elem
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
Holds xi^(1), xi^(2), and other data for a given mortar segment.
virtual SimpleRange< element_iterator > active_local_subdomain_elements_ptr_range(subdomain_id_type sid)=0
virtual Real volume() const
MeshBase & _mesh
Reference to the mesh stored in equation_systems.
std::vector< std::pair< SubdomainID, SubdomainID > > _primary_secondary_subdomain_id_pairs
A list of primary/secondary subdomain id pairs corresponding to each side of the mortar interface...
std::unique_ptr< MeshBase > _mortar_segment_mesh
1D Mesh of mortar segment elements which gets built by the call to build_mortar_segment_mesh().

◆ computeInactiveLMNodes()

void AutomaticMortarGeneration::computeInactiveLMNodes ( )

Get list of secondary nodes that don't contribute to interaction with any primary element.

Used to enforce zero values on inactive DoFs of nodal variables.

Definition at line 2000 of file AutomaticMortarGeneration.C.

Referenced by MortarInterfaceWarehouse::update().

2001 {
2003  {
2005  return;
2006  }
2007 
2008  std::unordered_map<processor_id_type, std::set<dof_id_type>> proc_to_active_nodes_set;
2009  const auto my_pid = _mesh.processor_id();
2010 
2011  // List of active nodes on local secondary elements
2012  std::unordered_set<dof_id_type> active_local_nodes;
2013 
2014  // Mark all active local nodes
2015  for (const auto msm_elem : _mortar_segment_mesh->active_local_element_ptr_range())
2016  {
2017  const MortarSegmentInfo & msinfo = _msm_elem_to_info.at(msm_elem);
2018  const Elem * secondary_elem = msinfo.secondary_elem;
2019 
2020  for (auto n : make_range(secondary_elem->n_nodes()))
2021  active_local_nodes.insert(secondary_elem->node_id(n));
2022  }
2023 
2024  // Assemble list of procs that nodes contribute to
2025  for (const auto & pr : _primary_secondary_subdomain_id_pairs)
2026  {
2027  const auto secondary_subd_id = pr.second;
2028 
2029  // Loop through all elements not on my processor
2030  for (const auto el : _mesh.active_subdomain_elements_ptr_range(secondary_subd_id))
2031  {
2032  // Get processor_id
2033  const auto pid = el->processor_id();
2034 
2035  // If element is in my subdomain, skip
2036  if (pid == my_pid)
2037  continue;
2038 
2039  // If element on proc pid shares any of my active nodes, mark to send
2040  for (const auto n : make_range(el->n_nodes()))
2041  {
2042  const auto node_id = el->node_id(n);
2043  if (active_local_nodes.find(node_id) != active_local_nodes.end())
2044  proc_to_active_nodes_set[pid].insert(node_id);
2045  }
2046  }
2047  }
2048 
2049  // Send list of active nodes
2050  {
2051  // Pack set into vector for sending (push_parallel_vector_data doesn't like sets)
2052  std::unordered_map<processor_id_type, std::vector<dof_id_type>> proc_to_active_nodes_vector;
2053  for (const auto & proc_set : proc_to_active_nodes_set)
2054  {
2055  proc_to_active_nodes_vector[proc_set.first].reserve(proc_to_active_nodes_set.size());
2056  for (const auto node_id : proc_set.second)
2057  proc_to_active_nodes_vector[proc_set.first].push_back(node_id);
2058  }
2059 
2060  // First push data
2061  auto action_functor = [this, &active_local_nodes](const processor_id_type pid,
2062  const std::vector<dof_id_type> & sent_data)
2063  {
2064  if (pid == _mesh.processor_id())
2065  mooseError("Should not be communicating with self.");
2066  active_local_nodes.insert(sent_data.begin(), sent_data.end());
2067  };
2068  TIMPI::push_parallel_vector_data(_mesh.comm(), proc_to_active_nodes_vector, action_functor);
2069  }
2070 
2071  // Every proc has correct list of active local nodes, now take complement (list of inactive nodes)
2072  // and store to use later to zero LM DoFs on inactive nodes
2073  _inactive_local_lm_nodes.clear();
2074  for (const auto & pr : _primary_secondary_subdomain_id_pairs)
2075  for (const auto el : _mesh.active_local_subdomain_elements_ptr_range(
2076  /*secondary_subd_id*/ pr.second))
2077  for (const auto n : make_range(el->n_nodes()))
2078  if (active_local_nodes.find(el->node_id(n)) == active_local_nodes.end())
2079  _inactive_local_lm_nodes.insert(el->node_ptr(n));
2080 }
void mooseError(Args &&... args)
Emit an error message with the given stringified, concatenated args and terminate the application...
Definition: MooseError.h:311
std::unordered_map< const Elem *, MortarSegmentInfo > _msm_elem_to_info
Map between Elems in the mortar segment mesh and their info structs.
const Parallel::Communicator & comm() const
void push_parallel_vector_data(const Communicator &comm, MapToVectors &&data, const ActionFunctor &act_on_data)
uint8_t processor_id_type
virtual unsigned int n_nodes() const=0
const bool _correct_edge_dropping
Flag to enable regressed treatment of edge dropping where all LM DoFs on edge dropping element are st...
void computeIncorrectEdgeDroppingInactiveLMNodes()
Computes inactive secondary nodes when incorrect edge dropping behavior is enabled (any node touching...
const Elem * secondary_elem
Holds xi^(1), xi^(2), and other data for a given mortar segment.
virtual SimpleRange< element_iterator > active_local_subdomain_elements_ptr_range(subdomain_id_type sid)=0
IntRange< T > make_range(T beg, T end)
unsigned int level ElemType type std::set< subdomain_id_type > ss processor_id_type pid unsigned int level std::set< subdomain_id_type > virtual ss SimpleRange< element_iterator > active_subdomain_elements_ptr_range(subdomain_id_type sid)=0
MeshBase & _mesh
Reference to the mesh stored in equation_systems.
std::vector< std::pair< SubdomainID, SubdomainID > > _primary_secondary_subdomain_id_pairs
A list of primary/secondary subdomain id pairs corresponding to each side of the mortar interface...
std::unordered_set< const Node * > _inactive_local_lm_nodes
processor_id_type processor_id() const
processor_id_type processor_id() const
dof_id_type node_id(const unsigned int i) const
std::unique_ptr< MeshBase > _mortar_segment_mesh
1D Mesh of mortar segment elements which gets built by the call to build_mortar_segment_mesh().

◆ computeIncorrectEdgeDroppingInactiveLMNodes()

void AutomaticMortarGeneration::computeIncorrectEdgeDroppingInactiveLMNodes ( )

Computes inactive secondary nodes when incorrect edge dropping behavior is enabled (any node touching a partially or fully dropped element is dropped)

Definition at line 1906 of file AutomaticMortarGeneration.C.

Referenced by computeInactiveLMNodes().

1907 {
1908  using std::abs;
1909 
1910  // Note that in 3D our trick to check whether an element has edge dropping needs loose tolerances
1911  // since the mortar segments are on the linearized element and comparing the volume of the
1912  // linearized element does not have the same volume as the warped element
1913  const Real tol = (dim() == 3) ? 0.1 : TOLERANCE;
1914 
1915  std::unordered_map<processor_id_type, std::set<dof_id_type>> proc_to_inactive_nodes_set;
1916  const auto my_pid = _mesh.processor_id();
1917 
1918  // List of inactive nodes on local secondary elements
1919  std::unordered_set<dof_id_type> inactive_node_ids;
1920 
1921  std::unordered_map<const Elem *, Real> active_volume{};
1922 
1923  for (const auto & pr : _primary_secondary_subdomain_id_pairs)
1924  for (const auto el : _mesh.active_subdomain_elements_ptr_range(pr.second))
1925  active_volume[el] = 0.;
1926 
1927  // Compute fraction of elements with corresponding primary elements
1928  for (const auto msm_elem : _mortar_segment_mesh->active_local_element_ptr_range())
1929  {
1930  const MortarSegmentInfo & msinfo = _msm_elem_to_info.at(msm_elem);
1931  const Elem * secondary_elem = msinfo.secondary_elem;
1932 
1933  active_volume[secondary_elem] += msm_elem->volume();
1934  }
1935 
1936  // Mark all inactive local nodes
1937  for (const auto & pr : _primary_secondary_subdomain_id_pairs)
1938  // Loop through all elements on my processor
1939  for (const auto el : _mesh.active_local_subdomain_elements_ptr_range(pr.second))
1940  // If elem fully or partially dropped
1941  if (abs(active_volume[el] / el->volume() - 1.0) > tol)
1942  {
1943  // Add all nodes to list of inactive
1944  for (auto n : make_range(el->n_nodes()))
1945  inactive_node_ids.insert(el->node_id(n));
1946  }
1947 
1948  // Assemble list of procs that nodes contribute to
1949  for (const auto & pr : _primary_secondary_subdomain_id_pairs)
1950  {
1951  const auto secondary_subd_id = pr.second;
1952 
1953  // Loop through all elements not on my processor
1954  for (const auto el : _mesh.active_subdomain_elements_ptr_range(secondary_subd_id))
1955  {
1956  // Get processor_id
1957  const auto pid = el->processor_id();
1958 
1959  // If element is in my subdomain, skip
1960  if (pid == my_pid)
1961  continue;
1962 
1963  // If element on proc pid shares any of my inactive nodes, mark to send
1964  for (const auto n : make_range(el->n_nodes()))
1965  {
1966  const auto node_id = el->node_id(n);
1967  if (inactive_node_ids.find(node_id) != inactive_node_ids.end())
1968  proc_to_inactive_nodes_set[pid].insert(node_id);
1969  }
1970  }
1971  }
1972 
1973  // Send list of inactive nodes
1974  {
1975  // Pack set into vector for sending (push_parallel_vector_data doesn't like sets)
1976  std::unordered_map<processor_id_type, std::vector<dof_id_type>> proc_to_inactive_nodes_vector;
1977  for (const auto & proc_set : proc_to_inactive_nodes_set)
1978  proc_to_inactive_nodes_vector[proc_set.first].insert(
1979  proc_to_inactive_nodes_vector[proc_set.first].end(),
1980  proc_set.second.begin(),
1981  proc_set.second.end());
1982 
1983  // First push data
1984  auto action_functor = [this, &inactive_node_ids](const processor_id_type pid,
1985  const std::vector<dof_id_type> & sent_data)
1986  {
1987  if (pid == _mesh.processor_id())
1988  mooseError("Should not be communicating with self.");
1989  for (const auto pr : sent_data)
1990  inactive_node_ids.insert(pr);
1991  };
1992  TIMPI::push_parallel_vector_data(_mesh.comm(), proc_to_inactive_nodes_vector, action_functor);
1993  }
1994  _inactive_local_lm_nodes.clear();
1995  for (const auto node_id : inactive_node_ids)
1996  _inactive_local_lm_nodes.insert(_mesh.node_ptr(node_id));
1997 }
MetaPhysicL::DualNumber< V, D, asd > abs(const MetaPhysicL::DualNumber< V, D, asd > &a)
Definition: EigenADReal.h:50
void mooseError(Args &&... args)
Emit an error message with the given stringified, concatenated args and terminate the application...
Definition: MooseError.h:311
std::unordered_map< const Elem *, MortarSegmentInfo > _msm_elem_to_info
Map between Elems in the mortar segment mesh and their info structs.
const Parallel::Communicator & comm() const
void push_parallel_vector_data(const Communicator &comm, MapToVectors &&data, const ActionFunctor &act_on_data)
uint8_t processor_id_type
const Elem * secondary_elem
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
Holds xi^(1), xi^(2), and other data for a given mortar segment.
virtual SimpleRange< element_iterator > active_local_subdomain_elements_ptr_range(subdomain_id_type sid)=0
virtual Real volume() const
IntRange< T > make_range(T beg, T end)
unsigned int level ElemType type std::set< subdomain_id_type > ss processor_id_type pid unsigned int level std::set< subdomain_id_type > virtual ss SimpleRange< element_iterator > active_subdomain_elements_ptr_range(subdomain_id_type sid)=0
MeshBase & _mesh
Reference to the mesh stored in equation_systems.
std::vector< std::pair< SubdomainID, SubdomainID > > _primary_secondary_subdomain_id_pairs
A list of primary/secondary subdomain id pairs corresponding to each side of the mortar interface...
virtual const Node * node_ptr(const dof_id_type i) const=0
std::unordered_set< const Node * > _inactive_local_lm_nodes
processor_id_type processor_id() const
std::unique_ptr< MeshBase > _mortar_segment_mesh
1D Mesh of mortar segment elements which gets built by the call to build_mortar_segment_mesh().

◆ computeMsmStatistics()

std::vector< AutomaticMortarGeneration::MsmSubdomainStats > AutomaticMortarGeneration::computeMsmStatistics ( )

Computes mortar segment mesh statistics and returns one entry per subdomain pair.

Must be called collectively on all ranks.

Definition at line 1784 of file AutomaticMortarGeneration.C.

Referenced by msmStatistics().

1785 {
1786  std::vector<MsmSubdomainStats> result;
1787  StatisticsVector<Real> primary;
1788  StatisticsVector<Real> secondary;
1790  std::unordered_map<dof_id_type, Real> primary_elems_to_volume;
1791 
1792  for (const auto & [primary_subd_id, secondary_subd_id] : _primary_secondary_subdomain_id_pairs)
1793  {
1794  for (const auto * const secondary_el :
1795  _mesh.active_local_subdomain_element_ptr_range(secondary_subd_id))
1796  {
1797  secondary.push_back(secondary_el->volume());
1798  // We may not have projected onto a primary face in which case we may not have created mortar
1799  // segments
1800  if (auto it = _secondary_elems_to_mortar_segments.find(secondary_el->id());
1802  for (const auto * const msm_elem : it->second)
1803  {
1804  msm.push_back(msm_elem->volume());
1805  const auto & msm_info = libmesh_map_find(_msm_elem_to_info, msm_elem);
1806  // Now it's also possible that we didn't project onto a primary face and we *did* create
1807  // mortar segments
1808  if (msm_info.primary_elem)
1809  {
1810  if (msm_info.primary_elem->subdomain_id() != primary_subd_id)
1811  mooseError("Unhandled primary-secondary pairing when computing mortar segment "
1812  "statistics. This could happen if you have the same secondary "
1813  "lower-dimensional subdomain ID paired with multiple lower-dimensional "
1814  "primary subdomain IDs. Contact a MOOSE developer for help.");
1815  if (const auto [it, inserted] =
1816  primary_elems_to_volume.emplace(msm_info.primary_elem->id(), Real{});
1817  inserted)
1818  it->second = msm_info.primary_elem->volume();
1819  else
1820  mooseAssert(
1821  MooseUtils::absoluteFuzzyEqual(it->second, msm_info.primary_elem->volume()),
1822  "Volumes should be consistent");
1823  }
1824  }
1825  }
1826 
1827  _mesh.comm().set_union(primary_elems_to_volume);
1828  _mesh.comm().allgather(static_cast<std::vector<Real> &>(secondary));
1829  _mesh.comm().allgather(static_cast<std::vector<Real> &>(msm));
1830  primary.reserve(primary_elems_to_volume.size());
1831  for (const auto [_, volume] : primary_elems_to_volume)
1832  primary.push_back(volume);
1833 
1834  MsmSubdomainStats stats;
1835  stats.primary_subd_id = primary_subd_id;
1836  stats.secondary_subd_id = secondary_subd_id;
1837  stats.secondary_lower_n_elems = secondary.size();
1838  stats.secondary_lower_max_volume = secondary.maximum();
1839  stats.secondary_lower_min_volume = secondary.minimum();
1840  stats.secondary_lower_median_volume = secondary.median();
1841  stats.primary_lower_n_elems = primary.size();
1842  stats.primary_lower_max_volume = primary.maximum();
1843  stats.primary_lower_min_volume = primary.minimum();
1844  stats.primary_lower_median_volume = primary.median();
1845  stats.msm_n_elems = msm.size();
1846  stats.msm_max_volume = msm.maximum();
1847  stats.msm_min_volume = msm.minimum();
1848  stats.msm_median_volume = msm.median();
1849  result.push_back(stats);
1850 
1851  primary.clear();
1852  secondary.clear();
1853  msm.clear();
1854  primary_elems_to_volume.clear();
1855  }
1856 
1857  return result;
1858 }
virtual T maximum() const
void allgather(const T &send_data, std::vector< T, A > &recv_data) const
void mooseError(Args &&... args)
Emit an error message with the given stringified, concatenated args and terminate the application...
Definition: MooseError.h:311
std::unordered_map< dof_id_type, std::set< Elem *, CompareDofObjectsByID > > _secondary_elems_to_mortar_segments
We maintain a mapping from lower-dimensional secondary elements in the original mesh to (sets of) ele...
std::unordered_map< const Elem *, MortarSegmentInfo > _msm_elem_to_info
Map between Elems in the mortar segment mesh and their info structs.
const Parallel::Communicator & comm() const
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
virtual T minimum() const
MeshBase & _mesh
Reference to the mesh stored in equation_systems.
std::vector< std::pair< SubdomainID, SubdomainID > > _primary_secondary_subdomain_id_pairs
A list of primary/secondary subdomain id pairs corresponding to each side of the mortar interface...
void set_union(T &data, const unsigned int root_id) const

◆ computeNodalGeometry()

void AutomaticMortarGeneration::computeNodalGeometry ( )

Computes and stores the nodal normal/tangent vectors in a local data structure instead of using the ExplicitSystem/NumericVector approach.

This design was triggered by the way that the GhostingFunctor operates, but I think it is a better/more efficient way to do it anyway.

The _periodic flag tells us whether we want to inward vs outward facing normals

Definition at line 2132 of file AutomaticMortarGeneration.C.

Referenced by MortarInterfaceWarehouse::update().

2133 {
2134  // The dimension according to Mesh::mesh_dimension().
2135  const auto dim = _mesh.mesh_dimension();
2136 
2137  mooseAssert(dim == 2 || dim == 3,
2138  "AutomaticMortarGeneration::computeNodalGeometry() is only valid for "
2139  "mortar constraints on 2D or 3D meshes.");
2140  // A nodal lower-dimensional nodal quadrature rule to be used on faces.
2141  QNodal qface(dim - 1);
2142 
2143  // A map from the node id to the attached elemental normals/weights evaluated at the node. Th
2144  // length of the vector will correspond to the number of elements attached to the node. If it is a
2145  // vertex node, for a 1D mortar mesh, the vector length will be two. If it is an interior node,
2146  // the vector will be length 1. The first member of the pair is that element's normal at the node.
2147  // The second member is that element's JxW at the node
2148  std::map<dof_id_type, std::vector<std::pair<Point, Real>>> node_to_normals_map;
2149 
2151  Real sign = _periodic ? -1 : 1;
2152 
2153  // First loop over lower-dimensional secondary side elements and compute/save the outward normal
2154  // for each one. We loop over all active elements currently, but this procedure could be
2155  // parallelized as well.
2156  for (MeshBase::const_element_iterator el = _mesh.active_elements_begin(),
2157  end_el = _mesh.active_elements_end();
2158  el != end_el;
2159  ++el)
2160  {
2161  const Elem * secondary_elem = *el;
2162 
2163  // If this is not one of the lower-dimensional secondary side elements, go on to the next one.
2164  if (!_secondary_boundary_subdomain_ids.count(secondary_elem->subdomain_id()))
2165  continue;
2166 
2167  // We will create an FE object and attach the nodal quadrature rule such that we can get out the
2168  // normals at the element nodes
2169  FEType nnx_fe_type(secondary_elem->default_order(), LAGRANGE);
2170  std::unique_ptr<FEBase> nnx_fe_face(FEBase::build(dim, nnx_fe_type));
2171  nnx_fe_face->attach_quadrature_rule(&qface);
2172  const auto & face_normals = nnx_fe_face->get_normals();
2173  const auto & face_points = nnx_fe_face->get_xyz();
2174 
2175  const auto & JxW = nnx_fe_face->get_JxW();
2176 
2177  // Which side of the parent are we? We need to know this to know
2178  // which side to reinit.
2179  const Elem * interior_parent = secondary_elem->interior_parent();
2180  mooseAssert(interior_parent,
2181  "No interior parent exists for element "
2182  << secondary_elem->id()
2183  << ". There may be a problem with your sideset set-up.");
2184 
2185  // Map to get lower dimensional element from interior parent on secondary surface
2186  // This map can be used to provide a handle to methods in this class that need to
2187  // operate on lower dimensional elements.
2188  _secondary_element_to_secondary_lowerd_element.emplace(interior_parent->id(), secondary_elem);
2189 
2190  // Look up which side of the interior parent secondary_elem is.
2191  auto s = interior_parent->which_side_am_i(secondary_elem);
2192 
2193  // Reinit the face FE object on side s.
2194  nnx_fe_face->reinit(interior_parent, s);
2195 
2196  // Match by physical location instead of assuming that parent-side nodal
2197  // quadrature ordering and lower-dimensional side-element node ordering are
2198  // identical.
2199  const auto qpoint_to_secondary_node =
2200  nodalQuadraturePointToSecondaryNodeMap(*secondary_elem, face_points);
2201 
2202  mooseAssert(face_normals.size() == face_points.size() && JxW.size() == face_points.size(),
2203  "Face nodal geometry vectors must have the same size.");
2204 
2205  for (const auto qp : make_range(face_points.size()))
2206  {
2207  const auto n = qpoint_to_secondary_node[qp];
2208  auto & normals_and_weights_vec = node_to_normals_map[secondary_elem->node_id(n)];
2209  normals_and_weights_vec.push_back(std::make_pair(sign * face_normals[qp], JxW[qp]));
2210  }
2211  }
2212 
2213  for (const auto & pr : node_to_normals_map)
2214  {
2215  // Compute normal vector
2216  const auto & node_id = pr.first;
2217  const auto & normals_and_weights_vec = pr.second;
2218 
2219  Point nodal_normal;
2220  for (const auto & norm_and_weight : normals_and_weights_vec)
2221  nodal_normal += norm_and_weight.first * norm_and_weight.second;
2222  nodal_normal = nodal_normal.unit();
2223 
2224  _secondary_node_to_nodal_normal[_mesh.node_ptr(node_id)] = nodal_normal;
2225 
2226  Point nodal_tangent_one;
2227  Point nodal_tangent_two;
2228  householderOrthogolization(nodal_normal, nodal_tangent_one, nodal_tangent_two);
2229 
2230  _secondary_node_to_hh_nodal_tangents[_mesh.node_ptr(node_id)][0] = nodal_tangent_one;
2231  _secondary_node_to_hh_nodal_tangents[_mesh.node_ptr(node_id)][1] = nodal_tangent_two;
2232  }
2233 }
LAGRANGE
std::unique_ptr< FEGenericBase< Real > > build(const unsigned int dim, const FEType &fet)
const Elem * interior_parent() const
std::set< SubdomainID > _secondary_boundary_subdomain_ids
The secondary/primary lower-dimensional boundary subdomain ids are the secondary/primary boundary ids...
const bool _periodic
Whether this object will be generating a mortar segment mesh for periodic constraints.
unsigned int which_side_am_i(const Elem *e) const
TypeVector< Real > unit() const
dof_id_type id() const
std::unordered_map< dof_id_type, const Elem * > _secondary_element_to_secondary_lowerd_element
Map from full dimensional secondary element id to lower dimensional secondary element.
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
subdomain_id_type subdomain_id() const
std::unordered_map< const Node *, std::array< Point, 2 > > _secondary_node_to_hh_nodal_tangents
Container for storing the nodal tangent/binormal vectors associated with each secondary node (Househo...
IntRange< T > make_range(T beg, T end)
KOKKOS_INLINE_FUNCTION T sign(T x)
Returns the sign of a value.
Definition: KokkosUtils.h:26
unsigned int mesh_dimension() const
MeshBase & _mesh
Reference to the mesh stored in equation_systems.
std::unordered_map< const Node *, Point > _secondary_node_to_nodal_normal
Container for storing the nodal normal vector associated with each secondary node.
virtual const Node * node_ptr(const dof_id_type i) const=0
virtual Order default_order() const=0
void householderOrthogolization(const Point &normal, Point &tangent_one, Point &tangent_two) const
Householder orthogonalization procedure to obtain proper basis for tangent and binormal vectors...
dof_id_type node_id(const unsigned int i) const

◆ dim()

int AutomaticMortarGeneration::dim ( ) const
inline

◆ getInactiveLMElems()

const std::unordered_set<const Elem *>& AutomaticMortarGeneration::getInactiveLMElems ( ) const
inline
Returns
The list of secondary elems on which mortar constraint is not active

Definition at line 350 of file AutomaticMortarGeneration.h.

Referenced by ComputeMortarFunctor::operator()().

351  {
353  }
std::unordered_set< const Elem * > _inactive_local_lm_elems
List of inactive lagrange multiplier nodes (for elemental variables)

◆ getInactiveLMNodes()

const std::unordered_set<const Node *>& AutomaticMortarGeneration::getInactiveLMNodes ( ) const
inline
Returns
The set of nodes on which mortar constraints are not active

Definition at line 342 of file AutomaticMortarGeneration.h.

Referenced by ComputeMortarFunctor::operator()().

343  {
345  }
std::unordered_set< const Node * > _inactive_local_lm_nodes

◆ getNodalNormals()

std::vector< Point > AutomaticMortarGeneration::getNodalNormals ( const Elem secondary_elem) const
Returns
The nodal normals associated with the provided secondary_elem

Definition at line 543 of file AutomaticMortarGeneration.C.

Referenced by buildMortarSegmentMesh3d(), getNormals(), and MortarConsumerInterface::setNormals().

544 {
545  std::vector<Point> nodal_normals(secondary_elem.n_nodes());
546  for (const auto n : make_range(secondary_elem.n_nodes()))
547  nodal_normals[n] = _secondary_node_to_nodal_normal.at(secondary_elem.node_ptr(n));
548 
549  return nodal_normals;
550 }
virtual unsigned int n_nodes() const=0
const Node * node_ptr(const unsigned int i) const
IntRange< T > make_range(T beg, T end)
std::unordered_map< const Node *, Point > _secondary_node_to_nodal_normal
Container for storing the nodal normal vector associated with each secondary node.

◆ getNodalTangents()

std::array< MooseUtils::SemidynamicVector< Point, 9 >, 2 > AutomaticMortarGeneration::getNodalTangents ( const Elem secondary_elem) const

Compute the two nodal tangents, which are built on-the-fly.

Returns
The nodal tangents associated with the provided secondary_elem

Definition at line 596 of file AutomaticMortarGeneration.C.

597 {
598  // MetaPhysicL will check if we ran out of allocated space.
599  MooseUtils::SemidynamicVector<Point, 9> nodal_tangents_one(0);
600  MooseUtils::SemidynamicVector<Point, 9> nodal_tangents_two(0);
601 
602  for (const auto n : make_range(secondary_elem.n_nodes()))
603  {
604  const auto & tangent_vectors =
605  libmesh_map_find(_secondary_node_to_hh_nodal_tangents, secondary_elem.node_ptr(n));
606  nodal_tangents_one.push_back(tangent_vectors[0]);
607  nodal_tangents_two.push_back(tangent_vectors[1]);
608  }
609 
610  return {{nodal_tangents_one, nodal_tangents_two}};
611 }
virtual unsigned int n_nodes() const=0
const Node * node_ptr(const unsigned int i) const
std::unordered_map< const Node *, std::array< Point, 2 > > _secondary_node_to_hh_nodal_tangents
Container for storing the nodal tangent/binormal vectors associated with each secondary node (Househo...
IntRange< T > make_range(T beg, T end)

◆ getNormals() [1/2]

std::vector< Point > AutomaticMortarGeneration::getNormals ( const Elem secondary_elem,
const std::vector< Point > &  xi1_pts 
) const

Compute the normals at given reference points on a secondary element.

Parameters
secondary_elemThe secondary element used to query for associated nodal normals
xi1_ptsThe reference points on the secondary element to evaluate the normals at. The points should only be non-zero in the zeroth entry because right now our mortar mesh elements are always 1D
Returns
The normals

Definition at line 625 of file AutomaticMortarGeneration.C.

Referenced by buildMortarSegmentMesh(), getNormals(), and MortarConsumerInterface::setNormals().

627 {
628  const auto mortar_dim = _mesh.mesh_dimension() - 1;
629  const auto num_qps = xi1_pts.size();
630  const auto nodal_normals = getNodalNormals(secondary_elem);
631  std::vector<Point> normals(num_qps);
632 
633  for (const auto n : make_range(secondary_elem.n_nodes()))
634  for (const auto qp : make_range(num_qps))
635  {
636  const auto phi =
637  (mortar_dim == 1)
638  ? Moose::fe_lagrange_1D_shape(secondary_elem.default_order(), n, xi1_pts[qp](0))
639  : Moose::fe_lagrange_2D_shape(secondary_elem.type(),
640  secondary_elem.default_order(),
641  n,
642  static_cast<const TypeVector<Real> &>(xi1_pts[qp]));
643  normals[qp] += phi * nodal_normals[n];
644  }
645 
646  if (_periodic)
647  for (auto & normal : normals)
648  normal *= -1;
649 
650  return normals;
651 }
T fe_lagrange_2D_shape(const libMesh::ElemType type, const Order order, const unsigned int i, const VectorType< T > &p)
const bool _periodic
Whether this object will be generating a mortar segment mesh for periodic constraints.
std::vector< Point > getNodalNormals(const Elem &secondary_elem) const
virtual unsigned int n_nodes() const=0
IntRange< T > make_range(T beg, T end)
unsigned int mesh_dimension() const
T fe_lagrange_1D_shape(const Order order, const unsigned int i, const T &xi)
MeshBase & _mesh
Reference to the mesh stored in equation_systems.
virtual Order default_order() const=0
virtual ElemType type() const=0

◆ getNormals() [2/2]

std::vector< Point > AutomaticMortarGeneration::getNormals ( const Elem secondary_elem,
const std::vector< Real > &  oned_xi1_pts 
) const

Compute the normals at given reference points on a secondary element.

Parameters
secondary_elemThe secondary element used to query for associated nodal normals
1d_xi1_ptsThe reference points on the secondary element to evaluate the normals at. The "points" are single reals corresponding to xi because right now our mortar mesh elements are always 1D
Returns
The normals

Definition at line 614 of file AutomaticMortarGeneration.C.

616 {
617  std::vector<Point> xi1_pts(oned_xi1_pts.size());
618  for (const auto qp : index_range(oned_xi1_pts))
619  xi1_pts[qp] = oned_xi1_pts[qp];
620 
621  return getNormals(secondary_elem, xi1_pts);
622 }
std::vector< Point > getNormals(const Elem &secondary_elem, const std::vector< Point > &xi1_pts) const
Compute the normals at given reference points on a secondary element.
auto index_range(const T &sizable)

◆ getPrimaryIpToLowerElementMap()

std::map< unsigned int, unsigned int > AutomaticMortarGeneration::getPrimaryIpToLowerElementMap ( const Elem primary_elem,
const Elem primary_elem_ip,
const Elem lower_secondary_elem 
) const

Compute on-the-fly mapping from primary interior parent nodes to its corresponding lower dimensional nodes.

Returns
The map from primary interior parent nodes to its corresponding lower dimensional nodes

Definition at line 579 of file AutomaticMortarGeneration.C.

583 {
584  std::map<unsigned int, unsigned int> primary_ip_i_to_lower_primary_i;
585 
586  for (const auto i : make_range(lower_primary_elem.n_nodes()))
587  {
588  const auto & nd = lower_primary_elem.node_ref(i);
589  primary_ip_i_to_lower_primary_i[primary_elem.get_node_index(&nd)] = i;
590  }
591 
592  return primary_ip_i_to_lower_primary_i;
593 }
unsigned int get_node_index(const Node *node_ptr) const
IntRange< T > make_range(T beg, T end)

◆ getSecondaryIpToLowerElementMap()

std::map< unsigned int, unsigned int > AutomaticMortarGeneration::getSecondaryIpToLowerElementMap ( const Elem lower_secondary_elem) const

Compute on-the-fly mapping from secondary interior parent nodes to lower dimensional nodes.

Returns
The map from secondary interior parent nodes to lower dimensional nodes

Definition at line 563 of file AutomaticMortarGeneration.C.

564 {
565  std::map<unsigned int, unsigned int> secondary_ip_i_to_lower_secondary_i;
566  const Elem * const secondary_ip = lower_secondary_elem.interior_parent();
567  mooseAssert(secondary_ip, "This should be non-null");
568 
569  for (const auto i : make_range(lower_secondary_elem.n_nodes()))
570  {
571  const auto & nd = lower_secondary_elem.node_ref(i);
572  secondary_ip_i_to_lower_secondary_i[secondary_ip->get_node_index(&nd)] = i;
573  }
574 
575  return secondary_ip_i_to_lower_secondary_i;
576 }
unsigned int get_node_index(const Node *node_ptr) const
const Elem * interior_parent() const
const Node & node_ref(const unsigned int i) const
virtual unsigned int n_nodes() const=0
IntRange< T > make_range(T beg, T end)

◆ getSecondaryLowerdElemFromSecondaryElem()

const Elem * AutomaticMortarGeneration::getSecondaryLowerdElemFromSecondaryElem ( dof_id_type  secondary_elem_id) const

Return lower dimensional secondary element given its interior parent.

Helpful outside the mortar generation to locate mortar-related quantities.

Parameters
secondary_elem_idThe secondary interior parent element id used to query for associated lower dimensional element
Returns
The corresponding lower dimensional secondary element

Definition at line 553 of file AutomaticMortarGeneration.C.

555 {
556  mooseAssert(_secondary_element_to_secondary_lowerd_element.count(secondary_elem_id),
557  "Map should locate secondary element");
558 
559  return _secondary_element_to_secondary_lowerd_element.at(secondary_elem_id);
560 }
std::unordered_map< dof_id_type, const Elem * > _secondary_element_to_secondary_lowerd_element
Map from full dimensional secondary element id to lower dimensional secondary element.

◆ householderOrthogolization()

void AutomaticMortarGeneration::householderOrthogolization ( const Point normal,
Point tangent_one,
Point tangent_two 
) const
private

Householder orthogonalization procedure to obtain proper basis for tangent and binormal vectors.

Definition at line 2236 of file AutomaticMortarGeneration.C.

Referenced by computeNodalGeometry().

2239 {
2240  using std::abs;
2241 
2242  mooseAssert(MooseUtils::absoluteFuzzyEqual(nodal_normal.norm(), 1),
2243  "The input nodal normal should have unity norm");
2244 
2245  const Real nx = nodal_normal(0);
2246  const Real ny = nodal_normal(1);
2247  const Real nz = nodal_normal(2);
2248 
2249  // See Lopes DS, Silva MT, Ambrosio JA. Tangent vectors to a 3-D surface normal: A geometric tool
2250  // to find orthogonal vectors based on the Householder transformation. Computer-Aided Design. 2013
2251  // Mar 1;45(3):683-94. We choose one definition of h_vector and deal with special case.
2252  const Point h_vector(nx + 1.0, ny, nz);
2253 
2254  // Avoid singularity of the equations at the end of routine by providing the solution to
2255  // (nx,ny,nz)=(-1,0,0) Normal/tangent fields can be visualized by outputting nodal geometry mesh
2256  // on a spherical problem.
2257  if (abs(h_vector(0)) < TOLERANCE)
2258  {
2259  nodal_tangent_one(0) = 0;
2260  nodal_tangent_one(1) = 1;
2261  nodal_tangent_one(2) = 0;
2262 
2263  nodal_tangent_two(0) = 0;
2264  nodal_tangent_two(1) = 0;
2265  nodal_tangent_two(2) = -1;
2266 
2267  return;
2268  }
2269 
2270  const Real h = h_vector.norm();
2271 
2272  nodal_tangent_one(0) = -2.0 * h_vector(0) * h_vector(1) / (h * h);
2273  nodal_tangent_one(1) = 1.0 - 2.0 * h_vector(1) * h_vector(1) / (h * h);
2274  nodal_tangent_one(2) = -2.0 * h_vector(1) * h_vector(2) / (h * h);
2275 
2276  nodal_tangent_two(0) = -2.0 * h_vector(0) * h_vector(2) / (h * h);
2277  nodal_tangent_two(1) = -2.0 * h_vector(1) * h_vector(2) / (h * h);
2278  nodal_tangent_two(2) = 1.0 - 2.0 * h_vector(2) * h_vector(2) / (h * h);
2279 }
MetaPhysicL::DualNumber< V, D, asd > abs(const MetaPhysicL::DualNumber< V, D, asd > &a)
Definition: EigenADReal.h:50
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real

◆ incorrectEdgeDropping()

bool AutomaticMortarGeneration::incorrectEdgeDropping ( ) const
inline

Definition at line 355 of file AutomaticMortarGeneration.h.

Referenced by ComputeMortarFunctor::operator()().

355 { return !_correct_edge_dropping; }
const bool _correct_edge_dropping
Flag to enable regressed treatment of edge dropping where all LM DoFs on edge dropping element are st...

◆ initOutput()

void AutomaticMortarGeneration::initOutput ( )

initialize mortar-mesh based output

Definition at line 450 of file AutomaticMortarGeneration.C.

451 {
452  if (!_debug)
453  return;
454 
455  _output_params = std::make_unique<InputParameters>(MortarNodalGeometryOutput::validParams());
456  _output_params->set<AutomaticMortarGeneration *>("_amg") = this;
457  _output_params->set<FEProblemBase *>("_fe_problem_base") = &_app.feProblem();
459  _output_params->set<std::string>(MooseBase::name_param) =
460  "mortar_nodal_geometry_" + mortarInterfaceName();
461  _output_params->finalize("MortarNodalGeometryOutput");
462  _app.getOutputWarehouse().addOutput(std::make_shared<MortarNodalGeometryOutput>(*_output_params));
463 }
static const std::string name_param
The name of the parameter that contains the object name.
Definition: MooseBase.h:55
static const std::string app_param
The name of the parameter that contains the MooseApp.
Definition: MooseBase.h:59
Base class for MOOSE-based applications.
Definition: MooseApp.h:109
Specialization of SubProblem for solving nonlinear equations plus auxiliary equations.
FEProblemBase & feProblem() const
Definition: MooseApp.C:1852
This class is a container/interface for the objects involved in automatic generation of mortar spaces...
const bool _debug
Whether to print debug output.
std::unique_ptr< InputParameters > _output_params
Storage for the input parameters used by the mortar nodal geometry output.
static InputParameters validParams()
void addOutput(std::shared_ptr< Output > output)
Adds an existing output object to the warehouse.
MooseApp & _app
The Moose app.
OutputWarehouse & getOutputWarehouse()
Get the OutputWarehouse objects.
Definition: MooseApp.C:2409

◆ meshChanged()

void AutomaticMortarGeneration::meshChanged ( )
inline

Invalidates the cached MSM node/element ID starting offset so that the next call to buildMortarSegmentMesh3d() recomputes it via allgather.

Call this when mesh topology changes.

Definition at line 217 of file AutomaticMortarGeneration.h.

217 { _msm_node_id_start = std::nullopt; }
std::optional< dof_id_type > _msm_node_id_start
Cached per-rank starting ID for 3D MSM nodes/elements.

◆ mortar3DQpMapping()

Mortar3DQuadraturePointMapping AutomaticMortarGeneration::mortar3DQpMapping ( ) const
inline

Return the 3D mortar quadrature-point mapping method.

Definition at line 333 of file AutomaticMortarGeneration.h.

Referenced by Moose::Mortar::loopOverMortarSegments().

333 { return _mortar_3d_qp_mapping; }
const Mortar3DQuadraturePointMapping _mortar_3d_qp_mapping
Method used to map 3D mortar segment quadrature points to primary and secondary faces.

◆ mortarInterfaceCoupling()

const std::unordered_map<dof_id_type, std::unordered_set<dof_id_type> >& AutomaticMortarGeneration::mortarInterfaceCoupling ( ) const
inline
Returns
The mortar interface coupling

Definition at line 307 of file AutomaticMortarGeneration.h.

308  {
310  }
std::unordered_map< dof_id_type, std::unordered_set< dof_id_type > > _mortar_interface_coupling
Used by the AugmentSparsityOnInterface functor to determine whether a given Elem is coupled to any ot...

◆ mortarInterfaceName()

std::string AutomaticMortarGeneration::mortarInterfaceName ( ) const
private
Returns
A string uniquely identifying this mortar interface

Definition at line 436 of file AutomaticMortarGeneration.C.

Referenced by initOutput(), and outputMortarMesh().

437 {
438  std::vector<std::string> string_vec(_primary_secondary_boundary_id_pairs.size() * 2 + 1);
440  {
441  const auto [primary_bnd_id, secondary_bnd_id] = _primary_secondary_boundary_id_pairs[i];
442  string_vec[2 * i] = std::to_string(primary_bnd_id);
443  string_vec[2 * i + 1] = std::to_string(secondary_bnd_id);
444  }
445  string_vec.back() = _on_displaced ? "displaced" : "undisplaced";
446  return MooseUtils::join(string_vec, "_");
447 }
std::vector< std::pair< BoundaryID, BoundaryID > > _primary_secondary_boundary_id_pairs
A list of primary/secondary boundary id pairs corresponding to each side of the mortar interface...
const bool _on_displaced
Whether this object is on the displaced mesh.
auto index_range(const T &sizable)

◆ mortarSegmentMesh()

const MeshBase& AutomaticMortarGeneration::mortarSegmentMesh ( ) const
inline
Returns
The mortar segment mesh

Definition at line 320 of file AutomaticMortarGeneration.h.

Referenced by Moose::Mortar::loopOverMortarSegments().

320 { return *_mortar_segment_mesh; }
std::unique_ptr< MeshBase > _mortar_segment_mesh
1D Mesh of mortar segment elements which gets built by the call to build_mortar_segment_mesh().

◆ mortarSegmentMeshElemToInfo()

const std::unordered_map<const Elem *, MortarSegmentInfo>& AutomaticMortarGeneration::mortarSegmentMeshElemToInfo ( ) const
inline
Returns
The mortar segment element to corresponding information

Definition at line 325 of file AutomaticMortarGeneration.h.

Referenced by Moose::Mortar::loopOverMortarSegments().

326  {
327  return _msm_elem_to_info;
328  }
std::unordered_map< const Elem *, MortarSegmentInfo > _msm_elem_to_info
Map between Elems in the mortar segment mesh and their info structs.

◆ mortarSegmentReferencePoints()

const MortarSegmentReferencePoints & AutomaticMortarGeneration::mortarSegmentReferencePoints ( const Elem mortar_segment_elem) const

Return the parent-face reference coordinates for a mortar segment.

Definition at line 488 of file AutomaticMortarGeneration.C.

Referenced by Moose::Mortar::loopOverMortarSegments().

489 {
490  if (_mortar_3d_qp_mapping != Mortar3DQuadraturePointMapping::REFERENCE_INTERPOLATION)
491  mooseError("Mortar segment reference points were requested for mortar segment element ",
492  mortar_segment_elem.id(),
493  ", but the reference-interpolation mapping mode is not enabled.");
494 
495  const auto reference_points_it = _msm_elem_to_reference_points.find(&mortar_segment_elem);
496  if (reference_points_it == _msm_elem_to_reference_points.end())
497  mooseError("No reference-point record was found for mortar segment element ",
498  mortar_segment_elem.id(),
499  ". The mortar segment info and reference-point maps are not aligned.");
500 
501  return reference_points_it->second;
502 }
void mooseError(Args &&... args)
Emit an error message with the given stringified, concatenated args and terminate the application...
Definition: MooseError.h:311
dof_id_type id() const
std::unordered_map< const Elem *, MortarSegmentReferencePoints > _msm_elem_to_reference_points
Reference-coordinate data used only by the reference-interpolation mapping mode.
const Mortar3DQuadraturePointMapping _mortar_3d_qp_mapping
Method used to map 3D mortar segment quadrature points to primary and secondary faces.

◆ msmStatistics()

void AutomaticMortarGeneration::msmStatistics ( )

Prints mortar segment mesh statistics to console (calls computeMsmStatistics internally)

Definition at line 1861 of file AutomaticMortarGeneration.C.

Referenced by buildMortarSegmentMesh3d().

1862 {
1863  const auto all_stats = computeMsmStatistics();
1864 
1865  if (_mesh.processor_id() != 0)
1866  return;
1867 
1868  Moose::out << "Mortar Interface Statistics:" << std::endl;
1869  for (const auto & stats : all_stats)
1870  {
1871  std::vector<std::string> col_names = {"mesh", "n_elems", "max", "min", "median"};
1872  std::vector<std::string> subds = {"secondary_lower", "primary_lower", "mortar_segment"};
1873  std::vector<size_t> n_elems = {
1874  stats.secondary_lower_n_elems, stats.primary_lower_n_elems, stats.msm_n_elems};
1875  std::vector<Real> maxs = {
1876  stats.secondary_lower_max_volume, stats.primary_lower_max_volume, stats.msm_max_volume};
1877  std::vector<Real> mins = {
1878  stats.secondary_lower_min_volume, stats.primary_lower_min_volume, stats.msm_min_volume};
1879  std::vector<Real> medians = {stats.secondary_lower_median_volume,
1880  stats.primary_lower_median_volume,
1881  stats.msm_median_volume};
1882 
1883  FormattedTable table;
1884  table.clear();
1885  for (auto i : index_range(subds))
1886  {
1887  table.addRow(i);
1888  table.addData<std::string>(col_names[0], subds[i]);
1889  table.addData<size_t>(col_names[1], n_elems[i]);
1890  table.addData<Real>(col_names[2], maxs[i]);
1891  table.addData<Real>(col_names[3], mins[i]);
1892  table.addData<Real>(col_names[4], medians[i]);
1893  }
1894 
1895  Moose::out << "secondary subdomain: " << stats.secondary_subd_id
1896  << " \tprimary subdomain: " << stats.primary_subd_id << std::endl;
1897  table.printTable(Moose::out, subds.size());
1898  }
1899 }
void addData(const std::string &name, const T &value)
Method for adding data to the output table.
This class is used for building, formatting, and outputting tables of numbers.
void printTable(std::ostream &out, unsigned int last_n_entries=0)
Methods for dumping the table to the stream - either by filename or by stream handle.
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
std::vector< MsmSubdomainStats > computeMsmStatistics()
Computes mortar segment mesh statistics and returns one entry per subdomain pair. ...
MeshBase & _mesh
Reference to the mesh stored in equation_systems.
void addRow(Real time)
Force a new row in the table with the passed in time.
processor_id_type processor_id() const
auto index_range(const T &sizable)

◆ nodesToSecondaryElem()

const std::unordered_map<dof_id_type, std::vector<const Elem *> >& AutomaticMortarGeneration::nodesToSecondaryElem ( ) const
inline
Returns
Map from node id to secondary lower-d element pointer

Definition at line 391 of file AutomaticMortarGeneration.h.

392  {
394  }
std::unordered_map< dof_id_type, std::vector< const Elem * > > _nodes_to_secondary_elem_map
Map from nodes to connected lower-dimensional elements on the secondary/primary subdomains.

◆ onDisplaced()

bool AutomaticMortarGeneration::onDisplaced ( ) const
inline

returns whether this object is on the displaced mesh

Definition at line 222 of file AutomaticMortarGeneration.h.

222 { return _on_displaced; }
const bool _on_displaced
Whether this object is on the displaced mesh.

◆ outputMortarMesh()

void AutomaticMortarGeneration::outputMortarMesh ( )
private

Write the mortar segment mesh to exodus.

Definition at line 1146 of file AutomaticMortarGeneration.C.

Referenced by buildMortarSegmentMesh(), and buildMortarSegmentMesh3d().

1147 {
1148  ExodusII_IO mortar_segment_mesh_writer(*_mortar_segment_mesh);
1149 
1150  // Default to non-HDF5 output for wider compatibility
1151  mortar_segment_mesh_writer.set_hdf5_writing(false);
1152 
1153  std::array<std::string, 3> file_pieces = {
1154  _app.getOutputFileBase(/*for_non_moose_build_output=*/true),
1156  "mortar_segment_mesh.e"};
1157  mortar_segment_mesh_writer.write(MooseUtils::join(file_pieces, "_"));
1158 }
std::string getOutputFileBase(bool for_non_moose_build_output=false) const
Get the output file base name.
Definition: MooseApp.C:1532
MooseApp & _app
The Moose app.
std::unique_ptr< MeshBase > _mortar_segment_mesh
1D Mesh of mortar segment elements which gets built by the call to build_mortar_segment_mesh().

◆ primaryIPSubIDs()

const std::set<SubdomainID>& AutomaticMortarGeneration::primaryIPSubIDs ( ) const
inline
Returns
All the primary interior parent subdomain IDs associated with the mortar mesh

Definition at line 386 of file AutomaticMortarGeneration.h.

Referenced by Moose::Mortar::setupMortarMaterials().

386 { return _primary_ip_sub_ids; }
std::set< SubdomainID > _primary_ip_sub_ids
All the primary interior parent subdomain IDs associated with the mortar mesh.

◆ primarySecondaryBoundaryIDPair()

const std::pair< BoundaryID, BoundaryID > & AutomaticMortarGeneration::primarySecondaryBoundaryIDPair ( ) const
inline
Returns
The primary-secondary boundary ID pair

Definition at line 622 of file AutomaticMortarGeneration.h.

Referenced by Moose::Mortar::loopOverMortarSegments(), and Moose::Mortar::setupMortarMaterials().

623 {
624  mooseAssert(_primary_secondary_boundary_id_pairs.size() == 1,
625  "We currently only support a single boundary pair per mortar generation object");
626 
628 }
std::vector< std::pair< BoundaryID, BoundaryID > > _primary_secondary_boundary_id_pairs
A list of primary/secondary boundary id pairs corresponding to each side of the mortar interface...

◆ processAlignedNodes()

bool AutomaticMortarGeneration::processAlignedNodes ( const Node secondary_node,
const Node primary_node,
const std::vector< const Elem *> *  secondary_node_neighbors,
const std::vector< const Elem *> *  primary_node_neighbors,
const VectorValue< Real > &  nodal_normal,
const Elem candidate_element,
std::set< const Elem *> &  rejected_element_candidates 
)
private

Process aligned nodes.

Returns
whether mortar segment(s) were created

Definition at line 2293 of file AutomaticMortarGeneration.C.

Referenced by projectPrimaryNodesSinglePair(), and projectSecondaryNodesSinglePair().

2301 {
2302  if (!secondary_node_neighbors)
2303  secondary_node_neighbors = &libmesh_map_find(_nodes_to_secondary_elem_map, secondary_node.id());
2304  if (!primary_node_neighbors)
2305  primary_node_neighbors = &libmesh_map_find(_nodes_to_primary_elem_map, primary_node.id());
2306 
2307  std::vector<bool> primary_elems_mapped(primary_node_neighbors->size(), false);
2308 
2309  // Add entries to secondary_node_and_elem_to_xi2_primary_elem container.
2310  //
2311  // First, determine "on left" vs. "on right" orientation of the nodal neighbors.
2312  // There can be a max of 2 nodal neighbors, and we want to make sure that the
2313  // secondary nodal neighbor on the "left" is associated with the primary nodal
2314  // neighbor on the "left" and similarly for the "right". We use cross products to determine
2315  // alignment. In the below diagram, 'x' denotes a node, and connected '|' are lower dimensional
2316  // elements.
2317  // x
2318  // x |
2319  // | |
2320  // secondary x ----> x primary
2321  // | |
2322  // | x
2323  // x
2324  //
2325  // Looking at the aligned nodes, the secondary node first, if we pick the top secondary lower
2326  // dimensional element, then the cross product as written a few lines below points out of the
2327  // screen towards you. (Point in the direction of the secondary nodal normal, and then curl your
2328  // hand towards the secondary element's opposite node, then the thumb points in the direction of
2329  // the cross product). Doing the same with the aligned primary node, if we pick the top primary
2330  // element, then the cross product also points out of the screen. Because the cross products
2331  // point in the same direction (positive dot product), then we know to associate the
2332  // secondary-primary element pair. If we had picked the bottom primary element whose cross
2333  // product points into the screen, then clearly the cross products point in the opposite
2334  // direction and we don't have a match
2335  std::array<Real, 2> secondary_node_neighbor_cps, primary_node_neighbor_cps;
2336 
2337  for (const auto nn : index_range(*secondary_node_neighbors))
2338  {
2339  const Elem * const secondary_neigh = (*secondary_node_neighbors)[nn];
2340  const Point opposite = (secondary_neigh->node_ptr(0) == &secondary_node)
2341  ? secondary_neigh->point(1)
2342  : secondary_neigh->point(0);
2343  const Point cp = nodal_normal.cross(opposite - secondary_node);
2344  secondary_node_neighbor_cps[nn] = cp(2);
2345  }
2346 
2347  for (const auto nn : index_range(*primary_node_neighbors))
2348  {
2349  const Elem * const primary_neigh = (*primary_node_neighbors)[nn];
2350  const Point opposite = (primary_neigh->node_ptr(0) == &primary_node) ? primary_neigh->point(1)
2351  : primary_neigh->point(0);
2352  const Point cp = nodal_normal.cross(opposite - primary_node);
2353  primary_node_neighbor_cps[nn] = cp(2);
2354  }
2355 
2356  // Associate secondary/primary elems on matching sides.
2357  bool found_match = false;
2358  for (const auto snn : index_range(*secondary_node_neighbors))
2359  for (const auto mnn : index_range(*primary_node_neighbors))
2360  if (secondary_node_neighbor_cps[snn] * primary_node_neighbor_cps[mnn] > 0)
2361  {
2362  found_match = true;
2363  if (primary_elems_mapped[mnn])
2364  continue;
2365  primary_elems_mapped[mnn] = true;
2366 
2367  // Figure out xi^(2) value by looking at which node primary_node is
2368  // of the current primary node neighbor.
2369  const Real xi2 = (&primary_node == (*primary_node_neighbors)[mnn]->node_ptr(0)) ? -1 : +1;
2370  const auto secondary_key =
2371  std::make_pair(&secondary_node, (*secondary_node_neighbors)[snn]);
2372  const auto primary_val = std::make_pair(xi2, (*primary_node_neighbors)[mnn]);
2373  _secondary_node_and_elem_to_xi2_primary_elem.emplace(secondary_key, primary_val);
2374 
2375  // Also map in the other direction.
2376  const Real xi1 =
2377  (&secondary_node == (*secondary_node_neighbors)[snn]->node_ptr(0)) ? -1 : +1;
2378 
2379  const auto primary_key =
2380  std::make_tuple(primary_node.id(), &primary_node, (*primary_node_neighbors)[mnn]);
2381  const auto secondary_val = std::make_pair(xi1, (*secondary_node_neighbors)[snn]);
2382  _primary_node_and_elem_to_xi1_secondary_elem.emplace(primary_key, secondary_val);
2383  }
2384 
2385  if (!found_match)
2386  {
2387  // There could be coincident nodes and this might be a bad primary candidate (see
2388  // issue #21680). Instead of giving up, let's try continuing
2389  rejected_elem_candidates.insert(&candidate_element);
2390  return false;
2391  }
2392 
2393  // We need to handle the case where we've exactly projected a secondary node onto a
2394  // primary node, but our secondary node is at one of the secondary boundary face endpoints and
2395  // our primary node is not.
2396  if (secondary_node_neighbors->size() == 1 && primary_node_neighbors->size() == 2)
2397  for (const auto i : index_range(primary_elems_mapped))
2398  if (!primary_elems_mapped[i])
2399  {
2401  std::make_tuple(primary_node.id(), &primary_node, (*primary_node_neighbors)[i]),
2402  std::make_pair(1, nullptr));
2403  }
2404 
2405  return found_match;
2406 }
dof_id_type id() const
std::map< std::tuple< dof_id_type, const Node *, const Elem * >, std::pair< Real, const Elem * > > _primary_node_and_elem_to_xi1_secondary_elem
Same type of container, but for mapping (Primary Node ID, Primary Node, Primary Elem) -> (xi^(1)...
std::unordered_map< dof_id_type, std::vector< const Elem * > > _nodes_to_primary_elem_map
TypeVector< typename CompareTypes< T, T2 >::supertype > cross(const TypeVector< T2 > &v) const
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
std::unordered_map< std::pair< const Node *, const Elem * >, std::pair< Real, const Elem * > > _secondary_node_and_elem_to_xi2_primary_elem
Similar to the map above, but associates a (Secondary Node, Secondary Elem) pair to a (xi^(2)...
const Node * node_ptr(const unsigned int i) const
const Point & point(const unsigned int i) const
auto index_range(const T &sizable)
std::unordered_map< dof_id_type, std::vector< const Elem * > > _nodes_to_secondary_elem_map
Map from nodes to connected lower-dimensional elements on the secondary/primary subdomains.

◆ projectPrimaryNodes()

void AutomaticMortarGeneration::projectPrimaryNodes ( )

(Inverse) project primary nodes to the points on the secondary surface where they would have come from (find (xi^(1) values)).

Inputs:

  • The nodal normals values
  • mesh
  • nodes_to_secondary_elem_map

Outputs:

  • primary_node_and_elem_to_xi1_secondary_elem

Defined in the file project_primary_nodes.C.

Definition at line 2654 of file AutomaticMortarGeneration.C.

Referenced by MortarInterfaceWarehouse::update().

2655 {
2656  // For each primary/secondary boundary id pair, call the
2657  // project_primary_nodes_single_pair() helper function.
2658  for (const auto & pr : _primary_secondary_subdomain_id_pairs)
2659  projectPrimaryNodesSinglePair(pr.first, pr.second);
2660 }
void projectPrimaryNodesSinglePair(SubdomainID lower_dimensional_primary_subdomain_id, SubdomainID lower_dimensional_secondary_subdomain_id)
Helper function used internally by AutomaticMortarGeneration::project_primary_nodes().
std::vector< std::pair< SubdomainID, SubdomainID > > _primary_secondary_subdomain_id_pairs
A list of primary/secondary subdomain id pairs corresponding to each side of the mortar interface...

◆ projectPrimaryNodesSinglePair()

void AutomaticMortarGeneration::projectPrimaryNodesSinglePair ( SubdomainID  lower_dimensional_primary_subdomain_id,
SubdomainID  lower_dimensional_secondary_subdomain_id 
)
private

Helper function used internally by AutomaticMortarGeneration::project_primary_nodes().

Definition at line 2663 of file AutomaticMortarGeneration.C.

Referenced by projectPrimaryNodes().

2666 {
2667  using std::abs;
2668 
2669  // Build a Nanoflann object on the lower-dimensional secondary elements of the Mesh.
2670  NanoflannMeshSubdomainAdaptor<3> mesh_adaptor(_mesh, lower_dimensional_secondary_subdomain_id);
2671  subdomain_kd_tree_t kd_tree(
2672  3, mesh_adaptor, nanoflann::KDTreeSingleIndexAdaptorParams(/*max leaf=*/10));
2673 
2674  // Construct the KD tree for lower-dimensional elements in the volume mesh.
2675  kd_tree.buildIndex();
2676 
2677  std::unordered_set<dof_id_type> primary_nodes_visited;
2678 
2679  for (const auto & primary_side_elem : _mesh.active_element_ptr_range())
2680  {
2681  // If this is not one of the lower-dimensional primary side elements, go on to the next one.
2682  if (primary_side_elem->subdomain_id() != lower_dimensional_primary_subdomain_id)
2683  continue;
2684 
2685  // For each node on this side, find the nearest node on the secondary side using the KDTree,
2686  // then search in nearby elements for where it projects along the nodal normal direction.
2687  for (MooseIndex(primary_side_elem->n_vertices()) n = 0; n < primary_side_elem->n_vertices();
2688  ++n)
2689  {
2690  // Get a pointer to this node.
2691  const Node * primary_node = primary_side_elem->node_ptr(n);
2692 
2693  // Get the nodal neighbors connected to this primary node.
2694  const std::vector<const Elem *> & primary_node_neighbors =
2695  _nodes_to_primary_elem_map.at(primary_node->id());
2696 
2697  // Check whether we have already successfully inverse mapped this primary node (whether during
2698  // secondary node projection or now during primary node projection) or we have already failed
2699  // to inverse map this primary node (now during primary node projection), and then skip if
2700  // either of those things is true
2701  auto primary_key =
2702  std::make_tuple(primary_node->id(), primary_node, primary_node_neighbors[0]);
2703  if (!primary_nodes_visited.insert(primary_node->id()).second ||
2705  continue;
2706 
2707  // Data structure for performing Nanoflann searches.
2708  Real query_pt[3] = {(*primary_node)(0), (*primary_node)(1), (*primary_node)(2)};
2709 
2710  // The number of results we want to get. We'll look for a
2711  // "few" nearest nodes, hopefully that is enough to let us
2712  // figure out which lower-dimensional Elem on the secondary side
2713  // we are across from.
2714  const size_t num_results = 3;
2715 
2716  // Initialize result_set and do the search.
2717  std::vector<size_t> ret_index(num_results);
2718  std::vector<Real> out_dist_sqr(num_results);
2719  nanoflann::KNNResultSet<Real> result_set(num_results);
2720  result_set.init(&ret_index[0], &out_dist_sqr[0]);
2721  kd_tree.findNeighbors(result_set, &query_pt[0], nanoflann::SearchParameters());
2722 
2723  // If this flag gets set in the loop below, we can break out of the outer r-loop as well.
2724  bool projection_succeeded = false;
2725 
2726  // Once we've rejected a candidate for a given
2727  // primary_node, there's no reason to check it
2728  // again.
2729  std::set<const Elem *> rejected_secondary_elem_candidates;
2730 
2731  // Loop over the closest nodes, check whether the secondary node successfully projects into
2732  // either of the closest neighbors, stop when the projection succeeds.
2733  for (MooseIndex(result_set) r = 0; r < result_set.size(); ++r)
2734  {
2735  // Verify that the squared distance we compute is the same as nanoflann's
2736  mooseAssert(abs((_mesh.point(ret_index[r]) - *primary_node).norm_sq() - out_dist_sqr[r]) <=
2737  TOLERANCE,
2738  "Lower-dimensional element squared distance verification failed.");
2739 
2740  // Get a reference to the vector of lower dimensional elements from the
2741  // nodes_to_secondary_elem_map.
2742  const std::vector<const Elem *> & secondary_elem_candidates =
2743  _nodes_to_secondary_elem_map.at(static_cast<dof_id_type>(ret_index[r]));
2744 
2745  // Print the Elems connected to this node on the secondary mesh side.
2746  for (MooseIndex(secondary_elem_candidates) e = 0; e < secondary_elem_candidates.size(); ++e)
2747  {
2748  const Elem * secondary_elem_candidate = secondary_elem_candidates[e];
2749 
2750  // If we've already rejected this candidate, we don't need to check it again.
2751  if (rejected_secondary_elem_candidates.count(secondary_elem_candidate))
2752  continue;
2753 
2754  std::vector<Point> nodal_normals(secondary_elem_candidate->n_nodes());
2755  for (const auto n : make_range(secondary_elem_candidate->n_nodes()))
2756  nodal_normals[n] =
2757  _secondary_node_to_nodal_normal.at(secondary_elem_candidate->node_ptr(n));
2758 
2759  // Use equation 2.4.6 from Bin Yang's dissertation to try and solve for
2760  // the position on the secondary element where this primary came from. This
2761  // requires a Newton iteration in general.
2762  DualNumber<Real> xi1_dn{0, 1}; // initial guess
2763  auto && order = secondary_elem_candidate->default_order();
2764  unsigned int current_iterate = 0, max_iterates = 10;
2765 
2766  VectorValue<DualNumber<Real>> normals(0);
2767 
2768  // Newton iteration loop - this to converge in 1 iteration when it
2769  // succeeds, and possibly two iterations when it converges to a
2770  // xi outside the reference element. I don't know any reason why it should
2771  // only take 1 iteration -- the Jacobian is not constant in general...
2772  do
2773  {
2775  for (MooseIndex(secondary_elem_candidate->n_nodes()) n = 0;
2776  n < secondary_elem_candidate->n_nodes();
2777  ++n)
2778  {
2779  const auto phi = Moose::fe_lagrange_1D_shape(order, n, xi1_dn);
2780  x1 += phi * secondary_elem_candidate->point(n);
2781  normals += phi * nodal_normals[n];
2782  }
2783 
2784  const auto u = x1 - (*primary_node);
2785 
2786  const auto F = u(0) * normals(1) - u(1) * normals(0);
2787 
2788  if (abs(F) < _newton_tolerance)
2789  break;
2790 
2791  // Unlike for projection of nodal normals onto primary surfaces, we should never have a
2792  // case where the nodal normal is completely orthogonal to the secondary surface, so we
2793  // do not have to guard against F.derivatives() == 0 here
2794  Real dxi1 = -F.value() / F.derivatives();
2795 
2796  xi1_dn += dxi1;
2797 
2798  normals = 0;
2799  } while (++current_iterate < max_iterates);
2800 
2801  Real xi1 = xi1_dn.value();
2802 
2803  // Check for convergence to a valid solution... The last condition checks for obliqueness
2804  // of the projection
2805  if ((current_iterate < max_iterates) && (abs(xi1) <= 1. + _xi_tolerance) &&
2806  (abs((primary_side_elem->point(0) - primary_side_elem->point(1)).unit() *
2807  MetaPhysicL::raw_value(normals).unit()) <
2809  {
2810  if (abs(abs(xi1) - 1.) < _xi_tolerance)
2811  {
2812  // Special case: xi1=+/-1.
2813  // It is unlikely that we get here, because this primary node should already
2814  // have been mapped during the project_secondary_nodes() routine, but
2815  // there is still a chance since the tolerances are applied to
2816  // the xi coordinate and that value may be different on a primary element and a
2817  // secondary element since they may have different sizes. It's also possible that we
2818  // may reach this point if the solve has yielded a non-physical configuration such as
2819  // one block being pushed way out into space
2820  const Node & secondary_node = (xi1 < 0) ? secondary_elem_candidate->node_ref(0)
2821  : secondary_elem_candidate->node_ref(1);
2822  bool created_mortar_segment = false;
2823 
2824  // If we have failed to project this secondary node, let's try again now
2825  if (_failed_secondary_node_projections.count(secondary_node.id()))
2826  created_mortar_segment = processAlignedNodes(secondary_node,
2827  *primary_node,
2828  nullptr,
2829  &primary_node_neighbors,
2830  MetaPhysicL::raw_value(normals),
2831  *secondary_elem_candidate,
2832  rejected_secondary_elem_candidates);
2833  else
2834  rejected_secondary_elem_candidates.insert(secondary_elem_candidate);
2835 
2836  if (!created_mortar_segment)
2837  // We used to throw an exception in this scope but now that we support processing
2838  // aligned nodes within this primary node projection method, I don't see any harm in
2839  // simply rejecting the secondary element candidate in the case of failure and
2840  // continuing just as we do when projecting secondary nodes
2841  continue;
2842  }
2843  else // somewhere in the middle of the Elem
2844  {
2845  // Add entry to primary_node_and_elem_to_xi1_secondary_elem
2846  //
2847  // Note: we originally duplicated the map values for the keys (node, left_neighbor)
2848  // and (node, right_neighbor) but I don't think that should be necessary. Instead we
2849  // just do it for neighbor 0, but really maybe we don't even need to do that since
2850  // we can always look up the neighbors later given the Node... keeping it like this
2851  // helps to maintain the "symmetry" of the two containers.
2852  const Elem * neigh = primary_node_neighbors[0];
2853  for (MooseIndex(neigh->n_vertices()) nid = 0; nid < neigh->n_vertices(); ++nid)
2854  {
2855  const Node * neigh_node = neigh->node_ptr(nid);
2856  if (primary_node == neigh_node)
2857  {
2858  auto key = std::make_tuple(neigh_node->id(), neigh_node, neigh);
2859  auto val = std::make_pair(xi1, secondary_elem_candidate);
2861  }
2862  }
2863  }
2864 
2865  projection_succeeded = true;
2866  break; // out of e-loop
2867  }
2868  else
2869  {
2870  // The current primary_point is not in this Elem, so keep track of the rejects.
2871  rejected_secondary_elem_candidates.insert(secondary_elem_candidate);
2872  }
2873  } // end e-loop over candidate elems
2874 
2875  if (projection_succeeded)
2876  break; // out of r-loop
2877  } // r-loop
2878 
2879  if (!projection_succeeded && _debug)
2880  {
2881  _console << "\nFailed to find point from which primary node "
2882  << static_cast<const Point &>(*primary_node) << " was projected." << std::endl
2883  << std::endl;
2884  }
2885  } // loop over side nodes
2886  } // end loop over elements for finding where primary points would have projected from.
2887 }
MetaPhysicL::DualNumber< V, D, asd > abs(const MetaPhysicL::DualNumber< V, D, asd > &a)
Definition: EigenADReal.h:50
auto norm_sq(const T &a)
Real _newton_tolerance
Newton solve tolerance for node projections.
auto raw_value(const Eigen::Map< T > &in)
Definition: EigenADReal.h:100
Special adaptor that works with subdomains of the Mesh.
const bool _debug
Whether to print debug output.
CTSub CT_OPERATOR_BINARY CTMul CTCompareLess CTCompareGreater CTCompareEqual _arg template cos(_arg) *_arg.template D< dtag >()) CT_SIMPLE_UNARY_FUNCTION(cos
const Node & node_ref(const unsigned int i) const
dof_id_type id() const
virtual unsigned int n_nodes() const=0
std::map< std::tuple< dof_id_type, const Node *, const Elem * >, std::pair< Real, const Elem * > > _primary_node_and_elem_to_xi1_secondary_elem
Same type of container, but for mapping (Primary Node ID, Primary Node, Primary Elem) -> (xi^(1)...
std::unordered_map< dof_id_type, std::vector< const Elem * > > _nodes_to_primary_elem_map
bool processAlignedNodes(const Node &secondary_node, const Node &primary_node, const std::vector< const Elem *> *secondary_node_neighbors, const std::vector< const Elem *> *primary_node_neighbors, const VectorValue< Real > &nodal_normal, const Elem &candidate_element, std::set< const Elem *> &rejected_element_candidates)
Process aligned nodes.
nanoflann::KDTreeSingleIndexAdaptor< subdomain_adatper_t, NanoflannMeshSubdomainAdaptor< 3 >, 3 > subdomain_kd_tree_t
Real _xi_tolerance
Tolerance for checking projection xi values.
std::unordered_set< dof_id_type > _failed_secondary_node_projections
Secondary nodes that failed to project.
virtual unsigned int n_vertices() const=0
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
const Node * node_ptr(const unsigned int i) const
IntRange< T > make_range(T beg, T end)
T fe_lagrange_1D_shape(const Order order, const unsigned int i, const T &xi)
const Real _minimum_projection_angle
Parameter to control which angle (in degrees) is admissible for the creation of mortar segments...
MeshBase & _mesh
Reference to the mesh stored in equation_systems.
virtual const Point & point(const dof_id_type i) const=0
const ConsoleStream _console
An instance of helper class to write streams to the Console objects.
std::unordered_map< const Node *, Point > _secondary_node_to_nodal_normal
Container for storing the nodal normal vector associated with each secondary node.
virtual Order default_order() const=0
SearchParams SearchParameters
const Point & point(const unsigned int i) const
std::unordered_map< dof_id_type, std::vector< const Elem * > > _nodes_to_secondary_elem_map
Map from nodes to connected lower-dimensional elements on the secondary/primary subdomains.
const Real pi

◆ projectSecondaryNodes()

void AutomaticMortarGeneration::projectSecondaryNodes ( )

Project secondary nodes (find xi^(2) values) to the closest points on the primary surface.

Inputs:

  • The nodal normals values
  • mesh
  • nodes_to_primary_elem_map

Outputs:

  • secondary_node_and_elem_to_xi2_primary_elem

Defined in the file project_secondary_nodes.C.

Definition at line 2284 of file AutomaticMortarGeneration.C.

Referenced by MortarInterfaceWarehouse::update().

2285 {
2286  // For each primary/secondary boundary id pair, call the
2287  // project_secondary_nodes_single_pair() helper function.
2288  for (const auto & pr : _primary_secondary_subdomain_id_pairs)
2289  projectSecondaryNodesSinglePair(pr.first, pr.second);
2290 }
void projectSecondaryNodesSinglePair(SubdomainID lower_dimensional_primary_subdomain_id, SubdomainID lower_dimensional_secondary_subdomain_id)
Helper function responsible for projecting secondary nodes onto primary elements for a single primary...
std::vector< std::pair< SubdomainID, SubdomainID > > _primary_secondary_subdomain_id_pairs
A list of primary/secondary subdomain id pairs corresponding to each side of the mortar interface...

◆ projectSecondaryNodesSinglePair()

void AutomaticMortarGeneration::projectSecondaryNodesSinglePair ( SubdomainID  lower_dimensional_primary_subdomain_id,
SubdomainID  lower_dimensional_secondary_subdomain_id 
)
private

Helper function responsible for projecting secondary nodes onto primary elements for a single primary/secondary pair.

Called by the class member AutomaticMortarGeneration::project_secondary_nodes().

Definition at line 2409 of file AutomaticMortarGeneration.C.

Referenced by projectSecondaryNodes().

2412 {
2413  using std::abs;
2414 
2415  // Build the "subdomain" adaptor based KD Tree.
2416  NanoflannMeshSubdomainAdaptor<3> mesh_adaptor(_mesh, lower_dimensional_primary_subdomain_id);
2417  subdomain_kd_tree_t kd_tree(
2418  3, mesh_adaptor, nanoflann::KDTreeSingleIndexAdaptorParams(/*max leaf=*/10));
2419 
2420  // Construct the KD tree.
2421  kd_tree.buildIndex();
2422 
2423  for (MeshBase::const_element_iterator el = _mesh.active_elements_begin(),
2424  end_el = _mesh.active_elements_end();
2425  el != end_el;
2426  ++el)
2427  {
2428  const Elem * secondary_side_elem = *el;
2429 
2430  // If this Elem is not in the current secondary subdomain, go on to the next one.
2431  if (secondary_side_elem->subdomain_id() != lower_dimensional_secondary_subdomain_id)
2432  continue;
2433 
2434  // For each node on the lower-dimensional element, find the nearest
2435  // node on the primary side using the KDTree, then
2436  // search in nearby elements for where it projects
2437  // along the nodal normal direction.
2438  for (MooseIndex(secondary_side_elem->n_vertices()) n = 0; n < secondary_side_elem->n_vertices();
2439  ++n)
2440  {
2441  const Node * secondary_node = secondary_side_elem->node_ptr(n);
2442 
2443  // Get the nodal neighbors for secondary_node, so we can check whether we've
2444  // already successfully projected it.
2445  const std::vector<const Elem *> & secondary_node_neighbors =
2446  this->_nodes_to_secondary_elem_map.at(secondary_node->id());
2447 
2448  // Check whether we've already mapped this secondary node
2449  // successfully for all of its nodal neighbors.
2450  bool is_mapped = true;
2451  for (MooseIndex(secondary_node_neighbors) snn = 0; snn < secondary_node_neighbors.size();
2452  ++snn)
2453  {
2454  auto secondary_key = std::make_pair(secondary_node, secondary_node_neighbors[snn]);
2455  if (!_secondary_node_and_elem_to_xi2_primary_elem.count(secondary_key))
2456  {
2457  is_mapped = false;
2458  break;
2459  }
2460  }
2461 
2462  // Go to the next node if this one has already been mapped.
2463  if (is_mapped)
2464  continue;
2465 
2466  // Look up the new nodal normal value in the local storage, error if not found.
2467  Point nodal_normal = _secondary_node_to_nodal_normal.at(secondary_node);
2468 
2469  // Data structure for performing Nanoflann searches.
2470  std::array<Real, 3> query_pt = {
2471  {(*secondary_node)(0), (*secondary_node)(1), (*secondary_node)(2)}};
2472 
2473  // The number of results we want to get. We'll look for a
2474  // "few" nearest nodes, hopefully that is enough to let us
2475  // figure out which lower-dimensional Elem on the primary
2476  // side we are across from.
2477  const std::size_t num_results = 3;
2478 
2479  // Initialize result_set and do the search.
2480  std::vector<size_t> ret_index(num_results);
2481  std::vector<Real> out_dist_sqr(num_results);
2482  nanoflann::KNNResultSet<Real> result_set(num_results);
2483  result_set.init(&ret_index[0], &out_dist_sqr[0]);
2484  kd_tree.findNeighbors(result_set, &query_pt[0], nanoflann::SearchParameters());
2485 
2486  // If this flag gets set in the loop below, we can break out of the outer r-loop as well.
2487  bool projection_succeeded = false;
2488 
2489  // Once we've rejected a candidate for a given secondary_node,
2490  // there's no reason to check it again.
2491  std::set<const Elem *> rejected_primary_elem_candidates;
2492 
2493  // Loop over the closest nodes, check whether
2494  // the secondary node successfully projects into
2495  // either of the closest neighbors, stop when
2496  // the projection succeeds.
2497  for (MooseIndex(result_set) r = 0; r < result_set.size(); ++r)
2498  {
2499  // Verify that the squared distance we compute is the same as nanoflann'sFss
2500  mooseAssert(abs((_mesh.point(ret_index[r]) - *secondary_node).norm_sq() -
2501  out_dist_sqr[r]) <= TOLERANCE,
2502  "Lower-dimensional element squared distance verification failed.");
2503 
2504  // Get a reference to the vector of lower dimensional elements from the
2505  // nodes_to_primary_elem_map.
2506  std::vector<const Elem *> & primary_elem_candidates =
2507  this->_nodes_to_primary_elem_map.at(static_cast<dof_id_type>(ret_index[r]));
2508 
2509  // Search the Elems connected to this node on the primary mesh side.
2510  for (MooseIndex(primary_elem_candidates) e = 0; e < primary_elem_candidates.size(); ++e)
2511  {
2512  const Elem * primary_elem_candidate = primary_elem_candidates[e];
2513 
2514  // If we've already rejected this candidate, we don't need to check it again.
2515  if (rejected_primary_elem_candidates.count(primary_elem_candidate))
2516  continue;
2517 
2518  // Now generically solve for xi2
2519  const auto order = primary_elem_candidate->default_order();
2520  DualNumber<Real> xi2_dn{0, 1};
2521  unsigned int current_iterate = 0, max_iterates = 10;
2522 
2523  // Newton loop
2524  do
2525  {
2527  for (MooseIndex(primary_elem_candidate->n_nodes()) n = 0;
2528  n < primary_elem_candidate->n_nodes();
2529  ++n)
2530  x2 +=
2531  Moose::fe_lagrange_1D_shape(order, n, xi2_dn) * primary_elem_candidate->point(n);
2532  const auto u = x2 - (*secondary_node);
2533  const auto F = u(0) * nodal_normal(1) - u(1) * nodal_normal(0);
2534 
2535  if (abs(F) < _newton_tolerance)
2536  break;
2537 
2538  if (F.derivatives())
2539  {
2540  Real dxi2 = -F.value() / F.derivatives();
2541 
2542  xi2_dn += dxi2;
2543  }
2544  else
2545  // It's possible that the secondary surface nodal normal is completely orthogonal to
2546  // the primary surface normal, in which case the derivative is 0. We know in this case
2547  // that the projection should be a failure
2548  current_iterate = max_iterates;
2549  } while (++current_iterate < max_iterates);
2550 
2551  Real xi2 = xi2_dn.value();
2552 
2553  // Check whether the projection worked. The last condition checks for obliqueness of the
2554  // projection
2555  //
2556  // We are projecting on one side first and the other side second. If we make the
2557  // tolerance bigger and remove the (5) factor we are going to continue to miss the
2558  // second projection and fall into the exception message in
2559  // projectPrimaryNodesSinglePair. What makes this modification to not fall in the
2560  // exception is that we are projecting on one side more xi than in the other. There
2561  // should be a better way of doing this by using actual distances and not parametric
2562  // coordinates. But I believe making the tolerance uniformly larger or smaller won't do
2563  // the trick here.
2564  if ((current_iterate < max_iterates) && (std::abs(xi2) <= 1. + 5 * _xi_tolerance) &&
2565  (abs((primary_elem_candidate->point(0) - primary_elem_candidate->point(1)).unit() *
2566  nodal_normal) < std::cos(_minimum_projection_angle * libMesh::pi / 180)))
2567  {
2568  // If xi2 == +1 or -1 then this secondary node mapped directly to a node on the primary
2569  // surface. This isn't as unlikely as you might think, it will happen if the meshes
2570  // on the interface start off being perfectly aligned. In this situation, we need to
2571  // associate the secondary node with two different elements (and two corresponding
2572  // xi^(2) values.
2573  if (abs(abs(xi2) - 1.) <= _xi_tolerance * 5.0)
2574  {
2575  const Node * primary_node = (xi2 < 0) ? primary_elem_candidate->node_ptr(0)
2576  : primary_elem_candidate->node_ptr(1);
2577  const bool created_mortar_segment =
2578  processAlignedNodes(*secondary_node,
2579  *primary_node,
2580  &secondary_node_neighbors,
2581  nullptr,
2582  nodal_normal,
2583  *primary_elem_candidate,
2584  rejected_primary_elem_candidates);
2585 
2586  if (!created_mortar_segment)
2587  continue;
2588  }
2589  else // Point falls somewhere in the middle of the Elem.
2590  {
2591  // Add two entries to secondary_node_and_elem_to_xi2_primary_elem.
2592  for (MooseIndex(secondary_node_neighbors) nn = 0;
2593  nn < secondary_node_neighbors.size();
2594  ++nn)
2595  {
2596  const Elem * neigh = secondary_node_neighbors[nn];
2597  for (MooseIndex(neigh->n_vertices()) nid = 0; nid < neigh->n_vertices(); ++nid)
2598  {
2599  const Node * neigh_node = neigh->node_ptr(nid);
2600  if (secondary_node == neigh_node)
2601  {
2602  auto key = std::make_pair(neigh_node, neigh);
2603  auto val = std::make_pair(xi2, primary_elem_candidate);
2605  }
2606  }
2607  }
2608  }
2609 
2610  projection_succeeded = true;
2611  break; // out of e-loop
2612  }
2613  else
2614  // The current secondary_node is not in this Elem, so keep track of the rejects.
2615  rejected_primary_elem_candidates.insert(primary_elem_candidate);
2616  }
2617 
2618  if (projection_succeeded)
2619  break; // out of r-loop
2620  } // r-loop
2621 
2622  if (!projection_succeeded)
2623  {
2624  _failed_secondary_node_projections.insert(secondary_node->id());
2625  if (_debug)
2626  _console << "Failed to find primary Elem into which secondary node "
2627  << static_cast<const Point &>(*secondary_node) << ", id '"
2628  << secondary_node->id() << "', projects onto\n"
2629  << std::endl;
2630  }
2631  else if (_debug)
2632  _projected_secondary_nodes.insert(secondary_node->id());
2633  } // loop over side nodes
2634  } // end loop over lower-dimensional elements
2635 
2636  if (_distributed)
2637  {
2638  if (_debug)
2641  }
2642 
2643  if (_debug)
2644  _console << "\n"
2645  << _projected_secondary_nodes.size() << " out of "
2647  << " secondary nodes were successfully projected\n"
2648  << std::endl;
2649 }
MetaPhysicL::DualNumber< V, D, asd > abs(const MetaPhysicL::DualNumber< V, D, asd > &a)
Definition: EigenADReal.h:50
auto norm_sq(const T &a)
Real _newton_tolerance
Newton solve tolerance for node projections.
const bool _distributed
Whether the mortar segment mesh is distributed.
const Parallel::Communicator & comm() const
Special adaptor that works with subdomains of the Mesh.
std::unordered_set< dof_id_type > _projected_secondary_nodes
Debugging container for printing information about fraction of successful projections for secondary n...
const bool _debug
Whether to print debug output.
CTSub CT_OPERATOR_BINARY CTMul CTCompareLess CTCompareGreater CTCompareEqual _arg template cos(_arg) *_arg.template D< dtag >()) CT_SIMPLE_UNARY_FUNCTION(cos
dof_id_type id() const
virtual unsigned int n_nodes() const=0
std::unordered_map< dof_id_type, std::vector< const Elem * > > _nodes_to_primary_elem_map
bool processAlignedNodes(const Node &secondary_node, const Node &primary_node, const std::vector< const Elem *> *secondary_node_neighbors, const std::vector< const Elem *> *primary_node_neighbors, const VectorValue< Real > &nodal_normal, const Elem &candidate_element, std::set< const Elem *> &rejected_element_candidates)
Process aligned nodes.
nanoflann::KDTreeSingleIndexAdaptor< subdomain_adatper_t, NanoflannMeshSubdomainAdaptor< 3 >, 3 > subdomain_kd_tree_t
Real _xi_tolerance
Tolerance for checking projection xi values.
std::unordered_set< dof_id_type > _failed_secondary_node_projections
Secondary nodes that failed to project.
virtual unsigned int n_vertices() const=0
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
std::unordered_map< std::pair< const Node *, const Elem * >, std::pair< Real, const Elem * > > _secondary_node_and_elem_to_xi2_primary_elem
Similar to the map above, but associates a (Secondary Node, Secondary Elem) pair to a (xi^(2)...
subdomain_id_type subdomain_id() const
const Node * node_ptr(const unsigned int i) const
T fe_lagrange_1D_shape(const Order order, const unsigned int i, const T &xi)
const Real _minimum_projection_angle
Parameter to control which angle (in degrees) is admissible for the creation of mortar segments...
MeshBase & _mesh
Reference to the mesh stored in equation_systems.
virtual const Point & point(const dof_id_type i) const=0
const ConsoleStream _console
An instance of helper class to write streams to the Console objects.
std::unordered_map< const Node *, Point > _secondary_node_to_nodal_normal
Container for storing the nodal normal vector associated with each secondary node.
virtual Order default_order() const=0
SearchParams SearchParameters
const Point & point(const unsigned int i) const
std::unordered_map< dof_id_type, std::vector< const Elem * > > _nodes_to_secondary_elem_map
Map from nodes to connected lower-dimensional elements on the secondary/primary subdomains.
const Real pi
void set_union(T &data, const unsigned int root_id) const

◆ secondariesToMortarSegments() [1/2]

std::vector< AutomaticMortarGeneration::MortarFilterIter > AutomaticMortarGeneration::secondariesToMortarSegments ( const Node node) const
Returns
A vector of iterators that point to the lower dimensional secondary elements and their associated mortar segment elements that would have nonzero values for a Lagrange shape function associated with the provided node. This method may return an empty container if the node is away from the mortar mesh

Definition at line 2890 of file AutomaticMortarGeneration.C.

Referenced by MortarUserObjectThread::operator()(), and ComputeMortarFunctor::operator()().

2891 {
2892  auto secondary_it = _nodes_to_secondary_elem_map.find(node.id());
2893  if (secondary_it == _nodes_to_secondary_elem_map.end())
2894  return {};
2895 
2896  const auto & secondary_elems = secondary_it->second;
2897  std::vector<MortarFilterIter> ret;
2898  ret.reserve(secondary_elems.size());
2899 
2900  for (const auto i : index_range(secondary_elems))
2901  {
2902  auto * const secondary_elem = secondary_elems[i];
2903  auto msm_it = _secondary_elems_to_mortar_segments.find(secondary_elem->id());
2904  if (msm_it == _secondary_elems_to_mortar_segments.end())
2905  // We may have removed this element key from this map
2906  continue;
2907 
2908  mooseAssert(secondary_elem->active(),
2909  "We loop over active elements when building the mortar segment mesh, so we golly "
2910  "well hope this is active.");
2911  mooseAssert(!msm_it->second.empty(),
2912  "We should have removed all secondaries from this map if they do not have any "
2913  "mortar segments associated with them.");
2914  ret.push_back(msm_it);
2915  }
2916 
2917  return ret;
2918 }
std::unordered_map< dof_id_type, std::set< Elem *, CompareDofObjectsByID > > _secondary_elems_to_mortar_segments
We maintain a mapping from lower-dimensional secondary elements in the original mesh to (sets of) ele...
dof_id_type id() const
auto index_range(const T &sizable)
std::unordered_map< dof_id_type, std::vector< const Elem * > > _nodes_to_secondary_elem_map
Map from nodes to connected lower-dimensional elements on the secondary/primary subdomains.

◆ secondariesToMortarSegments() [2/2]

const std::unordered_map<dof_id_type, std::set<Elem *, CompareDofObjectsByID> >& AutomaticMortarGeneration::secondariesToMortarSegments ( ) const
inline
Returns
the lower dimensional secondary element ids and their associated mortar segment elements

Definition at line 373 of file AutomaticMortarGeneration.h.

374  {
376  }
std::unordered_map< dof_id_type, std::set< Elem *, CompareDofObjectsByID > > _secondary_elems_to_mortar_segments
We maintain a mapping from lower-dimensional secondary elements in the original mesh to (sets of) ele...

◆ secondaryIPSubIDs()

const std::set<SubdomainID>& AutomaticMortarGeneration::secondaryIPSubIDs ( ) const
inline
Returns
All the secondary interior parent subdomain IDs associated with the mortar mesh

Definition at line 381 of file AutomaticMortarGeneration.h.

Referenced by Moose::Mortar::setupMortarMaterials().

381 { return _secondary_ip_sub_ids; }
std::set< SubdomainID > _secondary_ip_sub_ids
All the secondary interior parent subdomain IDs associated with the mortar mesh.

Friends And Related Function Documentation

◆ AugmentSparsityOnInterface

friend class AugmentSparsityOnInterface
friend

Definition at line 618 of file AutomaticMortarGeneration.h.

◆ MortarNodalGeometryOutput

friend class MortarNodalGeometryOutput
friend

Definition at line 617 of file AutomaticMortarGeneration.h.

Member Data Documentation

◆ _app

MooseApp& AutomaticMortarGeneration::_app
private

The Moose app.

Definition at line 418 of file AutomaticMortarGeneration.h.

Referenced by initOutput(), and outputMortarMesh().

◆ _console

const ConsoleStream ConsoleStreamInterface::_console
inherited

An instance of helper class to write streams to the Console objects.

Definition at line 31 of file ConsoleStreamInterface.h.

Referenced by IterationAdaptiveDT::acceptStep(), MeshOnlyAction::act(), SetupDebugAction::act(), MaterialOutputAction::act(), Adaptivity::adaptMesh(), FEProblemBase::adaptMesh(), PerfGraph::addToExecutionList(), SimplePredictor::apply(), SystemBase::applyScalingFactors(), MultiApp::backup(), FEProblemBase::backupMultiApps(), CoarsenedPiecewiseLinear::buildCoarsenedGrid(), DefaultSteadyStateConvergence::checkConvergence(), MeshDiagnosticsGenerator::checkElementOverlap(), MeshDiagnosticsGenerator::checkElementTypes(), MeshDiagnosticsGenerator::checkElementVolumes(), FEProblemBase::checkExceptionAndStopSolve(), SolverSystem::checkInvalidSolution(), MeshDiagnosticsGenerator::checkLocalJacobians(), MeshDiagnosticsGenerator::checkNonConformalMesh(), MeshDiagnosticsGenerator::checkNonConformalMeshFromAdaptivity(), MeshDiagnosticsGenerator::checkNonMatchingEdges(), MeshDiagnosticsGenerator::checkNonPlanarSides(), MeshDiagnosticsGenerator::checkPolygons(), FEProblemBase::checkProblemIntegrity(), ReferenceResidualConvergence::checkResidualConvergence(), MeshDiagnosticsGenerator::checkSidesetsOrientation(), MeshDiagnosticsGenerator::checkWatertightNodesets(), MeshDiagnosticsGenerator::checkWatertightSidesets(), IterationAdaptiveDT::computeAdaptiveDT(), TransientBase::computeConstrainedDT(), DefaultMultiAppFixedPointConvergence::computeCustomConvergencePostprocessor(), NonlinearSystemBase::computeDamping(), FixedPointIterationAdaptiveDT::computeDT(), IterationAdaptiveDT::computeDT(), IterationAdaptiveDT::computeFailedDT(), IterationAdaptiveDT::computeInitialDT(), IterationAdaptiveDT::computeInterpolationDT(), LinearSystem::computeLinearSystemTags(), FEProblemBase::computeLinearSystemTags(), NonlinearSystemBase::computeScaling(), Problem::console(), IterationAdaptiveDT::constrainStep(), TimeStepper::constrainStep(), MultiApp::createApp(), FEProblemBase::execMultiApps(), FEProblemBase::execMultiAppTransfers(), MFEMSteady::execute(), MessageFromInput::execute(), SteadyBase::execute(), Eigenvalue::execute(), ActionWarehouse::executeActionsWithAction(), ActionWarehouse::executeAllActions(), MeshGeneratorSystem::executeMeshGenerators(), ElementQualityChecker::finalize(), SidesetAroundSubdomainUpdater::finalize(), FEProblemBase::finishMultiAppStep(), MeshRepairGenerator::fixOverlappingNodes(), SurfaceSubdomainsDelaunayRemesher::General2DDelaunay(), OrientSurfaceMeshGenerator::generate(), CoarsenBlockGenerator::generate(), PolyLineMeshFollowingNodeSetGenerator::generate(), MeshGenerator::generateInternal(), VariableCondensationPreconditioner::getDofToCondense(), NonlinearEigen::init(), InversePowerMethod::init(), FEProblemBase::initialAdaptMesh(), DefaultMultiAppFixedPointConvergence::initialize(), SubProblem::initialSetup(), EigenExecutionerBase::inversePowerIteration(), FEProblemBase::joinAndFinalize(), TransientBase::keepGoing(), IterationAdaptiveDT::limitDTByFunction(), IterationAdaptiveDT::limitDTToPostprocessorValue(), FEProblemBase::logAdd(), EigenExecutionerBase::makeBXConsistent(), Console::meshChanged(), SurfaceDelaunayGeneratorBase::meshNormalDeviation2D(), MooseBase::mooseDeprecated(), MooseBase::mooseDeprecatedNoTrace(), MooseBase::mooseInfo(), MooseBase::mooseWarning(), MooseBase::mooseWarningNonPrefixed(), ReferenceResidualConvergence::nonlinearConvergenceSetup(), ReporterDebugOutput::output(), PerfGraphOutput::output(), SolutionInvalidityOutput::output(), MaterialPropertyDebugOutput::output(), DOFMapOutput::output(), VariableResidualNormsDebugOutput::output(), Console::output(), ControlOutput::outputActiveObjects(), ControlOutput::outputChangedControls(), ControlOutput::outputControls(), Console::outputInput(), WebServerControl::outputMessage(), Console::outputPostprocessors(), PseudoTimestep::outputPseudoTimestep(), Console::outputReporters(), DefaultMultiAppFixedPointConvergence::outputResidualNorm(), Console::outputScalarVariables(), Console::outputSystemInformation(), FEProblemBase::possiblyRebuildGeomSearchPatches(), EigenExecutionerBase::postExecute(), AB2PredictorCorrector::postSolve(), ActionWarehouse::printActionDependencySets(), BlockRestrictionDebugOutput::printBlockRestrictionGroups(), BlockRestrictionDebugOutput::printBlockRestrictionMap(), BlockRestrictionDebugOutput::printBoundaryRestrictionGroups(), SolutionInvalidity::printDebug(), EigenExecutionerBase::printEigenvalue(), SecantSolve::printFixedPointConvergenceHistory(), SteffensenSolve::printFixedPointConvergenceHistory(), PicardSolve::printFixedPointConvergenceHistory(), FixedPointSolve::printFixedPointConvergenceReason(), PerfGraphLivePrint::printLiveMessage(), MaterialPropertyDebugOutput::printMaterialMap(), PerfGraphLivePrint::printStats(), NEML2Action::printSummary(), projectPrimaryNodesSinglePair(), projectSecondaryNodesSinglePair(), CoarsenBlockGenerator::recursiveCoarsen(), SolutionTimeAdaptiveDT::rejectStep(), MultiApp::restore(), FEProblemBase::restoreMultiApps(), FEProblemBase::restoreSolutions(), NonlinearSystemBase::setInitialSolution(), MooseApp::setupOptions(), Checkpoint::shouldOutput(), SubProblem::showFunctorRequestors(), SubProblem::showFunctors(), FullSolveMultiApp::showStatusMessage(), EigenProblem::solve(), FEProblemSolve::solve(), NonlinearSystem::solve(), FixedPointSolve::solve(), LStableDirk2::solve(), LStableDirk3::solve(), ImplicitMidpoint::solve(), LinearSystem::solve(), ExplicitTVDRK2::solve(), LStableDirk4::solve(), AStableDirk4::solve(), ExplicitRK2::solve(), TransientMultiApp::solveStep(), FixedPointSolve::solveStep(), MeshRepairGenerator::splitNonConvexPolygons(), PerfGraphLivePrint::start(), WebServerControl::startServer(), AB2PredictorCorrector::step(), NonlinearEigen::takeStep(), MFEMTransient::takeStep(), TransientBase::takeStep(), TerminateChainControl::terminate(), SubProblem::timestepSetup(), FEProblemBase::updateMeshXFEM(), Convergence::verboseOutput(), Console::writeTimestepInformation(), Console::writeVariableNorms(), and FEProblemBase::~FEProblemBase().

◆ _correct_edge_dropping

const bool AutomaticMortarGeneration::_correct_edge_dropping
private

Flag to enable regressed treatment of edge dropping where all LM DoFs on edge dropping element are strongly set to 0.

Definition at line 579 of file AutomaticMortarGeneration.h.

Referenced by computeInactiveLMElems(), computeInactiveLMNodes(), and incorrectEdgeDropping().

◆ _debug

const bool AutomaticMortarGeneration::_debug
private

◆ _distributed

const bool AutomaticMortarGeneration::_distributed
private

Whether the mortar segment mesh is distributed.

Definition at line 567 of file AutomaticMortarGeneration.h.

Referenced by AutomaticMortarGeneration(), and projectSecondaryNodesSinglePair().

◆ _failed_secondary_node_projections

std::unordered_set<dof_id_type> AutomaticMortarGeneration::_failed_secondary_node_projections
private

Secondary nodes that failed to project.

Definition at line 615 of file AutomaticMortarGeneration.h.

Referenced by clear(), projectPrimaryNodesSinglePair(), and projectSecondaryNodesSinglePair().

◆ _inactive_local_lm_elems

std::unordered_set<const Elem *> AutomaticMortarGeneration::_inactive_local_lm_elems
private

List of inactive lagrange multiplier nodes (for elemental variables)

Definition at line 511 of file AutomaticMortarGeneration.h.

Referenced by computeInactiveLMElems(), and getInactiveLMElems().

◆ _inactive_local_lm_nodes

std::unordered_set<const Node *> AutomaticMortarGeneration::_inactive_local_lm_nodes
private

◆ _lower_elem_to_side_id

std::unordered_map<const Elem *, unsigned int> AutomaticMortarGeneration::_lower_elem_to_side_id
private

Keeps track of the mapping between lower-dimensional elements and the side_id of the interior_parent which they are.

Definition at line 476 of file AutomaticMortarGeneration.h.

Referenced by clear().

◆ _mesh

MeshBase& AutomaticMortarGeneration::_mesh
private

◆ _minimum_projection_angle

const Real AutomaticMortarGeneration::_minimum_projection_angle
private

Parameter to control which angle (in degrees) is admissible for the creation of mortar segments.

If set to a value close to zero, very oblique projections are allowed, which can result in mortar segments solving physics not meaningfully and overprojection of primary nodes onto the mortar segment mesh in extreme cases. This parameter is mostly intended for mortar mesh debugging purposes in 2D.

Definition at line 585 of file AutomaticMortarGeneration.h.

Referenced by buildMortarSegmentMesh3d(), projectPrimaryNodesSinglePair(), and projectSecondaryNodesSinglePair().

◆ _mortar_3d_qp_mapping

const Mortar3DQuadraturePointMapping AutomaticMortarGeneration::_mortar_3d_qp_mapping
private

Method used to map 3D mortar segment quadrature points to primary and secondary faces.

Definition at line 597 of file AutomaticMortarGeneration.h.

Referenced by buildMortarSegmentMesh3d(), mortar3DQpMapping(), and mortarSegmentReferencePoints().

◆ _mortar_3d_subpatch_plane

const Mortar3DSubpatchPlane AutomaticMortarGeneration::_mortar_3d_subpatch_plane
private

Method used to define the local projection planes for 3D secondary subpatches.

Definition at line 588 of file AutomaticMortarGeneration.h.

Referenced by buildMortarSegmentMesh3d().

◆ _mortar_interface_coupling

std::unordered_map<dof_id_type, std::unordered_set<dof_id_type> > AutomaticMortarGeneration::_mortar_interface_coupling
private

Used by the AugmentSparsityOnInterface functor to determine whether a given Elem is coupled to any others across the gap, and to explicitly set up the dependence between interior_parent() elements on the secondary side and their lower-dimensional sides which are on the interface.

This latter type of coupling must be explicitly declared when there is no primary_elem for a given mortar segment and you are using e.g. a P^1-P^0 discretization which does not induce the coupling automatically.

Definition at line 495 of file AutomaticMortarGeneration.h.

Referenced by buildCouplingInformation(), clear(), and mortarInterfaceCoupling().

◆ _mortar_segment_mesh

std::unique_ptr<MeshBase> AutomaticMortarGeneration::_mortar_segment_mesh
private

1D Mesh of mortar segment elements which gets built by the call to build_mortar_segment_mesh().

Definition at line 467 of file AutomaticMortarGeneration.h.

Referenced by AutomaticMortarGeneration(), buildMortarSegmentMesh(), buildMortarSegmentMesh3d(), clear(), computeInactiveLMElems(), computeInactiveLMNodes(), computeIncorrectEdgeDroppingInactiveLMNodes(), mortarSegmentMesh(), and outputMortarMesh().

◆ _msm_elem_to_info

std::unordered_map<const Elem *, MortarSegmentInfo> AutomaticMortarGeneration::_msm_elem_to_info
private

Map between Elems in the mortar segment mesh and their info structs.

This gets filled in by the call to build_mortar_segment_mesh().

Definition at line 472 of file AutomaticMortarGeneration.h.

Referenced by buildCouplingInformation(), buildMortarSegmentMesh(), buildMortarSegmentMesh3d(), clear(), computeInactiveLMElems(), computeInactiveLMNodes(), computeIncorrectEdgeDroppingInactiveLMNodes(), computeMsmStatistics(), and mortarSegmentMeshElemToInfo().

◆ _msm_elem_to_reference_points

std::unordered_map<const Elem *, MortarSegmentReferencePoints> AutomaticMortarGeneration::_msm_elem_to_reference_points
private

Reference-coordinate data used only by the reference-interpolation mapping mode.

Definition at line 600 of file AutomaticMortarGeneration.h.

Referenced by buildMortarSegmentMesh3d(), clear(), and mortarSegmentReferencePoints().

◆ _msm_node_id_start

std::optional<dof_id_type> AutomaticMortarGeneration::_msm_node_id_start
private

Cached per-rank starting ID for 3D MSM nodes/elements.

nullopt forces recomputation on next buildMortarSegmentMesh3d() call. Reset by meshChanged() on topology change so the allgather is skipped for displaced-mesh residual updates that only move nodes.

Definition at line 608 of file AutomaticMortarGeneration.h.

Referenced by buildMortarSegmentMesh3d(), and meshChanged().

◆ _newton_tolerance

Real AutomaticMortarGeneration::_newton_tolerance = 1e-12
private

Newton solve tolerance for node projections.

Definition at line 570 of file AutomaticMortarGeneration.h.

Referenced by projectPrimaryNodesSinglePair(), and projectSecondaryNodesSinglePair().

◆ _nodes_to_primary_elem_map

std::unordered_map<dof_id_type, std::vector<const Elem *> > AutomaticMortarGeneration::_nodes_to_primary_elem_map
private

◆ _nodes_to_secondary_elem_map

std::unordered_map<dof_id_type, std::vector<const Elem *> > AutomaticMortarGeneration::_nodes_to_secondary_elem_map
private

Map from nodes to connected lower-dimensional elements on the secondary/primary subdomains.

Definition at line 434 of file AutomaticMortarGeneration.h.

Referenced by buildNodeToElemMaps(), clear(), nodesToSecondaryElem(), processAlignedNodes(), projectPrimaryNodesSinglePair(), projectSecondaryNodesSinglePair(), and secondariesToMortarSegments().

◆ _on_displaced

const bool AutomaticMortarGeneration::_on_displaced
private

Whether this object is on the displaced mesh.

Definition at line 561 of file AutomaticMortarGeneration.h.

Referenced by mortarInterfaceName(), and onDisplaced().

◆ _output_params

std::unique_ptr<InputParameters> AutomaticMortarGeneration::_output_params
private

Storage for the input parameters used by the mortar nodal geometry output.

Definition at line 603 of file AutomaticMortarGeneration.h.

Referenced by initOutput().

◆ _periodic

const bool AutomaticMortarGeneration::_periodic
private

Whether this object will be generating a mortar segment mesh for periodic constraints.

Definition at line 564 of file AutomaticMortarGeneration.h.

Referenced by computeNodalGeometry(), and getNormals().

◆ _primary_boundary_subdomain_ids

std::set<SubdomainID> AutomaticMortarGeneration::_primary_boundary_subdomain_ids
private

Definition at line 485 of file AutomaticMortarGeneration.h.

Referenced by AutomaticMortarGeneration(), and buildNodeToElemMaps().

◆ _primary_ip_sub_ids

std::set<SubdomainID> AutomaticMortarGeneration::_primary_ip_sub_ids
private

All the primary interior parent subdomain IDs associated with the mortar mesh.

Definition at line 523 of file AutomaticMortarGeneration.h.

Referenced by buildMortarSegmentMesh(), buildMortarSegmentMesh3d(), clear(), and primaryIPSubIDs().

◆ _primary_node_and_elem_to_xi1_secondary_elem

std::map<std::tuple<dof_id_type, const Node *, const Elem *>, std::pair<Real, const Elem *> > AutomaticMortarGeneration::_primary_node_and_elem_to_xi1_secondary_elem
private

Same type of container, but for mapping (Primary Node ID, Primary Node, Primary Elem) -> (xi^(1), Secondary Elem) where they are inverse-projected along the nodal normal direction.

Note that the first item of the key, the primary node ID, is important for storing the key-value pairs in a consistent order across processes, e.g. this container has to be ordered!

Definition at line 463 of file AutomaticMortarGeneration.h.

Referenced by buildMortarSegmentMesh(), clear(), processAlignedNodes(), and projectPrimaryNodesSinglePair().

◆ _primary_requested_boundary_ids

std::set<BoundaryID> AutomaticMortarGeneration::_primary_requested_boundary_ids
private

The boundary ids corresponding to all the primary surfaces.

Definition at line 427 of file AutomaticMortarGeneration.h.

Referenced by AutomaticMortarGeneration(), and buildNodeToElemMaps().

◆ _primary_secondary_boundary_id_pairs

std::vector<std::pair<BoundaryID, BoundaryID> > AutomaticMortarGeneration::_primary_secondary_boundary_id_pairs
private

A list of primary/secondary boundary id pairs corresponding to each side of the mortar interface.

Definition at line 431 of file AutomaticMortarGeneration.h.

Referenced by AutomaticMortarGeneration(), buildMortarSegmentMesh(), mortarInterfaceName(), and primarySecondaryBoundaryIDPair().

◆ _primary_secondary_subdomain_id_pairs

std::vector<std::pair<SubdomainID, SubdomainID> > AutomaticMortarGeneration::_primary_secondary_subdomain_id_pairs
private

A list of primary/secondary subdomain id pairs corresponding to each side of the mortar interface.

Definition at line 480 of file AutomaticMortarGeneration.h.

Referenced by AutomaticMortarGeneration(), buildMortarSegmentMesh3d(), computeInactiveLMElems(), computeInactiveLMNodes(), computeIncorrectEdgeDroppingInactiveLMNodes(), computeMsmStatistics(), projectPrimaryNodes(), and projectSecondaryNodes().

◆ _projected_secondary_nodes

std::unordered_set<dof_id_type> AutomaticMortarGeneration::_projected_secondary_nodes
private

Debugging container for printing information about fraction of successful projections for secondary nodes.

If !_debug then this should always be empty

Definition at line 612 of file AutomaticMortarGeneration.h.

Referenced by clear(), and projectSecondaryNodesSinglePair().

◆ _secondary_boundary_subdomain_ids

std::set<SubdomainID> AutomaticMortarGeneration::_secondary_boundary_subdomain_ids
private

The secondary/primary lower-dimensional boundary subdomain ids are the secondary/primary boundary ids.

Definition at line 484 of file AutomaticMortarGeneration.h.

Referenced by AutomaticMortarGeneration(), buildMortarSegmentMesh(), buildNodeToElemMaps(), and computeNodalGeometry().

◆ _secondary_element_to_secondary_lowerd_element

std::unordered_map<dof_id_type, const Elem *> AutomaticMortarGeneration::_secondary_element_to_secondary_lowerd_element
private

Map from full dimensional secondary element id to lower dimensional secondary element.

Definition at line 505 of file AutomaticMortarGeneration.h.

Referenced by clear(), computeNodalGeometry(), and getSecondaryLowerdElemFromSecondaryElem().

◆ _secondary_elems_to_mortar_segments

std::unordered_map<dof_id_type, std::set<Elem *, CompareDofObjectsByID> > AutomaticMortarGeneration::_secondary_elems_to_mortar_segments
private

We maintain a mapping from lower-dimensional secondary elements in the original mesh to (sets of) elements in mortar_segment_mesh.

This allows us to quickly determine which elements need to be split.

Definition at line 517 of file AutomaticMortarGeneration.h.

Referenced by buildMortarSegmentMesh(), buildMortarSegmentMesh3d(), clear(), computeMsmStatistics(), and secondariesToMortarSegments().

◆ _secondary_ip_sub_ids

std::set<SubdomainID> AutomaticMortarGeneration::_secondary_ip_sub_ids
private

All the secondary interior parent subdomain IDs associated with the mortar mesh.

Definition at line 520 of file AutomaticMortarGeneration.h.

Referenced by buildMortarSegmentMesh(), buildMortarSegmentMesh3d(), clear(), and secondaryIPSubIDs().

◆ _secondary_node_and_elem_to_xi2_primary_elem

std::unordered_map<std::pair<const Node *, const Elem *>, std::pair<Real, const Elem *> > AutomaticMortarGeneration::_secondary_node_and_elem_to_xi2_primary_elem
private

Similar to the map above, but associates a (Secondary Node, Secondary Elem) pair to a (xi^(2), primary Elem) pair.

This allows a single secondary node, which is potentially connected to two elements on the secondary side, to be associated with multiple primary Elem/xi^(2) values to handle the case where the primary and secondary nodes are "matching". In this configuration:

A B o--—o--—o (secondary orientation ->) | v ---—x---— (primary orientation <-) C D

The entries in the map should be: (Elem A, Node 1) -> (Elem C, xi^(2)=-1) (Elem B, Node 0) -> (Elem D, xi^(2)=+1)

Definition at line 455 of file AutomaticMortarGeneration.h.

Referenced by buildMortarSegmentMesh(), clear(), processAlignedNodes(), and projectSecondaryNodesSinglePair().

◆ _secondary_node_to_hh_nodal_tangents

std::unordered_map<const Node *, std::array<Point, 2> > AutomaticMortarGeneration::_secondary_node_to_hh_nodal_tangents
private

Container for storing the nodal tangent/binormal vectors associated with each secondary node (Householder approach).

Definition at line 502 of file AutomaticMortarGeneration.h.

Referenced by clear(), computeNodalGeometry(), and getNodalTangents().

◆ _secondary_node_to_nodal_normal

std::unordered_map<const Node *, Point> AutomaticMortarGeneration::_secondary_node_to_nodal_normal
private

Container for storing the nodal normal vector associated with each secondary node.

Definition at line 498 of file AutomaticMortarGeneration.h.

Referenced by clear(), computeNodalGeometry(), getNodalNormals(), projectPrimaryNodesSinglePair(), and projectSecondaryNodesSinglePair().

◆ _secondary_requested_boundary_ids

std::set<BoundaryID> AutomaticMortarGeneration::_secondary_requested_boundary_ids
private

The boundary ids corresponding to all the secondary surfaces.

Definition at line 424 of file AutomaticMortarGeneration.h.

Referenced by AutomaticMortarGeneration(), and buildNodeToElemMaps().

◆ _triangulate_triangles

const bool AutomaticMortarGeneration::_triangulate_triangles
private

Whether already-triangular clipped polygons should still be centroid-subdivided.

Definition at line 594 of file AutomaticMortarGeneration.h.

Referenced by buildMortarSegmentMesh3d().

◆ _triangulation_mode

const MortarSegmentTriangulationMode AutomaticMortarGeneration::_triangulation_mode
private

Triangulation mode used for clipped 3D mortar polygons.

Definition at line 591 of file AutomaticMortarGeneration.h.

Referenced by buildMortarSegmentMesh3d().

◆ _xi_tolerance

Real AutomaticMortarGeneration::_xi_tolerance = 1e-6
private

Tolerance for checking projection xi values.

Usually we are checking whether we projected onto a certain element (in which case -1 <= xi <= 1) or whether we should have already projected a primary node (in which case we error if abs(xi) is sufficiently close to 1)

Definition at line 575 of file AutomaticMortarGeneration.h.

Referenced by buildMortarSegmentMesh(), projectPrimaryNodesSinglePair(), and projectSecondaryNodesSinglePair().

◆ system_name

const std::string AutomaticMortarGeneration::system_name
static

The name of the nodal normals system.

We store this in one place so it's easy to change later.

Definition at line 75 of file AutomaticMortarGeneration.h.


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