22 #include "libmesh/mesh_tools.h" 23 #include "libmesh/explicit_system.h" 24 #include "libmesh/numeric_vector.h" 25 #include "libmesh/elem.h" 26 #include "libmesh/node.h" 27 #include "libmesh/dof_map.h" 28 #include "libmesh/edge_edge2.h" 29 #include "libmesh/edge_edge3.h" 30 #include "libmesh/face_tri3.h" 31 #include "libmesh/face_tri6.h" 32 #include "libmesh/face_tri7.h" 33 #include "libmesh/face_quad4.h" 34 #include "libmesh/face_quad8.h" 35 #include "libmesh/face_quad9.h" 36 #include "libmesh/exodusII_io.h" 37 #include "libmesh/quadrature_gauss.h" 38 #include "libmesh/quadrature_nodal.h" 39 #include "libmesh/distributed_mesh.h" 40 #include "libmesh/replicated_mesh.h" 41 #include "libmesh/enum_to_string.h" 42 #include "libmesh/statistics.h" 43 #include "libmesh/equation_systems.h" 45 #include "metaphysicl/dualnumber.h" 58 #if NANOFLANN_VERSION < 0x150 72 params.addPrivateParam<
MooseApp *>(
"_moose_app",
nullptr);
73 params.set<std::string>(
"_type") =
"MortarNodalGeometryOutput";
85 if (_amg._secondary_node_to_nodal_normal.empty() ||
86 _amg._secondary_node_to_hh_nodal_tangents.empty())
87 mooseError(
"No entries found in the secondary node -> nodal geometry map.");
89 auto & problem = _app.feProblem();
90 auto & subproblem = _amg._on_displaced
91 ?
static_cast<SubProblem &
>(*problem.getDisplacedProblem())
92 : static_cast<SubProblem &>(problem);
93 auto & nodal_normals_es = subproblem.
es();
95 const std::string nodal_normals_sys_name =
"nodal_normals";
97 if (!_nodal_normals_system)
99 for (
const auto s :
make_range(nodal_normals_es.n_systems()))
100 if (!nodal_normals_es.get_system(s).is_initialized())
105 _nodal_normals_system =
106 &nodal_normals_es.template add_system<ExplicitSystem>(nodal_normals_sys_name);
107 _nnx_var_num = _nodal_normals_system->add_variable(
"nodal_normal_x",
FEType(
FIRST,
LAGRANGE)),
108 _nny_var_num = _nodal_normals_system->add_variable(
"nodal_normal_y",
FEType(
FIRST,
LAGRANGE));
109 _nnz_var_num = _nodal_normals_system->add_variable(
"nodal_normal_z",
FEType(
FIRST,
LAGRANGE));
124 nodal_normals_es.
reinit();
127 const DofMap & dof_map = _nodal_normals_system->get_dof_map();
128 std::vector<dof_id_type> dof_indices_nnx, dof_indices_nny, dof_indices_nnz;
129 std::vector<dof_id_type> dof_indices_t1x, dof_indices_t1y, dof_indices_t1z;
130 std::vector<dof_id_type> dof_indices_t2x, dof_indices_t2y, dof_indices_t2z;
132 for (MeshBase::const_element_iterator el = _amg._mesh.elements_begin(),
133 end_el = _amg._mesh.elements_end();
137 const Elem * elem = *el;
140 dof_map.
dof_indices(elem, dof_indices_nnx, _nnx_var_num);
141 dof_map.
dof_indices(elem, dof_indices_nny, _nny_var_num);
142 dof_map.
dof_indices(elem, dof_indices_nnz, _nnz_var_num);
144 dof_map.
dof_indices(elem, dof_indices_t1x, _t1x_var_num);
145 dof_map.
dof_indices(elem, dof_indices_t1y, _t1y_var_num);
146 dof_map.
dof_indices(elem, dof_indices_t1z, _t1z_var_num);
148 dof_map.
dof_indices(elem, dof_indices_t2x, _t2x_var_num);
149 dof_map.
dof_indices(elem, dof_indices_t2y, _t2y_var_num);
150 dof_map.
dof_indices(elem, dof_indices_t2z, _t2z_var_num);
158 auto it = _amg._secondary_node_to_nodal_normal.find(elem->
node_ptr(n));
159 if (it != _amg._secondary_node_to_nodal_normal.end())
161 _nodal_normals_system->solution->set(dof_indices_nnx[n], it->second(0));
162 _nodal_normals_system->solution->set(dof_indices_nny[n], it->second(1));
163 _nodal_normals_system->solution->set(dof_indices_nnz[n], it->second(2));
166 auto it_tangent = _amg._secondary_node_to_hh_nodal_tangents.find(elem->
node_ptr(n));
167 if (it_tangent != _amg._secondary_node_to_hh_nodal_tangents.end())
169 _nodal_normals_system->solution->set(dof_indices_t1x[n], it_tangent->second[0](0));
170 _nodal_normals_system->solution->set(dof_indices_t1y[n], it_tangent->second[0](1));
171 _nodal_normals_system->solution->set(dof_indices_t1z[n], it_tangent->second[0](2));
173 _nodal_normals_system->solution->set(dof_indices_t2x[n], it_tangent->second[1](0));
174 _nodal_normals_system->solution->set(dof_indices_t2y[n], it_tangent->second[1](1));
175 _nodal_normals_system->solution->set(dof_indices_t2z[n], it_tangent->second[1](2));
182 _nodal_normals_system->solution->close();
184 std::set<std::string> sys_names = {nodal_normals_sys_name};
193 "nodal_geometry_only.e", nodal_normals_es, &sys_names);
220 const std::pair<BoundaryID, BoundaryID> & boundary_key,
221 const std::pair<SubdomainID, SubdomainID> & subdomain_key,
225 const bool correct_edge_dropping,
226 const Real minimum_projection_angle)
231 _on_displaced(on_displaced),
233 _distributed(!_mesh.is_replicated()),
234 _correct_edge_dropping(correct_edge_dropping),
235 _minimum_projection_angle(minimum_projection_angle)
261 "mortar_nodal_geometry_" +
293 "Must specify secondary and primary boundary ids before building node-to-elem maps.");
296 for (
const auto & secondary_elem :
303 for (
const auto & nd : secondary_elem->node_ref_range())
306 vec.push_back(secondary_elem);
311 for (
const auto & primary_elem :
318 for (
const auto & nd : primary_elem->node_ref_range())
321 vec.push_back(primary_elem);
329 std::vector<Point> nodal_normals(secondary_elem.
n_nodes());
333 return nodal_normals;
341 "Map should locate secondary element");
346 std::map<unsigned int, unsigned int>
349 std::map<unsigned int, unsigned int> secondary_ip_i_to_lower_secondary_i;
351 mooseAssert(secondary_ip,
"This should be non-null");
355 const auto & nd = lower_secondary_elem.
node_ref(i);
356 secondary_ip_i_to_lower_secondary_i[secondary_ip->
get_node_index(&nd)] = i;
359 return secondary_ip_i_to_lower_secondary_i;
362 std::map<unsigned int, unsigned int>
364 const Elem & lower_primary_elem,
365 const Elem & primary_elem,
368 std::map<unsigned int, unsigned int> primary_ip_i_to_lower_primary_i;
372 const auto & nd = lower_primary_elem.
node_ref(i);
373 primary_ip_i_to_lower_primary_i[primary_elem.
get_node_index(&nd)] = i;
376 return primary_ip_i_to_lower_primary_i;
379 std::array<MooseUtils::SemidynamicVector<Point, 9>, 2>
388 const auto & tangent_vectors =
390 nodal_tangents_one.
push_back(tangent_vectors[0]);
391 nodal_tangents_two.
push_back(tangent_vectors[1]);
394 return {{nodal_tangents_one, nodal_tangents_two}};
399 const std::vector<Real> & oned_xi1_pts)
const 401 std::vector<Point> xi1_pts(oned_xi1_pts.size());
403 xi1_pts[qp] = oned_xi1_pts[qp];
410 const std::vector<Point> & xi1_pts)
const 413 const auto num_qps = xi1_pts.size();
415 std::vector<Point> normals(num_qps);
427 normals[qp] += phi * nodal_normals[n];
431 for (
auto & normal : normals)
441 std::size_t node_unique_id_offset = 0;
449 const auto primary_bnd_id = pr.first;
450 const auto secondary_bnd_id = pr.second;
451 const auto num_primary_nodes =
452 std::distance(
_mesh.bid_nodes_begin(primary_bnd_id),
_mesh.bid_nodes_end(primary_bnd_id));
453 const auto num_secondary_nodes = std::distance(
_mesh.bid_nodes_begin(secondary_bnd_id),
454 _mesh.bid_nodes_end(secondary_bnd_id));
455 mooseAssert(num_primary_nodes,
456 "There are no primary nodes on boundary ID " 457 << primary_bnd_id <<
". Does that bondary ID even exist on the mesh?");
458 mooseAssert(num_secondary_nodes,
459 "There are no secondary nodes on boundary ID " 460 << secondary_bnd_id <<
". Does that bondary ID even exist on the mesh?");
462 node_unique_id_offset += num_primary_nodes + 2 * num_secondary_nodes;
466 for (MeshBase::const_element_iterator el =
_mesh.active_elements_begin(),
467 end_el =
_mesh.active_elements_end();
471 const Elem * secondary_elem = *el;
477 std::vector<Node *> new_nodes;
478 for (MooseIndex(secondary_elem->
n_nodes()) n = 0; n < secondary_elem->
n_nodes(); ++n)
482 Node *
const new_node = new_nodes.back();
486 std::unique_ptr<Elem> new_elem;
488 new_elem = std::make_unique<Edge3>();
490 new_elem = std::make_unique<Edge2>();
492 new_elem->processor_id() = secondary_elem->
processor_id();
493 new_elem->subdomain_id() = secondary_elem->
subdomain_id();
494 new_elem->set_id(local_id_index++);
495 new_elem->set_unique_id(new_elem->id());
497 for (MooseIndex(new_elem->n_nodes()) n = 0; n < new_elem->n_nodes(); ++n)
498 new_elem->set_node(n, new_nodes[n]);
509 std::make_pair(secondary_elem->
node_ptr(0), secondary_elem)),
511 std::make_pair(secondary_elem->
node_ptr(1), secondary_elem));
513 bool new_container_node0_found =
515 new_container_node1_found =
518 const Elem * node0_primary_candidate =
nullptr;
519 const Elem * node1_primary_candidate =
nullptr;
521 if (new_container_node0_found)
523 const auto & xi2_primary_elem_pair = new_container_it0->second;
524 msinfo.
xi2_a = xi2_primary_elem_pair.first;
525 node0_primary_candidate = xi2_primary_elem_pair.second;
528 if (new_container_node1_found)
530 const auto & xi2_primary_elem_pair = new_container_it1->second;
531 msinfo.
xi2_b = xi2_primary_elem_pair.first;
532 node1_primary_candidate = xi2_primary_elem_pair.second;
539 if (node0_primary_candidate == node1_primary_candidate)
554 auto val = pr.second;
556 const Node * primary_node = std::get<1>(key);
557 Real xi1 = val.first;
558 const Elem * secondary_elem = val.second;
568 for (MooseIndex(secondary_elem->
n_nodes()) n = 0; n < secondary_elem->
n_nodes(); ++n)
573 Elem * current_mortar_segment =
nullptr;
576 for (
const auto & mortar_segment_candidate : mortar_segment_set)
582 catch (std::out_of_range &)
584 mooseError(
"MortarSegmentInfo not found for the mortar segment candidate");
586 if (
info->xi1_a <= xi1 && xi1 <= info->xi1_b)
588 current_mortar_segment = mortar_segment_candidate;
594 if (current_mortar_segment ==
nullptr)
595 mooseError(
"Unable to find appropriate mortar segment during linear search!");
602 if (
info->xi1_a == xi1 || xi1 ==
info->xi1_b)
607 "new_id must be the same on all processes");
608 Node *
const new_node =
610 new_node->set_unique_id(new_id + node_unique_id_offset);
615 const Point normal =
getNormals(*secondary_elem, std::vector<Real>({xi1}))[0];
620 mooseError(
"We should already have built this primary node to elem pair!");
621 const std::vector<const Elem *> & primary_node_neighbors =
625 if (primary_node_neighbors.size() == 0 || primary_node_neighbors.size() > 2)
626 mooseError(
"We must have either 1 or 2 primary side nodal neighbors, but we had ",
627 primary_node_neighbors.size());
634 const Elem * left_primary_elem = primary_node_neighbors[0];
635 const Elem * right_primary_elem =
636 (primary_node_neighbors.size() == 2) ? primary_node_neighbors[1] :
nullptr;
642 std::array<Real, 2> secondary_node_cps;
643 std::vector<Real> primary_node_cps(primary_node_neighbors.size());
646 for (
unsigned int nid = 0; nid < 2; ++nid)
647 secondary_node_cps[nid] = normal.
cross(secondary_elem->
point(nid) - new_pt)(2);
649 for (MooseIndex(primary_node_neighbors) mnn = 0; mnn < primary_node_neighbors.size(); ++mnn)
651 const Elem * primary_neigh = primary_node_neighbors[mnn];
652 Point opposite = (primary_neigh->
node_ptr(0) == primary_node) ? primary_neigh->
point(1)
653 : primary_neigh->
point(0);
655 primary_node_cps[mnn] = cp(2);
659 bool orientation1_valid =
false, orientation2_valid =
false;
661 if (primary_node_neighbors.size() == 2)
664 orientation1_valid = (secondary_node_cps[0] * primary_node_cps[0] > 0.) &&
665 (secondary_node_cps[1] * primary_node_cps[1] > 0.);
667 orientation2_valid = (secondary_node_cps[0] * primary_node_cps[1] > 0.) &&
668 (secondary_node_cps[1] * primary_node_cps[0] > 0.);
670 else if (primary_node_neighbors.size() == 1)
673 orientation1_valid = (secondary_node_cps[0] * primary_node_cps[0] > 0.);
674 orientation2_valid = (secondary_node_cps[1] * primary_node_cps[0] > 0.);
677 mooseError(
"Invalid primary node neighbors size ", primary_node_neighbors.size());
683 if (orientation1_valid && orientation2_valid)
685 "AutomaticMortarGeneration: Both orientations cannot simultaneously be valid.");
691 if (!orientation1_valid && !orientation2_valid)
694 "AutomaticMortarGeneration: Unable to determine valid secondary-primary orientation. " 695 "Consequently we will consider projection of the primary node invalid and not split the " 697 "This situation can indicate there are very oblique projections between primary (mortar) " 698 "and secondary (non-mortar) surfaces for a good problem set up. It can also mean your " 699 "time step is too large. This message is only printed once."));
704 std::unique_ptr<Elem> new_elem_left;
706 new_elem_left = std::make_unique<Edge3>();
708 new_elem_left = std::make_unique<Edge2>();
710 new_elem_left->processor_id() = current_mortar_segment->
processor_id();
711 new_elem_left->subdomain_id() = current_mortar_segment->
subdomain_id();
712 new_elem_left->set_id(local_id_index++);
713 new_elem_left->set_unique_id(new_elem_left->id());
714 new_elem_left->set_node(0, current_mortar_segment->
node_ptr(0));
715 new_elem_left->set_node(1, new_node);
718 std::unique_ptr<Elem> new_elem_right;
720 new_elem_right = std::make_unique<Edge3>();
722 new_elem_right = std::make_unique<Edge2>();
724 new_elem_right->processor_id() = current_mortar_segment->
processor_id();
725 new_elem_right->subdomain_id() = current_mortar_segment->
subdomain_id();
726 new_elem_right->set_id(local_id_index++);
727 new_elem_right->set_unique_id(new_elem_right->id());
728 new_elem_right->set_node(0, new_node);
729 new_elem_right->set_node(1, current_mortar_segment->
node_ptr(1));
734 Point left_interior_point(0);
735 Real left_interior_xi = (xi1 +
info->xi1_a) / 2;
738 Real current_left_interior_eta =
739 (2. * left_interior_xi -
info->xi1_a -
info->xi1_b) / (
info->xi1_b -
info->xi1_a);
741 for (MooseIndex(current_mortar_segment->
n_nodes()) n = 0;
742 n < current_mortar_segment->
n_nodes();
745 current_mortar_segment->
point(n);
749 "new_id must be the same on all processes");
751 left_interior_point, new_interior_left_id, new_elem_left->processor_id());
752 new_elem_left->set_node(2, new_interior_node_left);
753 new_interior_node_left->set_unique_id(new_interior_left_id + node_unique_id_offset);
756 Point right_interior_point(0);
757 Real right_interior_xi = (xi1 +
info->xi1_b) / 2;
759 Real current_right_interior_eta =
760 (2. * right_interior_xi -
info->xi1_a -
info->xi1_b) / (
info->xi1_b -
info->xi1_a);
762 for (MooseIndex(current_mortar_segment->
n_nodes()) n = 0;
763 n < current_mortar_segment->
n_nodes();
766 current_mortar_segment->
point(n);
770 "new_id must be the same on all processes");
772 right_interior_point, new_interior_id_right, new_elem_right->processor_id());
773 new_elem_right->set_node(2, new_interior_node_right);
774 new_interior_node_right->set_unique_id(new_interior_id_right + node_unique_id_offset);
778 if (orientation2_valid)
779 std::swap(left_primary_elem, right_primary_elem);
783 if (left_primary_elem)
784 left_xi2 = (primary_node == left_primary_elem->
node_ptr(0)) ? -1 : +1;
785 if (right_primary_elem)
786 right_xi2 = (primary_node == right_primary_elem->
node_ptr(0)) ? -1 : +1;
795 mooseError(
"MortarSegmentInfo not found for current_mortar_segment.");
811 new_msinfo_left.
xi1_b = xi1;
812 new_msinfo_left.
xi2_b = left_xi2;
820 mortar_segment_set.insert(msm_new_elem);
830 new_msinfo_right.
xi1_b = current_msinfo.
xi1_b;
831 new_msinfo_right.
xi2_b = current_msinfo.
xi2_b;
833 new_msinfo_right.
xi1_a = xi1;
834 new_msinfo_right.
xi2_a = right_xi2;
839 mortar_segment_set.insert(msm_new_elem);
847 mortar_segment_set.erase(current_mortar_segment);
864 if (primary_elem ==
nullptr ||
std::abs(msinfo.
xi2_a) > 1.0 + TOLERANCE ||
870 "We should have found the element");
871 auto & msm_set = it->second;
872 msm_set.erase(msm_elem);
895 std::unordered_set<Node *> msm_connected_nodes;
900 for (
auto & n : element->node_ref_range())
901 msm_connected_nodes.insert(&n);
904 if (!msm_connected_nodes.count(node))
913 "All mortar segment elements should have valid " 928 mortar_segment_mesh_writer.
write(
"mortar_segment_mesh.e");
947 const auto primary_subd_id = pr.first;
948 const auto secondary_subd_id = pr.second;
953 3, mesh_adaptor, nanoflann::KDTreeSingleIndexAdaptorParams(10));
956 kd_tree.buildIndex();
959 auto get_sub_elem_nodes = [](
const ElemType type,
960 const unsigned int sub_elem) -> std::vector<unsigned int>
967 return {{0, 1, 2, 3}};
981 mooseError(
"get_sub_elem_nodes: Invalid sub_elem: ", sub_elem);
995 return {{4, 5, 6, 7}};
997 mooseError(
"get_sub_elem_nodes: Invalid sub_elem: ", sub_elem);
1003 return {{0, 4, 8, 7}};
1005 return {{4, 1, 5, 8}};
1007 return {{8, 5, 2, 6}};
1009 return {{7, 8, 6, 3}};
1011 mooseError(
"get_sub_elem_nodes: Invalid sub_elem: ", sub_elem);
1014 mooseError(
"get_sub_elem_inds: Face element type: ",
1015 libMesh::Utility::enum_to_string<ElemType>(type),
1016 " invalid for 3D mortar");
1023 for (MeshBase::const_element_iterator el =
_mesh.active_local_elements_begin(),
1024 end_el =
_mesh.active_local_elements_end();
1028 const Elem * secondary_side_elem = *el;
1030 const Real secondary_volume = secondary_side_elem->
volume();
1033 if (secondary_side_elem->
subdomain_id() != secondary_subd_id)
1036 auto [secondary_elem_to_msm_map_it, insertion_happened] =
1038 std::set<Elem *, CompareDofObjectsByID>{});
1040 auto & secondary_to_msm_element_set = secondary_elem_to_msm_map_it->second;
1042 std::vector<std::unique_ptr<MortarSegmentHelper>> mortar_segment_helper(
1059 auto sub_elem_nodes = get_sub_elem_nodes(secondary_side_elem->
type(), sel);
1064 std::vector<Point> nodes(sub_elem_nodes.size());
1067 for (
auto iv :
make_range(sub_elem_nodes.size()))
1069 const auto n = sub_elem_nodes[iv];
1070 nodes[iv] = secondary_side_elem->
point(n);
1071 center += secondary_side_elem->
point(n);
1072 normal += nodal_normals[n];
1074 center /= sub_elem_nodes.size();
1075 normal = normal.
unit();
1078 mortar_segment_helper[sel] = std::make_unique<MortarSegmentHelper>(nodes, center, normal);
1089 std::array<Real, 3> query_pt;
1091 switch (secondary_side_elem->
type())
1095 center_point = mortar_segment_helper[0]->center();
1096 query_pt = {{center_point(0), center_point(1), center_point(2)}};
1100 center_point = mortar_segment_helper[1]->center();
1101 query_pt = {{center_point(0), center_point(1), center_point(2)}};
1104 center_point = mortar_segment_helper[4]->center();
1105 query_pt = {{center_point(0), center_point(1), center_point(2)}};
1108 center_point = secondary_side_elem->
point(8);
1109 query_pt = {{center_point(0), center_point(1), center_point(2)}};
1113 "Face element type: ", secondary_side_elem->
type(),
"not supported for 3D mortar");
1119 const std::size_t num_results = 3;
1122 std::vector<size_t> ret_index(num_results);
1123 std::vector<Real> out_dist_sqr(num_results);
1124 nanoflann::KNNResultSet<Real> result_set(num_results);
1125 result_set.init(&ret_index[0], &out_dist_sqr[0]);
1129 std::set<const Elem *, CompareDofObjectsByID> processed_primary_elems;
1133 bool primary_elem_found =
false;
1134 std::set<const Elem *, CompareDofObjectsByID> primary_elem_candidates;
1141 out_dist_sqr[r]) <= TOLERANCE,
1142 "Lower-dimensional element squared distance verification failed.");
1145 std::vector<const Elem *> & node_elems =
1149 for (
auto elem : node_elems)
1150 primary_elem_candidates.insert(elem);
1160 while (!primary_elem_candidates.empty())
1162 const Elem * primary_elem_candidate = *primary_elem_candidates.begin();
1165 if (processed_primary_elems.count(primary_elem_candidate))
1169 std::vector<Point> nodal_points;
1172 std::vector<std::vector<unsigned int>> elem_to_node_map;
1175 std::vector<std::pair<unsigned int, unsigned int>> sub_elem_map;
1184 auto sub_elem_nodes = get_sub_elem_nodes(primary_elem_candidate->
type(), p_el);
1187 std::vector<Point> primary_sub_elem(sub_elem_nodes.size());
1188 for (
auto iv :
make_range(sub_elem_nodes.size()))
1190 const auto n = sub_elem_nodes[iv];
1191 primary_sub_elem[iv] = primary_elem_candidate->
point(n);
1204 mortar_segment_helper[s_el]->getMortarSegments(
1205 primary_sub_elem, nodal_points, elem_to_node_map);
1208 for (
auto i = sub_elem_map.size(); i < elem_to_node_map.size(); ++i)
1209 sub_elem_map.push_back(std::make_pair(s_el, p_el));
1214 processed_primary_elems.insert(primary_elem_candidate);
1215 primary_elem_candidates.erase(primary_elem_candidate);
1218 if (!elem_to_node_map.empty())
1223 if (!primary_elem_found)
1225 primary_elem_found =
true;
1226 primary_elem_candidates.clear();
1233 if (neighbor ==
nullptr || neighbor->subdomain_id() != primary_subd_id)
1236 if (processed_primary_elems.count(neighbor))
1239 primary_elem_candidates.insert(neighbor);
1245 std::vector<Node *> new_nodes;
1246 for (
auto pt : nodal_points)
1254 std::unique_ptr<Elem> new_elem;
1255 if (elem_to_node_map[el].size() == 3)
1256 new_elem = std::make_unique<Tri3>();
1258 mooseError(
"Active mortar segments only supports TRI elements, 3 nodes expected " 1260 elem_to_node_map[el].size(),
1263 new_elem->processor_id() = secondary_side_elem->
processor_id();
1264 new_elem->subdomain_id() = secondary_side_elem->
subdomain_id();
1265 new_elem->set_id(local_id_index++);
1269 new_elem->set_node(i, new_nodes[elem_to_node_map[el][i]]);
1272 if (new_elem->volume() / secondary_volume < TOLERANCE)
1290 secondary_to_msm_element_set.insert(msm_new_elem);
1304 if (mortar_segment_helper[sel]->remainder() == 1.0)
1306 mooseWarning(
"Some secondary elements on mortar interface were unable to identify" 1307 " a corresponding primary element; this may be expected depending on" 1308 " problem geometry but may indicate a failure of the element search" 1312 if (secondary_to_msm_element_set.empty())
1325 if (msm_el->type() !=
TRI3)
1326 msm_el->subdomain_id()++;
1333 mortar_segment_mesh_writer.
write(
"mortar_segment_mesh.e");
1337 if (msm_el->type() !=
TRI3)
1338 msm_el->subdomain_id()--;
1349 mooseWarning(
"Mortar segment mesh statistics intended for debugging purposes in serial only, " 1350 "parallel will only provide statistics for local mortar segment mesh.");
1357 std::unordered_map<processor_id_type, std::vector<std::pair<dof_id_type, dof_id_type>>>
1368 const Elem * secondary_elem = pr.second.secondary_elem;
1369 const Elem * primary_elem = pr.second.primary_elem;
1372 coupling_info[secondary_elem->
processor_id()].emplace_back(
1381 coupling_info[secondary_elem->
processor_id()].emplace_back(
1390 coupling_info[secondary_elem->
processor_id()].emplace_back(secondary_elem->
id(),
1391 primary_elem->
id());
1415 auto action_functor =
1417 const std::vector<std::pair<dof_id_type, dof_id_type>> & coupling_info)
1419 for (
auto [i, j] : coupling_info)
1429 Moose::out <<
"Mortar Interface Statistics:" << std::endl;
1434 const auto primary_subd_id = pr.first;
1435 const auto secondary_subd_id = pr.second;
1442 for (
auto * el :
_mesh.active_element_ptr_range())
1445 if (el->subdomain_id() == secondary_subd_id)
1446 secondary.push_back(el->volume());
1447 else if (el->subdomain_id() == primary_subd_id)
1448 primary.push_back(el->volume());
1456 msm.push_back(msm_elem->volume());
1460 std::vector<std::string> col_names = {
"mesh",
"n_elems",
"max",
"min",
"median"};
1461 std::vector<std::string> subds = {
"secondary_lower",
"primary_lower",
"mortar_segment"};
1462 std::vector<size_t> n_elems = {secondary.size(), primary.size(), msm.size()};
1472 table.
addData<std::string>(col_names[0], subds[i]);
1473 table.
addData<
size_t>(col_names[1], n_elems[i]);
1479 Moose::out <<
"secondary subdomain: " << secondary_subd_id
1480 <<
" \tprimary subdomain: " << primary_subd_id << std::endl;
1495 const Real tol = (
dim() == 3) ? 0.1 : TOLERANCE;
1497 std::unordered_map<processor_id_type, std::set<dof_id_type>> proc_to_inactive_nodes_set;
1501 std::unordered_set<dof_id_type> inactive_node_ids;
1503 std::unordered_map<const Elem *, Real> active_volume{};
1507 active_volume[el] = 0.;
1515 active_volume[secondary_elem] += msm_elem->
volume();
1523 if (
std::abs(active_volume[el] / el->volume() - 1.0) > tol)
1527 inactive_node_ids.insert(el->node_id(n));
1533 const auto secondary_subd_id = pr.second;
1539 const auto pid = el->processor_id();
1546 for (
const auto n :
make_range(el->n_nodes()))
1548 const auto node_id = el->node_id(n);
1549 if (inactive_node_ids.find(node_id) != inactive_node_ids.end())
1550 proc_to_inactive_nodes_set[pid].insert(node_id);
1558 std::unordered_map<processor_id_type, std::vector<dof_id_type>> proc_to_inactive_nodes_vector;
1559 for (
const auto & proc_set : proc_to_inactive_nodes_set)
1560 proc_to_inactive_nodes_vector[proc_set.first].insert(
1561 proc_to_inactive_nodes_vector[proc_set.first].end(),
1562 proc_set.second.begin(),
1563 proc_set.second.end());
1567 const std::vector<dof_id_type> & sent_data)
1570 mooseError(
"Should not be communicating with self.");
1571 for (
const auto pr : sent_data)
1572 inactive_node_ids.insert(pr);
1577 for (
const auto node_id : inactive_node_ids)
1590 std::unordered_map<processor_id_type, std::set<dof_id_type>> proc_to_active_nodes_set;
1594 std::unordered_set<dof_id_type> active_local_nodes;
1603 active_local_nodes.insert(secondary_elem->
node_id(n));
1609 const auto secondary_subd_id = pr.second;
1622 for (
const auto n :
make_range(el->n_nodes()))
1624 const auto node_id = el->node_id(n);
1625 if (active_local_nodes.find(node_id) != active_local_nodes.end())
1626 proc_to_active_nodes_set[pid].insert(node_id);
1634 std::unordered_map<processor_id_type, std::vector<dof_id_type>> proc_to_active_nodes_vector;
1635 for (
const auto & proc_set : proc_to_active_nodes_set)
1637 proc_to_active_nodes_vector[proc_set.first].reserve(proc_to_active_nodes_set.size());
1638 for (
const auto node_id : proc_set.second)
1639 proc_to_active_nodes_vector[proc_set.first].push_back(node_id);
1644 const std::vector<dof_id_type> & sent_data)
1647 mooseError(
"Should not be communicating with self.");
1648 active_local_nodes.insert(sent_data.begin(), sent_data.end());
1659 for (
const auto n :
make_range(el->n_nodes()))
1660 if (active_local_nodes.find(el->node_id(n)) == active_local_nodes.end())
1669 std::unordered_set<const Elem *> active_local_elems;
1675 const Real tol = (
dim() == 3) ? 0.1 : TOLERANCE;
1677 std::unordered_map<const Elem *, Real> active_volume;
1686 active_volume[secondary_elem] += msm_elem->
volume();
1697 if (
std::abs(active_volume[secondary_elem] / secondary_elem->
volume() - 1.0) > tol)
1701 active_local_elems.insert(secondary_elem);
1709 if (active_local_elems.find(el) == active_local_elems.end())
1727 std::map<dof_id_type, std::vector<std::pair<Point, Real>>> node_to_normals_map;
1735 for (MeshBase::const_element_iterator el =
_mesh.active_elements_begin(),
1736 end_el =
_mesh.active_elements_end();
1740 const Elem * secondary_elem = *el;
1750 nnx_fe_face->attach_quadrature_rule(&qface);
1751 const std::vector<Point> & face_normals = nnx_fe_face->get_normals();
1753 const auto & JxW = nnx_fe_face->get_JxW();
1758 mooseAssert(interior_parent,
1759 "No interior parent exists for element " 1760 << secondary_elem->
id()
1761 <<
". There may be a problem with your sideset set-up.");
1772 nnx_fe_face->reinit(interior_parent, s);
1774 for (MooseIndex(secondary_elem->
n_nodes()) n = 0; n < secondary_elem->
n_nodes(); ++n)
1776 auto & normals_and_weights_vec = node_to_normals_map[secondary_elem->
node_id(n)];
1777 normals_and_weights_vec.push_back(std::make_pair(
sign * face_normals[n], JxW[n]));
1784 for (
const auto & pr : node_to_normals_map)
1787 const auto & node_id = pr.first;
1788 const auto & normals_and_weights_vec = pr.second;
1791 for (
const auto & norm_and_weight : normals_and_weights_vec)
1792 nodal_normal += norm_and_weight.first * norm_and_weight.second;
1793 nodal_normal = nodal_normal.
unit();
1797 Point nodal_tangent_one;
1798 Point nodal_tangent_two;
1808 Point & nodal_tangent_one,
1809 Point & nodal_tangent_two)
const 1812 "The input nodal normal should have unity norm");
1814 const Real nx = nodal_normal(0);
1815 const Real ny = nodal_normal(1);
1816 const Real nz = nodal_normal(2);
1821 const Point h_vector(nx + 1.0, ny, nz);
1826 if (
std::abs(h_vector(0)) < TOLERANCE)
1828 nodal_tangent_one(0) = 0;
1829 nodal_tangent_one(1) = 1;
1830 nodal_tangent_one(2) = 0;
1832 nodal_tangent_two(0) = 0;
1833 nodal_tangent_two(1) = 0;
1834 nodal_tangent_two(2) = -1;
1841 nodal_tangent_one(0) = -2.0 * h_vector(0) * h_vector(1) / (h * h);
1842 nodal_tangent_one(1) = 1.0 - 2.0 * h_vector(1) * h_vector(1) / (h * h);
1843 nodal_tangent_one(2) = -2.0 * h_vector(1) * h_vector(2) / (h * h);
1845 nodal_tangent_two(0) = -2.0 * h_vector(0) * h_vector(2) / (h * h);
1846 nodal_tangent_two(1) = -2.0 * h_vector(1) * h_vector(2) / (h * h);
1847 nodal_tangent_two(2) = 1.0 - 2.0 * h_vector(2) * h_vector(2) / (h * h);
1863 SubdomainID lower_dimensional_primary_subdomain_id,
1864 SubdomainID lower_dimensional_secondary_subdomain_id)
1869 3, mesh_adaptor, nanoflann::KDTreeSingleIndexAdaptorParams(10));
1872 kd_tree.buildIndex();
1874 for (MeshBase::const_element_iterator el =
_mesh.active_elements_begin(),
1875 end_el =
_mesh.active_elements_end();
1879 const Elem * secondary_side_elem = *el;
1882 if (secondary_side_elem->
subdomain_id() != lower_dimensional_secondary_subdomain_id)
1889 for (MooseIndex(secondary_side_elem->
n_vertices()) n = 0; n < secondary_side_elem->
n_vertices();
1892 const Node * secondary_node = secondary_side_elem->
node_ptr(n);
1896 const std::vector<const Elem *> & secondary_node_neighbors =
1901 bool is_mapped =
true;
1902 for (MooseIndex(secondary_node_neighbors) snn = 0; snn < secondary_node_neighbors.size();
1905 auto secondary_key = std::make_pair(secondary_node, secondary_node_neighbors[snn]);
1921 std::array<Real, 3> query_pt = {
1922 {(*secondary_node)(0), (*secondary_node)(1), (*secondary_node)(2)}};
1928 const std::size_t num_results = 3;
1931 std::vector<size_t> ret_index(num_results);
1932 std::vector<Real> out_dist_sqr(num_results);
1933 nanoflann::KNNResultSet<Real> result_set(num_results);
1934 result_set.init(&ret_index[0], &out_dist_sqr[0]);
1938 bool projection_succeeded =
false;
1942 std::set<const Elem *> rejected_primary_elem_candidates;
1948 for (MooseIndex(result_set) r = 0; r < result_set.size(); ++r)
1952 out_dist_sqr[r]) <= TOLERANCE,
1953 "Lower-dimensional element squared distance verification failed.");
1957 std::vector<const Elem *> & primary_elem_candidates =
1961 for (MooseIndex(primary_elem_candidates) e = 0; e < primary_elem_candidates.size(); ++e)
1963 const Elem * primary_elem_candidate = primary_elem_candidates[e];
1966 if (rejected_primary_elem_candidates.count(primary_elem_candidate))
1972 unsigned int current_iterate = 0, max_iterates = 10;
1978 for (MooseIndex(primary_elem_candidate->
n_nodes()) n = 0;
1979 n < primary_elem_candidate->
n_nodes();
1983 const auto u = x2 - (*secondary_node);
1984 const auto F = u(0) * nodal_normal(1) - u(1) * nodal_normal(0);
1989 if (F.derivatives())
1991 Real dxi2 = -F.value() / F.derivatives();
1999 current_iterate = max_iterates;
2000 }
while (++current_iterate < max_iterates);
2002 Real xi2 = xi2_dn.value();
2008 (primary_elem_candidate->
point(0) - primary_elem_candidate->
point(1)).unit() *
2027 const Node * primary_node = (xi2 < 0) ? primary_elem_candidate->
node_ptr(0)
2028 : primary_elem_candidate->
node_ptr(1);
2030 const std::vector<const Elem *> & primary_node_neighbors =
2033 std::vector<bool> primary_elems_mapped(primary_node_neighbors.size(),
false);
2041 std::vector<Real> secondary_node_neighbor_cps(2), primary_node_neighbor_cps(2);
2045 for (MooseIndex(secondary_node_neighbors) nn = 0;
2046 nn < secondary_node_neighbors.size();
2049 const Elem * secondary_neigh = secondary_node_neighbors[nn];
2050 Point opposite = (secondary_neigh->
node_ptr(0) == secondary_node)
2051 ? secondary_neigh->
point(1)
2052 : secondary_neigh->
point(0);
2053 Point cp = nodal_normal.
cross(opposite - *secondary_node);
2054 secondary_node_neighbor_cps[nn] = cp(2);
2059 for (MooseIndex(primary_node_neighbors) nn = 0; nn < primary_node_neighbors.size();
2062 const Elem * primary_neigh = primary_node_neighbors[nn];
2063 Point opposite = (primary_neigh->
node_ptr(0) == primary_node)
2064 ? primary_neigh->
point(1)
2065 : primary_neigh->
point(0);
2066 Point cp = nodal_normal.
cross(opposite - *secondary_node);
2067 primary_node_neighbor_cps[nn] = cp(2);
2071 bool found_match =
false;
2072 for (MooseIndex(secondary_node_neighbors) snn = 0;
2073 snn < secondary_node_neighbors.size();
2075 for (MooseIndex(primary_node_neighbors) mnn = 0;
2076 mnn < primary_node_neighbors.size();
2078 if (secondary_node_neighbor_cps[snn] * primary_node_neighbor_cps[mnn] > 0)
2081 primary_elems_mapped[mnn] =
true;
2085 Real xi2 = (primary_node == primary_node_neighbors[mnn]->node_ptr(0)) ? -1 : +1;
2086 auto secondary_key =
2087 std::make_pair(secondary_node, secondary_node_neighbors[snn]);
2088 auto primary_val = std::make_pair(xi2, primary_node_neighbors[mnn]);
2094 (secondary_node == secondary_node_neighbors[snn]->node_ptr(0)) ? -1 : +1;
2096 auto primary_key = std::make_tuple(
2097 primary_node->
id(), primary_node, primary_node_neighbors[mnn]);
2098 auto secondary_val = std::make_pair(xi1, secondary_node_neighbors[snn]);
2107 rejected_primary_elem_candidates.insert(primary_elem_candidate);
2114 if (secondary_node_neighbors.size() == 1 && primary_node_neighbors.size() == 2)
2115 for (
auto it = primary_elems_mapped.begin(); it != primary_elems_mapped.end(); ++it)
2118 auto index = std::distance(primary_elems_mapped.begin(), it);
2121 primary_node->
id(), primary_node, primary_node_neighbors[index]),
2122 std::make_pair(1,
nullptr));
2128 for (MooseIndex(secondary_node_neighbors) nn = 0;
2129 nn < secondary_node_neighbors.size();
2132 const Elem * neigh = secondary_node_neighbors[nn];
2136 if (secondary_node == neigh_node)
2138 auto key = std::make_pair(neigh_node, neigh);
2139 auto val = std::make_pair(xi2, primary_elem_candidate);
2146 projection_succeeded =
true;
2151 rejected_primary_elem_candidates.insert(primary_elem_candidate);
2154 if (projection_succeeded)
2158 if (!projection_succeeded &&
_debug)
2159 _console <<
"Failed to find primary Elem into which secondary node " 2160 <<
static_cast<const Point &
>(*secondary_node) <<
" was projected." << std::endl
2179 SubdomainID lower_dimensional_primary_subdomain_id,
2180 SubdomainID lower_dimensional_secondary_subdomain_id)
2185 3, mesh_adaptor, nanoflann::KDTreeSingleIndexAdaptorParams(10));
2188 kd_tree.buildIndex();
2190 std::unordered_set<dof_id_type> primary_nodes_visited;
2192 for (
const auto & primary_side_elem :
_mesh.active_element_ptr_range())
2195 if (primary_side_elem->subdomain_id() != lower_dimensional_primary_subdomain_id)
2200 for (MooseIndex(primary_side_elem->n_vertices()) n = 0; n < primary_side_elem->n_vertices();
2204 const Node * primary_node = primary_side_elem->node_ptr(n);
2207 const std::vector<const Elem *> & primary_node_neighbors =
2215 std::make_tuple(primary_node->
id(), primary_node, primary_node_neighbors[0]);
2216 if (!primary_nodes_visited.insert(primary_node->
id()).second ||
2221 Real query_pt[3] = {(*primary_node)(0), (*primary_node)(1), (*primary_node)(2)};
2227 const size_t num_results = 3;
2230 std::vector<size_t> ret_index(num_results);
2231 std::vector<Real> out_dist_sqr(num_results);
2232 nanoflann::KNNResultSet<Real> result_set(num_results);
2233 result_set.init(&ret_index[0], &out_dist_sqr[0]);
2237 bool projection_succeeded =
false;
2242 std::set<const Elem *> rejected_secondary_elem_candidates;
2246 for (MooseIndex(result_set) r = 0; r < result_set.size(); ++r)
2250 out_dist_sqr[r]) <= TOLERANCE,
2251 "Lower-dimensional element squared distance verification failed.");
2255 const std::vector<const Elem *> & secondary_elem_candidates =
2259 for (MooseIndex(secondary_elem_candidates) e = 0; e < secondary_elem_candidates.size(); ++e)
2261 const Elem * secondary_elem_candidate = secondary_elem_candidates[e];
2264 if (rejected_secondary_elem_candidates.count(secondary_elem_candidate))
2267 std::vector<Point> nodal_normals(secondary_elem_candidate->
n_nodes());
2277 unsigned int current_iterate = 0, max_iterates = 10;
2288 for (MooseIndex(secondary_elem_candidate->
n_nodes()) n = 0;
2289 n < secondary_elem_candidate->
n_nodes();
2293 x1 += phi * secondary_elem_candidate->
point(n);
2294 normals += phi * nodal_normals[n];
2297 const auto u = x1 - (*primary_node);
2299 const auto F = u(0) * normals(1) - u(1) * normals(0);
2307 Real dxi1 = -F.value() / F.derivatives();
2312 }
while (++current_iterate < max_iterates);
2314 Real xi1 = xi1_dn.value();
2319 (
std::abs((primary_side_elem->point(0) - primary_side_elem->point(1)).unit() *
2331 throw MooseException(
"Nodes on primary and secondary surfaces are aligned. This is " 2332 "causing trouble when identifying projections from secondary " 2333 "nodes when performing primary node projections.");
2344 const Elem * neigh = primary_node_neighbors[0];
2348 if (primary_node == neigh_node)
2350 auto key = std::make_tuple(neigh_node->
id(), neigh_node, neigh);
2351 auto val = std::make_pair(xi1, secondary_elem_candidate);
2357 projection_succeeded =
true;
2363 rejected_secondary_elem_candidates.insert(secondary_elem_candidate);
2367 if (projection_succeeded)
2371 if (!projection_succeeded &&
_debug)
2373 _console <<
"Failed to find point from which primary node " 2374 <<
static_cast<const Point &
>(*primary_node) <<
" was projected." << std::endl
2381 std::vector<AutomaticMortarGeneration::MortarFilterIter>
2388 const auto & secondary_elems = secondary_it->second;
2389 std::vector<MortarFilterIter> ret;
2390 ret.reserve(secondary_elems.size());
2394 auto *
const secondary_elem = secondary_elems[i];
2400 mooseAssert(secondary_elem->active(),
2401 "We loop over active elements when building the mortar segment mesh, so we golly " 2402 "well hope this is active.");
2403 mooseAssert(!msm_it->second.empty(),
2404 "We should have removed all secondaries from this map if they do not have any " 2405 "mortar segments associated with them.");
2406 ret.push_back(msm_it);
virtual T maximum() const
std::set< SubdomainID > _primary_boundary_subdomain_ids
MetaPhysicL::DualNumber< V, D, asd > abs(const MetaPhysicL::DualNumber< V, D, asd > &a)
const Elem * getSecondaryLowerdElemFromSecondaryElem(dof_id_type secondary_elem_id) const
Return lower dimensional secondary element given its interior parent.
T fe_lagrange_2D_shape(const libMesh::ElemType type, const Order order, const unsigned int i, const VectorType< T > &p)
std::unique_ptr< FEGenericBase< Real > > build(const unsigned int dim, const FEType &fet)
unique_id_type & set_unique_id()
std::unordered_set< const Elem * > _inactive_local_lm_elems
List of inactive lagrange multiplier nodes (for elemental variables)
auto norm() const -> decltype(std::norm(Real()))
bool absoluteFuzzyEqual(const T &var1, const T2 &var2, const T3 &tol=libMesh::TOLERANCE *libMesh::TOLERANCE)
Function to check whether two variables are equal within an absolute tolerance.
unsigned int get_node_index(const Node *node_ptr) const
void computeInactiveLMNodes()
Get list of secondary nodes that don't contribute to interaction with any primary element...
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
unsigned int _t2x_var_num
void clear()
Clears the mortar segment mesh and accompanying data structures.
void dof_indices(const Elem *const elem, std::vector< dof_id_type > &di) const
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...
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::map< unsigned int, unsigned int > 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 ...
unsigned int _t2z_var_num
const bool _periodic
Whether this object will be generating a mortar segment mesh for periodic constraints.
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.
unsigned int which_side_am_i(const Elem *e) const
Real _newton_tolerance
Newton solve tolerance for node projections.
void mooseWarning(Args &&... args)
Emit a warning message with the given stringified, concatenated args.
void output() override
Overload this function with the desired output activities.
MortarNodalGeometryOutput(const InputParameters ¶ms)
void push_back(const T &v)
std::unordered_map< const Elem *, MortarSegmentInfo > _msm_elem_to_info
Map between Elems in the mortar segment mesh and their info structs.
const bool _distributed
Whether the mortar segment mesh is distributed.
Base class for MOOSE-based applications.
void buildMortarSegmentMesh()
Builds the mortar segment mesh once the secondary and primary node projections have been completed...
const Parallel::Communicator & comm() const
void msmStatistics()
Outputs mesh statistics for mortar segment mesh.
unsigned int _nnx_var_num
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.
The following methods are specializations for using the libMesh::Parallel::packed_range_* routines fo...
static const Real invalid_xi
unsigned int _nny_var_num
Specialization of SubProblem for solving nonlinear equations plus auxiliary equations.
FEProblemBase & feProblem() const
unsigned int _nnz_var_num
void projectSecondaryNodes()
Project secondary nodes (find xi^(2) values) to the closest points on the primary surface...
const std::unordered_map< dof_id_type, std::set< Elem *, CompareDofObjectsByID > > & secondariesToMortarSegments() const
Utility class template for a semidynamic vector with a maximum size N and a chosen dynamic size...
This class is a container/interface for the objects involved in automatic generation of mortar spaces...
const bool _debug
Whether to print debug output.
Based class for output objects.
const Elem * primary_elem
void push_parallel_vector_data(const Communicator &comm, MapToVectors &&data, const ActionFunctor &act_on_data)
uint8_t processor_id_type
std::set< SubdomainID > _secondary_ip_sub_ids
All the secondary interior parent subdomain IDs associated with the mortar mesh.
virtual libMesh::EquationSystems & es()=0
processor_id_type n_processors() const
TypeVector< Real > unit() const
void libmesh_ignore(const Args &...)
void projectPrimaryNodes()
(Inverse) project primary nodes to the points on the secondary surface where they would have come fro...
unsigned int _t1z_var_num
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
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)
Must be constructed with a reference to the Mesh we are generating mortar spaces for.
virtual unsigned int n_nodes() const=0
unsigned int _t2y_var_num
void buildNodeToElemMaps()
Once the secondary_requested_boundary_ids and primary_requested_boundary_ids containers have been fil...
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.
virtual void write_equation_systems(const std::string &fname, const EquationSystems &es, const std::set< std::string > *system_names=nullptr) override
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)...
void projectPrimaryNodesSinglePair(SubdomainID lower_dimensional_primary_subdomain_id, SubdomainID lower_dimensional_secondary_subdomain_id)
Helper function used internally by AutomaticMortarGeneration::project_primary_nodes().
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)
const bool _correct_edge_dropping
Flag to enable regressed treatment of edge dropping where all LM DoFs on edge dropping element are st...
An inteface for the _console for outputting to the Console object.
unsigned int _t1y_var_num
std::array< MooseUtils::SemidynamicVector< Point, 9 >, 2 > getNodalTangents(const Elem &secondary_elem) const
Compute the two nodal tangents, which are built on-the-fly.
std::unique_ptr< InputParameters > _output_params
Storage for the input parameters used by the mortar nodal geometry output.
std::set< BoundaryID > _primary_requested_boundary_ids
The boundary ids corresponding to all the primary surfaces.
std::map< unsigned int, unsigned int > getSecondaryIpToLowerElementMap(const Elem &lower_secondary_elem) const
Compute on-the-fly mapping from secondary interior parent nodes to lower dimensional nodes...
TypeVector< typename CompareTypes< Real, T2 >::supertype > cross(const TypeVector< T2 > &v) const
nanoflann::KDTreeSingleIndexAdaptor< subdomain_adatper_t, NanoflannMeshSubdomainAdaptor< 3 >, 3 > subdomain_kd_tree_t
void computeNodalGeometry()
Computes and stores the nodal normal/tangent vectors in a local data structure instead of using the E...
Real _xi_tolerance
Tolerance for checking projection xi values.
static InputParameters validParams()
void computeIncorrectEdgeDroppingInactiveLMNodes()
Computes inactive secondary nodes when incorrect edge dropping behavior is enabled (any node touching...
void buildMortarSegmentMesh3d()
Builds the mortar segment mesh once the secondary and primary node projections have been completed...
const Elem * secondary_elem
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...
Provides a way for users to bail out of the current solve.
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)...
Holds xi^(1), xi^(2), and other data for a given mortar segment.
Generic class for solving transient nonlinear problems.
virtual SimpleRange< element_iterator > active_local_subdomain_elements_ptr_range(subdomain_id_type sid)=0
subdomain_id_type subdomain_id() const
const bool _on_displaced
Whether this object is on the displaced mesh.
const Node * node_ptr(const unsigned int i) const
virtual void write(const std::string &fname) override
void initOutput()
initialize mortar-mesh based output
void addOutput(std::shared_ptr< Output > output)
Adds an existing output object to the warehouse.
virtual Real volume() 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)
unsigned int _t1x_var_num
unsigned int mesh_dimension() const
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
virtual T minimum() 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
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...
const ConsoleStream _console
An instance of helper class to write streams to the Console objects.
SimpleRange< NeighborPtrIter > neighbor_ptr_range()
virtual unsigned int n_sub_elem() const=0
std::set< BoundaryID > _secondary_requested_boundary_ids
The boundary ids corresponding to all the secondary surfaces.
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
std::unordered_set< const Node * > _inactive_local_lm_nodes
processor_id_type processor_id() const
virtual Order default_order() const=0
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.
SearchParams SearchParameters
AutomaticMortarGeneration & _amg
The mortar generation object that we will query for nodal normal and tangent information.
void set_hdf5_writing(bool write_hdf5)
processor_id_type processor_id() const
void householderOrthogolization(const Point &normal, Point &tangent_one, Point &tangent_two) const
Householder orthogonalization procedure to obtain proper basis for tangent and binormal vectors...
virtual ElemType type() const=0
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.
static InputParameters validParams()
MooseApp & _app
The Moose app.
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)
OutputWarehouse & getOutputWarehouse()
Get the OutputWarehouse objects.
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 ...
void computeInactiveLMElems()
Get list of secondary elems without any corresponding primary elements.