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 "SideSetsFromBoundingBoxGenerator.h"
      11             : #include "Conversion.h"
      12             : #include "MooseTypes.h"
      13             : #include "MooseMeshUtils.h"
      14             : #include "CastUniquePointer.h"
      15             : #include "MooseUtils.h"
      16             : 
      17             : #include "libmesh/distributed_mesh.h"
      18             : #include "libmesh/elem.h"
      19             : 
      20             : #include <typeinfo>
      21             : 
      22             : registerMooseObject("MooseApp", SideSetsFromBoundingBoxGenerator);
      23             : 
      24             : InputParameters
      25        3521 : SideSetsFromBoundingBoxGenerator::validParams()
      26             : {
      27        3521 :   InputParameters params = SideSetsGeneratorBase::validParams();
      28             : 
      29        7042 :   params.addClassDescription("Defines new sidesets using currently-defined sideset IDs inside or "
      30             :                              "outside of a bounding box.");
      31             : 
      32       14084 :   MooseEnum location("INSIDE OUTSIDE", "INSIDE");
      33             : 
      34       14084 :   params.addRequiredParam<RealVectorValue>(
      35             :       "bottom_left", "The bottom left point (in x,y,z with spaces in-between).");
      36       14084 :   params.addRequiredParam<RealVectorValue>(
      37             :       "top_right", "The bottom left point (in x,y,z with spaces in-between).");
      38       21126 :   params.addDeprecatedParam<subdomain_id_type>(
      39             :       "block_id",
      40             :       "Subdomain id to set for inside/outside the bounding box",
      41             :       "The parameter 'block_id' is not used.");
      42       14084 :   params.addRequiredParam<BoundaryName>(
      43             :       "boundary_new", "Boundary on specified block within the bounding box to assign");
      44       10563 :   params.addParam<bool>("boundary_id_overlap",
      45        7042 :                         false,
      46             :                         "Set to true if boundaries need to overlap on sideset to be detected.");
      47       14084 :   params.addParam<MooseEnum>(
      48             :       "location", location, "Control of where the subdomain id is to be set");
      49             : 
      50             :   // TODO: Implement each of these in the generate() routine using utilities in SidesetGeneratorBase
      51        7042 :   params.suppressParameter<bool>("fixed_normal");
      52        3521 :   params.suppressParameter<std::vector<BoundaryName>>("new_boundary");
      53             : 
      54        7042 :   return params;
      55        3521 : }
      56             : 
      57         230 : SideSetsFromBoundingBoxGenerator::SideSetsFromBoundingBoxGenerator(
      58         230 :     const InputParameters & parameters)
      59             :   : SideSetsGeneratorBase(parameters),
      60         230 :     _location(parameters.get<MooseEnum>("location")),
      61         230 :     _bounding_box(MooseUtils::buildBoundingBox(parameters.get<RealVectorValue>("bottom_left"),
      62         230 :                                                parameters.get<RealVectorValue>("top_right"))),
      63         690 :     _boundary_id_overlap(parameters.get<bool>("boundary_id_overlap"))
      64             : {
      65         230 :   if (_boundary_id_overlap)
      66             :   {
      67             :     const std::vector<std::string> incompatible_params = {"normal",
      68             :                                                           "replace",
      69             :                                                           "include_only_external_sides",
      70             :                                                           "included_subdomains",
      71         126 :                                                           "included_neighbors"};
      72         378 :     for (const auto & param_name : incompatible_params)
      73         315 :       if (isParamSetByUser(param_name))
      74           0 :         paramError(param_name, "Parameter should not be used with boundary_id_overlap = true.");
      75          63 :   }
      76             : 
      77         230 :   _boundary_names.push_back(parameters.get<BoundaryName>("boundary_new"));
      78         230 : }
      79             : 
      80             : std::unique_ptr<MeshBase>
      81         227 : SideSetsFromBoundingBoxGenerator::generate()
      82             : {
      83         227 :   std::unique_ptr<MeshBase> mesh = std::move(_input);
      84         227 :   if (!mesh->is_replicated())
      85           0 :     mooseError("SideSetsFromBoundingBoxGenerator is not implemented for distributed meshes");
      86             : 
      87             :   // construct the FE object so we can compute normals of faces
      88         227 :   setup(*mesh);
      89             : 
      90             :   // Get a reference to our BoundaryInfo object for later use
      91         227 :   BoundaryInfo & boundary_info = mesh->get_boundary_info();
      92         227 :   boundary_info.build_node_list_from_side_list();
      93             : 
      94         227 :   bool found_element = false;
      95         227 :   bool found_side_sets = false;
      96         227 :   const bool inside = (_location == "INSIDE");
      97             : 
      98             :   // Attempt to get the new boundary id from the name
      99         227 :   auto boundary_id_new = MooseMeshUtils::getBoundaryID(_boundary_names[0], *mesh);
     100             : 
     101             :   // If the new boundary id is not valid, make it instead
     102         227 :   if (boundary_id_new == Moose::INVALID_BOUNDARY_ID)
     103             :   {
     104           8 :     boundary_id_new = MooseMeshUtils::getNextFreeBoundaryID(*mesh);
     105             : 
     106             :     // Write the name alias of the boundary id to the mesh boundary info
     107           8 :     boundary_info.sideset_name(boundary_id_new) = _boundary_names[0];
     108           8 :     boundary_info.nodeset_name(boundary_id_new) = _boundary_names[0];
     109             :   }
     110             : 
     111         227 :   if (!_boundary_id_overlap)
     112             :   {
     113             :     // Request to compute normal vectors
     114         167 :     const std::vector<Point> & face_normals = _fe_face->get_normals();
     115             : 
     116             :     // Loop over the elements
     117       80111 :     for (const auto & elem : mesh->active_element_ptr_range())
     118             :     {
     119             :       // boolean if element centroid is in bounding box
     120       39972 :       bool contains = _bounding_box.contains_point(elem->vertex_average());
     121             : 
     122             :       // check if active elements are found either in or out of the bounding box, apropos "inside"
     123       39972 :       if (contains == inside)
     124             :       {
     125        3537 :         found_element = true;
     126             :         // loop over sides of elements within bounding box
     127       20805 :         for (const auto & side : make_range(elem->n_sides()))
     128             :         {
     129       17268 :           _fe_face->reinit(elem, side);
     130             :           // We'll just use the normal of the first qp
     131       17268 :           const Point face_normal = face_normals[0];
     132             : 
     133       17268 :           if (elemSideSatisfiesRequirements(elem, side, *mesh, _normal, face_normal))
     134             :           {
     135             :             // assign new boundary value to boundary which meets meshmodifier criteria
     136         932 :             if (_replace)
     137           0 :               boundary_info.remove_side(elem, side);
     138         932 :             boundary_info.add_side(elem, side, boundary_id_new);
     139         932 :             found_side_sets = true;
     140             :           }
     141             :         }
     142             :       }
     143         167 :     }
     144         167 :     if (!found_element && inside)
     145           3 :       mooseError("No elements found within the bounding box");
     146             : 
     147         164 :     if (!found_element && !inside)
     148           0 :       mooseError("No elements found outside the bounding box");
     149             : 
     150         164 :     if (!found_side_sets)
     151           3 :       mooseError("No side sets found on active elements within the bounding box");
     152             :   }
     153             : 
     154          60 :   else if (_boundary_id_overlap)
     155             :   {
     156          60 :     if (_included_boundary_ids.size() < 2)
     157           6 :       mooseError("boundary_id_old out of bounds: ",
     158           3 :                  _included_boundary_ids.size(),
     159             :                  " Must be 2 boundary inputs or more.");
     160             : 
     161          57 :     bool found_node = false;
     162             : 
     163             :     // Loop over the elements and assign node set id to nodes within the bounding box
     164      144531 :     for (auto node = mesh->active_nodes_begin(); node != mesh->active_nodes_end(); ++node)
     165             :     {
     166             :       // check if nodes are inside of bounding box
     167       72237 :       if (_bounding_box.contains_point(**node) == inside)
     168             :       {
     169             :         // read out boundary ids for nodes
     170       71970 :         std::vector<boundary_id_type> node_boundary_ids;
     171       71970 :         boundary_info.boundary_ids(*node, node_boundary_ids);
     172             : 
     173             :         // sort boundary ids on node and sort boundary ids provided in input file
     174       71970 :         std::sort(node_boundary_ids.begin(), node_boundary_ids.end());
     175       71970 :         std::sort(_included_boundary_ids.begin(), _included_boundary_ids.end());
     176             : 
     177             :         // check if input boundary ids are all contained in the node
     178             :         // if true, write new boundary id on respective node
     179       71970 :         if (std::includes(node_boundary_ids.begin(),
     180             :                           node_boundary_ids.end(),
     181             :                           _included_boundary_ids.begin(),
     182             :                           _included_boundary_ids.end()))
     183             :         {
     184         391 :           boundary_info.add_node(*node, boundary_id_new);
     185         391 :           found_node = true;
     186             :         }
     187       71970 :       }
     188          57 :     }
     189             : 
     190          57 :     if (!found_node)
     191           6 :       mooseError("No nodes found within the bounding box");
     192             :   }
     193             : 
     194         212 :   mesh->unset_is_prepared();
     195         424 :   return dynamic_pointer_cast<MeshBase>(mesh);
     196         212 : }