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 "SideSetsFromPointsGenerator.h"
      11             : #include "Parser.h"
      12             : #include "InputParameters.h"
      13             : #include "MeshTraversingUtils.h"
      14             : #include "MooseMeshUtils.h"
      15             : #include "CastUniquePointer.h"
      16             : 
      17             : #include "libmesh/mesh.h"
      18             : #include "libmesh/mesh_generation.h"
      19             : #include "libmesh/mesh_serializer.h"
      20             : #include "libmesh/string_to_enum.h"
      21             : #include "libmesh/quadrature_gauss.h"
      22             : #include "libmesh/point_locator_base.h"
      23             : #include "libmesh/enum_point_locator_type.h"
      24             : #include "libmesh/distributed_mesh.h"
      25             : #include "libmesh/elem.h"
      26             : #include "libmesh/fe_base.h"
      27             : #include "libmesh/mesh_tools.h"
      28             : 
      29             : #include <typeinfo>
      30             : 
      31             : registerMooseObject("MooseApp", SideSetsFromPointsGenerator);
      32             : 
      33             : InputParameters
      34        3101 : SideSetsFromPointsGenerator::validParams()
      35             : {
      36        3101 :   InputParameters params = SideSetsGeneratorBase::validParams();
      37             : 
      38        6202 :   params.addClassDescription("Adds a new sideset starting at the specified point containing all "
      39             :                              "connected element faces with the same normal.");
      40       12404 :   params.addRequiredParam<std::vector<Point>>(
      41             :       "points", "A list of points from which to start painting sidesets");
      42             : 
      43        6202 :   params.suppressParameter<Point>("normal");
      44        3101 :   params.suppressParameter<Real>("normal_tol");
      45             : 
      46             :   // It doesn't make sense to allow internal sides for this side set generator.
      47        9303 :   params.setParameters("include_only_external_sides", true);
      48        3101 :   params.suppressParameter<bool>("include_only_external_sides");
      49             : 
      50        3101 :   return params;
      51           0 : }
      52             : 
      53          20 : SideSetsFromPointsGenerator::SideSetsFromPointsGenerator(const InputParameters & parameters)
      54          40 :   : SideSetsGeneratorBase(parameters), _points(getParam<std::vector<Point>>("points"))
      55             : {
      56          20 :   if (_points.size() != _boundary_names.size())
      57           0 :     mooseError("point list and boundary list are not the same length");
      58          20 :   _using_normal = true;
      59          20 : }
      60             : 
      61             : std::unique_ptr<MeshBase>
      62          19 : SideSetsFromPointsGenerator::generate()
      63             : {
      64          19 :   std::unique_ptr<MeshBase> mesh = std::move(_input);
      65             : 
      66             :   // Our flood fill doesn't do any communication, so it requires a
      67             :   // serialized mesh
      68          19 :   libMesh::MeshSerializer serial(*mesh);
      69             : 
      70             :   // Get the BoundaryIDs from the mesh
      71             :   std::vector<BoundaryID> boundary_ids =
      72          19 :       MooseMeshUtils::getBoundaryIDs(*mesh, _boundary_names, true);
      73             : 
      74          19 :   setup(*mesh);
      75             : 
      76          19 :   _visited.clear();
      77             : 
      78             :   std::unique_ptr<libMesh::PointLocatorBase> pl =
      79          19 :       libMesh::PointLocatorBase::build(libMesh::TREE, *mesh);
      80             : 
      81          75 :   for (const auto i : index_range(boundary_ids))
      82             :   {
      83          59 :     std::set<const Elem *> candidate_elements;
      84          59 :     (*pl)(_points[i], candidate_elements);
      85             : 
      86          59 :     const Elem * elem_to_flood = nullptr;
      87          59 :     Point normal_to_flood;
      88             : 
      89         155 :     for (const Elem * elem : candidate_elements)
      90         550 :       for (const auto side : make_range(elem->n_sides()))
      91             :       {
      92         454 :         if (elem->neighbor_ptr(side))
      93         352 :           continue;
      94             : 
      95             :         // See if this point is on this side
      96         102 :         std::unique_ptr<const Elem> elem_side = elem->side_ptr(side);
      97             : 
      98         102 :         if (elem_side->contains_point(_points[i]))
      99             :         {
     100             :           // This is a good side to paint our sideset with.
     101             :           // Get the normal.
     102         102 :           const std::vector<Point> & normals = _fe_face->get_normals();
     103         102 :           _fe_face->reinit(elem, side);
     104             : 
     105             :           // If we *already* found a good but different side to paint
     106             :           // our sideset with, we've got an ambiguity here.
     107         145 :           if (elem_to_flood &&
     108          43 :               (!MeshTraversingUtils::normalsWithinTol(normal_to_flood, normals[0], _normal_tol) ||
     109          40 :                elem_to_flood->which_neighbor_am_i(elem) == libMesh::invalid_uint))
     110           3 :             mooseError("Two ambiguous potential sideset sources found for boundary `",
     111           3 :                        _boundary_names[i],
     112             :                        "' at ",
     113           3 :                        _points[i],
     114             :                        ":\nElement ",
     115           3 :                        elem_to_flood->id(),
     116             :                        " normal ",
     117             :                        normal_to_flood,
     118             :                        " and\n",
     119             :                        ":\nElement ",
     120           3 :                        elem->id(),
     121             :                        " normal ",
     122           3 :                        normals[0]);
     123             : 
     124          99 :           elem_to_flood = elem;
     125          99 :           normal_to_flood = normals[0];
     126             :           // Don't just flood here; keep looking for possible ambiguity
     127             :         }
     128          99 :       }
     129             : 
     130          56 :     flood(elem_to_flood, normal_to_flood, boundary_ids[i], *mesh);
     131          56 :   }
     132             : 
     133             : #ifdef DEBUG
     134             :   MeshTools::libmesh_assert_valid_boundary_ids(*mesh);
     135             : #endif
     136             : 
     137          16 :   finalize();
     138             : 
     139          72 :   for (const auto i : index_range(boundary_ids))
     140          56 :     mesh->get_boundary_info().sideset_name(boundary_ids[i]) = _boundary_names[i];
     141             : 
     142          16 :   mesh->unset_is_prepared();
     143          32 :   return dynamic_pointer_cast<MeshBase>(mesh);
     144          16 : }