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