Line data Source code
1 : //* This file is part of the MOOSE framework
2 : //* https://mooseframework.inl.gov
3 : //*
4 : //* All rights reserved, see COPYRIGHT for full restrictions
5 : //* https://github.com/idaholab/moose/blob/master/COPYRIGHT
6 : //*
7 : //* Licensed under LGPL 2.1, please see LICENSE for details
8 : //* https://www.gnu.org/licenses/lgpl-2.1.html
9 :
10 : #pragma once
11 :
12 : #include "MortarSegmentInfo.h"
13 : #include "MooseHashing.h"
14 : #include "ConsoleStreamInterface.h"
15 : #include "MooseError.h"
16 : #include "MooseUtils.h"
17 :
18 : // libMesh includes
19 : #include "libmesh/id_types.h"
20 : #include "libmesh/equation_systems.h"
21 : #include "libmesh/elem.h"
22 : #include "libmesh/int_range.h"
23 :
24 : // C++ includes
25 : #include <set>
26 : #include <memory>
27 : #include <vector>
28 : #include <unordered_map>
29 :
30 : // Forward declarations
31 : namespace libMesh
32 : {
33 : class MeshBase;
34 : class System;
35 : }
36 : class GetPot;
37 :
38 : // Using statements
39 : using libMesh::boundary_id_type;
40 : using libMesh::CompareDofObjectsByID;
41 : using libMesh::dof_id_type;
42 : using libMesh::Elem;
43 : using libMesh::MeshBase;
44 : using libMesh::Node;
45 : using libMesh::Point;
46 : using libMesh::Real;
47 : using libMesh::subdomain_id_type;
48 :
49 : typedef boundary_id_type BoundaryID;
50 : typedef subdomain_id_type SubdomainID;
51 :
52 : /**
53 : * This class is a container/interface for the objects involved in
54 : * automatic generation of mortar spaces.
55 : */
56 : class AutomaticMortarGeneration : public ConsoleStreamInterface
57 : {
58 : public:
59 : /**
60 : * The name of the nodal normals system. We store this in one place
61 : * so it's easy to change later.
62 : */
63 : const static std::string system_name;
64 :
65 : /**
66 : * Must be constructed with a reference to the Mesh we are
67 : * generating mortar spaces for.
68 : */
69 : AutomaticMortarGeneration(MooseApp & app,
70 : MeshBase & mesh_in,
71 : const std::pair<BoundaryID, BoundaryID> & boundary_key,
72 : const std::pair<SubdomainID, SubdomainID> & subdomain_key,
73 : bool on_displaced,
74 : bool periodic,
75 : const bool debug,
76 : const bool correct_edge_dropping,
77 : const Real minimum_projection_angle);
78 :
79 : /**
80 : * Once the secondary_requested_boundary_ids and
81 : * primary_requested_boundary_ids containers have been filled in,
82 : * call this function to build node-to-Elem maps for the
83 : * lower-dimensional elements.
84 : */
85 : void buildNodeToElemMaps();
86 :
87 : /**
88 : * Computes and stores the nodal normal/tangent vectors in a local data
89 : * structure instead of using the ExplicitSystem/NumericVector
90 : * approach. This design was triggered by the way that the
91 : * GhostingFunctor operates, but I think it is a better/more
92 : * efficient way to do it anyway.
93 : */
94 : void computeNodalGeometry();
95 :
96 : /**
97 : * Project secondary nodes (find xi^(2) values) to the closest points on
98 : * the primary surface.
99 : * Inputs:
100 : * - The nodal normals values
101 : * - mesh
102 : * - nodes_to_primary_elem_map
103 : *
104 : * Outputs:
105 : * - secondary_node_and_elem_to_xi2_primary_elem
106 : *
107 : * Defined in the file project_secondary_nodes.C.
108 : */
109 : void projectSecondaryNodes();
110 :
111 : /**
112 : * (Inverse) project primary nodes to the points on the secondary surface
113 : * where they would have come from (find (xi^(1) values)).
114 : *
115 : * Inputs:
116 : * - The nodal normals values
117 : * - mesh
118 : * - nodes_to_secondary_elem_map
119 : *
120 : * Outputs:
121 : * - primary_node_and_elem_to_xi1_secondary_elem
122 : *
123 : * Defined in the file project_primary_nodes.C.
124 : */
125 : void projectPrimaryNodes();
126 :
127 : /**
128 : * Builds the mortar segment mesh once the secondary and primary node
129 : * projections have been completed.
130 : *
131 : * Inputs:
132 : * - mesh
133 : * - primary_node_and_elem_to_xi1_secondary_elem
134 : * - secondary_node_and_elem_to_xi2_primary_elem
135 : * - nodes_to_primary_elem_map
136 : *
137 : * Outputs:
138 : * - mortar_segment_mesh
139 : * - msm_elem_to_info
140 : *
141 : * Defined in the file build_mortar_segment_mesh.C.
142 : */
143 : void buildMortarSegmentMesh();
144 :
145 : /**
146 : * Builds the mortar segment mesh once the secondary and primary node
147 : * projections have been completed.
148 : *
149 : * Inputs:
150 : * - mesh
151 : *
152 : * Outputs:
153 : * - mortar_segment_mesh
154 : * - msm_elem_to_info
155 : */
156 : void buildMortarSegmentMesh3d();
157 :
158 : /**
159 : * Outputs mesh statistics for mortar segment mesh
160 : */
161 : void msmStatistics();
162 :
163 : /**
164 : * Clears the mortar segment mesh and accompanying data structures
165 : */
166 : void clear();
167 :
168 : /**
169 : * returns whether this object is on the displaced mesh
170 : */
171 : bool onDisplaced() const { return _on_displaced; }
172 :
173 : /**
174 : * @return The nodal normals associated with the provided \p secondary_elem
175 : */
176 : std::vector<Point> getNodalNormals(const Elem & secondary_elem) const;
177 :
178 : /**
179 : * Compute the two nodal tangents, which are built on-the-fly.
180 : * @return The nodal tangents associated with the provided \p secondary_elem
181 : */
182 : std::array<MooseUtils::SemidynamicVector<Point, 9>, 2>
183 : getNodalTangents(const Elem & secondary_elem) const;
184 :
185 : /**
186 : * Compute on-the-fly mapping from secondary interior parent nodes to lower dimensional nodes
187 : * @return The map from secondary interior parent nodes to lower dimensional nodes
188 : */
189 : std::map<unsigned int, unsigned int>
190 : getSecondaryIpToLowerElementMap(const Elem & lower_secondary_elem) const;
191 :
192 : /**
193 : * Compute on-the-fly mapping from primary interior parent nodes to its corresponding lower
194 : * dimensional nodes
195 : * @return The map from primary interior parent nodes to its corresponding lower dimensional
196 : * nodes
197 : */
198 : std::map<unsigned int, unsigned int>
199 : getPrimaryIpToLowerElementMap(const Elem & primary_elem,
200 : const Elem & primary_elem_ip,
201 : const Elem & lower_secondary_elem) const;
202 :
203 : /**
204 : * Compute the normals at given reference points on a secondary element
205 : * @param secondary_elem The secondary element used to query for associated nodal normals
206 : * @param xi1_pts The reference points on the secondary element to evaluate the normals at. The
207 : * points should only be non-zero in the zeroth entry because right now our mortar mesh elements
208 : * are always 1D
209 : * @return The normals
210 : */
211 : std::vector<Point> getNormals(const Elem & secondary_elem,
212 : const std::vector<Point> & xi1_pts) const;
213 :
214 : /**
215 : * Compute the normals at given reference points on a secondary element
216 : * @param secondary_elem The secondary element used to query for associated nodal normals
217 : * @param 1d_xi1_pts The reference points on the secondary element to evaluate the normals at. The
218 : * "points" are single reals corresponding to xi because right now our mortar mesh elements are
219 : * always 1D
220 : * @return The normals
221 : */
222 : std::vector<Point> getNormals(const Elem & secondary_elem,
223 : const std::vector<Real> & oned_xi1_pts) const;
224 :
225 : /**
226 : * Return lower dimensional secondary element given its interior parent. Helpful outside the
227 : * mortar generation to locate mortar-related quantities.
228 : * @param secondary_elem_id The secondary interior parent element id used to query for associated
229 : * lower dimensional element
230 : * @return The corresponding lower dimensional secondary element
231 : */
232 : const Elem * getSecondaryLowerdElemFromSecondaryElem(dof_id_type secondary_elem_id) const;
233 :
234 : /**
235 : * Get list of secondary nodes that don't contribute to interaction with any primary element.
236 : * Used to enforce zero values on inactive DoFs of nodal variables.
237 : */
238 : void computeInactiveLMNodes();
239 :
240 : /**
241 : * Computes inactive secondary nodes when incorrect edge dropping behavior is enabled
242 : * (any node touching a partially or fully dropped element is dropped)
243 : */
244 : void computeIncorrectEdgeDroppingInactiveLMNodes();
245 :
246 : /**
247 : * Get list of secondary elems without any corresponding primary elements.
248 : * Used to enforce zero values on inactive DoFs of elemental variables.
249 : */
250 : void computeInactiveLMElems();
251 :
252 : /**
253 : * @return The mortar interface coupling
254 : */
255 : const std::unordered_map<dof_id_type, std::unordered_set<dof_id_type>> &
256 349378 : mortarInterfaceCoupling() const
257 : {
258 349378 : return _mortar_interface_coupling;
259 : }
260 :
261 : /**
262 : * @return The primary-secondary boundary ID pair
263 : */
264 : const std::pair<BoundaryID, BoundaryID> & primarySecondaryBoundaryIDPair() const;
265 :
266 : /**
267 : * @return The mortar segment mesh
268 : */
269 1488 : const MeshBase & mortarSegmentMesh() const { return *_mortar_segment_mesh; }
270 :
271 : /**
272 : * @return The mortar segment element to corresponding information
273 : */
274 526172 : const std::unordered_map<const Elem *, MortarSegmentInfo> & mortarSegmentMeshElemToInfo() const
275 : {
276 526172 : return _msm_elem_to_info;
277 : }
278 :
279 23185 : int dim() const { return _mesh.mesh_dimension(); }
280 :
281 : /**
282 : * @return The set of nodes on which mortar constraints are not active
283 : */
284 22641 : const std::unordered_set<const Node *> & getInactiveLMNodes() const
285 : {
286 22641 : return _inactive_local_lm_nodes;
287 : }
288 :
289 : /**
290 : * @return The list of secondary elems on which mortar constraint is not active
291 : */
292 15741 : const std::unordered_set<const Elem *> & getInactiveLMElems() const
293 : {
294 15741 : return _inactive_local_lm_elems;
295 : }
296 :
297 15741 : bool incorrectEdgeDropping() const { return !_correct_edge_dropping; }
298 :
299 : using MortarFilterIter =
300 : std::unordered_map<dof_id_type, std::set<Elem *, CompareDofObjectsByID>>::const_iterator;
301 :
302 : /**
303 : * @return A vector of iterators that point to the lower dimensional secondary elements and their
304 : * associated mortar segment elements that would have nonzero values for a Lagrange shape function
305 : * associated with the provided node. This method may return an empty container if the node is
306 : * away from the mortar mesh
307 : */
308 : std::vector<MortarFilterIter> secondariesToMortarSegments(const Node & node) const;
309 :
310 : /**
311 : * @return the lower dimensional secondary element ids and their associated mortar segment
312 : * elements
313 : */
314 : const std::unordered_map<dof_id_type, std::set<Elem *, CompareDofObjectsByID>> &
315 13084 : secondariesToMortarSegments() const
316 : {
317 13084 : return _secondary_elems_to_mortar_segments;
318 : }
319 :
320 : /**
321 : * @return All the secondary interior parent subdomain IDs associated with the mortar mesh
322 : */
323 997 : const std::set<SubdomainID> & secondaryIPSubIDs() const { return _secondary_ip_sub_ids; }
324 :
325 : /**
326 : * @return All the primary interior parent subdomain IDs associated with the mortar mesh
327 : */
328 997 : const std::set<SubdomainID> & primaryIPSubIDs() const { return _primary_ip_sub_ids; }
329 :
330 : /**
331 : * @return Map from node id to secondary lower-d element pointer
332 : */
333 : const std::unordered_map<dof_id_type, std::vector<const Elem *>> & nodesToSecondaryElem() const
334 : {
335 : return _nodes_to_secondary_elem_map;
336 : }
337 :
338 : /**
339 : * initialize mortar-mesh based output
340 : */
341 : void initOutput();
342 :
343 : private:
344 : /**
345 : * build the \p _mortar_interface_coupling data
346 : */
347 : void buildCouplingInformation();
348 :
349 : /// The Moose app
350 : MooseApp & _app;
351 :
352 : /// Reference to the mesh stored in equation_systems.
353 : MeshBase & _mesh;
354 :
355 : /// The boundary ids corresponding to all the secondary surfaces.
356 : std::set<BoundaryID> _secondary_requested_boundary_ids;
357 :
358 : /// The boundary ids corresponding to all the primary surfaces.
359 : std::set<BoundaryID> _primary_requested_boundary_ids;
360 :
361 : /// A list of primary/secondary boundary id pairs corresponding to each
362 : /// side of the mortar interface.
363 : std::vector<std::pair<BoundaryID, BoundaryID>> _primary_secondary_boundary_id_pairs;
364 :
365 : /// Map from nodes to connected lower-dimensional elements on the secondary/primary subdomains.
366 : std::unordered_map<dof_id_type, std::vector<const Elem *>> _nodes_to_secondary_elem_map;
367 : std::unordered_map<dof_id_type, std::vector<const Elem *>> _nodes_to_primary_elem_map;
368 :
369 : /// Similar to the map above, but associates a (Secondary Node, Secondary Elem)
370 : /// pair to a (xi^(2), primary Elem) pair. This allows a single secondary node, which is
371 : /// potentially connected to two elements on the secondary side, to be associated with
372 : /// multiple primary Elem/xi^(2) values to handle the case where the primary and secondary
373 : /// nodes are "matching".
374 : /// In this configuration:
375 : ///
376 : /// A B
377 : /// o-----o-----o (secondary orientation ->)
378 : /// |
379 : /// v
380 : /// ------x------ (primary orientation <-)
381 : /// C D
382 : ///
383 : /// The entries in the map should be:
384 : /// (Elem A, Node 1) -> (Elem C, xi^(2)=-1)
385 : /// (Elem B, Node 0) -> (Elem D, xi^(2)=+1)
386 : std::unordered_map<std::pair<const Node *, const Elem *>, std::pair<Real, const Elem *>>
387 : _secondary_node_and_elem_to_xi2_primary_elem;
388 :
389 : /// Same type of container, but for mapping (Primary Node ID, Primary Node,
390 : /// Primary Elem) -> (xi^(1), Secondary Elem) where they are inverse-projected along
391 : /// the nodal normal direction. Note that the first item of the key, the primary
392 : /// node ID, is important for storing the key-value pairs in a consistent order
393 : /// across processes, e.g. this container has to be ordered!
394 : std::map<std::tuple<dof_id_type, const Node *, const Elem *>, std::pair<Real, const Elem *>>
395 : _primary_node_and_elem_to_xi1_secondary_elem;
396 :
397 : /// 1D Mesh of mortar segment elements which gets built by the call
398 : /// to build_mortar_segment_mesh().
399 : std::unique_ptr<MeshBase> _mortar_segment_mesh;
400 :
401 : /// Map between Elems in the mortar segment mesh and their info
402 : /// structs. This gets filled in by the call to
403 : /// build_mortar_segment_mesh().
404 : std::unordered_map<const Elem *, MortarSegmentInfo> _msm_elem_to_info;
405 :
406 : /// Keeps track of the mapping between lower-dimensional elements and
407 : /// the side_id of the interior_parent which they are.
408 : std::unordered_map<const Elem *, unsigned int> _lower_elem_to_side_id;
409 :
410 : /// A list of primary/secondary subdomain id pairs corresponding to each
411 : /// side of the mortar interface.
412 : std::vector<std::pair<SubdomainID, SubdomainID>> _primary_secondary_subdomain_id_pairs;
413 :
414 : /// The secondary/primary lower-dimensional boundary subdomain ids are the
415 : /// secondary/primary *boundary* ids
416 : std::set<SubdomainID> _secondary_boundary_subdomain_ids;
417 : std::set<SubdomainID> _primary_boundary_subdomain_ids;
418 :
419 : /// Used by the AugmentSparsityOnInterface functor to determine
420 : /// whether a given Elem is coupled to any others across the gap, and
421 : /// to explicitly set up the dependence between interior_parent()
422 : /// elements on the secondary side and their lower-dimensional sides
423 : /// which are on the interface. This latter type of coupling must be
424 : /// explicitly declared when there is no primary_elem for a given
425 : /// mortar segment and you are using e.g. a P^1-P^0 discretization
426 : /// which does not induce the coupling automatically.
427 : std::unordered_map<dof_id_type, std::unordered_set<dof_id_type>> _mortar_interface_coupling;
428 :
429 : /// Container for storing the nodal normal vector associated with each secondary node.
430 : std::unordered_map<const Node *, Point> _secondary_node_to_nodal_normal;
431 :
432 : /// Container for storing the nodal tangent/binormal vectors associated with each secondary node
433 : /// (Householder approach).
434 : std::unordered_map<const Node *, std::array<Point, 2>> _secondary_node_to_hh_nodal_tangents;
435 :
436 : /// Map from full dimensional secondary element id to lower dimensional secondary element
437 : std::unordered_map<dof_id_type, const Elem *> _secondary_element_to_secondary_lowerd_element;
438 :
439 : // List of inactive lagrange multiplier nodes (for nodal variables)
440 : std::unordered_set<const Node *> _inactive_local_lm_nodes;
441 :
442 : /// List of inactive lagrange multiplier nodes (for elemental variables)
443 : std::unordered_set<const Elem *> _inactive_local_lm_elems;
444 :
445 : /// We maintain a mapping from lower-dimensional secondary elements in the original mesh to (sets
446 : /// of) elements in mortar_segment_mesh. This allows us to quickly determine which elements need
447 : /// to be split.
448 : std::unordered_map<dof_id_type, std::set<Elem *, CompareDofObjectsByID>>
449 : _secondary_elems_to_mortar_segments;
450 :
451 : /// All the secondary interior parent subdomain IDs associated with the mortar mesh
452 : std::set<SubdomainID> _secondary_ip_sub_ids;
453 :
454 : /// All the primary interior parent subdomain IDs associated with the mortar mesh
455 : std::set<SubdomainID> _primary_ip_sub_ids;
456 :
457 : /**
458 : * Helper function responsible for projecting secondary nodes
459 : * onto primary elements for a single primary/secondary pair. Called by the class member
460 : * AutomaticMortarGeneration::project_secondary_nodes().
461 : */
462 : void projectSecondaryNodesSinglePair(SubdomainID lower_dimensional_primary_subdomain_id,
463 : SubdomainID lower_dimensional_secondary_subdomain_id);
464 :
465 : /**
466 : * Helper function used internally by AutomaticMortarGeneration::project_primary_nodes().
467 : */
468 : void projectPrimaryNodesSinglePair(SubdomainID lower_dimensional_primary_subdomain_id,
469 : SubdomainID lower_dimensional_secondary_subdomain_id);
470 :
471 : /**
472 : * Householder orthogonalization procedure to obtain proper basis for tangent and binormal vectors
473 : */
474 : void
475 : householderOrthogolization(const Point & normal, Point & tangent_one, Point & tangent_two) const;
476 :
477 : /// Whether to print debug output
478 : const bool _debug;
479 :
480 : /// Whether this object is on the displaced mesh
481 : const bool _on_displaced;
482 :
483 : /// Whether this object will be generating a mortar segment mesh for periodic constraints
484 : const bool _periodic;
485 :
486 : /// Whether the mortar segment mesh is distributed
487 : const bool _distributed;
488 :
489 : /// Newton solve tolerance for node projections
490 : Real _newton_tolerance = 1e-12;
491 :
492 : /// Tolerance for checking projection xi values. Usually we are checking whether we projected onto
493 : /// a certain element (in which case -1 <= xi <= 1) or whether we should have *already* projected
494 : /// a primary node (in which case we error if abs(xi) is sufficiently close to 1)
495 : Real _xi_tolerance = 1e-6;
496 :
497 : /// Flag to enable regressed treatment of edge dropping where all LM DoFs on edge dropping element
498 : /// are strongly set to 0.
499 : const bool _correct_edge_dropping;
500 :
501 : /// Parameter to control which angle (in degrees) is admissible for the creation of mortar segments.
502 : /// If set to a value close to zero, very oblique projections are allowed, which can result in mortar
503 : /// segments solving physics not meaningfully and overprojection of primary nodes onto the mortar
504 : /// segment mesh in extreme cases. This parameter is mostly intended for mortar mesh debugging purposes in 2D.
505 : const Real _minimum_projection_angle;
506 :
507 : /// Storage for the input parameters used by the mortar nodal geometry output
508 : std::unique_ptr<InputParameters> _output_params;
509 :
510 : friend class MortarNodalGeometryOutput;
511 : friend class AugmentSparsityOnInterface;
512 : };
513 :
514 : inline const std::pair<BoundaryID, BoundaryID> &
515 14676 : AutomaticMortarGeneration::primarySecondaryBoundaryIDPair() const
516 : {
517 : mooseAssert(_primary_secondary_boundary_id_pairs.size() == 1,
518 : "We currently only support a single boundary pair per mortar generation object");
519 :
520 14676 : return _primary_secondary_boundary_id_pairs.front();
521 : }
|
