https://mooseframework.inl.gov
SolutionUserObjectBase.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 "SolutionUserObjectBase.h"
11 
12 // MOOSE includes
13 #include "ConsoleUtils.h"
14 #include "MooseError.h"
15 #include "MooseMesh.h"
16 #include "MooseUtils.h"
17 #include "MooseVariableFE.h"
18 #include "RotationMatrix.h"
19 #include "Function.h"
20 
21 // libMesh includes
22 #include "libmesh/equation_systems.h"
23 #include "libmesh/mesh_function.h"
24 #include "libmesh/numeric_vector.h"
25 #include "libmesh/nonlinear_implicit_system.h"
26 #include "libmesh/transient_system.h"
27 #include "libmesh/parallel_mesh.h"
28 #include "libmesh/serial_mesh.h"
29 #include "libmesh/exodusII_io.h"
30 #include "libmesh/exodusII_io_helper.h"
31 #include "libmesh/enum_xdr_mode.h"
32 #include "libmesh/string_to_enum.h"
33 
36 {
37  // Get the input parameters from the parent class
39 
40  // Add required parameters
41  params.addRequiredParam<MeshFileName>(
42  "mesh", "The name of the mesh file (must be xda/xdr or exodusII file).");
43  params.addParam<std::vector<std::string>>(
44  "system_variables",
45  std::vector<std::string>(),
46  "The name of the nodal and elemental variables from the file you want to use for values");
47 
48  // When using XDA/XDR files the following must be defined
49  params.addParam<FileName>(
50  "es",
51  "<not supplied>",
52  "The name of the file holding the equation system info in xda/xdr format (xda/xdr only).");
53  params.addParam<std::string>(
54  "system",
55  "nl0",
56  "The name of the system to pull values out of (xda/xdr only). The default name for the "
57  "nonlinear system is 'nl0', auxiliary system is 'aux0'");
58 
59  // When using ExodusII a specific time is extracted
60  params.addParam<std::string>("timestep",
61  "Index of the single timestep used or \"LATEST\" for "
62  "the last timestep (exodusII only). If not supplied, "
63  "time interpolation will occur.");
64 
65  // Choose the interpolation order for incoming nodal variables
66  params.addParam<MooseEnum>("nodal_variable_order",
67  MooseEnum("FIRST SECOND", "FIRST"),
68  "Specifies the order of the nodal solution data.");
69 
70  // Add ability to perform coordinate transformation: scale, factor
71  params.addParam<std::vector<Real>>(
72  "scale", std::vector<Real>(LIBMESH_DIM, 1), "Scale factor for points in the simulation");
73  params.addParam<std::vector<Real>>("scale_multiplier",
74  std::vector<Real>(LIBMESH_DIM, 1),
75  "Scale multiplying factor for points in the simulation");
76  params.addParam<std::vector<Real>>("translation",
77  std::vector<Real>(LIBMESH_DIM, 0),
78  "Translation factors for x,y,z coordinates of the simulation");
79  params.addParam<RealVectorValue>("rotation0_vector",
80  RealVectorValue(0, 0, 1),
81  "Vector about which to rotate points of the simulation.");
82  params.addParam<Real>(
83  "rotation0_angle",
84  0.0,
85  "Anticlockwise rotation angle (in degrees) to use for rotation about rotation0_vector.");
86  params.addParam<RealVectorValue>("rotation1_vector",
87  RealVectorValue(0, 0, 1),
88  "Vector about which to rotate points of the simulation.");
89  params.addParam<Real>(
90  "rotation1_angle",
91  0.0,
92  "Anticlockwise rotation angle (in degrees) to use for rotation about rotation1_vector.");
93 
94  // following lines build the default_transformation_order
95  MultiMooseEnum default_transformation_order(
96  "rotation0 translation scale rotation1 scale_multiplier", "translation scale");
97  params.addParam<MultiMooseEnum>(
98  "transformation_order",
99  default_transformation_order,
100  "The order to perform the operations in. Define R0 to be the rotation matrix encoded by "
101  "rotation0_vector and rotation0_angle. Similarly for R1. Denote the scale by s, the "
102  "scale_multiplier by m, and the translation by t. Then, given a point x in the simulation, "
103  "if transformation_order = 'rotation0 scale_multiplier translation scale rotation1' then "
104  "form p = R1*(R0*x*m - t)/s. Then the values provided by the SolutionUserObjectBase at "
105  "point x "
106  "in the simulation are the variable values at point p in the mesh.");
107  params.addParamNamesToGroup("scale scale_multiplier translation rotation0_vector rotation0_angle "
108  "rotation1_angle transformation_order",
109  "Coordinate system transformation");
110  params.addClassDescription("Reads a variable from a mesh in one simulation to another");
111  return params;
112 }
113 
114 // Static mutex definition
116 
118  : GeneralUserObject(parameters),
119  _file_type(MooseEnum("xda=0 exodusII=1 xdr=2")),
120  _mesh_file(getParam<MeshFileName>("mesh")),
121  _es_file(getParam<FileName>("es")),
122  _system_name(getParam<std::string>("system")),
123  _system_variables(getParam<std::vector<std::string>>("system_variables")),
124  _exodus_time_index(-1),
125  _interpolate_times(false),
126  _system(nullptr),
127  _system2(nullptr),
128  _interpolation_time(0.0),
129  _interpolation_factor(0.0),
130  _exodus_times(nullptr),
131  _exodus_index1(-1),
132  _exodus_index2(-1),
133  _nodal_variable_order(getParam<MooseEnum>("nodal_variable_order")),
134  _scale(getParam<std::vector<Real>>("scale")),
135  _scale_multiplier(getParam<std::vector<Real>>("scale_multiplier")),
136  _translation(getParam<std::vector<Real>>("translation")),
137  _rotation0_vector(getParam<RealVectorValue>("rotation0_vector")),
138  _rotation0_angle(getParam<Real>("rotation0_angle")),
139  _r0(RealTensorValue()),
140  _rotation1_vector(getParam<RealVectorValue>("rotation1_vector")),
141  _rotation1_angle(getParam<Real>("rotation1_angle")),
142  _r1(RealTensorValue()),
143  _transformation_order(getParam<MultiMooseEnum>("transformation_order")),
144  _initialized(false)
145 {
146  // form rotation matrices with the specified angles
147  Real halfPi = libMesh::pi / 2.0;
148  Real a;
149  Real b;
150 
151  a = std::cos(halfPi * -_rotation0_angle / 90.0);
152  b = std::sin(halfPi * -_rotation0_angle / 90.0);
153  // the following is an anticlockwise rotation about z
154  RealTensorValue rot0_z(a, -b, 0, b, a, 0, 0, 0, 1);
155  // form the rotation matrix that will take rotation0_vector to the z axis
157  // _r0 is then: rotate points so vec0 lies along z; then rotate about angle0; then rotate points
158  // back
159  _r0 = vec0_to_z.transpose() * (rot0_z * vec0_to_z);
160 
161  a = std::cos(halfPi * -_rotation1_angle / 90.0);
162  b = std::sin(halfPi * -_rotation1_angle / 90.0);
163  // the following is an anticlockwise rotation about z
164  RealTensorValue rot1_z(a, -b, 0, b, a, 0, 0, 0, 1);
165  // form the rotation matrix that will take rotation1_vector to the z axis
167  // _r1 is then: rotate points so vec1 lies along z; then rotate about angle1; then rotate points
168  // back
169  _r1 = vec1_to_z.transpose() * (rot1_z * vec1_to_z);
170 
171  if (isParamValid("timestep") && getParam<std::string>("timestep") == "-1")
172  mooseError("A \"timestep\" of -1 is no longer supported for interpolation. Instead simply "
173  "remove this parameter altogether for interpolation");
174 
175  // check for possible inconsistencies between mesh and input data
177  mooseInfo(
178  "The mesh has second-order elements, be sure to set 'nodal_variable_order' if needed.");
179 }
180 
181 void
183 {
184  if (!isParamSetByUser("es"))
185  paramError("es", "Equation system file (.xda or .xdr) should have been specified");
186 
187  // Check that the required files exist
190 
191  // Read the libmesh::mesh from the xda file
192  _mesh->read(_mesh_file);
193 
194  // Create the libmesh::EquationSystems
195  _es = std::make_unique<EquationSystems>(*_mesh);
196 
197  // Use new read syntax (binary)
198  if (_file_type == "xdr")
199  _es->read(_es_file,
201  EquationSystems::READ_HEADER | EquationSystems::READ_DATA |
202  EquationSystems::READ_ADDITIONAL_DATA);
203 
204  // Use new read syntax
205  else if (_file_type == "xda")
206  _es->read(_es_file,
208  EquationSystems::READ_HEADER | EquationSystems::READ_DATA |
209  EquationSystems::READ_ADDITIONAL_DATA);
210 
211  // This should never occur, just in case produce an error
212  else
213  mooseError("Failed to determine proper read method for XDA/XDR equation system file: ",
214  _es_file);
215 
216  // Update and store the EquationSystems name locally
217  _es->update();
218  _system = &_es->get_system(_system_name);
219 }
220 
221 void
223 {
224  // Define a default system name
225  if (_system_name == "")
226  _system_name = "SolutionUserObjectSystem";
227 
228  // Read the Exodus file
229  _exodusII_io = std::make_unique<libMesh::ExodusII_IO>(*_mesh);
230  _exodusII_io->read(_mesh_file);
232  _exodus_times = &_exodusII_io->get_time_steps();
233 
234  if (isParamValid("timestep"))
235  {
236  std::string s_timestep = getParam<std::string>("timestep");
237  int n_steps = _exodusII_io->get_num_time_steps();
238  if (s_timestep == "LATEST")
239  _exodus_time_index = n_steps;
240  else
241  {
242  std::istringstream ss(s_timestep);
243  if (!((ss >> _exodus_time_index) && ss.eof()) || _exodus_time_index > n_steps)
244  mooseError("Invalid value passed as \"timestep\". Expected \"LATEST\" or a valid integer "
245  "less than ",
246  n_steps,
247  ", received ",
248  s_timestep);
249  }
250  }
251  else
252  // Interpolate between times rather than using values from a set timestep
253  _interpolate_times = true;
254 
255  // Check that the number of time steps is valid
256  int num_exo_times = _exodus_times->size();
257  if (num_exo_times == 0)
258  mooseError("In SolutionUserObjectBase, exodus file contains no timesteps.");
259 
260  // Account for parallel mesh
261  if (dynamic_cast<DistributedMesh *>(_mesh.get()))
262  {
263  _mesh->allow_renumbering(true);
264  _mesh->prepare_for_use();
265  }
266  else
267  {
268  _mesh->allow_renumbering(false);
269  _mesh->prepare_for_use();
270  }
271 
272  // Create EquationSystems object for solution
273  _es = std::make_unique<EquationSystems>(*_mesh);
275  _system = &_es->get_system(_system_name);
276 
277  // Get the variable name lists as set; these need to be sets to perform set_intersection
278  const std::vector<std::string> & all_nodal(_exodusII_io->get_nodal_var_names());
279  const std::vector<std::string> & all_elemental(_exodusII_io->get_elem_var_names());
280  const std::vector<std::string> & all_scalar(_exodusII_io->get_global_var_names());
281 
282  // Build nodal/elemental variable lists, limit to variables listed in 'system_variables', if
283  // provided
284  // This function could be called more than once, so clear the member variables so we don't keep
285  // adding to the vectors
286  _nodal_variables.clear();
287  _elemental_variables.clear();
288  _scalar_variables.clear();
289  if (!_system_variables.empty())
290  {
291  for (const auto & var_name : _system_variables)
292  {
293  if (std::find(all_nodal.begin(), all_nodal.end(), var_name) != all_nodal.end())
294  _nodal_variables.push_back(var_name);
295  if (std::find(all_elemental.begin(), all_elemental.end(), var_name) != all_elemental.end())
296  _elemental_variables.push_back(var_name);
297  if (std::find(all_scalar.begin(), all_scalar.end(), var_name) != all_scalar.end())
298  // Check if the scalar matches any field variables, and ignore the var if it does. This
299  // means its a Postprocessor.
300  if (std::find(begin(_nodal_variables), end(_nodal_variables), var_name) ==
301  _nodal_variables.end() &&
302  std::find(begin(_elemental_variables), end(_elemental_variables), var_name) ==
303  _elemental_variables.end())
304  _scalar_variables.push_back(var_name);
305  }
306  }
307  else
308  {
309  _nodal_variables = all_nodal;
310  _elemental_variables = all_elemental;
311 
312  for (auto var_name : all_scalar)
313  // Check if the scalar matches any field variables, and ignore the var if it does. This means
314  // its a Postprocessor.
315  if (std::find(begin(_nodal_variables), end(_nodal_variables), var_name) ==
316  _nodal_variables.end() &&
317  std::find(begin(_elemental_variables), end(_elemental_variables), var_name) ==
318  _elemental_variables.end())
319  _scalar_variables.push_back(var_name);
320  }
321 
322  // Add the variables to the system
323  for (const auto & var_name : _nodal_variables)
324  _system->add_variable(var_name, Utility::string_to_enum<Order>(_nodal_variable_order));
325 
326  for (const auto & var_name : _elemental_variables)
327  _system->add_variable(var_name, CONSTANT, MONOMIAL);
328 
329  for (const auto & var_name : _scalar_variables)
330  _system->add_variable(var_name, FIRST, SCALAR);
331 
332  // Initialize the equations systems
333  _es->init();
334 
335  // Interpolate times
336  if (_interpolate_times)
337  {
338  // Create a second equation system
339  _es2 = std::make_unique<EquationSystems>(*_mesh);
341  _system2 = &_es2->get_system(_system_name);
342 
343  // Add the variables to the system
344  for (const auto & var_name : _nodal_variables)
345  _system2->add_variable(var_name, Utility::string_to_enum<Order>(_nodal_variable_order));
346 
347  for (const auto & var_name : _elemental_variables)
348  _system2->add_variable(var_name, CONSTANT, MONOMIAL);
349 
350  for (const auto & var_name : _scalar_variables)
351  _system2->add_variable(var_name, FIRST, SCALAR);
352 
353  // Initialize
354  _es2->init();
355 
356  // Update the times for interpolation (initially start at 0)
358 
359  // Copy the solutions from the first system
360  for (const auto & var_name : _nodal_variables)
361  {
362  _exodusII_io->copy_nodal_solution(*_system, var_name, var_name, _exodus_index1 + 1);
363  _exodusII_io->copy_nodal_solution(*_system2, var_name, var_name, _exodus_index2 + 1);
364  }
365 
366  for (const auto & var_name : _elemental_variables)
367  {
368  _exodusII_io->copy_elemental_solution(*_system, var_name, var_name, _exodus_index1 + 1);
369  _exodusII_io->copy_elemental_solution(*_system2, var_name, var_name, _exodus_index2 + 1);
370  }
371 
372  if (_scalar_variables.size() > 0)
373  {
374  _exodusII_io->copy_scalar_solution(
376  _exodusII_io->copy_scalar_solution(
378  }
379 
380  // Update the systems
381  _system->update();
382  _es->update();
383  _system2->update();
384  _es2->update();
385  }
386 
387  // Non-interpolated times
388  else
389  {
390  if (_exodus_time_index > num_exo_times)
391  mooseError("In SolutionUserObjectBase, timestep = ",
393  ", but there are only ",
394  num_exo_times,
395  " time steps.");
396 
397  // Copy the values from the ExodusII file
398  for (const auto & var_name : _nodal_variables)
399  _exodusII_io->copy_nodal_solution(*_system, var_name, var_name, _exodus_time_index);
400 
401  for (const auto & var_name : _elemental_variables)
402  _exodusII_io->copy_elemental_solution(*_system, var_name, var_name, _exodus_time_index);
403 
404  if (!_scalar_variables.empty())
405  _exodusII_io->copy_scalar_solution(
407 
408  // Update the equations systems
409  _system->update();
410  _es->update();
411  }
412 }
413 
414 Real
415 SolutionUserObjectBase::directValue(const Node * node, const std::string & var_name) const
416 {
417  // Get the libmesh variable and system numbers
418  unsigned int var_num = _system->variable_number(var_name);
419  unsigned int sys_num = _system->number();
420 
421  // Get the node id and associated dof
422  dof_id_type node_id = node->id();
423  const Node & sys_node = _system->get_mesh().node_ref(node_id);
424  mooseAssert(sys_node.n_dofs(sys_num, var_num) > 0,
425  "Variable " << var_name << " has no DoFs on node " << sys_node.id());
426  dof_id_type dof_id = sys_node.dof_number(sys_num, var_num, 0);
427 
428  // Return the desired value for the dof
429  return directValue(dof_id);
430 }
431 
432 Real
433 SolutionUserObjectBase::directValue(const Elem * elem, const std::string & var_name) const
434 {
435  // Get the libmesh variable and system numbers
436  unsigned int var_num = _system->variable_number(var_name);
437  unsigned int sys_num = _system->number();
438 
439  // Get the element id and associated dof
440  dof_id_type elem_id = elem->id();
441  const Elem & sys_elem = _system->get_mesh().elem_ref(elem_id);
442  mooseAssert(sys_elem.n_dofs(sys_num, var_num) > 0,
443  "Variable " << var_name << " has no DoFs on element " << sys_elem.id());
444  dof_id_type dof_id = sys_elem.dof_number(sys_num, var_num, 0);
445 
446  // Return the desired value
447  return directValue(dof_id);
448 }
449 
450 void
452 {
453 }
454 
455 void
457 {
458 }
459 
460 void
462 {
463  // Update time interpolation for ExodusII solution
464  if (_file_type == 1 && _interpolate_times)
466 }
467 
468 void
470 {
471 }
472 
473 void
475 {
476 
477  // Make sure this only happens once
478  if (_initialized)
479  return;
480 
481  // Create a libmesh::Mesh object for storing the loaded data.
482  // Several aspects of SolutionUserObjectBase won't work with a DistributedMesh:
483  // .) ExodusII_IO::copy_nodal_solution() doesn't work in parallel.
484  // .) We don't know if directValue will be used, which may request
485  // a value on a Node we don't have.
486  // We force the Mesh used here to be a ReplicatedMesh.
487  _mesh = std::make_unique<ReplicatedMesh>(_communicator);
488 
489  // ExodusII mesh file supplied
490  if (MooseUtils::hasExtension(_mesh_file, "e", /*strip_exodus_ext =*/true) ||
491  MooseUtils::hasExtension(_mesh_file, "exo", true))
492  {
493  _file_type = "exodusII";
494  readExodusII();
495  }
496 
497  // XDA mesh file supplied
498  else if (MooseUtils::hasExtension(_mesh_file, "xda"))
499  {
500  _file_type = "xda";
501  readXda();
502  }
503 
504  // XDR mesh file supplied
505  else if (MooseUtils::hasExtension(_mesh_file, "xdr"))
506  {
507  _file_type = "xdr";
508  readXda();
509  }
510 
511  // Produce an error for an unknown file type
512  else
513  mooseError(
514  "In SolutionUserObjectBase, invalid file type: only .xda, .xdr, .exo and .e supported");
515 
516  // Initialize the serial solution vector
519 
520  // Pull down a full copy of this vector on every processor so we can get values in parallel
522 
523  // Vector of variable numbers to apply the MeshFunction to
524  std::vector<unsigned int> var_nums;
525 
526  // If no variables were given, use all of them
527  if (_system_variables.empty())
528  {
530  for (const auto & var_num : var_nums)
531  _system_variables.push_back(_system->variable_name(var_num));
532  }
533 
534  // Otherwise, gather the numbers for the variables given
535  else
536  {
537  for (const auto & var_name : _system_variables)
538  var_nums.push_back(_system->variable_number(var_name));
539  }
540 
541  // Create the MeshFunction for working with the solution data
542  _mesh_function = std::make_unique<libMesh::MeshFunction>(
543  *_es, *_serialized_solution, _system->get_dof_map(), var_nums);
544  _mesh_function->init();
545 
546  // Tell the MeshFunctions that we might be querying them outside the
547  // mesh, so we can handle any errors at the MOOSE rather than at the
548  // libMesh level.
549  DenseVector<Number> default_values;
550  _mesh_function->enable_out_of_mesh_mode(default_values);
551 
552  // Build second MeshFunction for interpolation
553  if (_interpolate_times)
554  {
555  // Need to pull down a full copy of this vector on every processor so we can get values in
556  // parallel
560 
561  // Create the MeshFunction for the second copy of the data
562  _mesh_function2 = std::make_unique<libMesh::MeshFunction>(
564  _mesh_function2->init();
565  _mesh_function2->enable_out_of_mesh_mode(default_values);
566  }
567 
568  // Populate the MeshFunction variable index
569  for (unsigned int i = 0; i < _system_variables.size(); ++i)
570  {
571  std::string name = _system_variables[i];
573  }
574 
575  // Set initialization flag
576  _initialized = true;
577 }
578 
579 MooseEnum
581 {
582  return _file_type;
583 }
584 
585 void
587 {
588  if (_t != _interpolation_time)
589  {
591  {
592 
593  for (const auto & var_name : _nodal_variables)
594  _exodusII_io->copy_nodal_solution(*_system, var_name, var_name, _exodus_index1 + 1);
595  for (const auto & var_name : _elemental_variables)
596  _exodusII_io->copy_elemental_solution(*_system, var_name, var_name, _exodus_index1 + 1);
597  if (_scalar_variables.size() > 0)
598  _exodusII_io->copy_scalar_solution(
600 
601  _system->update();
602  _es->update();
604 
605  for (const auto & var_name : _nodal_variables)
606  _exodusII_io->copy_nodal_solution(*_system2, var_name, var_name, _exodus_index2 + 1);
607  for (const auto & var_name : _elemental_variables)
608  _exodusII_io->copy_elemental_solution(*_system2, var_name, var_name, _exodus_index2 + 1);
609  if (_scalar_variables.size() > 0)
610  _exodusII_io->copy_scalar_solution(
612 
613  _system2->update();
614  _es2->update();
616  }
618  }
619 }
620 
621 bool
623 {
624  if (_file_type != 1)
625  mooseError("In SolutionUserObjectBase, getTimeInterpolationData only applicable for exodusII "
626  "file type");
627 
628  int old_index1 = _exodus_index1;
629  int old_index2 = _exodus_index2;
630 
631  int num_exo_times = _exodus_times->size();
632 
633  const auto solution_time = solutionSampleTime();
634 
635  if (solution_time < (*_exodus_times)[0])
636  {
637  _exodus_index1 = 0;
638  _exodus_index2 = 0;
639  _interpolation_factor = 0.0;
640  }
641  else
642  {
643  for (int i = 0; i < num_exo_times - 1; ++i)
644  {
645  if (solution_time <= (*_exodus_times)[i + 1])
646  {
647  _exodus_index1 = i;
648  _exodus_index2 = i + 1;
650  (solution_time - (*_exodus_times)[i]) / ((*_exodus_times)[i + 1] - (*_exodus_times)[i]);
651  break;
652  }
653  else if (i == num_exo_times - 2)
654  {
655  _exodus_index1 = num_exo_times - 1;
656  _exodus_index2 = num_exo_times - 1;
657  _interpolation_factor = 1.0;
658  break;
659  }
660  }
661  }
662 
663  bool indices_modified(false);
664 
665  if (_exodus_index1 != old_index1 || _exodus_index2 != old_index2)
666  indices_modified = true;
667 
668  return indices_modified;
669 }
670 
671 unsigned int
672 SolutionUserObjectBase::getLocalVarIndex(const std::string & var_name) const
673 {
674  // Extract the variable index for the MeshFunction(s)
675  std::map<std::string, unsigned int>::const_iterator it = _local_variable_index.find(var_name);
676  if (it == _local_variable_index.end())
677  mooseError("Value requested for nonexistent variable '",
678  var_name,
679  "' in the '",
680  name(),
681  "' SolutionUserObjectBase.\nSystem selected: ",
682  _system_name,
683  "\nAvailable variables:\n",
685  return it->second;
686 }
687 
688 Real
690  const Point & p,
691  const std::string & var_name,
692  const MooseEnum & weighting_type,
693  const std::set<subdomain_id_type> * subdomain_ids) const
694 {
695  // first check if the FE type is continuous because in that case the value is
696  // unique and we can take a short cut, the default weighting_type found_first also
697  // shortcuts out
698  auto family =
700  .getVariable(
702  .feType()
703  .family;
704 
705  if (weighting_type == 1 ||
706  (family != L2_LAGRANGE && family != MONOMIAL && family != L2_HIERARCHIC))
707  return pointValue(t, p, var_name, subdomain_ids);
708 
709  // the shape function is discontinuous so we need to compute a suitable unique value
710  std::map<const Elem *, Real> values = discontinuousPointValue(t, p, var_name);
711  switch (weighting_type)
712  {
713  case 2:
714  {
715  Real average = 0.0;
716  for (auto & v : values)
717  average += v.second;
718  return average / Real(values.size());
719  }
720  case 4:
721  {
722  Real smallest_elem_id_value = std::numeric_limits<Real>::max();
723  dof_id_type smallest_elem_id = _fe_problem.mesh().maxElemId();
724  for (auto & v : values)
725  if (v.first->id() < smallest_elem_id)
726  {
727  smallest_elem_id = v.first->id();
728  smallest_elem_id_value = v.second;
729  }
730  return smallest_elem_id_value;
731  }
732  case 8:
733  {
734  Real largest_elem_id_value = std::numeric_limits<Real>::lowest();
735  dof_id_type largest_elem_id = 0;
736  for (auto & v : values)
737  if (v.first->id() > largest_elem_id)
738  {
739  largest_elem_id = v.first->id();
740  largest_elem_id_value = v.second;
741  }
742  return largest_elem_id_value;
743  }
744  }
745 
746  mooseError("SolutionUserObjectBase::pointValueWrapper reaches line that it should not be able to "
747  "reach.");
748  return 0.0;
749 }
750 
751 Real
753  const Point & p,
754  const std::string & var_name,
755  const std::set<subdomain_id_type> * subdomain_ids) const
756 {
757  const unsigned int local_var_index = getLocalVarIndex(var_name);
758  return pointValue(t, p, local_var_index, subdomain_ids);
759 }
760 
761 Real
762 SolutionUserObjectBase::pointValue(Real libmesh_dbg_var(t),
763  const Point & p,
764  const unsigned int local_var_index,
765  const std::set<subdomain_id_type> * subdomain_ids) const
766 {
767  // Create copy of point
768  Point pt(p);
769 
770  // do the transformations
771  for (unsigned int trans_num = 0; trans_num < _transformation_order.size(); ++trans_num)
772  {
773  if (_transformation_order[trans_num] == "rotation0")
774  pt = _r0 * pt;
775  else if (_transformation_order[trans_num] == "translation")
776  for (const auto i : make_range(Moose::dim))
777  pt(i) -= _translation[i];
778  else if (_transformation_order[trans_num] == "scale")
779  for (const auto i : make_range(Moose::dim))
780  pt(i) /= _scale[i];
781  else if (_transformation_order[trans_num] == "scale_multiplier")
782  for (const auto i : make_range(Moose::dim))
783  pt(i) *= _scale_multiplier[i];
784  else if (_transformation_order[trans_num] == "rotation1")
785  pt = _r1 * pt;
786  }
787 
788  // Extract the value at the current point
789  Real val = evalMeshFunction(pt, local_var_index, 1, subdomain_ids);
790 
791  // Interpolate
792  if (_file_type == 1 && _interpolate_times)
793  {
794  mooseAssert(t == _interpolation_time,
795  "Time passed into value() must match time at last call to timestepSetup()");
796  Real val2 = evalMeshFunction(pt, local_var_index, 2, subdomain_ids);
797  val = val + (val2 - val) * _interpolation_factor;
798  }
799 
800  return val;
801 }
802 
803 std::map<const Elem *, Real>
805  Real t,
806  const Point & p,
807  const std::string & var_name,
808  const std::set<subdomain_id_type> * subdomain_ids) const
809 {
810  const unsigned int local_var_index = getLocalVarIndex(var_name);
811  return discontinuousPointValue(t, p, local_var_index, subdomain_ids);
812 }
813 
814 std::map<const Elem *, Real>
816  Real libmesh_dbg_var(t),
817  Point pt,
818  const unsigned int local_var_index,
819  const std::set<subdomain_id_type> * subdomain_ids) const
820 {
821  // do the transformations
822  for (unsigned int trans_num = 0; trans_num < _transformation_order.size(); ++trans_num)
823  {
824  if (_transformation_order[trans_num] == "rotation0")
825  pt = _r0 * pt;
826  else if (_transformation_order[trans_num] == "translation")
827  for (const auto i : make_range(Moose::dim))
828  pt(i) -= _translation[i];
829  else if (_transformation_order[trans_num] == "scale")
830  for (const auto i : make_range(Moose::dim))
831  pt(i) /= _scale[i];
832  else if (_transformation_order[trans_num] == "scale_multiplier")
833  for (const auto i : make_range(Moose::dim))
834  pt(i) *= _scale_multiplier[i];
835  else if (_transformation_order[trans_num] == "rotation1")
836  pt = _r1 * pt;
837  }
838 
839  // Extract the value at the current point
840  std::map<const Elem *, Real> map =
841  evalMultiValuedMeshFunction(pt, local_var_index, 1, subdomain_ids);
842 
843  // Interpolate
844  if (_file_type == 1 && _interpolate_times)
845  {
846  mooseAssert(t == _interpolation_time,
847  "Time passed into value() must match time at last call to timestepSetup()");
848  std::map<const Elem *, Real> map2 = evalMultiValuedMeshFunction(pt, local_var_index, 2);
849 
850  if (map.size() != map2.size())
851  mooseError(
852  "In SolutionUserObjectBase::discontinuousPointValue map and map2 have different size");
853 
854  // construct the interpolated map
855  for (auto & k : map)
856  {
857  if (map2.find(k.first) == map2.end())
858  mooseError(
859  "In SolutionUserObjectBase::discontinuousPointValue map and map2 have differing keys");
860  Real val = k.second;
861  Real val2 = map2[k.first];
862  map[k.first] = val + (val2 - val) * _interpolation_factor;
863  }
864  }
865 
866  return map;
867 }
868 
871  Real t,
872  const Point & p,
873  const std::string & var_name,
874  const MooseEnum & weighting_type,
875  const std::set<subdomain_id_type> * subdomain_ids) const
876 {
877  // the default weighting_type found_first shortcuts out
878  if (weighting_type == 1)
879  return pointValueGradient(t, p, var_name, subdomain_ids);
880 
881  // the shape function is discontinuous so we need to compute a suitable unique value
882  std::map<const Elem *, RealGradient> values =
883  discontinuousPointValueGradient(t, p, var_name, subdomain_ids);
884  switch (weighting_type)
885  {
886  case 2:
887  {
888  RealGradient average = RealGradient(0.0, 0.0, 0.0);
889  for (auto & v : values)
890  average += v.second;
891  return average / Real(values.size());
892  }
893  case 4:
894  {
895  RealGradient smallest_elem_id_value;
896  dof_id_type smallest_elem_id = _fe_problem.mesh().maxElemId();
897  for (auto & v : values)
898  if (v.first->id() < smallest_elem_id)
899  {
900  smallest_elem_id = v.first->id();
901  smallest_elem_id_value = v.second;
902  }
903  return smallest_elem_id_value;
904  }
905  case 8:
906  {
907  RealGradient largest_elem_id_value;
908  dof_id_type largest_elem_id = 0;
909  for (auto & v : values)
910  if (v.first->id() > largest_elem_id)
911  {
912  largest_elem_id = v.first->id();
913  largest_elem_id_value = v.second;
914  }
915  return largest_elem_id_value;
916  }
917  }
918 
919  mooseError("SolutionUserObjectBase::pointValueGradientWrapper reaches line that it should not be "
920  "able to reach.");
921  return RealGradient(0.0, 0.0, 0.0);
922 }
923 
926  const Point & p,
927  const std::string & var_name,
928  const std::set<subdomain_id_type> * subdomain_ids) const
929 {
930  const unsigned int local_var_index = getLocalVarIndex(var_name);
931  return pointValueGradient(t, p, local_var_index, subdomain_ids);
932 }
933 
935 SolutionUserObjectBase::pointValueGradient(Real libmesh_dbg_var(t),
936  Point pt,
937  const unsigned int local_var_index,
938  const std::set<subdomain_id_type> * subdomain_ids) const
939 {
940  // do the transformations
941  for (unsigned int trans_num = 0; trans_num < _transformation_order.size(); ++trans_num)
942  {
943  if (_transformation_order[trans_num] == "rotation0")
944  pt = _r0 * pt;
945  else if (_transformation_order[trans_num] == "translation")
946  for (const auto i : make_range(Moose::dim))
947  pt(i) -= _translation[i];
948  else if (_transformation_order[trans_num] == "scale")
949  for (const auto i : make_range(Moose::dim))
950  pt(i) /= _scale[i];
951  else if (_transformation_order[trans_num] == "scale_multiplier")
952  for (const auto i : make_range(Moose::dim))
953  pt(i) *= _scale_multiplier[i];
954  else if (_transformation_order[trans_num] == "rotation1")
955  pt = _r1 * pt;
956  }
957 
958  // Extract the value at the current point
959  RealGradient val = evalMeshFunctionGradient(pt, local_var_index, 1, subdomain_ids);
960 
961  // Interpolate
962  if (_file_type == 1 && _interpolate_times)
963  {
964  mooseAssert(t == _interpolation_time,
965  "Time passed into value() must match time at last call to timestepSetup()");
966  RealGradient val2 = evalMeshFunctionGradient(pt, local_var_index, 2, subdomain_ids);
967  val = val + (val2 - val) * _interpolation_factor;
968  }
969 
970  return val;
971 }
972 
973 std::map<const Elem *, RealGradient>
975  Real t,
976  const Point & p,
977  const std::string & var_name,
978  const std::set<subdomain_id_type> * subdomain_ids) const
979 {
980  const unsigned int local_var_index = getLocalVarIndex(var_name);
981  return discontinuousPointValueGradient(t, p, local_var_index, subdomain_ids);
982 }
983 
984 std::map<const Elem *, RealGradient>
986  Real libmesh_dbg_var(t),
987  Point pt,
988  const unsigned int local_var_index,
989  const std::set<subdomain_id_type> * subdomain_ids) const
990 {
991  // do the transformations
992  for (unsigned int trans_num = 0; trans_num < _transformation_order.size(); ++trans_num)
993  {
994  if (_transformation_order[trans_num] == "rotation0")
995  pt = _r0 * pt;
996  else if (_transformation_order[trans_num] == "translation")
997  for (const auto i : make_range(Moose::dim))
998  pt(i) -= _translation[i];
999  else if (_transformation_order[trans_num] == "scale")
1000  for (const auto i : make_range(Moose::dim))
1001  pt(i) /= _scale[i];
1002  else if (_transformation_order[trans_num] == "scale_multiplier")
1003  for (const auto i : make_range(Moose::dim))
1004  pt(i) *= _scale_multiplier[i];
1005  else if (_transformation_order[trans_num] == "rotation1")
1006  pt = _r1 * pt;
1007  }
1008 
1009  // Extract the value at the current point
1010  std::map<const Elem *, RealGradient> map =
1011  evalMultiValuedMeshFunctionGradient(pt, local_var_index, 1, subdomain_ids);
1012 
1013  // Interpolate
1014  if (_file_type == 1 && _interpolate_times)
1015  {
1016  mooseAssert(t == _interpolation_time,
1017  "Time passed into value() must match time at last call to timestepSetup()");
1018  std::map<const Elem *, RealGradient> map2 =
1019  evalMultiValuedMeshFunctionGradient(pt, local_var_index, 1, subdomain_ids);
1020 
1021  if (map.size() != map2.size())
1022  mooseError(
1023  "In SolutionUserObjectBase::discontinuousPointValue map and map2 have different size");
1024 
1025  // construct the interpolated map
1026  for (auto & k : map)
1027  {
1028  if (map2.find(k.first) == map2.end())
1029  mooseError(
1030  "In SolutionUserObjectBase::discontinuousPointValue map and map2 have differing keys");
1031  RealGradient val = k.second;
1032  RealGradient val2 = map2[k.first];
1033  map[k.first] = val + (val2 - val) * _interpolation_factor;
1034  }
1035  }
1036 
1037  return map;
1038 }
1039 
1040 Real
1042 {
1043  Real val = (*_serialized_solution)(dof_index);
1044  if (_file_type == 1 && _interpolate_times)
1045  {
1046  Real val2 = (*_serialized_solution2)(dof_index);
1047  val = val + (val2 - val) * _interpolation_factor;
1048  }
1049  return val;
1050 }
1051 
1052 Real
1054  const unsigned int local_var_index,
1055  unsigned int func_num,
1056  const std::set<subdomain_id_type> * subdomain_ids) const
1057 {
1058  // Storage for mesh function output
1059  DenseVector<Number> output;
1060 
1061  // Extract a value from the _mesh_function
1062  {
1063  Threads::spin_mutex::scoped_lock lock(_solution_user_object_mutex);
1064  if (func_num == 1)
1065  (*_mesh_function)(p, 0.0, output, subdomain_ids);
1066 
1067  // Extract a value from _mesh_function2
1068  else if (func_num == 2)
1069  (*_mesh_function2)(p, 0.0, output, subdomain_ids);
1070 
1071  else
1072  mooseError("The func_num must be 1 or 2");
1073  }
1074 
1075  // Error if the data is out-of-range, which will be the case if the mesh functions are evaluated
1076  // outside the domain
1077  if (output.size() == 0)
1078  {
1079  std::ostringstream oss;
1080  p.print(oss);
1081  mooseError("Failed to access the data for variable '",
1082  _system_variables[local_var_index],
1083  "' at point ",
1084  oss.str(),
1085  " in the '",
1086  name(),
1087  "' SolutionUserObjectBase");
1088  }
1089  return output(local_var_index);
1090 }
1091 
1092 std::map<const Elem *, Real>
1094  const Point & p,
1095  const unsigned int local_var_index,
1096  unsigned int func_num,
1097  const std::set<subdomain_id_type> * subdomain_ids) const
1098 {
1099  // Storage for mesh function output
1100  std::map<const Elem *, DenseVector<Number>> temporary_output;
1101 
1102  // Extract a value from the _mesh_function
1103  {
1104  Threads::spin_mutex::scoped_lock lock(_solution_user_object_mutex);
1105  if (func_num == 1)
1106  _mesh_function->discontinuous_value(p, 0.0, temporary_output, subdomain_ids);
1107 
1108  // Extract a value from _mesh_function2
1109  else if (func_num == 2)
1110  _mesh_function2->discontinuous_value(p, 0.0, temporary_output, subdomain_ids);
1111 
1112  else
1113  mooseError("The func_num must be 1 or 2");
1114  }
1115 
1116  // Error if the data is out-of-range, which will be the case if the mesh functions are evaluated
1117  // outside the domain
1118  if (temporary_output.size() == 0)
1119  {
1120  std::ostringstream oss;
1121  p.print(oss);
1122  mooseError("Failed to access the data for variable '",
1123  _system_variables[local_var_index],
1124  "' at point ",
1125  oss.str(),
1126  " in the '",
1127  name(),
1128  "' SolutionUserObjectBase");
1129  }
1130 
1131  // Fill the actual map that is returned
1132  std::map<const Elem *, Real> output;
1133  for (auto & k : temporary_output)
1134  {
1135  mooseAssert(
1136  k.second.size() > local_var_index,
1137  "In SolutionUserObjectBase::evalMultiValuedMeshFunction variable with local_var_index "
1138  << local_var_index << " does not exist");
1139  output[k.first] = k.second(local_var_index);
1140  }
1141 
1142  return output;
1143 }
1144 
1147  const Point & p,
1148  const unsigned int local_var_index,
1149  unsigned int func_num,
1150  const std::set<subdomain_id_type> * subdomain_ids) const
1151 {
1152  // Storage for mesh function output
1153  std::vector<Gradient> output;
1154 
1155  // Extract a value from the _mesh_function
1156  {
1157  Threads::spin_mutex::scoped_lock lock(_solution_user_object_mutex);
1158  if (func_num == 1)
1159  _mesh_function->gradient(p, 0.0, output, subdomain_ids);
1160 
1161  // Extract a value from _mesh_function2
1162  else if (func_num == 2)
1163  _mesh_function2->gradient(p, 0.0, output, subdomain_ids);
1164 
1165  else
1166  mooseError("The func_num must be 1 or 2");
1167  }
1168 
1169  // Error if the data is out-of-range, which will be the case if the mesh functions are evaluated
1170  // outside the domain
1171  if (output.size() == 0)
1172  {
1173  std::ostringstream oss;
1174  p.print(oss);
1175  mooseError("Failed to access the data for variable '",
1176  _system_variables[local_var_index],
1177  "' at point ",
1178  oss.str(),
1179  " in the '",
1180  name(),
1181  "' SolutionUserObjectBase");
1182  }
1183  return output[local_var_index];
1184 }
1185 
1186 std::map<const Elem *, RealGradient>
1188  const Point & p,
1189  const unsigned int local_var_index,
1190  unsigned int func_num,
1191  const std::set<subdomain_id_type> * subdomain_ids) const
1192 {
1193  // Storage for mesh function output
1194  std::map<const Elem *, std::vector<Gradient>> temporary_output;
1195 
1196  // Extract a value from the _mesh_function
1197  {
1198  Threads::spin_mutex::scoped_lock lock(_solution_user_object_mutex);
1199  if (func_num == 1)
1200  _mesh_function->discontinuous_gradient(p, 0.0, temporary_output, subdomain_ids);
1201 
1202  // Extract a value from _mesh_function2
1203  else if (func_num == 2)
1204  _mesh_function2->discontinuous_gradient(p, 0.0, temporary_output, subdomain_ids);
1205 
1206  else
1207  mooseError("The func_num must be 1 or 2");
1208  }
1209 
1210  // Error if the data is out-of-range, which will be the case if the mesh functions are evaluated
1211  // outside the domain
1212  if (temporary_output.size() == 0)
1213  {
1214  std::ostringstream oss;
1215  p.print(oss);
1216  mooseError("Failed to access the data for variable '",
1217  _system_variables[local_var_index],
1218  "' at point ",
1219  oss.str(),
1220  " in the '",
1221  name(),
1222  "' SolutionUserObjectBase");
1223  }
1224 
1225  // Fill the actual map that is returned
1226  std::map<const Elem *, RealGradient> output;
1227  for (auto & k : temporary_output)
1228  {
1229  mooseAssert(
1230  k.second.size() > local_var_index,
1231  "In SolutionUserObjectBase::evalMultiValuedMeshFunction variable with local_var_index "
1232  << local_var_index << " does not exist");
1233  output[k.first] = k.second[local_var_index];
1234  }
1235 
1236  return output;
1237 }
1238 
1239 const std::vector<std::string> &
1241 {
1242  return _system_variables;
1243 }
1244 
1245 bool
1246 SolutionUserObjectBase::isVariableNodal(const std::string & var_name) const
1247 {
1248  return std::find(_nodal_variables.begin(), _nodal_variables.end(), var_name) !=
1249  _nodal_variables.end();
1250 }
1251 
1252 Real
1253 SolutionUserObjectBase::scalarValue(Real /*t*/, const std::string & var_name) const
1254 {
1255  unsigned int var_num = _system->variable_number(var_name);
1256  const DofMap & dof_map = _system->get_dof_map();
1257  std::vector<dof_id_type> dofs;
1258  dof_map.SCALAR_dof_indices(dofs, var_num);
1259  // We can handle only FIRST order scalar variables
1260  return directValue(dofs[0]);
1261 }
1262 
1263 void
1265 {
1266 #ifdef LIBMESH_HAVE_EXODUS_API
1267  libMesh::ExodusII_IO_Helper & exio_helper = _exodusII_io->get_exio_helper();
1268  const auto & id_to_block = exio_helper.id_to_block_names;
1269  _block_name_to_id.clear();
1270  _block_id_to_name.clear();
1271  for (const auto & it : id_to_block)
1272  {
1273  _block_name_to_id[it.second] = it.first;
1274  _block_id_to_name[it.first] = it.second;
1275  }
1276 #endif
1277 }
L2_HIERARCHIC
MultiMooseEnum _transformation_order
transformations (rotations, translation, scales) are performed in this order
std::vector< std::string > _system_variables
A list of variables to extract from the read system.
static InputParameters validParams()
SolutionUserObjectBase(const InputParameters &parameters)
dof_id_type dof_number(const unsigned int s, const unsigned int var, const unsigned int comp) const
std::unique_ptr< libMesh::MeshFunction > _mesh_function2
Pointer to second libMesh::MeshFuntion, used for interpolation.
const libMesh::FEType & feType() const
Get the type of finite element object.
std::unique_ptr< libMesh::EquationSystems > _es
Pointer to the libMesh::EquationSystems object.
CTSub CT_OPERATOR_BINARY CTMul CTCompareLess CTCompareGreater CTCompareEqual _arg template * sin(_arg) *_arg.template D< dtag >()) CT_SIMPLE_UNARY_FUNCTION(tan
RealVectorValue _rotation0_vector
vector about which to rotate
const std::vector< Real > * _exodus_times
The times available in the ExodusII file.
virtual void initialSetup() override
Initialize the System and Mesh objects for the solution being read.
virtual Real solutionSampleTime()=0
Get the time at which to sample the solution.
SCALAR
static InputParameters validParams()
void readXda()
Method for reading XDA mesh and equation systems file(s) This method is called by the constructor whe...
static Threads::spin_mutex _solution_user_object_mutex
std::map< const Elem *, libMesh::RealGradient > discontinuousPointValueGradient(Real t, const Point &p, const std::string &var_name, const std::set< subdomain_id_type > *subdomain_ids=nullptr) const
Returns the gradient at a specific location and variable for cases where the gradient is multivalued ...
void mooseInfo(Args &&... args) const
FIRST
TypeTensor< T > transpose() const
std::unique_ptr< NumericVector< Number > > _serialized_solution2
Pointer to second serial solution, used for interpolation.
std::vector< Real > _scale
Scale parameter.
std::vector< Real > _translation
Translation.
Real _rotation1_angle
angle (in degrees) which to rotate through about vector _rotation1_vector
unsigned int size() const
Return the number of active items in the MultiMooseEnum.
The main MOOSE class responsible for handling user-defined parameters in almost every MOOSE system...
virtual dof_id_type maxElemId() const
Definition: MooseMesh.C:3088
Real _rotation0_angle
angle (in degrees) which to rotate through about vector _rotation0_vector
static constexpr std::size_t dim
This is the dimension of all vector and tensor datastructures used in MOOSE.
Definition: Moose.h:154
MooseEnum getSolutionFileType() const
Get the type of file that was read.
void get_all_variable_numbers(std::vector< unsigned int > &all_variable_numbers) const
int _exodus_time_index
Current ExodusII time index.
std::unique_ptr< libMesh::ExodusII_IO > _exodusII_io
Pointer to the libMesh::ExodusII used to read the files.
Real _interpolation_time
Time in the current simulation at which the solution interpolation was last updated.
std::unique_ptr< libMesh::EquationSystems > _es2
Pointer to second libMesh::EquationSystems object, used for interpolation.
const Parallel::Communicator & _communicator
std::map< const Elem *, Real > evalMultiValuedMeshFunction(const Point &p, const unsigned int local_var_index, unsigned int func_num, const std::set< subdomain_id_type > *subdomain_ids=nullptr) const
A wrapper method for calling the various MeshFunctions that calls the mesh function functionality for...
libMesh::RealGradient evalMeshFunctionGradient(const Point &p, const unsigned int local_var_index, unsigned int func_num, const std::set< subdomain_id_type > *subdomain_ids=nullptr) const
A wrapper method interfacing with the libMesh mesh function for evaluating the gradient.
std::unique_ptr< NumericVector< Number > > _serialized_solution
Pointer to the serial solution vector.
bool _initialized
True if initial_setup has executed.
const MeshBase & get_mesh() const
virtual const std::string & name() const
Get the name of the class.
Definition: MooseBase.h:57
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...
virtual const MooseVariableFieldBase & getVariable(const THREAD_ID tid, const std::string &var_name, Moose::VarKindType expected_var_type=Moose::VarKindType::VAR_ANY, Moose::VarFieldType expected_var_field_type=Moose::VarFieldType::VAR_FIELD_ANY) const override
Returns the variable reference for requested variable which must be of the expected_var_type (Nonline...
dof_id_type n_dofs() const
std::map< SubdomainName, SubdomainID > _block_name_to_id
Map from block ID to block names. Read from the ExodusII file.
void SCALAR_dof_indices(std::vector< dof_id_type > &di, const unsigned int vn, const bool old_dofs=false) const
auto max(const L &left, const R &right)
unsigned int variable_number(std::string_view var) const
bool isParamValid(const std::string &name) const
Test if the supplied parameter is valid.
virtual void timestepSetup() override
When reading ExodusII files, this will update the interpolation times.
virtual void execute() override
Execute method.
std::unique_ptr< libMesh::MeshBase > _mesh
Pointer the libMesh::mesh object.
libMesh::System * _system
Pointer libMesh::System class storing the read solution.
bool updateExodusBracketingTimeIndices()
Updates the time indices to interpolate between for ExodusII data.
CONSTANT
unsigned int number() const
RealVectorValue _rotation1_vector
vector about which to rotate
CTSub CT_OPERATOR_BINARY CTMul CTCompareLess CTCompareGreater CTCompareEqual _arg template cos(_arg) *_arg.template D< dtag >()) CT_SIMPLE_UNARY_FUNCTION(cos
std::string formatString(std::string message, const std::string &prefix)
Add new lines and prefixes to a string for pretty display in output NOTE: This makes a copy of the st...
Definition: ConsoleUtils.C:582
bool checkFileReadable(const std::string &filename, bool check_line_endings=false, bool throw_on_unreadable=true, bool check_for_git_lfs_pointer=true)
Checks to see if a file is readable (exists and permissions)
Definition: MooseUtils.C:250
std::string _mesh_file
The XDA or ExodusII file that is being read.
std::unique_ptr< NumericVector< Number > > solution
std::vector< std::string > _scalar_variables
Stores names of scalar variables.
This is a "smart" enum class intended to replace many of the shortcomings in the C++ enum type It sho...
Definition: MooseEnum.h:33
RealTensorValue _r1
Rotation matrix that performs the "_rotation1_angle about rotation1_vector".
bool hasExtension(const std::string &filename, std::string ext, bool strip_exodus_ext=false)
Function tests if the supplied filename as the desired extension.
Definition: MooseUtils.C:383
unsigned int add_variable(std::string_view var, const FEType &type, const std::set< subdomain_id_type > *const active_subdomains=nullptr)
const std::string & variable_name(const unsigned int i) const
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 ...
void readExodusII()
Method for reading an ExodusII file, which is called when &#39;file_type = exodusII is set in the input f...
std::string stringify(const T &t)
conversion to string
Definition: Conversion.h:64
bool hasSecondOrderElements()
check if the mesh has SECOND order elements
Definition: MooseMesh.C:3706
virtual void initialize() override
Called before execute() is ever called so that data can be cleared.
Real evalMeshFunction(const Point &p, const unsigned int local_var_index, unsigned int func_num, const std::set< subdomain_id_type > *subdomain_ids=nullptr) const
A wrapper method for calling the various MeshFunctions used for reading the data. ...
MONOMIAL
void updateExodusTimeInterpolation()
Updates the times for interpolating ExodusII data.
Real scalarValue(Real t, const std::string &var_name) const
Returns a value of a global variable.
const MooseEnum _nodal_variable_order
Nodal variable order, used when reading in solution data.
libMesh::System * _system2
Pointer to a second libMesh::System object, used for interpolation.
GenericRealTensorValue< is_ad > rotVecToZ(GenericRealVectorValue< is_ad > vec)
provides a rotation matrix that will rotate the vector vec to the z axis (the "2" direction) ...
bool isParamSetByUser(const std::string &nm) const
Test if the supplied parameter is set by a user, as opposed to not set or set to default.
virtual void update()
Real directValue(const Node *node, const std::string &var_name) const
Return a value directly from a Node.
void readBlockIdMapFromExodusII()
Read block ID map from the ExodusII file.
const std::vector< std::string > & variableNames() const
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
FEProblemBase & _fe_problem
Reference to the FEProblemBase for this user object.
Definition: UserObject.h:211
bool isVariableNodal(const std::string &var_name) const
std::map< std::string, unsigned int > _local_variable_index
Stores the local index need by MeshFunction.
std::map< const Elem *, libMesh::RealGradient > evalMultiValuedMeshFunctionGradient(const Point &p, const unsigned int local_var_index, unsigned int func_num, const std::set< subdomain_id_type > *subdomain_ids=nullptr) const
A wrapper method interfacing with the libMesh mesh function that calls the gradient functionality for...
Real pointValueWrapper(Real t, const Point &p, const std::string &var_name, const MooseEnum &weighting_type=weightingType(), const std::set< subdomain_id_type > *subdomain_ids=nullptr) const
Returns a value at a specific location and variable checking for multiple values and weighting these ...
int _exodus_index1
Time index 1, used for interpolation.
virtual const Elem & elem_ref(const dof_id_type i) const
std::vector< std::string > _nodal_variables
Stores names of nodal variables.
std::map< int, std::string > id_to_block_names
std::string _es_file
The XDA file that contians the EquationSystems data (xda only)
const THREAD_ID _tid
Thread ID of this postprocessor.
Definition: UserObject.h:218
IntRange< T > make_range(T beg, T end)
virtual MooseMesh & mesh() override
static std::unique_ptr< NumericVector< Number > > build(const Parallel::Communicator &comm, SolverPackage solver_package=libMesh::default_solver_package(), ParallelType parallel_type=AUTOMATIC)
RealTensorValue _r0
Rotation matrix that performs the "_rotation0_angle about rotation0_vector".
void mooseError(Args &&... args) const
Emits an error prefixed with object name and type.
virtual unsigned int size() const override final
L2_LAGRANGE
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...
std::vector< std::string > _elemental_variables
Stores names of elemental variables.
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...
virtual const Node & node_ref(const dof_id_type i) const
This is a "smart" enum class intended to replace many of the shortcomings in the C++ enum type...
MooseEnum _file_type
File type to read (0 = xda; 1 = ExodusII)
libMesh::RealGradient pointValueGradientWrapper(Real t, const Point &p, const std::string &var_name, const MooseEnum &weighting_type=weightingType(), const std::set< subdomain_id_type > *subdomain_ids=nullptr) const
Returns the gradient at a specific location and variable checking for multiple values and weighting t...
std::vector< Real > _scale_multiplier
scale_multiplier parameter
std::unique_ptr< libMesh::MeshFunction > _mesh_function
Pointer the libMesh::MeshFunction object that the read data is stored.
bool _interpolate_times
Flag for triggering interpolation of ExodusII data.
std::string _system_name
The system name to extract from the XDA file (xda only)
const DofMap & get_dof_map() const
std::map< const Elem *, Real > discontinuousPointValue(Real t, Point pt, const unsigned int local_var_index, const std::set< subdomain_id_type > *subdomain_ids=nullptr) const
Returns a value at a specific location and variable for cases where the solution is multivalued at el...
Real pointValue(Real t, const Point &p, const unsigned int local_var_index, const std::set< subdomain_id_type > *subdomain_ids=nullptr) const
Returns a value at a specific location and variable (see SolutionFunction)
unsigned int getLocalVarIndex(const std::string &var_name) const
Returns the local index for a given variable name.
int _exodus_index2
Time index 2, used for interpolation.
Real _interpolation_factor
Interpolation weight factor.
std::map< SubdomainID, SubdomainName > _block_id_to_name
Map from block names to block IDs. Read from the ExodusII file.
virtual void finalize() override
Finalize.
uint8_t dof_id_type
const Real pi
void addParamNamesToGroup(const std::string &space_delim_names, const std::string group_name)
This method takes a space delimited list of parameter names and adds them to the specified group name...
libMesh::RealGradient pointValueGradient(Real t, const Point &p, const std::string &var_name, const std::set< subdomain_id_type > *subdomain_ids=nullptr) const
Returns the gradient at a specific location and variable (see SolutionFunction)