https://mooseframework.inl.gov
MeshInfo.C
Go to the documentation of this file.
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 "MeshInfo.h"
11 #include "SubProblem.h"
12 #include "libmesh/system.h"
13 #include "libmesh/equation_systems.h"
14 #include "libmesh/parallel_sync.h"
15 
16 registerMooseObject("MooseApp", MeshInfo);
17 
20 {
22  params.addClassDescription(
23  "Report mesh information, such as the number of elements, nodes, and degrees of freedom.");
24 
25  MultiMooseEnum items(
26  "num_dofs num_dofs_nonlinear num_dofs_auxiliary num_elements num_nodes num_local_dofs "
27  "num_local_dofs_nonlinear num_local_dofs_auxiliary num_local_elements num_local_nodes "
28  "local_sidesets local_sideset_elems sidesets sideset_elems local_subdomains "
29  "local_subdomain_elems subdomains subdomain_elems");
30  params.addParam<MultiMooseEnum>(
31  "items",
32  items,
33  "The iteration information to output, if nothing is provided everything will be output.");
34 
35  return params;
36 }
37 
39  : GeneralReporter(parameters),
40  _items(getParam<MultiMooseEnum>("items")),
41  _num_dofs(declareHelper<unsigned int>("num_dofs", REPORTER_MODE_REPLICATED)),
42  _num_dofs_nl(declareHelper<unsigned int>("num_dofs_nonlinear", REPORTER_MODE_REPLICATED)),
43  _num_dofs_aux(declareHelper<unsigned int>("num_dofs_auxiliary", REPORTER_MODE_REPLICATED)),
44  _num_dofs_constrained(
45  declareHelper<unsigned int>("num_dofs_constrained", REPORTER_MODE_REPLICATED)),
46  _num_elem(declareHelper<unsigned int>("num_elements", REPORTER_MODE_REPLICATED)),
47  _num_node(declareHelper<unsigned int>("num_nodes", REPORTER_MODE_REPLICATED)),
48  _num_local_dofs(declareHelper<unsigned int>("num_local_dofs", REPORTER_MODE_DISTRIBUTED)),
49  _num_local_dofs_nl(
50  declareHelper<unsigned int>("num_dofs_local_nonlinear", REPORTER_MODE_DISTRIBUTED)),
51  _num_local_dofs_aux(
52  declareHelper<unsigned int>("num_dofs_local_auxiliary", REPORTER_MODE_DISTRIBUTED)),
53  _num_local_elem(declareHelper<unsigned int>("num_local_elements", REPORTER_MODE_DISTRIBUTED)),
54  _num_local_node(declareHelper<unsigned int>("num_local_nodes", REPORTER_MODE_DISTRIBUTED)),
55 
56  _local_sidesets(declareHelper<std::map<BoundaryID, SidesetInfo>>("local_sidesets",
58  _local_sideset_elems(declareHelper<std::map<BoundaryID, SidesetInfo>>(
59  "local_sideset_elems", REPORTER_MODE_DISTRIBUTED)),
60  _sidesets(declareHelper<std::map<BoundaryID, SidesetInfo>>("sidesets", REPORTER_MODE_ROOT)),
61  _sideset_elems(
62  declareHelper<std::map<BoundaryID, SidesetInfo>>("sideset_elems", REPORTER_MODE_ROOT)),
63 
64  _local_subdomains(declareHelper<std::map<SubdomainID, SubdomainInfo>>(
65  "local_subdomains", REPORTER_MODE_DISTRIBUTED)),
66  _local_subdomain_elems(declareHelper<std::map<SubdomainID, SubdomainInfo>>(
67  "local_subdomain_elems", REPORTER_MODE_DISTRIBUTED)),
68  _subdomains(
69  declareHelper<std::map<SubdomainID, SubdomainInfo>>("subdomains", REPORTER_MODE_ROOT)),
70  _subdomain_elems(
71  declareHelper<std::map<SubdomainID, SubdomainInfo>>("subdomain_elems", REPORTER_MODE_ROOT)),
72 
73  _equation_systems(_fe_problem.es()),
74  _nonlinear_system(_fe_problem.es().get_system("nl0")),
75  _aux_system(_fe_problem.es().get_system("aux0")),
76  _mesh(_fe_problem.mesh().getMesh())
77 {
78 }
79 
80 void
82 {
87  for (auto s : make_range(_equation_systems.n_systems()))
88  _num_dofs_constrained += _equation_systems.get_system(s).n_constrained_dofs();
89 
97 
100 }
101 
102 void
104 {
105  // Helper for adding the sideset names to a given map of sidesets
106  auto add_sideset_names = [&](std::map<BoundaryID, SidesetInfo> & sidesets)
107  {
108  for (auto & pair : sidesets)
109  pair.second.name = _mesh.get_boundary_info().get_sideset_name(pair.second.id);
110  };
111 
112  // Helper for sorting all of the sides in each sideset
113  auto sort_sides = [](std::map<BoundaryID, SidesetInfo> & sidesets)
114  {
115  for (auto & pair : sidesets)
116  std::sort(pair.second.sides.begin(), pair.second.sides.end());
117  };
118 
119  const bool include_all = !_items.isValid();
120 
121  if (include_all || _items.isValueSet("local_sidesets") ||
122  _items.isValueSet("local_sideset_elems") || _items.isValueSet("sideset_elems"))
123  {
124  _local_sidesets.clear();
125  _local_sideset_elems.clear();
126  _sideset_elems.clear();
127 
128  // Fill the local sideset information; all cases need it
129  std::map<BoundaryID, SidesetInfo> sidesets;
130  for (const auto & bnd_elem :
132  if (bnd_elem->_elem->processor_id() == processor_id())
133  {
134  auto & entry = sidesets[bnd_elem->_bnd_id];
135  entry.id = bnd_elem->_bnd_id;
136  entry.sides.emplace_back(bnd_elem->_elem->id(), bnd_elem->_side);
137  }
138 
139  // For local sidesets: copy over the local info, remove the sides, and add the names
140  if (include_all || _items.isValueSet("local_sidesets"))
141  {
142  // Copy over the local sideset info, remove the sides, and add the names
143  _local_sidesets = sidesets;
144  for (auto & pair : _local_sidesets)
145  pair.second.sides.clear();
146  add_sideset_names(_local_sidesets);
147  }
148 
149  // For local sideset elems: copy over the local info, and add the names
150  if (include_all || _items.isValueSet("local_sideset_elems"))
151  {
152  _local_sideset_elems = sidesets;
153  sort_sides(_local_sideset_elems);
154  add_sideset_names(_local_sideset_elems);
155  }
156 
157  // For the global sideset elems, we need to communicate all of the elems
158  if (include_all || _items.isValueSet("sideset_elems"))
159  {
160  // Set up a structure for sending each (id, elem id, side) tuple to root
161  std::map<processor_id_type,
162  std::vector<std::tuple<boundary_id_type, dof_id_type, unsigned int>>>
163  send_info;
164  // Avoid empty sends
165  if (!sidesets.empty())
166  {
167  auto & root_info = send_info[0];
168  for (const auto & pair : sidesets)
169  for (const auto & side : pair.second.sides)
170  root_info.emplace_back(pair.second.id, side.first, side.second);
171  }
172 
173  // Take the received information and insert it into _sideset_elems
174  auto accumulate_info =
175  [this](processor_id_type,
176  const std::vector<std::tuple<boundary_id_type, dof_id_type, unsigned int>> & info)
177  {
178  for (const auto & tuple : info)
179  {
180  const auto id = std::get<0>(tuple);
181  auto & entry = _sideset_elems[id];
182  entry.id = id;
183  entry.sides.emplace_back(std::get<1>(tuple), std::get<2>(tuple));
184  }
185  };
186 
187  // Push the information and insert it into _sideset_elems on root
188  Parallel::push_parallel_vector_data(comm(), send_info, accumulate_info);
189 
190  sort_sides(_sideset_elems);
191  add_sideset_names(_sideset_elems);
192  }
193  }
194 
195  // For global sideset information without elements, we can simplify communication.
196  // All we need are the boundary IDs from libMesh (may not be reduced, so take the union)
197  // and then add the names (global)
198  if (include_all || _items.isValueSet("sidesets"))
199  {
200  _sidesets.clear();
201 
202  auto boundary_ids = _mesh.get_boundary_info().get_boundary_ids();
203  comm().set_union(boundary_ids, 0);
204  if (processor_id() == 0)
205  {
206  for (const auto id : boundary_ids)
207  _sidesets[id].id = id;
208  add_sideset_names(_sidesets);
209  }
210  }
211 }
212 
213 void
214 to_json(nlohmann::json & json, const std::map<BoundaryID, MeshInfo::SidesetInfo> & sidesets)
215 {
216  for (const auto & pair : sidesets)
217  {
218  const MeshInfo::SidesetInfo & sideset_info = pair.second;
219 
220  nlohmann::json sideset_json;
221  sideset_json["id"] = sideset_info.id;
222  if (sideset_info.name.size())
223  sideset_json["name"] = sideset_info.name;
224  if (sideset_info.sides.size())
225  {
226  auto & sides_json = sideset_json["sides"];
227 
228  for (const std::pair<dof_id_type, unsigned int> & pair : sideset_info.sides)
229  {
230  nlohmann::json side_json;
231  side_json["elem_id"] = pair.first;
232  side_json["side"] = pair.second;
233  sides_json.push_back(side_json);
234  }
235  }
236 
237  json.push_back(sideset_json);
238  }
239 }
240 
241 void
242 dataStore(std::ostream & stream, MeshInfo::SidesetInfo & sideset_info, void * context)
243 {
244  storeHelper(stream, sideset_info.id, context);
245  storeHelper(stream, sideset_info.name, context);
246  storeHelper(stream, sideset_info.sides, context);
247 }
248 
249 void
250 dataLoad(std::istream & stream, MeshInfo::SidesetInfo & sideset_info, void * context)
251 {
252  loadHelper(stream, sideset_info.id, context);
253  loadHelper(stream, sideset_info.name, context);
254  loadHelper(stream, sideset_info.sides, context);
255 }
256 
257 void
259 {
260  // Helper for adding the subdomain names to a given map of subdomains
261  auto add_subdomain_names = [&](std::map<SubdomainID, SubdomainInfo> & subdomains)
262  {
263  for (auto & pair : subdomains)
264  pair.second.name = _mesh.subdomain_name(pair.second.id);
265  };
266 
267  // Helper for sorting all of the elems in each subdomain
268  auto sort_elems = [](std::map<SubdomainID, SubdomainInfo> & subdomains)
269  {
270  for (auto & pair : subdomains)
271  std::sort(pair.second.elems.begin(), pair.second.elems.end());
272  };
273 
274  const bool include_all = !_items.isValid();
275 
276  if (include_all || _items.isValueSet("local_subdomains") ||
277  _items.isValueSet("local_subdomain_elems") || _items.isValueSet("subdomain_elems"))
278  {
279  _local_subdomains.clear();
280  _local_subdomain_elems.clear();
281  _subdomain_elems.clear();
282 
283  // Fill the local subdomain information; all cases need it
284  std::map<SubdomainID, SubdomainInfo> subdomains;
285  for (const auto & elem : *_fe_problem.mesh().getActiveLocalElementRange())
286  {
287  auto & entry = subdomains[elem->subdomain_id()];
288  entry.id = elem->subdomain_id();
289  entry.elems.push_back(elem->id());
290  }
291 
292  // For local subdomains: copy over the local info, remove the elems, and add the names
293  if (include_all || _items.isValueSet("local_subdomains"))
294  {
295  _local_subdomains = subdomains;
296  for (auto & pair : _local_subdomains)
297  pair.second.elems.clear();
298  add_subdomain_names(_local_subdomains);
299  }
300 
301  // For local subdomain elems: copy over the local info, and add the names
302  if (include_all || _items.isValueSet("local_subdomain_elems"))
303  {
304  _local_subdomain_elems = subdomains;
305  sort_elems(_local_subdomain_elems);
306  add_subdomain_names(_local_subdomain_elems);
307  }
308 
309  // For the global subdomain elems, we need to communicate all of the elems
310  if (include_all || _items.isValueSet("subdomain_elems"))
311  {
312  // Set up a structure for sending each (id, elem id) to root
313  std::map<processor_id_type, std::vector<std::pair<subdomain_id_type, dof_id_type>>> send_info;
314  // Avoid creating empty entry
315  if (!subdomains.empty())
316  {
317  auto & root_info = send_info[0];
318  for (const auto & pair : subdomains)
319  for (const auto elem_id : pair.second.elems)
320  root_info.emplace_back(pair.second.id, elem_id);
321  }
322 
323  // Take the received information and insert it into _subdomain_elems
324  auto accumulate_info =
325  [this](processor_id_type,
326  const std::vector<std::pair<subdomain_id_type, dof_id_type>> & info)
327  {
328  for (const auto & subdomain_elem_pair : info)
329  {
330  auto & entry = _subdomain_elems[subdomain_elem_pair.first];
331  entry.id = subdomain_elem_pair.first;
332  entry.elems.emplace_back(subdomain_elem_pair.second);
333  }
334  };
335 
336  // Push the information and insert it into _subdomain_elems on root
337  Parallel::push_parallel_vector_data(comm(), send_info, accumulate_info);
338 
339  if (processor_id() == 0)
340  {
341  sort_elems(_subdomain_elems);
342  add_subdomain_names(_subdomain_elems);
343  }
344  }
345  }
346 
347  // For global subdomain information without elements, we can simplify communication.
348  // All we need are the subdomain IDs from libMesh and then add the names (global)
349  if (include_all || _items.isValueSet("subdomains"))
350  {
351  _subdomains.clear();
352 
353  std::set<subdomain_id_type> subdomain_ids;
354  _mesh.subdomain_ids(subdomain_ids);
355 
356  if (processor_id() == 0)
357  {
358  for (const auto id : subdomain_ids)
359  _subdomains[id].id = id;
360  add_subdomain_names(_subdomains);
361  }
362  }
363 }
364 
365 void
366 to_json(nlohmann::json & json, const std::map<SubdomainID, MeshInfo::SubdomainInfo> & subdomains)
367 {
368  for (const auto & pair : subdomains)
369  {
370  const MeshInfo::SubdomainInfo & subdomain_info = pair.second;
371 
372  nlohmann::json subdomain_json;
373  subdomain_json["id"] = subdomain_info.id;
374  if (subdomain_info.name.size())
375  subdomain_json["name"] = subdomain_info.name;
376  if (subdomain_info.elems.size())
377  {
378  auto & sides_json = subdomain_json["elems"];
379  for (const auto & id : subdomain_info.elems)
380  sides_json.push_back(id);
381  }
382 
383  json.push_back(subdomain_json);
384  }
385 }
386 
387 void
388 dataStore(std::ostream & stream, MeshInfo::SubdomainInfo & subdomain_info, void * context)
389 {
390  storeHelper(stream, subdomain_info.id, context);
391  storeHelper(stream, subdomain_info.name, context);
392  storeHelper(stream, subdomain_info.elems, context);
393 }
394 
395 void
396 dataLoad(std::istream & stream, MeshInfo::SubdomainInfo & subdomain_info, void * context)
397 {
398  loadHelper(stream, subdomain_info.id, context);
399  loadHelper(stream, subdomain_info.name, context);
400  loadHelper(stream, subdomain_info.elems, context);
401 }
virtual bnd_elem_iterator bndElemsEnd()
Definition: MooseMesh.C:1607
Helper struct for defining information about a single sideset.
Definition: MeshInfo.h:36
std::vector< dof_id_type > elems
Definition: MeshInfo.h:50
libMesh::ConstElemRange * getActiveLocalElementRange()
Return pointers to range objects for various types of ranges (local nodes, boundary elems...
Definition: MooseMesh.C:1276
std::string name
Definition: MeshInfo.h:39
std::string name
Definition: MeshInfo.h:49
unsigned int n_systems() const
static InputParameters validParams()
Definition: MeshInfo.C:19
unsigned int & _num_local_elem
Definition: MeshInfo.h:66
void possiblyAddSidesetInfo()
Possibly add to _local_sidesets, _local_sideset_elems, _sidesets, and _sideset_elems.
Definition: MeshInfo.C:103
const libMesh::MeshBase & _mesh
Definition: MeshInfo.h:90
unsigned int & _num_dofs_constrained
Definition: MeshInfo.h:60
MPI_Info info
Reporter object that has a single execution of the "execute" method for for each execute flag...
std::size_t n_dofs() const
const ReporterMode REPORTER_MODE_ROOT
virtual void execute() override
Execute method.
Definition: MeshInfo.C:81
virtual bool isValid() const override
IsValid.
std::map< SubdomainID, SubdomainInfo > & _subdomains
Definition: MeshInfo.h:74
MeshBase & mesh
unsigned int & _num_dofs_aux
Definition: MeshInfo.h:59
dof_id_type n_local_nodes() const
The main MOOSE class responsible for handling user-defined parameters in almost every MOOSE system...
std::map< BoundaryID, SidesetInfo > & _sideset_elems
Definition: MeshInfo.h:71
const Parallel::Communicator & comm() const
virtual bnd_elem_iterator bndElemsBegin()
Return iterators to the beginning/end of the boundary elements list.
Definition: MooseMesh.C:1599
static InputParameters validParams()
dof_id_type n_local_dofs() const
const BoundaryInfo & get_boundary_info() const
const T_sys & get_system(std::string_view name) const
unsigned int & _num_local_dofs_aux
Definition: MeshInfo.h:65
unsigned int & _num_local_node
Definition: MeshInfo.h:67
dof_id_type n_local_elem() const
void to_json(nlohmann::json &json, const std::map< BoundaryID, MeshInfo::SidesetInfo > &sidesets)
Definition: MeshInfo.C:214
dof_id_type n_dofs() const
MeshInfo(const InputParameters &parameters)
Definition: MeshInfo.C:38
unsigned int & _num_dofs_nl
Definition: MeshInfo.h:58
Report mesh information, such as the number of elements, nodes, and degrees of freedom.
Definition: MeshInfo.h:24
void storeHelper(std::ostream &stream, P &data, void *context)
Scalar helper routine.
Definition: DataIO.h:893
uint8_t processor_id_type
const libMesh::EquationSystems & _equation_systems
Definition: MeshInfo.h:87
boundary_id_type BoundaryID
SimpleRange< IndexType > as_range(const std::pair< IndexType, IndexType > &p)
std::vector< std::pair< dof_id_type, unsigned int > > sides
Definition: MeshInfo.h:40
void subdomain_ids(std::set< subdomain_id_type > &ids, const bool global=true) const
std::string & subdomain_name(subdomain_id_type id)
void possiblyAddSubdomainInfo()
Possibly add to _local_subdomains, _local_subdomain_elems, _subdomains, and _subdomain_elems.
Definition: MeshInfo.C:258
std::map< BoundaryID, SidesetInfo > & _local_sidesets
Definition: MeshInfo.h:68
bool isValueSet(const std::string &value) const
Methods for seeing if a value is set in the MultiMooseEnum.
const ReporterMode REPORTER_MODE_DISTRIBUTED
unsigned int & _num_node
Definition: MeshInfo.h:62
Helper struct for defining information about a single subdomain.
Definition: MeshInfo.h:46
const std::set< boundary_id_type > & get_boundary_ids() const
const libMesh::System & _aux_system
Definition: MeshInfo.h:89
FEProblemBase & _fe_problem
Reference to the FEProblemBase for this user object.
Definition: UserObject.h:211
std::map< BoundaryID, SidesetInfo > & _sidesets
Definition: MeshInfo.h:70
if(!dmm->_nl) SETERRQ(PETSC_COMM_WORLD
void dataLoad(std::istream &stream, MeshInfo::SidesetInfo &sideset_info, void *context)
Definition: MeshInfo.C:250
const std::string & get_sideset_name(boundary_id_type id) const
std::map< SubdomainID, SubdomainInfo > & _local_subdomain_elems
Definition: MeshInfo.h:73
IntRange< T > make_range(T beg, T end)
virtual MooseMesh & mesh() override
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 optional parameter and a documentation string to the InputParameters object...
std::map< SubdomainID, SubdomainInfo > & _local_subdomains
Definition: MeshInfo.h:72
const ReporterMode REPORTER_MODE_REPLICATED
unsigned int & _num_local_dofs
Definition: MeshInfo.h:63
unsigned int & _num_dofs
Definition: MeshInfo.h:57
std::map< SubdomainID, SubdomainInfo > & _subdomain_elems
Definition: MeshInfo.h:75
This is a "smart" enum class intended to replace many of the shortcomings in the C++ enum type...
const MultiMooseEnum & _items
Definition: MeshInfo.h:54
void dataStore(std::ostream &stream, MeshInfo::SidesetInfo &sideset_info, void *context)
Definition: MeshInfo.C:242
virtual dof_id_type n_elem() const=0
processor_id_type processor_id() const
unsigned int & _num_local_dofs_nl
Definition: MeshInfo.h:64
std::map< BoundaryID, SidesetInfo > & _local_sideset_elems
Definition: MeshInfo.h:69
void loadHelper(std::istream &stream, P &data, void *context)
Scalar helper routine.
Definition: DataIO.h:985
unsigned int & _num_elem
Definition: MeshInfo.h:61
const libMesh::System & _nonlinear_system
Definition: MeshInfo.h:88
void ErrorVector unsigned int
registerMooseObject("MooseApp", MeshInfo)
virtual dof_id_type n_nodes() const=0
void set_union(T &data, const unsigned int root_id) const