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 : // MOOSE Includes
11 : #include "MortarConsumerInterface.h"
12 : #include "InputParameters.h"
13 : #include "MooseObject.h"
14 : #include "FEProblemBase.h"
15 : #include "MooseMesh.h"
16 : #include "Mortar3DSubpatchPlane.h"
17 : #include "MortarInterfaceWarehouse.h"
18 : #include "Assembly.h"
19 : #include "AutomaticMortarGeneration.h"
20 : #include "libmesh/quadrature.h"
21 :
22 : #include <algorithm>
23 :
24 : InputParameters
25 58574 : MortarConsumerInterface::validParams()
26 : {
27 : // Create InputParameters object that will be appended to the parameters for the inheriting object
28 58574 : InputParameters params = emptyInputParameters();
29 : // On a displaced mesh this will geometrically and algebraically ghost the entire interface
30 175722 : params.addRelationshipManager(
31 : "AugmentSparsityOnInterface",
32 : Moose::RelationshipManagerType::GEOMETRIC | Moose::RelationshipManagerType::ALGEBRAIC,
33 58574 : [](const InputParameters & obj_params, InputParameters & rm_params)
34 : {
35 8940 : rm_params.set<bool>("use_displaced_mesh") = obj_params.get<bool>("use_displaced_mesh");
36 8940 : rm_params.set<BoundaryName>("secondary_boundary") =
37 8940 : obj_params.get<BoundaryName>("secondary_boundary");
38 8940 : rm_params.set<BoundaryName>("primary_boundary") =
39 8940 : obj_params.get<BoundaryName>("primary_boundary");
40 8940 : rm_params.set<SubdomainName>("secondary_subdomain") =
41 8940 : obj_params.get<SubdomainName>("secondary_subdomain");
42 8940 : rm_params.set<SubdomainName>("primary_subdomain") =
43 8940 : obj_params.get<SubdomainName>("primary_subdomain");
44 4470 : rm_params.set<bool>("ghost_point_neighbors") =
45 4470 : obj_params.get<bool>("ghost_point_neighbors");
46 4470 : rm_params.set<bool>("ghost_higher_d_neighbors") =
47 4470 : obj_params.get<bool>("ghost_higher_d_neighbors");
48 4470 : });
49 :
50 234296 : params.addRequiredParam<BoundaryName>("primary_boundary",
51 : "The name of the primary boundary sideset.");
52 234296 : params.addRequiredParam<BoundaryName>("secondary_boundary",
53 : "The name of the secondary boundary sideset.");
54 234296 : params.addRequiredParam<SubdomainName>("primary_subdomain", "The name of the primary subdomain.");
55 234296 : params.addRequiredParam<SubdomainName>("secondary_subdomain",
56 : "The name of the secondary subdomain.");
57 175722 : params.addParam<bool>(
58 : "periodic",
59 117148 : false,
60 : "Whether this constraint is going to be used to enforce a periodic condition. This has the "
61 : "effect of changing the normals vector for projection from outward to inward facing");
62 :
63 175722 : params.addParam<bool>(
64 : "debug_mesh",
65 117148 : false,
66 : "Whether this constraint is going to enable mortar segment mesh debug information. An exodus"
67 : "file will be generated if the user sets this flag to true");
68 :
69 175722 : params.addParam<bool>(
70 : "correct_edge_dropping",
71 117148 : false,
72 : "Whether to enable correct edge dropping treatment for mortar constraints. When disabled "
73 : "any Lagrange Multiplier degree of freedom on a secondary element without full primary "
74 : "contributions will be set (strongly) to 0.");
75 :
76 175722 : params.addParam<bool>(
77 : "interpolate_normals",
78 117148 : true,
79 : "Whether to interpolate the nodal normals (e.g. classic idea of evaluating field at "
80 : "quadrature points). If this is set to false, then non-interpolated nodal normals will be "
81 : "used, and then the _normals member should be indexed with _i instead of _qp");
82 :
83 175722 : params.addParam<bool>("ghost_point_neighbors",
84 117148 : false,
85 : "Whether we should ghost point neighbors of secondary face elements, their "
86 : "interior parents, and their mortar interface couples.");
87 175722 : params.addParam<Real>(
88 : "minimum_projection_angle",
89 117148 : 40.0,
90 : "Parameter to control which angle (in degrees) is admissible for the creation of mortar "
91 : "segments. If set to a value close to zero, very oblique projections are allowed, which "
92 : "can result in mortar segments solving physics not meaningfully, and overprojection of "
93 : "primary nodes onto the mortar segment mesh in extreme cases. In 3D with "
94 : "mortar_3d_subpatch_plane = GEOMETRIC_NORMAL, this parameter also controls which primary "
95 : "and secondary subpatch normal pairings are admissible before polygon clipping.");
96 117148 : params.addParam<MooseEnum>(
97 : "mortar_3d_subpatch_plane",
98 175722 : MooseEnum(getMortar3DSubpatchPlaneOptions(), "GEOMETRIC_NORMAL"),
99 : "Method used to construct the local 3D mortar subpatch planes used for projection and "
100 : "clipping. GEOMETRIC_NORMAL uses the geometric normal of each linear or bilinear secondary "
101 : "subpatch. AVERAGED_NODAL_NORMAL uses the averaged nodal normal on each secondary subpatch.");
102 58574 : params += MortarConsumerInterface::triangulationParams();
103 :
104 117148 : MooseEnum mortar_3d_qp_mapping(getMortar3DQuadraturePointMappingOptions(), "normal_projection");
105 234296 : mortar_3d_qp_mapping.addDocumentation(
106 : "normal_projection",
107 : "Project each mortar-segment quadrature point onto the linearized primary and secondary "
108 : "sub-elements along the mortar-segment normal.");
109 234296 : mortar_3d_qp_mapping.addDocumentation(
110 : "reference_interpolation",
111 : "Interpolate stored primary and secondary parent-face reference coordinates over each "
112 : "triangular mortar segment.");
113 234296 : params.addParam<MooseEnum>(
114 : "mortar_3d_qp_mapping",
115 : mortar_3d_qp_mapping,
116 : "Method used to map quadrature points on 3D mortar segments to primary and secondary "
117 : "parent-face reference coordinates. This parameter has no effect for 2D mortar.");
118 :
119 117148 : params.addParam<bool>(
120 : "ghost_higher_d_neighbors",
121 117148 : false,
122 : "Whether we should ghost higher-dimensional neighbors. This is necessary when we are doing "
123 : "second order mortar with finite volume primal variables, because in order for the method to "
124 : "be second order we must use cell gradients, which couples in the neighbor cells.");
125 :
126 117148 : return params;
127 58574 : }
128 :
129 : InputParameters
130 58574 : MortarConsumerInterface::triangulationParams()
131 : {
132 58574 : InputParameters params = emptyInputParameters();
133 :
134 : MooseEnum triangulation(
135 : #if defined(LIBMESH_HAVE_TRIANGLE) || defined(LIBMESH_HAVE_POLY2TRI)
136 : "vertex centroid ear_clipping delaunay",
137 : #else
138 : "vertex centroid ear_clipping",
139 : #endif
140 234296 : "centroid");
141 234296 : triangulation.addDocumentation(
142 : "vertex",
143 : "Triangulate clipped 3D mortar polygons by forming a fan from an existing polygon vertex.");
144 234296 : triangulation.addDocumentation(
145 : "centroid",
146 : "Triangulate clipped 3D mortar polygons by forming a fan from a polygon centroid.");
147 234296 : triangulation.addDocumentation(
148 : "ear_clipping", "Triangulate clipped 3D mortar polygons with an ear-clipping algorithm.");
149 : #if defined(LIBMESH_HAVE_TRIANGLE) || defined(LIBMESH_HAVE_POLY2TRI)
150 234296 : triangulation.addDocumentation(
151 : "delaunay",
152 : "Triangulate clipped 3D mortar polygons using libMesh's constrained-Delaunay PSLG "
153 : "triangulation backend while preserving polygon boundary edges.");
154 : #endif
155 234296 : params.addParam<MooseEnum>(
156 : "triangulation",
157 : triangulation,
158 : "Strategy used to triangulate clipped 3D mortar polygons into mortar segments. The default "
159 : "is 'centroid' to preserve the legacy 3D mortar segmentation behavior.");
160 117148 : params.addParam<bool>(
161 : "triangulate_triangles",
162 117148 : false,
163 : "Whether a clipped 3D mortar polygon that is already a triangle should still be subdivided "
164 : "during triangulation. When enabled, already-triangular polygons are subdivided with the "
165 : "centroid-based path because the vertex-fan, ear-clipping, and Delaunay backends cannot "
166 : "refine a triangle any further on their own.");
167 :
168 117148 : return params;
169 58574 : }
170 :
171 : // Standard constructor
172 1484 : MortarConsumerInterface::MortarConsumerInterface(const MooseObject * moose_object)
173 5936 : : _mci_fe_problem(*moose_object->getCheckedPointerParam<FEProblemBase *>("_fe_problem_base")),
174 5936 : _mci_subproblem(*moose_object->getCheckedPointerParam<SubProblem *>("_subproblem")),
175 2968 : _mci_tid(moose_object->getParam<THREAD_ID>("_tid")),
176 1484 : _mci_mesh(_mci_subproblem.mesh()),
177 : // all geometric assembly information should be correct for nl system number 0
178 1484 : _mci_assembly(_mci_subproblem.assembly(_mci_tid, 0)),
179 1484 : _mortar_data(_mci_fe_problem.mortarData()),
180 1484 : _secondary_id(
181 2968 : _mci_mesh.getBoundaryID(moose_object->getParam<BoundaryName>("secondary_boundary"))),
182 2968 : _primary_id(_mci_mesh.getBoundaryID(moose_object->getParam<BoundaryName>("primary_boundary"))),
183 1484 : _secondary_subdomain_id(
184 2968 : _mci_mesh.getSubdomainID(moose_object->getParam<SubdomainName>("secondary_subdomain"))),
185 1484 : _primary_subdomain_id(
186 2968 : _mci_mesh.getSubdomainID(moose_object->getParam<SubdomainName>("primary_subdomain"))),
187 2968 : _secondary_set({_secondary_id}),
188 2968 : _interpolate_normals(moose_object->getParam<bool>("interpolate_normals")),
189 1484 : _phys_points_secondary(_mci_assembly.qPointsFace()),
190 1484 : _phys_points_primary(_mci_assembly.qPointsFaceNeighbor()),
191 1484 : _qrule_msm(_mci_assembly.qRuleMortar()),
192 1484 : _qrule_face(_mci_assembly.qRuleFace()),
193 1484 : _lower_secondary_elem(_mci_assembly.lowerDElem()),
194 1484 : _lower_primary_elem(_mci_assembly.neighborLowerDElem()),
195 1484 : _JxW_msm(_mci_assembly.jxWMortar()),
196 2968 : _msm_elem(_mci_assembly.msmElem())
197 : {
198 2968 : const bool displaced = moose_object->isParamValid("use_displaced_mesh")
199 5936 : ? moose_object->getParam<bool>("use_displaced_mesh")
200 1484 : : false;
201 2968 : const auto minimum_projection_angle = moose_object->getParam<Real>("minimum_projection_angle");
202 : const auto mortar_3d_subpatch_plane =
203 4452 : moose_object->getParam<MooseEnum>("mortar_3d_subpatch_plane")
204 1484 : .getEnum<Mortar3DSubpatchPlane>();
205 :
206 : // Only geometric 3D mode interprets this value as an angle between subpatch normals. Averaged
207 : // normal and 2D interfaces continue to use the other projection-angle checks.
208 1902 : if (_mci_mesh.dimension() == 3 &&
209 2292 : mortar_3d_subpatch_plane == Mortar3DSubpatchPlane::GEOMETRIC_NORMAL &&
210 390 : (minimum_projection_angle < 0.0 || minimum_projection_angle > 90.0))
211 6 : moose_object->paramError("minimum_projection_angle",
212 : "Must be between 0 and 90 degrees when "
213 : "mortar_3d_subpatch_plane = GEOMETRIC_NORMAL.");
214 :
215 : // Create the mortar interface if it hasn't already been created
216 2962 : _mci_fe_problem.createMortarInterface(
217 1481 : std::make_pair(_primary_id, _secondary_id),
218 0 : std::make_pair(_primary_subdomain_id, _secondary_subdomain_id),
219 : displaced,
220 2956 : moose_object->getParam<bool>("periodic"),
221 2956 : moose_object->getParam<bool>("debug_mesh"),
222 5912 : moose_object->getParam<bool>("correct_edge_dropping"),
223 : minimum_projection_angle,
224 : mortar_3d_subpatch_plane,
225 : moose_object->getParam<MooseEnum>("triangulation"),
226 4437 : moose_object->getParam<bool>("triangulate_triangles"),
227 4437 : moose_object->getParam<MooseEnum>("mortar_3d_qp_mapping")
228 : .getEnum<Mortar3DQuadraturePointMapping>());
229 :
230 1475 : _amg = &_mci_fe_problem.getMortarInterface(
231 1475 : std::make_pair(_primary_id, _secondary_id),
232 1475 : std::make_pair(_primary_subdomain_id, _secondary_subdomain_id),
233 : displaced);
234 :
235 1475 : const auto & secondary_set = _mortar_data.getHigherDimSubdomainIDs(_secondary_subdomain_id);
236 1475 : const auto & primary_set = _mortar_data.getHigherDimSubdomainIDs(_primary_subdomain_id);
237 :
238 2950 : std::set_union(secondary_set.begin(),
239 : secondary_set.end(),
240 : primary_set.begin(),
241 : primary_set.end(),
242 1475 : std::inserter(_higher_dim_subdomain_ids, _higher_dim_subdomain_ids.begin()));
243 1475 : _boundary_ids = {_secondary_id, _primary_id};
244 1475 : }
245 :
246 : void
247 581232 : MortarConsumerInterface::setNormals()
248 : {
249 581232 : if (interpolateNormals())
250 580036 : _normals = amg().getNormals(*_lower_secondary_elem, _qrule_face->get_points());
251 : else
252 1196 : _normals = amg().getNodalNormals(*_lower_secondary_elem);
253 581232 : }
254 :
255 : void
256 196448 : MortarConsumerInterface::trimDerivative(const dof_id_type remove_derivative_index,
257 : ADReal & dual_number)
258 : {
259 196448 : auto md_it = dual_number.derivatives().nude_data().begin();
260 196448 : auto mi_it = dual_number.derivatives().nude_indices().begin();
261 :
262 196448 : auto d_it = dual_number.derivatives().nude_data().begin();
263 :
264 196448 : for (auto i_it = dual_number.derivatives().nude_indices().begin();
265 554252 : i_it != dual_number.derivatives().nude_indices().end();
266 357804 : ++i_it, ++d_it)
267 357804 : if (*i_it != remove_derivative_index)
268 : {
269 333806 : *mi_it = *i_it;
270 333806 : *md_it = *d_it;
271 333806 : ++mi_it;
272 333806 : ++md_it;
273 : }
274 :
275 196448 : std::size_t n_indices = md_it - dual_number.derivatives().nude_data().begin();
276 196448 : dual_number.derivatives().nude_indices().resize(n_indices);
277 196448 : dual_number.derivatives().nude_data().resize(n_indices);
278 196448 : }
|
