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 60037 : GeometricSearchData::GeometricSearchData(SubProblem & subproblem, MooseMesh & mesh)
24 60037 : : _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 60037 : const auto & nodeset_ids = _mesh.meshNodesetIds();
28 293257 : for (const auto id : nodeset_ids)
29 233220 : 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 60037 : }
33 :
34 55724 : GeometricSearchData::~GeometricSearchData()
35 : {
36 57400 : for (auto & it : _penetration_locators)
37 1676 : delete it.second;
38 :
39 57484 : for (auto & it : _nearest_node_locators)
40 1760 : delete it.second;
41 55724 : }
42 :
43 : void
44 265480 : GeometricSearchData::update(GeometricSearchType type)
45 : {
46 265480 : if (type == ALL || type == QUADRATURE || type == NEAREST_NODE)
47 : {
48 265480 : if (_first) // Only do this once
49 : {
50 57268 : _first = false;
51 :
52 57402 : for (const auto & it : _secondary_to_qsecondary)
53 134 : generateQuadratureNodes(it.first, it.second);
54 :
55 : // reinit on displaced mesh before update
56 57268 : 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 285344 : for (const auto & qbnd : _quadrature_boundaries)
65 19864 : updateQuadratureNodes(qbnd);
66 : }
67 :
68 265480 : if (type == ALL || type == NEAREST_NODE)
69 : {
70 437586 : for (const auto & nnl_it : _nearest_node_locators)
71 : {
72 172106 : NearestNodeLocator * nnl = nnl_it.second;
73 172106 : nnl->findNodes();
74 : }
75 : }
76 :
77 265480 : if (type == ALL || type == PENETRATION)
78 : {
79 378242 : for (const auto & pl_it : _penetration_locators)
80 : {
81 170034 : PenetrationLocator * pl = pl_it.second;
82 170034 : pl->detectPenetration();
83 : }
84 : }
85 :
86 265476 : if (type == ALL || type == PENETRATION)
87 : {
88 208208 : 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 265476 : }
95 :
96 : void
97 6777 : GeometricSearchData::reinit()
98 : {
99 6777 : _mesh.clearQuadratureNodes();
100 : // Update the position of quadrature nodes first
101 6832 : for (const auto & qbnd : _quadrature_boundaries)
102 55 : reinitQuadratureNodes(qbnd);
103 :
104 6909 : for (const auto & nnl_it : _nearest_node_locators)
105 : {
106 132 : NearestNodeLocator * nnl = nnl_it.second;
107 132 : nnl->reinit();
108 : }
109 :
110 6887 : for (const auto & pl_it : _penetration_locators)
111 : {
112 110 : PenetrationLocator * pl = pl_it.second;
113 110 : pl->reinit();
114 : }
115 :
116 6777 : for (const auto & epl_it : _element_pair_locators)
117 : {
118 0 : ElementPairLocator & epl = *(epl_it.second);
119 0 : epl.reinit();
120 : }
121 6777 : }
122 :
123 : void
124 858 : GeometricSearchData::clearNearestNodeLocators()
125 : {
126 1716 : for (const auto & nnl_it : _nearest_node_locators)
127 : {
128 858 : NearestNodeLocator * nnl = nnl_it.second;
129 858 : nnl->reinit();
130 : }
131 858 : }
132 :
133 : Real
134 330 : GeometricSearchData::maxPatchPercentage()
135 : {
136 330 : Real max = 0.0;
137 :
138 990 : for (const auto & nnl_it : _nearest_node_locators)
139 : {
140 660 : NearestNodeLocator * nnl = nnl_it.second;
141 :
142 660 : if (nnl->_max_patch_percentage > max)
143 319 : max = nnl->_max_patch_percentage;
144 : }
145 :
146 330 : return max;
147 : }
148 :
149 : PenetrationLocator &
150 12585 : GeometricSearchData::getPenetrationLocator(const BoundaryName & primary,
151 : const BoundaryName & secondary,
152 : Order order)
153 : {
154 12585 : auto primary_id = _mesh.getBoundaryID(primary);
155 12585 : auto secondary_id = _mesh.getBoundaryID(secondary);
156 :
157 12585 : _subproblem.addGhostedBoundary(primary_id);
158 12585 : _subproblem.addGhostedBoundary(secondary_id);
159 :
160 : PenetrationLocator * pl =
161 12585 : _penetration_locators[std::pair<BoundaryID, BoundaryID>(primary_id, secondary_id)];
162 :
163 12585 : if (!pl)
164 : {
165 0 : pl = new PenetrationLocator(_subproblem,
166 : *this,
167 : _mesh,
168 : primary_id,
169 : secondary_id,
170 : order,
171 1574 : getNearestNodeLocator(primary_id, secondary_id));
172 1574 : _penetration_locators[std::pair<BoundaryID, BoundaryID>(primary_id, secondary_id)] = pl;
173 : }
174 :
175 12585 : return *pl;
176 : }
177 :
178 : PenetrationLocator &
179 120 : GeometricSearchData::getQuadraturePenetrationLocator(const BoundaryName & primary,
180 : const BoundaryName & secondary,
181 : Order order)
182 : {
183 120 : const auto primary_id = _mesh.getBoundaryID(primary);
184 120 : const auto secondary_id = _mesh.getBoundaryID(secondary);
185 :
186 120 : _subproblem.addGhostedBoundary(primary_id);
187 120 : _subproblem.addGhostedBoundary(secondary_id);
188 :
189 : // Generate a new boundary id (we use the negative numbers and subtract 1 to disambiguate id=0)
190 120 : const auto qsecondary_id = -secondary_id - 1;
191 :
192 120 : _secondary_to_qsecondary[secondary_id] = qsecondary_id;
193 :
194 : PenetrationLocator * pl =
195 120 : _penetration_locators[std::pair<BoundaryID, BoundaryID>(primary_id, qsecondary_id)];
196 :
197 120 : if (!pl)
198 : {
199 0 : pl = new PenetrationLocator(_subproblem,
200 : *this,
201 : _mesh,
202 : primary_id,
203 : qsecondary_id,
204 : order,
205 110 : getQuadratureNearestNodeLocator(primary_id, secondary_id));
206 110 : _penetration_locators[std::pair<BoundaryID, BoundaryID>(primary_id, qsecondary_id)] = pl;
207 : }
208 :
209 120 : return *pl;
210 : }
211 :
212 : NearestNodeLocator &
213 78 : GeometricSearchData::getNearestNodeLocator(const BoundaryName & primary,
214 : const BoundaryName & secondary)
215 : {
216 78 : const auto primary_id = _mesh.getBoundaryID(primary);
217 78 : const auto secondary_id = _mesh.getBoundaryID(secondary);
218 :
219 78 : _subproblem.addGhostedBoundary(primary_id);
220 78 : _subproblem.addGhostedBoundary(secondary_id);
221 :
222 78 : return getNearestNodeLocator(primary_id, secondary_id);
223 : }
224 :
225 : NearestNodeLocator &
226 1788 : GeometricSearchData::getNearestNodeLocator(const BoundaryID primary_id,
227 : const BoundaryID secondary_id)
228 : {
229 : NearestNodeLocator * nnl =
230 1788 : _nearest_node_locators[std::pair<BoundaryID, BoundaryID>(primary_id, secondary_id)];
231 :
232 1788 : _subproblem.addGhostedBoundary(primary_id);
233 1788 : _subproblem.addGhostedBoundary(secondary_id);
234 :
235 1788 : if (!nnl)
236 : {
237 1768 : nnl = new NearestNodeLocator(_subproblem, _mesh, primary_id, secondary_id);
238 1768 : _nearest_node_locators[std::pair<BoundaryID, BoundaryID>(primary_id, secondary_id)] = nnl;
239 : }
240 :
241 1788 : return *nnl;
242 : }
243 :
244 : NearestNodeLocator &
245 26 : GeometricSearchData::getQuadratureNearestNodeLocator(const BoundaryName & primary,
246 : const BoundaryName & secondary)
247 : {
248 26 : const auto primary_id = _mesh.getBoundaryID(primary);
249 26 : const auto secondary_id = _mesh.getBoundaryID(secondary);
250 :
251 26 : _subproblem.addGhostedBoundary(primary_id);
252 26 : _subproblem.addGhostedBoundary(secondary_id);
253 :
254 26 : return getQuadratureNearestNodeLocator(primary_id, secondary_id);
255 : }
256 :
257 : NearestNodeLocator &
258 136 : 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 136 : const auto qsecondary_id = -secondary_id - 1;
263 :
264 136 : _secondary_to_qsecondary[secondary_id] = qsecondary_id;
265 136 : return getNearestNodeLocator(primary_id, qsecondary_id);
266 : }
267 :
268 : void
269 189 : 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 189 : if (_quadrature_boundaries.find(secondary_id) != _quadrature_boundaries.end())
275 : {
276 55 : if (!reiniting)
277 0 : return;
278 : }
279 : else
280 134 : _quadrature_boundaries.insert(secondary_id);
281 :
282 189 : const MooseArray<Point> & points_face = _subproblem.assembly(0, 0).qPointsFace();
283 :
284 189 : ConstBndElemRange & range = *_mesh.getBoundaryElementRange();
285 34364 : for (const auto & belem : range)
286 : {
287 34175 : const Elem * elem = belem->_elem;
288 34175 : const auto side = belem->_side;
289 34175 : const auto boundary_id = belem->_bnd_id;
290 :
291 34175 : if (elem->processor_id() == _subproblem.processor_id())
292 : {
293 25463 : 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 4796 : _subproblem.setCurrentSubdomainID(elem, 0);
298 4796 : _subproblem.assembly(0, 0).reinit(elem, side);
299 :
300 10151 : for (const auto qp : index_range(points_face))
301 5355 : _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 19864 : GeometricSearchData::updateQuadratureNodes(const BoundaryID secondary_id)
316 : {
317 19864 : const MooseArray<Point> & points_face = _subproblem.assembly(0, 0).qPointsFace();
318 :
319 19864 : ConstBndElemRange & range = *_mesh.getBoundaryElementRange();
320 3933239 : for (const auto & belem : range)
321 : {
322 3913375 : const Elem * elem = belem->_elem;
323 3913375 : const auto side = belem->_side;
324 3913375 : const auto boundary_id = belem->_bnd_id;
325 :
326 3913375 : if (elem->processor_id() == _subproblem.processor_id())
327 : {
328 3047406 : 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 40267 : _subproblem.setCurrentSubdomainID(elem, 0);
333 40267 : _subproblem.assembly(0, 0).reinit(elem, side);
334 :
335 112327 : for (const auto qp : index_range(points_face))
336 72060 : (*_mesh.getQuadratureNode(elem, side, qp)) = points_face[qp];
337 : }
338 : }
339 : }
340 19864 : }
341 :
342 : void
343 55 : GeometricSearchData::reinitQuadratureNodes(const BoundaryID /*secondary_id*/)
344 : {
345 : // Regenerate the quadrature nodes
346 110 : for (const auto & it : _secondary_to_qsecondary)
347 55 : generateQuadratureNodes(it.first, it.second, /*reiniting=*/true);
348 55 : }
349 :
350 : void
351 660 : GeometricSearchData::updateGhostedElems()
352 : {
353 1320 : for (const auto & nnl_it : _nearest_node_locators)
354 : {
355 660 : NearestNodeLocator * nnl = nnl_it.second;
356 660 : nnl->updateGhostedElems();
357 : }
358 660 : }
|
