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 "MeshRepairGenerator.h"
      11             : #include "CastUniquePointer.h"
      12             : #include "MooseMeshUtils.h"
      13             : 
      14             : #include "libmesh/mesh_tools.h"
      15             : #include "libmesh/mesh_modification.h"
      16             : 
      17             : registerMooseObject("MooseApp", MeshRepairGenerator);
      18             : 
      19             : InputParameters
      20       14389 : MeshRepairGenerator::validParams()
      21             : {
      22             : 
      23       14389 :   InputParameters params = MeshGenerator::validParams();
      24       14389 :   params.addClassDescription(
      25             :       "Mesh generator to perform various improvement / fixing operations on an input mesh");
      26       14389 :   params.addRequiredParam<MeshGeneratorName>("input",
      27             :                                              "Name of the mesh generator providing the mesh");
      28             : 
      29       14389 :   params.addParam<bool>("fix_node_overlap", false, "Whether to merge overlapping nodes");
      30       43167 :   params.addParam<Real>(
      31       28778 :       "node_overlap_tol", 1e-8, "Absolute tolerance for merging overlapping nodes");
      32             : 
      33       43167 :   params.addParam<bool>(
      34       28778 :       "fix_elements_orientation", false, "Whether to flip elements with negative volumes");
      35             : 
      36       43167 :   params.addParam<bool>("separate_blocks_by_element_types",
      37       28778 :                         false,
      38             :                         "Create new blocks if multiple element types are present in a block");
      39             : 
      40       43167 :   params.addParam<bool>("merge_boundary_ids_with_same_name",
      41       28778 :                         false,
      42             :                         "Merge boundaries if they have the same name but different boundary IDs");
      43             : 
      44       14389 :   return params;
      45           0 : }
      46             : 
      47          60 : MeshRepairGenerator::MeshRepairGenerator(const InputParameters & parameters)
      48             :   : MeshGenerator(parameters),
      49          60 :     _input(getMesh("input")),
      50          60 :     _fix_overlapping_nodes(getParam<bool>("fix_node_overlap")),
      51          60 :     _node_overlap_tol(getParam<Real>("node_overlap_tol")),
      52          60 :     _fix_element_orientation(getParam<bool>("fix_elements_orientation")),
      53          60 :     _elem_type_separation(getParam<bool>("separate_blocks_by_element_types")),
      54         120 :     _boundary_id_merge(getParam<bool>("merge_boundary_ids_with_same_name"))
      55             : {
      56          60 :   if (!_fix_overlapping_nodes && !_fix_element_orientation && !_elem_type_separation &&
      57          10 :       !_boundary_id_merge)
      58           0 :     mooseError("No specific item to fix. Are any of the parameters misspelled?");
      59          60 : }
      60             : 
      61             : std::unique_ptr<MeshBase>
      62          52 : MeshRepairGenerator::generate()
      63             : {
      64          52 :   std::unique_ptr<MeshBase> mesh = std::move(_input);
      65          52 :   mesh->prepare_for_use();
      66             : 
      67             :   // Blanket ban on distributed. This can be relaxed for some operations if needed
      68          52 :   if (!mesh->is_serial())
      69           0 :     mooseError("MeshRepairGenerator requires a serial mesh. The mesh should not be distributed.");
      70             : 
      71          52 :   if (_fix_overlapping_nodes)
      72          25 :     fixOverlappingNodes(mesh);
      73             : 
      74             :   // Flip orientation of elements to keep positive volumes
      75          52 :   if (_fix_element_orientation)
      76           9 :     MeshTools::Modification::orient_elements(*mesh);
      77             : 
      78             :   // Disambiguate any block that has elements of multiple types
      79          52 :   if (_elem_type_separation)
      80           9 :     separateSubdomainsByElementType(mesh);
      81             : 
      82             :   // Assign a single boundary ID to boundaries that have the same name
      83          52 :   if (_boundary_id_merge)
      84           9 :     MooseMeshUtils::mergeBoundaryIDsWithSameName(*mesh);
      85             : 
      86          52 :   mesh->set_isnt_prepared();
      87          52 :   return mesh;
      88           0 : }
      89             : 
      90             : void
      91          25 : MeshRepairGenerator::fixOverlappingNodes(std::unique_ptr<MeshBase> & mesh) const
      92             : {
      93          25 :   unsigned int num_fixed_nodes = 0;
      94          25 :   auto pl = mesh->sub_point_locator();
      95          25 :   pl->set_close_to_point_tol(_node_overlap_tol);
      96             : 
      97          25 :   std::unordered_set<dof_id_type> nodes_removed;
      98             :   // loop on nodes
      99         517 :   for (auto & node : mesh->node_ptr_range())
     100             :   {
     101             :     // this node has already been removed
     102         246 :     if (nodes_removed.count(node->id()))
     103          66 :       continue;
     104             : 
     105             :     // find all the elements around this node
     106         180 :     std::set<const Elem *> elements;
     107         180 :     (*pl)(*node, elements);
     108             : 
     109         426 :     for (auto & elem : elements)
     110             :     {
     111         246 :       bool found = false;
     112        1103 :       for (auto & elem_node : elem->node_ref_range())
     113             :       {
     114        1037 :         if (node->id() == elem_node.id())
     115             :         {
     116         180 :           found = true;
     117         180 :           break;
     118             :         }
     119             :       }
     120         246 :       if (!found)
     121             :       {
     122         440 :         for (auto & elem_node : elem->node_ref_range())
     123             :         {
     124         374 :           if (elem_node.id() == node->id())
     125           0 :             continue;
     126         374 :           const Real tol = _node_overlap_tol;
     127             :           // Compares the coordinates
     128         374 :           const auto x_node = (*node)(0);
     129         374 :           const auto x_elem_node = elem_node(0);
     130         374 :           const auto y_node = (*node)(1);
     131         374 :           const auto y_elem_node = elem_node(1);
     132         374 :           const auto z_node = (*node)(2);
     133         374 :           const auto z_elem_node = elem_node(2);
     134             : 
     135         374 :           if (MooseUtils::absoluteFuzzyEqual(x_node, x_elem_node, tol) &&
     136         472 :               MooseUtils::absoluteFuzzyEqual(y_node, y_elem_node, tol) &&
     137          98 :               MooseUtils::absoluteFuzzyEqual(z_node, z_elem_node, tol))
     138             :           {
     139             :             // Merging two nodes from the same element is almost never a good idea
     140          66 :             if (elem->get_node_index(node) != libMesh::invalid_uint)
     141             :             {
     142           0 :               _console << "Two overlapping nodes in element " << elem->id() << " right by "
     143           0 :                        << elem->vertex_average() << ".\n They will not be stitched" << std::endl;
     144           0 :               continue;
     145             :             }
     146             : 
     147             :             // Coordinates are the same but it's not the same node
     148             :             // Replace the node in the element
     149          66 :             const_cast<Elem *>(elem)->set_node(elem->get_node_index(&elem_node), node);
     150          66 :             nodes_removed.insert(elem_node.id());
     151             : 
     152          66 :             num_fixed_nodes++;
     153          66 :             if (num_fixed_nodes < 10)
     154          66 :               _console << "Stitching nodes " << *node << " and            " << elem_node
     155          66 :                        << std::endl;
     156           0 :             else if (num_fixed_nodes == 10)
     157           0 :               _console << "Node stitching will now proceed silently." << std::endl;
     158             :           }
     159             :         }
     160             :       }
     161             :     }
     162         205 :   }
     163          25 :   _console << "Number of overlapping nodes which got merged: " << num_fixed_nodes << std::endl;
     164          25 :   if (mesh->allow_renumbering())
     165           9 :     mesh->renumber_nodes_and_elements();
     166             :   else
     167             :   {
     168          16 :     mesh->remove_orphaned_nodes();
     169          16 :     mesh->update_parallel_id_counts();
     170             :   }
     171          25 : }
     172             : 
     173             : void
     174           9 : MeshRepairGenerator::separateSubdomainsByElementType(std::unique_ptr<MeshBase> & mesh) const
     175             : {
     176           9 :   std::set<subdomain_id_type> ids;
     177           9 :   mesh->subdomain_ids(ids);
     178             :   // loop on sub-domain
     179          18 :   for (const auto id : ids)
     180             :   {
     181             :     // Gather all the element types and blocks
     182             :     // ElemType defines an enum for geometric element types
     183           9 :     std::set<ElemType> types;
     184             :     // loop on elements within this sub-domain
     185          27 :     for (auto & elem : mesh->active_subdomain_elements_ptr_range(id))
     186          27 :       types.insert(elem->type());
     187             : 
     188             :     // This call must be performed on all processes
     189           9 :     auto next_block_id = MooseMeshUtils::getNextFreeSubdomainID(*mesh);
     190             : 
     191           9 :     if (types.size() > 1)
     192             :     {
     193           9 :       subdomain_id_type i = 0;
     194          27 :       for (const auto type : types)
     195             :       {
     196          18 :         auto new_id = next_block_id + i++;
     197             :         // Create blocks when a block has multiple element types
     198          18 :         mesh->subdomain_name(new_id) = mesh->subdomain_name(id) + "_" + Moose::stringify(type);
     199             : 
     200             :         // Re-assign elements to the new blocks
     201          72 :         for (auto elem : mesh->active_subdomain_elements_ptr_range(id))
     202          27 :           if (elem->type() == type)
     203          36 :             elem->subdomain_id() = new_id;
     204             :       }
     205             :     }
     206           9 :   }
     207           9 : }