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 : #include "GeometricSearchData.h"
11 : // Moose includes
12 : #include "NearestNodeLocator.h"
13 : #include "PenetrationLocator.h"
14 : #include "ElementPairLocator.h"
15 : #include "SubProblem.h"
16 : #include "MooseMesh.h"
17 : #include "Assembly.h"
18 :
19 : #include "libmesh/elem.h"
20 : #include "libmesh/node.h"
21 : #include "libmesh/int_range.h"
22 :
23 64731 : GeometricSearchData::GeometricSearchData(SubProblem & subproblem, MooseMesh & mesh)
24 64731 : : _subproblem(subproblem), _mesh(mesh), _first(true)
25 : {
26 : // do a check on the boundary IDs to see if any conflict will rise from computing QP node set IDs
27 64731 : const auto & nodeset_ids = _mesh.meshNodesetIds();
28 315895 : for (const auto id : nodeset_ids)
29 251164 : if (nodeset_ids.find(-id - 1) != nodeset_ids.end())
30 0 : mooseError("Your mesh contains nodesets with negative IDs that interfere with QP nodeset IDs "
31 : "potentially generated by the GeometricSearch system.");
32 64731 : }
33 :
34 60383 : GeometricSearchData::~GeometricSearchData()
35 : {
36 62201 : for (auto & it : _penetration_locators)
37 1818 : delete it.second;
38 :
39 62292 : for (auto & it : _nearest_node_locators)
40 1909 : delete it.second;
41 60383 : }
42 :
43 : void
44 287845 : GeometricSearchData::update(GeometricSearchType type)
45 : {
46 287845 : if (type == ALL || type == QUADRATURE || type == NEAREST_NODE)
47 : {
48 287845 : if (_first) // Only do this once
49 : {
50 61950 : _first = false;
51 :
52 62095 : for (const auto & it : _secondary_to_qsecondary)
53 145 : generateQuadratureNodes(it.first, it.second);
54 :
55 : // reinit on displaced mesh before update
56 61950 : for (const auto & epl_it : _element_pair_locators)
57 : {
58 0 : ElementPairLocator & epl = *(epl_it.second);
59 0 : epl.reinit();
60 : }
61 : }
62 :
63 : // Update the position of quadrature nodes first
64 309450 : for (const auto & qbnd : _quadrature_boundaries)
65 21605 : updateQuadratureNodes(qbnd);
66 : }
67 :
68 287845 : if (type == ALL || type == NEAREST_NODE)
69 : {
70 475976 : for (const auto & nnl_it : _nearest_node_locators)
71 : {
72 188131 : NearestNodeLocator * nnl = nnl_it.second;
73 188131 : nnl->findNodes();
74 : }
75 : }
76 :
77 287845 : if (type == ALL || type == PENETRATION)
78 : {
79 411777 : for (const auto & pl_it : _penetration_locators)
80 : {
81 185886 : PenetrationLocator * pl = pl_it.second;
82 185886 : pl->detectPenetration();
83 : }
84 : }
85 :
86 287841 : if (type == ALL || type == PENETRATION)
87 : {
88 225891 : for (auto & elem_pair_locator_pair : _element_pair_locators)
89 : {
90 0 : ElementPairLocator & epl = (*elem_pair_locator_pair.second);
91 0 : epl.update();
92 : }
93 : }
94 287841 : }
95 :
96 : void
97 7378 : GeometricSearchData::reinit()
98 : {
99 7378 : _mesh.clearQuadratureNodes();
100 : // Update the position of quadrature nodes first
101 7440 : for (const auto & qbnd : _quadrature_boundaries)
102 62 : reinitQuadratureNodes(qbnd);
103 :
104 7526 : for (const auto & nnl_it : _nearest_node_locators)
105 : {
106 148 : NearestNodeLocator * nnl = nnl_it.second;
107 148 : nnl->reinit();
108 : }
109 :
110 7502 : for (const auto & pl_it : _penetration_locators)
111 : {
112 124 : PenetrationLocator * pl = pl_it.second;
113 124 : pl->reinit();
114 : }
115 :
116 7378 : for (const auto & epl_it : _element_pair_locators)
117 : {
118 0 : ElementPairLocator & epl = *(epl_it.second);
119 0 : epl.reinit();
120 : }
121 7378 : }
122 :
123 : void
124 942 : GeometricSearchData::clearNearestNodeLocators()
125 : {
126 1884 : for (const auto & nnl_it : _nearest_node_locators)
127 : {
128 942 : NearestNodeLocator * nnl = nnl_it.second;
129 942 : nnl->reinit();
130 : }
131 942 : }
132 :
133 : Real
134 361 : GeometricSearchData::maxPatchPercentage()
135 : {
136 361 : Real max = 0.0;
137 :
138 1083 : for (const auto & nnl_it : _nearest_node_locators)
139 : {
140 722 : NearestNodeLocator * nnl = nnl_it.second;
141 :
142 722 : if (nnl->_max_patch_percentage > max)
143 349 : max = nnl->_max_patch_percentage;
144 : }
145 :
146 361 : return max;
147 : }
148 :
149 : PenetrationLocator &
150 13563 : GeometricSearchData::getPenetrationLocator(const BoundaryName & primary,
151 : const BoundaryName & secondary,
152 : Order order)
153 : {
154 13563 : auto primary_id = _mesh.getBoundaryID(primary);
155 13563 : auto secondary_id = _mesh.getBoundaryID(secondary);
156 :
157 13563 : _subproblem.addGhostedBoundary(primary_id);
158 13563 : _subproblem.addGhostedBoundary(secondary_id);
159 :
160 : PenetrationLocator * pl =
161 13563 : _penetration_locators[std::pair<BoundaryID, BoundaryID>(primary_id, secondary_id)];
162 :
163 13563 : if (!pl)
164 : {
165 0 : pl = new PenetrationLocator(_subproblem,
166 : *this,
167 : _mesh,
168 : primary_id,
169 : secondary_id,
170 : order,
171 1707 : getNearestNodeLocator(primary_id, secondary_id));
172 1707 : _penetration_locators[std::pair<BoundaryID, BoundaryID>(primary_id, secondary_id)] = pl;
173 : }
174 :
175 13563 : return *pl;
176 : }
177 :
178 : PenetrationLocator &
179 129 : GeometricSearchData::getQuadraturePenetrationLocator(const BoundaryName & primary,
180 : const BoundaryName & secondary,
181 : Order order)
182 : {
183 129 : const auto primary_id = _mesh.getBoundaryID(primary);
184 129 : const auto secondary_id = _mesh.getBoundaryID(secondary);
185 :
186 129 : _subproblem.addGhostedBoundary(primary_id);
187 129 : _subproblem.addGhostedBoundary(secondary_id);
188 :
189 : // Generate a new boundary id (we use the negative numbers and subtract 1 to disambiguate id=0)
190 129 : const auto qsecondary_id = -secondary_id - 1;
191 :
192 129 : _secondary_to_qsecondary[secondary_id] = qsecondary_id;
193 :
194 : PenetrationLocator * pl =
195 129 : _penetration_locators[std::pair<BoundaryID, BoundaryID>(primary_id, qsecondary_id)];
196 :
197 129 : if (!pl)
198 : {
199 0 : pl = new PenetrationLocator(_subproblem,
200 : *this,
201 : _mesh,
202 : primary_id,
203 : qsecondary_id,
204 : order,
205 119 : getQuadratureNearestNodeLocator(primary_id, secondary_id));
206 119 : _penetration_locators[std::pair<BoundaryID, BoundaryID>(primary_id, qsecondary_id)] = pl;
207 : }
208 :
209 129 : return *pl;
210 : }
211 :
212 : NearestNodeLocator &
213 84 : GeometricSearchData::getNearestNodeLocator(const BoundaryName & primary,
214 : const BoundaryName & secondary)
215 : {
216 84 : const auto primary_id = _mesh.getBoundaryID(primary);
217 84 : const auto secondary_id = _mesh.getBoundaryID(secondary);
218 :
219 84 : _subproblem.addGhostedBoundary(primary_id);
220 84 : _subproblem.addGhostedBoundary(secondary_id);
221 :
222 84 : return getNearestNodeLocator(primary_id, secondary_id);
223 : }
224 :
225 : NearestNodeLocator &
226 1938 : GeometricSearchData::getNearestNodeLocator(const BoundaryID primary_id,
227 : const BoundaryID secondary_id)
228 : {
229 : NearestNodeLocator * nnl =
230 1938 : _nearest_node_locators[std::pair<BoundaryID, BoundaryID>(primary_id, secondary_id)];
231 :
232 1938 : _subproblem.addGhostedBoundary(primary_id);
233 1938 : _subproblem.addGhostedBoundary(secondary_id);
234 :
235 1938 : if (!nnl)
236 : {
237 1917 : nnl = new NearestNodeLocator(_subproblem, _mesh, primary_id, secondary_id);
238 1917 : _nearest_node_locators[std::pair<BoundaryID, BoundaryID>(primary_id, secondary_id)] = nnl;
239 : }
240 :
241 1938 : return *nnl;
242 : }
243 :
244 : NearestNodeLocator &
245 28 : GeometricSearchData::getQuadratureNearestNodeLocator(const BoundaryName & primary,
246 : const BoundaryName & secondary)
247 : {
248 28 : const auto primary_id = _mesh.getBoundaryID(primary);
249 28 : const auto secondary_id = _mesh.getBoundaryID(secondary);
250 :
251 28 : _subproblem.addGhostedBoundary(primary_id);
252 28 : _subproblem.addGhostedBoundary(secondary_id);
253 :
254 28 : return getQuadratureNearestNodeLocator(primary_id, secondary_id);
255 : }
256 :
257 : NearestNodeLocator &
258 147 : GeometricSearchData::getQuadratureNearestNodeLocator(const BoundaryID primary_id,
259 : const BoundaryID secondary_id)
260 : {
261 : // Generate a new boundary id (we use the negative numbers and subtract 1 to disambiguate id=0)
262 147 : const auto qsecondary_id = -secondary_id - 1;
263 :
264 147 : _secondary_to_qsecondary[secondary_id] = qsecondary_id;
265 147 : return getNearestNodeLocator(primary_id, qsecondary_id);
266 : }
267 :
268 : void
269 207 : GeometricSearchData::generateQuadratureNodes(const BoundaryID secondary_id,
270 : const BoundaryID qsecondary_id,
271 : bool reiniting)
272 : {
273 : // Have we already generated quadrature nodes for this boundary id?
274 207 : if (_quadrature_boundaries.find(secondary_id) != _quadrature_boundaries.end())
275 : {
276 62 : if (!reiniting)
277 0 : return;
278 : }
279 : else
280 145 : _quadrature_boundaries.insert(secondary_id);
281 :
282 207 : const MooseArray<Point> & points_face = _subproblem.assembly(0, 0).qPointsFace();
283 :
284 207 : ConstBndElemRange & range = *_mesh.getBoundaryElementRange();
285 37763 : for (const auto & belem : range)
286 : {
287 37556 : const Elem * elem = belem->_elem;
288 37556 : const auto side = belem->_side;
289 37556 : const auto boundary_id = belem->_bnd_id;
290 :
291 37556 : if (elem->processor_id() == _subproblem.processor_id())
292 : {
293 28844 : if (boundary_id == (BoundaryID)secondary_id)
294 : {
295 : // All we should need to do here is reinit the underlying libMesh::FE object because that
296 : // will get us the correct points for the current element and side
297 5466 : _subproblem.setCurrentSubdomainID(elem, 0);
298 5466 : _subproblem.assembly(0, 0).reinit(elem, side);
299 :
300 11556 : for (const auto qp : index_range(points_face))
301 6090 : _mesh.addQuadratureNode(elem, side, qp, qsecondary_id, points_face[qp]);
302 : }
303 : }
304 : }
305 : }
306 :
307 : void
308 0 : GeometricSearchData::addElementPairLocator(const BoundaryID interface_id,
309 : std::shared_ptr<ElementPairLocator> epl)
310 : {
311 0 : _element_pair_locators[interface_id] = epl;
312 0 : }
313 :
314 : void
315 21605 : GeometricSearchData::updateQuadratureNodes(const BoundaryID secondary_id)
316 : {
317 21605 : const MooseArray<Point> & points_face = _subproblem.assembly(0, 0).qPointsFace();
318 :
319 21605 : ConstBndElemRange & range = *_mesh.getBoundaryElementRange();
320 4291442 : for (const auto & belem : range)
321 : {
322 4269837 : const Elem * elem = belem->_elem;
323 4269837 : const auto side = belem->_side;
324 4269837 : const auto boundary_id = belem->_bnd_id;
325 :
326 4269837 : if (elem->processor_id() == _subproblem.processor_id())
327 : {
328 3403868 : if (boundary_id == (BoundaryID)secondary_id)
329 : {
330 : // All we should need to do here is reinit the underlying libMesh::FE object because that
331 : // will get us the correct points for the current element and side
332 45178 : _subproblem.setCurrentSubdomainID(elem, 0);
333 45178 : _subproblem.assembly(0, 0).reinit(elem, side);
334 :
335 125850 : for (const auto qp : index_range(points_face))
336 80672 : (*_mesh.getQuadratureNode(elem, side, qp)) = points_face[qp];
337 : }
338 : }
339 : }
340 21605 : }
341 :
342 : void
343 62 : GeometricSearchData::reinitQuadratureNodes(const BoundaryID /*secondary_id*/)
344 : {
345 : // Regenerate the quadrature nodes
346 124 : for (const auto & it : _secondary_to_qsecondary)
347 62 : generateQuadratureNodes(it.first, it.second, /*reiniting=*/true);
348 62 : }
349 :
350 : void
351 722 : GeometricSearchData::updateGhostedElems()
352 : {
353 1444 : for (const auto & nnl_it : _nearest_node_locators)
354 : {
355 722 : NearestNodeLocator * nnl = nnl_it.second;
356 722 : nnl->updateGhostedElems();
357 : }
358 722 : }
|
