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 "MeshCut2DUserObjectBase.h"
11 : #include "MeshCut2DNucleationBase.h"
12 : #include "CrackFrontDefinition.h"
13 :
14 : #include "XFEMFuncs.h"
15 : #include "MooseError.h"
16 : #include "MooseMesh.h"
17 : #include "libmesh/edge_edge2.h"
18 : #include "libmesh/serial_mesh.h"
19 : #include "libmesh/mesh_tools.h"
20 :
21 : InputParameters
22 99 : MeshCut2DUserObjectBase::validParams()
23 : {
24 99 : InputParameters params = MeshCutUserObjectBase::validParams();
25 198 : params.addParam<UserObjectName>("nucleate_uo", "The MeshCutNucleation UO for nucleating cracks.");
26 198 : params.addParam<UserObjectName>("crack_front_definition",
27 : "crackFrontDefinition",
28 : "The CrackFrontDefinition user object name");
29 99 : params.addClassDescription("Creates a UserObject base class for a mesh cutter in 2D problems");
30 99 : return params;
31 0 : }
32 :
33 50 : MeshCut2DUserObjectBase::MeshCut2DUserObjectBase(const InputParameters & parameters)
34 : : MeshCutUserObjectBase(parameters),
35 100 : _mesh(_subproblem.mesh()),
36 50 : _nucleate_uo(isParamValid("nucleate_uo")
37 67 : ? &getUserObject<MeshCut2DNucleationBase>("nucleate_uo")
38 : : nullptr),
39 100 : _is_mesh_modified(false)
40 : {
41 150 : _depend_uo.insert(getParam<UserObjectName>("crack_front_definition"));
42 :
43 : // test element type; only line elements are allowed
44 318 : for (const auto & cut_elem : _cutter_mesh->element_ptr_range())
45 : {
46 109 : if (cut_elem->n_nodes() != 2)
47 0 : mooseError("The input cut mesh should include EDGE2 elements only!");
48 109 : if (cut_elem->dim() != 1)
49 0 : mooseError("The input cut mesh should have 1D elements (in a 2D space) only!");
50 50 : }
51 :
52 : // find node fronts of the original cutmesh. This is used to order EVERYTHING.
53 50 : findOriginalCrackFrontNodes();
54 50 : }
55 :
56 : void
57 49 : MeshCut2DUserObjectBase::initialSetup()
58 : {
59 49 : const auto uo_name = getParam<UserObjectName>("crack_front_definition");
60 49 : _crack_front_definition = &_fe_problem.getUserObject<CrackFrontDefinition>(uo_name);
61 49 : }
62 :
63 : bool
64 266083 : MeshCut2DUserObjectBase::cutElementByGeometry(const Elem * elem,
65 : std::vector<Xfem::CutEdge> & cut_edges,
66 : std::vector<Xfem::CutNode> & cut_nodes) const
67 : {
68 : // With the crack defined by a line, this method cuts a 2D elements by a line
69 : // Fixme lynn Copy and paste from InterfaceMeshCut2DUserObject::cutElementByGeometry
70 : mooseAssert(elem->dim() == 2, "Dimension of element to be cut must be 2");
71 :
72 : bool elem_cut = false;
73 :
74 3651334 : for (const auto & cut_elem : _cutter_mesh->element_ptr_range())
75 : {
76 1559584 : unsigned int n_sides = elem->n_sides();
77 :
78 7797920 : for (unsigned int i = 0; i < n_sides; ++i)
79 : {
80 6238336 : std::unique_ptr<const Elem> curr_side = elem->side_ptr(i);
81 :
82 : mooseAssert(curr_side->type() == EDGE2, "Element side type must be EDGE2.");
83 :
84 : const Node * node1 = curr_side->node_ptr(0);
85 : const Node * node2 = curr_side->node_ptr(1);
86 6238336 : Real seg_int_frac = 0.0;
87 :
88 6238336 : const std::pair<Point, Point> elem_endpoints(cut_elem->node_ref(0), cut_elem->node_ref(1));
89 :
90 6238336 : if (Xfem::intersectSegmentWithCutLine(*node1, *node2, elem_endpoints, 1.0, seg_int_frac))
91 : {
92 18680 : if (seg_int_frac > Xfem::tol && seg_int_frac < 1.0 - Xfem::tol)
93 : {
94 : elem_cut = true;
95 : Xfem::CutEdge mycut;
96 18680 : mycut._id1 = node1->id();
97 18680 : mycut._id2 = node2->id();
98 18680 : mycut._distance = seg_int_frac;
99 18680 : mycut._host_side_id = i;
100 18680 : cut_edges.push_back(mycut);
101 18680 : }
102 0 : else if (seg_int_frac < Xfem::tol)
103 : {
104 : elem_cut = true;
105 : Xfem::CutNode mycut;
106 0 : mycut._id = node1->id();
107 0 : mycut._host_id = i;
108 0 : cut_nodes.push_back(mycut);
109 : }
110 : }
111 6238336 : }
112 266083 : }
113 266083 : return elem_cut;
114 : }
115 :
116 : bool
117 0 : MeshCut2DUserObjectBase::cutElementByGeometry(const Elem * /*elem*/,
118 : std::vector<Xfem::CutFace> & /*cut_faces*/) const
119 : {
120 0 : mooseError("Invalid method for 2D mesh cutting.");
121 : return false;
122 : }
123 :
124 : bool
125 8949 : MeshCut2DUserObjectBase::cutFragmentByGeometry(std::vector<std::vector<Point>> & frag_edges,
126 : std::vector<Xfem::CutEdge> & cut_edges) const
127 : {
128 : bool cut_frag = false;
129 :
130 247488 : for (const auto & cut_elem : _cutter_mesh->element_ptr_range())
131 : {
132 114795 : const std::pair<Point, Point> elem_endpoints(cut_elem->node_ref(0), cut_elem->node_ref(1));
133 114795 : unsigned int n_sides = frag_edges.size();
134 583174 : for (unsigned int i = 0; i < n_sides; ++i)
135 : {
136 468379 : Real seg_int_frac = 0.0;
137 468379 : if (Xfem::intersectSegmentWithCutLine(
138 468379 : frag_edges[i][0], frag_edges[i][1], elem_endpoints, 1, seg_int_frac))
139 : {
140 : cut_frag = true;
141 : Xfem::CutEdge mycut;
142 16494 : mycut._id1 = i;
143 16494 : mycut._id2 = (i < (n_sides - 1) ? (i + 1) : 0);
144 16494 : mycut._distance = seg_int_frac;
145 16494 : mycut._host_side_id = i;
146 16494 : cut_edges.push_back(mycut);
147 : }
148 : }
149 8949 : }
150 8949 : return cut_frag;
151 : }
152 :
153 : bool
154 0 : MeshCut2DUserObjectBase::cutFragmentByGeometry(std::vector<std::vector<Point>> & /*frag_faces*/,
155 : std::vector<Xfem::CutFace> & /*cut_faces*/) const
156 : {
157 0 : mooseError("Invalid method for 2D mesh fragment cutting.");
158 : return false;
159 : }
160 :
161 : const std::vector<Point>
162 457 : MeshCut2DUserObjectBase::getCrackFrontPoints(unsigned int number_crack_front_points) const
163 : {
164 457 : std::vector<Point> crack_front_points(number_crack_front_points);
165 : // number_crack_front_points is updated via
166 : // _crack_front_definition->updateNumberOfCrackFrontPoints(_crack_front_points.size())
167 457 : if (number_crack_front_points != _original_and_current_front_node_ids.size())
168 0 : mooseError("Number of nodes in CrackFrontDefinition does not match the number of nodes in the "
169 0 : "cutter_mesh.\nCrackFrontDefinition nodes = " +
170 0 : Moose::stringify(number_crack_front_points) + "\ncutter_mesh nodes = " +
171 0 : Moose::stringify(_original_and_current_front_node_ids.size()));
172 :
173 1769 : for (unsigned int i = 0; i < number_crack_front_points; ++i)
174 : {
175 1312 : dof_id_type id = _original_and_current_front_node_ids[i].second;
176 1312 : Node * this_node = _cutter_mesh->node_ptr(id);
177 : mooseAssert(this_node, "Node is NULL");
178 : Point & this_point = *this_node;
179 1312 : crack_front_points[i] = this_point;
180 : }
181 457 : return crack_front_points;
182 0 : }
183 :
184 : const std::vector<RealVectorValue>
185 457 : MeshCut2DUserObjectBase::getCrackPlaneNormals(unsigned int number_crack_front_points) const
186 : {
187 457 : if (number_crack_front_points != _original_and_current_front_node_ids.size())
188 0 : mooseError("MeshCut2DFractureUserObject::getCrackPlaneNormals: number_crack_front_points=" +
189 0 : Moose::stringify(number_crack_front_points) +
190 : " does not match the number of nodes given in "
191 0 : "_original_and_current_front_node_ids=" +
192 : Moose::stringify(_original_and_current_front_node_ids.size()) +
193 : ". This will happen if a crack front exits the boundary because the number of "
194 : "points in the CrackFrontDefinition is never updated.");
195 :
196 : std::vector<std::pair<dof_id_type, RealVectorValue>> crack_plane_normals;
197 8248 : for (const auto & elem : _cutter_mesh->element_ptr_range())
198 : {
199 3667 : dof_id_type id0 = elem->node_id(0);
200 3667 : dof_id_type id1 = elem->node_id(1);
201 : dof_id_type id;
202 :
203 3667 : auto it0 = std::find_if(_original_and_current_front_node_ids.begin(),
204 : _original_and_current_front_node_ids.end(),
205 : [&id0](const std::pair<dof_id_type, dof_id_type> & element)
206 8816 : { return element.second == id0; });
207 3667 : auto it1 = std::find_if(_original_and_current_front_node_ids.begin(),
208 : _original_and_current_front_node_ids.end(),
209 : [&id1](const std::pair<dof_id_type, dof_id_type> & element)
210 8040 : { return element.second == id1; });
211 :
212 : bool found_it0 = (it0 != _original_and_current_front_node_ids.end());
213 : bool found_it1 = (it1 != _original_and_current_front_node_ids.end());
214 :
215 : // Newly nucleated crack elements can have one normal if they are on the edge OR
216 : // two normals if they are in the bulk.
217 3667 : if (found_it0)
218 : {
219 : Point end_pt, connecting_pt;
220 :
221 304 : end_pt = elem->node_ref(0);
222 304 : connecting_pt = elem->node_ref(1);
223 304 : id = it0->first; // sort by original crack front node ids
224 :
225 : Point fracture_dir = end_pt - connecting_pt;
226 : // The crack normal is orthogonal to the crack extension direction (fracture_dir),
227 : // and is defined in this implementation as the cross product of the direction of crack
228 : // extension with the tangent direction, which is always (0, 0, 1) in 2D.
229 304 : RealVectorValue normal_dir{fracture_dir(1), -fracture_dir(0), 0};
230 304 : normal_dir /= normal_dir.norm();
231 304 : crack_plane_normals.push_back(std::make_pair(id, normal_dir));
232 : }
233 :
234 3667 : if (found_it1)
235 : {
236 : Point end_pt, connecting_pt;
237 :
238 1008 : end_pt = elem->node_ref(1);
239 1008 : connecting_pt = elem->node_ref(0);
240 1008 : id = it1->first; // sort by original crack front node ids
241 :
242 : Point fracture_dir = end_pt - connecting_pt;
243 : // The crack normal is orthogonal to the crack extension direction (fracture_dir),
244 : // and is defined in this implementation as the cross product of the direction of crack
245 : // extension with the tangent direction, which is always (0, 0, 1) in 2D.
246 1008 : RealVectorValue normal_dir{fracture_dir(1), -fracture_dir(0), 0};
247 1008 : normal_dir /= normal_dir.norm();
248 1008 : crack_plane_normals.push_back(std::make_pair(id, normal_dir));
249 : }
250 457 : }
251 : mooseAssert(
252 : _original_and_current_front_node_ids.size() == crack_plane_normals.size(),
253 : "Boundary nodes are attached to more than one element. This should not happen for a 1D "
254 : "cutter mesh."
255 : "\n Number of _original_and_current_front_node_ids=" +
256 : Moose::stringify(_original_and_current_front_node_ids.size()) +
257 : "\n Number of crack_plane_normals=" + Moose::stringify(crack_plane_normals.size()));
258 :
259 : // the crack_plane_normals are now sorted by the ORIGINAL crack front ids
260 457 : std::sort(crack_plane_normals.begin(), crack_plane_normals.end());
261 : std::vector<RealVectorValue> sorted_crack_plane_normals;
262 1769 : for (auto & crack : crack_plane_normals)
263 1312 : sorted_crack_plane_normals.push_back(crack.second);
264 :
265 457 : return sorted_crack_plane_normals;
266 457 : }
267 :
268 : void
269 50 : MeshCut2DUserObjectBase::findOriginalCrackFrontNodes()
270 : {
271 50 : std::unique_ptr<PointLocatorBase> pl = _mesh.getPointLocator();
272 50 : pl->enable_out_of_mesh_mode();
273 50 : std::unordered_set boundary_nodes = MeshTools::find_boundary_nodes(*_cutter_mesh);
274 150 : for (const auto & node : boundary_nodes)
275 : {
276 100 : auto node_id = node;
277 100 : Node * this_node = _cutter_mesh->node_ptr(node_id);
278 : mooseAssert(this_node, "Node is NULL");
279 100 : Point & this_point = *this_node;
280 :
281 100 : bool is_on_bc = isPointOnEdgeBoundary(this_point);
282 :
283 100 : const Elem * elem = (*pl)(this_point);
284 100 : if (elem != NULL && !is_on_bc)
285 : {
286 58 : _original_and_current_front_node_ids.push_back(std::make_pair(node, node));
287 : }
288 : }
289 50 : std::sort(_original_and_current_front_node_ids.begin(),
290 : _original_and_current_front_node_ids.end());
291 50 : }
292 :
293 : bool
294 100 : MeshCut2DUserObjectBase::isPointOnEdgeBoundary(const Point & point, Real tolerance)
295 : {
296 28876 : for (const auto & bnd_elem : as_range(_mesh.bndElemsBegin(), _mesh.bndElemsEnd()))
297 : {
298 : // if (bnd_elem->_elem->processor_id() == processor_id())
299 : {
300 14296 : auto side = bnd_elem->_elem->side_ptr(bnd_elem->_side);
301 :
302 14296 : if (side->n_nodes() == 2)
303 : {
304 : // Create vectors for the endpoints of the edge
305 14296 : Point p1 = side->point(0);
306 14296 : Point p2 = side->point(1);
307 :
308 : // Compute the vector from p1 to p2 (edge vector) and from p1 to the point
309 : Point edge_vector = p2 - p1;
310 : Point point_vector = point - p1;
311 :
312 : // Compute the cross product magnitude to find if the point is collinear with the edge
313 14296 : Real cross_product_magnitude = std::sqrt(
314 14296 : std::pow(edge_vector(1) * point_vector(2) - edge_vector(2) * point_vector(1), 2) +
315 14296 : std::pow(edge_vector(2) * point_vector(0) - edge_vector(0) * point_vector(2), 2) +
316 14296 : std::pow(edge_vector(0) * point_vector(1) - edge_vector(1) * point_vector(0), 2));
317 :
318 : // Check if the point lies on the edge within a tolerance
319 14296 : if (cross_product_magnitude < tolerance)
320 : {
321 : // Check if the point lies within the edge segment bounds
322 : Real dot_product = point_vector * edge_vector;
323 : Real edge_length_squared = edge_vector * edge_vector;
324 :
325 16 : if (dot_product >= 0 && dot_product <= edge_length_squared)
326 : return true;
327 : }
328 : }
329 : else
330 0 : mooseError("MeshCut2DUserObjectBase currently only supports linear elements.\n");
331 14296 : }
332 100 : }
333 92 : return false;
334 : }
335 :
336 : void
337 374 : MeshCut2DUserObjectBase::growFront()
338 : {
339 : dof_id_type current_front_node_id;
340 1320 : for (std::size_t i = 0; i < _original_and_current_front_node_ids.size(); ++i)
341 : {
342 946 : current_front_node_id = _original_and_current_front_node_ids[i].second;
343 : // check if node front node id is active
344 : auto direction_iter =
345 946 : std::find_if(_active_front_node_growth_vectors.begin(),
346 : _active_front_node_growth_vectors.end(),
347 : [¤t_front_node_id](const std::pair<dof_id_type, Point> & element)
348 480 : { return element.first == current_front_node_id; });
349 : // only add an element for active node front ids
350 946 : if (direction_iter != _active_front_node_growth_vectors.end())
351 : {
352 260 : Node * this_node = _cutter_mesh->node_ptr(current_front_node_id);
353 : mooseAssert(this_node, "Node is NULL");
354 : Point & this_point = *this_node;
355 :
356 260 : Point new_node_offset = direction_iter->second;
357 : Point x = this_point + new_node_offset;
358 :
359 : // TODO: Should check if cut line segment created between "this_point" and "x" crosses
360 : // another line element in the cutter mesh or solid mesh boundary.
361 : // Crossing another line element would be a special case that still needs to be handled,
362 : // however, it doesnot cause an error, it will just ignore the other line segment and recut
363 : // the solid mesh element.
364 : // Crossing a solid mesh boundary would be for aesthetics reasons so
365 : // that element was trimmed close to the boundary but would have not effect on the simulation.
366 : // Crossing a solid mesh boundary should be handled by something like
367 : // MeshCut2DRankTwoTensorNucleation::lineLineIntersect2D
368 :
369 : // add node to front
370 520 : this_node = Node::build(x, _cutter_mesh->n_nodes()).release();
371 260 : _cutter_mesh->add_node(this_node);
372 260 : dof_id_type new_front_node_id = _cutter_mesh->n_nodes() - 1;
373 :
374 : // add element to front
375 : std::vector<dof_id_type> elem;
376 260 : elem.push_back(current_front_node_id);
377 260 : elem.push_back(new_front_node_id);
378 260 : Elem * new_elem = Elem::build(EDGE2).release();
379 780 : for (unsigned int i = 0; i < new_elem->n_nodes(); ++i)
380 : {
381 : mooseAssert(_cutter_mesh->node_ptr(elem[i]) != nullptr, "Node is NULL");
382 520 : new_elem->set_node(i, _cutter_mesh->node_ptr(elem[i]));
383 : }
384 260 : _cutter_mesh->add_elem(new_elem);
385 : // now push to the end of _original_and_current_front_node_ids for tracking and fracture
386 : // integrals
387 260 : _original_and_current_front_node_ids[i].second = new_front_node_id;
388 260 : _is_mesh_modified = true;
389 260 : }
390 : }
391 374 : _cutter_mesh->prepare_for_use();
392 374 : }
393 :
394 : void
395 374 : MeshCut2DUserObjectBase::addNucleatedCracksToMesh()
396 : {
397 374 : if (_nucleate_uo)
398 : {
399 : std::map<unsigned int, std::pair<RealVectorValue, RealVectorValue>> nucleated_elems_map =
400 : _nucleate_uo->getNucleatedElemsMap();
401 180 : const Real nucleationRadius = _nucleate_uo->getNucleationRadius();
402 :
403 180 : removeNucleatedCracksTooCloseToEachOther(nucleated_elems_map, nucleationRadius);
404 180 : removeNucleatedCracksTooCloseToExistingCracks(nucleated_elems_map, nucleationRadius);
405 :
406 180 : std::unique_ptr<PointLocatorBase> pl = _mesh.getPointLocator();
407 180 : pl->enable_out_of_mesh_mode();
408 252 : for (const auto & elem_nodes : nucleated_elems_map)
409 : {
410 72 : std::pair<RealVectorValue, RealVectorValue> nodes = elem_nodes.second;
411 : // add nodes for the elements that define the nucleated cracks
412 144 : Node * node_0 = Node::build(nodes.first, _cutter_mesh->n_nodes()).release();
413 72 : _cutter_mesh->add_node(node_0);
414 72 : dof_id_type node_id_0 = _cutter_mesh->n_nodes() - 1;
415 144 : Node * node_1 = Node::build(nodes.second, _cutter_mesh->n_nodes()).release();
416 72 : _cutter_mesh->add_node(node_1);
417 72 : dof_id_type node_id_1 = _cutter_mesh->n_nodes() - 1;
418 : // add elements that define nucleated cracks
419 : std::vector<dof_id_type> elem;
420 72 : elem.push_back(node_id_0);
421 72 : elem.push_back(node_id_1);
422 72 : Elem * new_elem = Elem::build(EDGE2).release();
423 216 : for (unsigned int i = 0; i < new_elem->n_nodes(); ++i)
424 : {
425 : mooseAssert(_cutter_mesh->node_ptr(elem[i]) != nullptr, "Node is NULL");
426 144 : new_elem->set_node(i, _cutter_mesh->node_ptr(elem[i]));
427 : }
428 72 : _cutter_mesh->add_elem(new_elem);
429 : // now add the nucleated nodes to the crack id data struct
430 : // edge nucleated cracks will add one node to _original_and_current_front_node_ids
431 : // bulk nucleated cracks will add two nodes to _original_and_current_front_node_ids
432 72 : Point & point_0 = *node_0;
433 72 : const Elem * crack_front_elem_0 = (*pl)(point_0);
434 72 : if (crack_front_elem_0 != NULL)
435 36 : _original_and_current_front_node_ids.push_back(std::make_pair(node_id_0, node_id_0));
436 :
437 72 : Point & point_1 = *node_1;
438 72 : const Elem * crack_front_elem_1 = (*pl)(point_1);
439 72 : if (crack_front_elem_1 != NULL)
440 60 : _original_and_current_front_node_ids.push_back(std::make_pair(node_id_1, node_id_1));
441 :
442 72 : _is_mesh_modified = true;
443 72 : }
444 180 : _cutter_mesh->prepare_for_use();
445 180 : }
446 374 : }
447 :
448 : void
449 180 : MeshCut2DUserObjectBase::removeNucleatedCracksTooCloseToEachOther(
450 : std::map<unsigned int, std::pair<RealVectorValue, RealVectorValue>> & nucleated_elems_map,
451 : const Real nucleationRadius)
452 : {
453 : // remove nucleated elements that are too close too each other. Lowest key wins
454 710 : for (auto it1 = nucleated_elems_map.begin(); it1 != nucleated_elems_map.end(); ++it1)
455 : {
456 530 : std::pair<RealVectorValue, RealVectorValue> nodes = it1->second;
457 : Point p2 = nodes.first;
458 : Point p1 = nodes.second;
459 : Point p = p1 + (p2 - p1) / 2;
460 4608 : for (auto it2 = nucleated_elems_map.begin(); it2 != nucleated_elems_map.end();)
461 : {
462 4078 : if (it1 == it2)
463 : {
464 : ++it2;
465 530 : continue;
466 : }
467 :
468 : nodes = it2->second;
469 : p2 = nodes.first;
470 : p1 = nodes.second;
471 : Point q = p1 + (p2 - p1) / 2;
472 : Point pq = q - p;
473 3548 : if (pq.norm() <= nucleationRadius)
474 478 : it2 = nucleated_elems_map.erase(it2);
475 : else
476 : ++it2;
477 : }
478 : }
479 180 : }
480 :
481 : void
482 180 : MeshCut2DUserObjectBase::removeNucleatedCracksTooCloseToExistingCracks(
483 : std::map<unsigned int, std::pair<RealVectorValue, RealVectorValue>> & nucleated_elems_map,
484 : const Real nucleationRadius)
485 : {
486 710 : for (auto it = nucleated_elems_map.begin(); it != nucleated_elems_map.end();)
487 : {
488 530 : std::pair<RealVectorValue, RealVectorValue> nodes = it->second;
489 : Point p2 = nodes.first;
490 : Point p1 = nodes.second;
491 : Point p = p1 + (p2 - p1) / 2;
492 : bool removeNucleatedElem = false;
493 : for (const auto & cutter_elem :
494 8312 : as_range(_cutter_mesh->active_elements_begin(), _cutter_mesh->active_elements_end()))
495 : {
496 3819 : const Node * const * cutter_elem_nodes = cutter_elem->get_nodes();
497 3819 : Point m = *cutter_elem_nodes[1] - *cutter_elem_nodes[0];
498 3819 : Real t = m * (p - *cutter_elem_nodes[0]) / m.norm_sq();
499 : Real d = std::numeric_limits<Real>::max();
500 3819 : if (t <= 0)
501 : {
502 : Point j = p - *cutter_elem_nodes[0];
503 1610 : d = j.norm();
504 : }
505 2209 : else if (t >= 0)
506 : {
507 : Point j = p - *cutter_elem_nodes[1];
508 2209 : d = j.norm();
509 : }
510 : else
511 : {
512 : Point j = p - (*cutter_elem_nodes[0] + t * m);
513 0 : d = j.norm();
514 : }
515 3819 : if (d <= nucleationRadius)
516 : {
517 : removeNucleatedElem = true;
518 : break;
519 : }
520 530 : }
521 530 : if (removeNucleatedElem)
522 458 : it = nucleated_elems_map.erase(it);
523 : else
524 : ++it;
525 : }
526 180 : }
|
