www.mooseframework.org
SideSetsAroundSubdomainGenerator.C
Go to the documentation of this file.
1 //* This file is part of the MOOSE framework
2 //* https://www.mooseframework.org
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 
11 #include "InputParameters.h"
12 #include "MooseTypes.h"
13 #include "MooseMeshUtils.h"
14 #include "CastUniquePointer.h"
15 
16 #include "libmesh/mesh.h"
17 #include "libmesh/remote_elem.h"
18 #include "libmesh/point.h"
19 #include "libmesh/fe_base.h"
20 
22 
25 {
27 
28  params.addRequiredParam<MeshGeneratorName>("input", "The mesh we want to modify");
29  params.addRequiredParam<std::vector<BoundaryName>>(
30  "new_boundary", "The list of boundary IDs to create on the supplied subdomain");
31  params.addRequiredParam<std::vector<SubdomainName>>("block",
32  "The blocks around which to create sidesets");
33  params.addParam<Point>("normal",
34  "If supplied, only faces with normal equal to this, up to "
35  "normal_tol, will be added to the sidesets specified");
36  params.addParam<bool>("external_only", false, "Only apply the sideset to external boundaries");
37  params.addRangeCheckedParam<Real>("normal_tol",
38  0.1,
39  "normal_tol>=0 & normal_tol<=2",
40  "If normal is supplied then faces are "
41  "only added if face_normal.normal_hat >= "
42  "1 - normal_tol, where normal_hat = "
43  "normal/|normal|");
44 
45  params.addClassDescription(
46  "Adds element faces that are on the exterior of the given block to the sidesets specified");
47 
48  return params;
49 }
50 
52  const InputParameters & parameters)
53  : SideSetsGeneratorBase(parameters),
54  _input(getMesh("input")),
55  _boundary_names(getParam<std::vector<BoundaryName>>("new_boundary")),
56  _using_normal(isParamValid("normal")),
57  _external_only(getParam<bool>("external_only")),
58  _normal_tol(getParam<Real>("normal_tol")),
59  _normal(_using_normal ? getParam<Point>("normal") : Point())
60 {
61  if (_using_normal)
62  {
63  // normalize
64  mooseAssert(_normal.norm() >= 1E-5, "Normal is zero");
65  _normal /= _normal.norm();
66  }
67 }
68 
69 std::unique_ptr<MeshBase>
71 {
72  std::unique_ptr<MeshBase> mesh = std::move(_input);
73 
74  std::vector<SubdomainName> block_names = getParam<std::vector<SubdomainName>>("block");
75  // check that the blocks exist in the mesh
76  for (const auto & name : block_names)
78  paramError("block", "The block '", name, "' was not found in the mesh");
79 
80  auto blocks = MooseMeshUtils::getSubdomainIDs(*mesh, block_names);
81  std::set<subdomain_id_type> block_ids(blocks.begin(), blocks.end());
82 
83  // Create the boundary IDs from the list of names provided (the true flag creates ids from unknown
84  // names)
85  std::vector<boundary_id_type> boundary_ids =
87 
88  // construct the FE object so we can compute normals of faces
89  setup(*mesh);
90  Point face_normal;
91  bool add_to_bdy = true;
92 
93  // Get a reference to our BoundaryInfo object for later use
94  BoundaryInfo & boundary_info = mesh->get_boundary_info();
95 
96  // Prepare to query about sides adjacent to remote elements if we're
97  // on a distributed mesh
98  const processor_id_type my_n_proc = mesh->n_processors();
99  const processor_id_type my_proc_id = mesh->processor_id();
100  typedef std::vector<std::pair<dof_id_type, unsigned int>> vec_type;
101  std::vector<vec_type> queries(my_n_proc);
102 
103  // Loop over the elements
104  for (const auto & elem : mesh->active_element_ptr_range())
105  {
106  subdomain_id_type curr_subdomain = elem->subdomain_id();
107 
108  // We only need to loop over elements in the source subdomain
109  if (block_ids.count(curr_subdomain) == 0)
110  continue;
111 
112  for (unsigned int side = 0; side < elem->n_sides(); ++side)
113  {
114  const Elem * neighbor = elem->neighbor_ptr(side);
115 
116  // On a replicated mesh, we add all subdomain sides ourselves.
117  // On a distributed mesh, we may have missed sides which
118  // neighbor remote elements. We should query any such cases.
119  if (neighbor == remote_elem)
120  {
121  queries[elem->processor_id()].push_back(std::make_pair(elem->id(), side));
122  }
123  else if (neighbor == nullptr || // element on boundary OR
124  (!_external_only && block_ids.count(neighbor->subdomain_id()) ==
125  0)) // neighboring element is on a different subdomain
126  {
127  if (_using_normal)
128  {
129  const std::vector<Point> & normals = _fe_face->get_normals();
130  _fe_face->reinit(elem, side);
131  face_normal = normals[0];
132  add_to_bdy = (_normal * face_normal >= 1.0 - _normal_tol);
133  }
134 
135  // Add the boundaries, if appropriate
136  if (add_to_bdy)
137  for (const auto & boundary_id : boundary_ids)
138  boundary_info.add_side(elem, side, boundary_id);
139  }
140  }
141  }
142 
143  if (!mesh->is_serial())
144  {
145  const auto queries_tag = mesh->comm().get_unique_tag(),
146  replies_tag = mesh->comm().get_unique_tag();
147 
148  std::vector<Parallel::Request> side_requests(my_n_proc - 1), reply_requests(my_n_proc - 1);
149 
150  // Make all requests
151  for (processor_id_type p = 0; p != my_n_proc; ++p)
152  {
153  if (p == my_proc_id)
154  continue;
155 
156  Parallel::Request & request = side_requests[p - (p > my_proc_id)];
157 
158  mesh->comm().send(p, queries[p], request, queries_tag);
159  }
160 
161  // Reply to all requests
162  std::vector<vec_type> responses(my_n_proc - 1);
163 
164  for (processor_id_type p = 1; p != my_n_proc; ++p)
165  {
166  vec_type query;
167 
168  Parallel::Status status(mesh->comm().probe(Parallel::any_source, queries_tag));
169  const processor_id_type source_pid = cast_int<processor_id_type>(status.source());
170 
171  mesh->comm().receive(source_pid, query, queries_tag);
172 
173  Parallel::Request & request = reply_requests[p - 1];
174 
175  for (const auto & q : query)
176  {
177  const Elem * elem = mesh->elem_ptr(q.first);
178  const unsigned int side = q.second;
179  const Elem * neighbor = elem->neighbor_ptr(side);
180 
181  if (neighbor == nullptr || // element on boundary OR
182  block_ids.count(neighbor->subdomain_id()) ==
183  0) // neighboring element is on a different subdomain
184  {
185  if (_using_normal)
186  {
187  const std::vector<Point> & normals = _fe_face->get_normals();
188  _fe_face->reinit(elem, side);
189  face_normal = normals[0];
190  add_to_bdy = (_normal * face_normal >= 1.0 - _normal_tol);
191  }
192 
193  // Add the boundaries, if appropriate
194  if (add_to_bdy)
195  responses[p - 1].push_back(std::make_pair(elem->id(), side));
196  }
197  }
198 
199  mesh->comm().send(source_pid, responses[p - 1], request, replies_tag);
200  }
201 
202  // Process all incoming replies
203  for (processor_id_type p = 1; p != my_n_proc; ++p)
204  {
205  Parallel::Status status(this->comm().probe(Parallel::any_source, replies_tag));
206  const processor_id_type source_pid = cast_int<processor_id_type>(status.source());
207 
208  vec_type response;
209 
210  this->comm().receive(source_pid, response, replies_tag);
211 
212  for (const auto & r : response)
213  {
214  const Elem * elem = mesh->elem_ptr(r.first);
215  const unsigned int side = r.second;
216 
217  for (const auto & boundary_id : boundary_ids)
218  boundary_info.add_side(elem, side, boundary_id);
219  }
220  }
221 
222  Parallel::wait(side_requests);
223  Parallel::wait(reply_requests);
224  }
225 
226  // Assign the supplied names to the newly created side sets
227  for (unsigned int i = 0; i < boundary_ids.size(); ++i)
228  boundary_info.sideset_name(boundary_ids[i]) = _boundary_names[i];
229 
230  mesh->set_isnt_prepared();
231  return dynamic_pointer_cast<MeshBase>(mesh);
232 }
SideSetsAroundSubdomainGenerator(const InputParameters &parameters)
MPI_Request request
std::vector< BoundaryName > _boundary_names
names of the sidesets to which the faces will be added
Point _normal
if specified, then faces are only added if their normal is close to this
MeshBase & mesh
The main MOOSE class responsible for handling user-defined parameters in almost every MOOSE system...
const Parallel::Communicator & comm() const
std::unique_ptr< T_DEST, T_DELETER > dynamic_pointer_cast(std::unique_ptr< T_SRC, T_DELETER > &src)
These are reworked from https://stackoverflow.com/a/11003103.
std::vector< subdomain_id_type > getSubdomainIDs(const libMesh::MeshBase &mesh, const std::vector< SubdomainName > &subdomain_name)
Get the associated subdomainIDs for the subdomain names that are passed in.
virtual const std::string & name() const
Get the name of the class.
Definition: MooseBase.h:56
void addRequiredParam(const std::string &name, const std::string &doc_string)
This method adds a parameter and documentation string to the InputParameters object that will be extr...
MPI_Status status
uint8_t processor_id_type
std::unique_ptr< MeshBase > generate() override
Generate / modify the mesh.
Real _normal_tol
if normal is specified, then faces are only added if face_normal.normal_hat <= 1 - normal_tol where n...
Status receive(const unsigned int dest_processor_id, T &buf, const MessageTag &tag=any_tag) const
std::unique_ptr< FEBase > _fe_face
bool hasSubdomainName(MeshBase &input_mesh, const SubdomainName &name)
Whether a particular subdomain name exists in the mesh.
std::vector< BoundaryID > getBoundaryIDs(const libMesh::MeshBase &mesh, const std::vector< BoundaryName > &boundary_name, bool generate_unknown, const std::set< BoundaryID > &mesh_boundary_ids)
Gets the boundary IDs with their names.
void paramError(const std::string &param, Args... args) const
Emits an error prefixed with the file and line number of the given param (from the input file) along ...
static InputParameters validParams()
query_obj query
const bool _external_only
Only apply the sideset to external boundaries.
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
registerMooseObject("MooseApp", SideSetsAroundSubdomainGenerator)
bool _using_normal
true if only faces close to "normal" will be added
void addClassDescription(const std::string &doc_string)
This method adds a description of the class that will be displayed in the input file syntax dump...
void addParam(const std::string &name, const S &value, const std::string &doc_string)
These methods add an option parameter and a documentation string to the InputParameters object...
void addRangeCheckedParam(const std::string &name, const T &value, const std::string &parsed_function, const std::string &doc_string)
Adds the faces on the boundary of given block to the sidesets specified by "boundary" Optionally...
void setup(MeshBase &mesh)
This method is used to construct the FE object so we can compute normals of faces.
const RemoteElem * remote_elem