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" 42 "mesh",
"The name of the mesh file (must be xda/xdr or exodusII file).");
43 params.
addParam<std::vector<std::string>>(
45 std::vector<std::string>(),
46 "The name of the nodal and elemental variables from the file you want to use for values");
52 "The name of the file holding the equation system info in xda/xdr format (xda/xdr only).");
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'");
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.");
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");
76 "Vector about which to rotate points of the simulation.");
80 "Anticlockwise rotation angle (in degrees) to use for rotation about rotation0_vector.");
83 "Vector about which to rotate points of the simulation.");
87 "Anticlockwise rotation angle (in degrees) to use for rotation about rotation1_vector.");
91 "rotation0 translation scale rotation1 scale_multiplier",
"translation scale");
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");
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),
122 _interpolation_time(0.0),
123 _interpolation_factor(0.0),
124 _exodus_times(nullptr),
127 _scale(getParam<
std::vector<
Real>>(
"scale")),
128 _scale_multiplier(getParam<
std::vector<
Real>>(
"scale_multiplier")),
129 _translation(getParam<
std::vector<
Real>>(
"translation")),
131 _rotation0_angle(getParam<
Real>(
"rotation0_angle")),
134 _rotation1_angle(getParam<
Real>(
"rotation1_angle")),
136 _transformation_order(getParam<
MultiMooseEnum>(
"transformation_order")),
140 Real halfPi = std::acos(0.0);
152 _r0 = vec0_to_z.transpose() * (rot0_z * vec0_to_z);
162 _r1 = vec1_to_z.transpose() * (rot1_z * vec1_to_z);
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");
175 paramError(
"es",
"Equation system file (.xda or .xdr) should have been specified");
185 _es = std::make_unique<EquationSystems>(*_mesh);
191 EquationSystems::READ_HEADER | EquationSystems::READ_DATA |
192 EquationSystems::READ_ADDITIONAL_DATA);
198 EquationSystems::READ_HEADER | EquationSystems::READ_DATA |
199 EquationSystems::READ_ADDITIONAL_DATA);
203 mooseError(
"Failed to determine proper read method for XDA/XDR equation system file: ",
226 std::string s_timestep = getParam<std::string>(
"timestep");
228 if (s_timestep ==
"LATEST")
232 std::istringstream ss(s_timestep);
234 mooseError(
"Invalid value passed as \"timestep\". Expected \"LATEST\" or a valid integer " 247 if (num_exo_times == 0)
248 mooseError(
"In SolutionUserObject, exodus file contains no timesteps.");
251 if (dynamic_cast<DistributedMesh *>(
_mesh.get()))
253 _mesh->allow_renumbering(
true);
254 _mesh->prepare_for_use();
258 _mesh->allow_renumbering(
false);
259 _mesh->prepare_for_use();
263 _es = std::make_unique<EquationSystems>(*_mesh);
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());
283 if (std::find(all_nodal.begin(), all_nodal.end(), var_name) != all_nodal.end())
285 if (std::find(all_elemental.begin(), all_elemental.end(), var_name) != all_elemental.end())
287 if (std::find(all_scalar.begin(), all_scalar.end(), var_name) != all_scalar.end())
302 for (
auto var_name : all_scalar)
329 _es2 = std::make_unique<EquationSystems>(*_mesh);
381 mooseError(
"In SolutionUserObject, timestep = ",
383 ", but there are only ",
408 unsigned int var_num =
_system->variable_number(var_name);
409 unsigned int sys_num =
_system->number();
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);
426 unsigned int var_num =
_system->variable_number(var_name);
427 unsigned int sys_num =
_system->number();
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);
503 mooseError(
"In SolutionUserObject, invalid file type: only .xda, .xdr, .exo and .e supported");
513 std::vector<unsigned int> var_nums;
518 _system->get_all_variable_numbers(var_nums);
519 for (
const auto & var_num : var_nums)
527 var_nums.push_back(
_system->variable_number(var_name));
615 "In SolutionUserObject, getTimeInterpolationData only applicable for exodusII file type");
630 for (
int i = 0; i < num_exo_times - 1; ++i)
637 (time - (*_exodus_times)[i]) / ((*
_exodus_times)[i + 1] - (*_exodus_times)[i]);
640 else if (i == num_exo_times - 2)
650 bool indices_modified(
false);
653 indices_modified =
true;
655 return indices_modified;
664 mooseError(
"Value requested for nonexistent variable '",
668 "' SolutionUserObject.\nSystem selected: ",
670 "\nAvailable variables:\n",
678 const std::string & var_name,
680 const std::set<subdomain_id_type> * subdomain_ids)
const 692 if (weighting_type == 1 ||
694 return pointValue(t, p, var_name, subdomain_ids);
698 switch (weighting_type)
703 for (
auto & v : values)
705 return average /
Real(values.size());
711 for (
auto & v : values)
712 if (v.first->id() < smallest_elem_id)
714 smallest_elem_id = v.first->id();
715 smallest_elem_id_value = v.second;
717 return smallest_elem_id_value;
721 Real largest_elem_id_value = std::numeric_limits<Real>::lowest();
723 for (
auto & v : values)
724 if (v.first->id() > largest_elem_id)
726 largest_elem_id = v.first->id();
727 largest_elem_id_value = v.second;
729 return largest_elem_id_value;
734 "SolutionUserObject::pointValueWrapper reaches line that it should not be able to reach.");
741 const std::string & var_name,
742 const std::set<subdomain_id_type> * subdomain_ids)
const 745 return pointValue(t, p, local_var_index, subdomain_ids);
751 const unsigned int local_var_index,
752 const std::set<subdomain_id_type> * subdomain_ids)
const 782 "Time passed into value() must match time at last call to timestepSetup()");
790 std::map<const Elem *, Real>
793 const std::string & var_name,
794 const std::set<subdomain_id_type> * subdomain_ids)
const 800 std::map<const Elem *, Real>
803 const unsigned int local_var_index,
804 const std::set<subdomain_id_type> * subdomain_ids)
const 825 std::map<const Elem *, Real> map =
832 "Time passed into value() must match time at last call to timestepSetup()");
835 if (map.size() != map2.size())
836 mooseError(
"In SolutionUserObject::discontinuousPointValue map and map2 have different size");
841 if (map2.find(k.first) == map2.end())
843 "In SolutionUserObject::discontinuousPointValue map and map2 have differing keys");
845 Real val2 = map2[k.first];
857 const std::string & var_name,
859 const std::set<subdomain_id_type> * subdomain_ids)
const 862 if (weighting_type == 1)
866 std::map<const Elem *, RealGradient> values =
868 switch (weighting_type)
873 for (
auto & v : values)
875 return average /
Real(values.size());
881 for (
auto & v : values)
882 if (v.first->id() < smallest_elem_id)
884 smallest_elem_id = v.first->id();
885 smallest_elem_id_value = v.second;
887 return smallest_elem_id_value;
893 for (
auto & v : values)
894 if (v.first->id() > largest_elem_id)
896 largest_elem_id = v.first->id();
897 largest_elem_id_value = v.second;
899 return largest_elem_id_value;
903 mooseError(
"SolutionUserObject::pointValueGradientWrapper reaches line that it should not be " 911 const std::string & var_name,
912 const std::set<subdomain_id_type> * subdomain_ids)
const 921 const unsigned int local_var_index,
922 const std::set<subdomain_id_type> * subdomain_ids)
const 949 "Time passed into value() must match time at last call to timestepSetup()");
957 std::map<const Elem *, RealGradient>
961 const std::string & var_name,
962 const std::set<subdomain_id_type> * subdomain_ids)
const 968 std::map<const Elem *, RealGradient>
970 Real libmesh_dbg_var(t),
972 const unsigned int local_var_index,
973 const std::set<subdomain_id_type> * subdomain_ids)
const 994 std::map<const Elem *, RealGradient> map =
1001 "Time passed into value() must match time at last call to timestepSetup()");
1002 std::map<const Elem *, RealGradient> map2 =
1005 if (map.size() != map2.size())
1006 mooseError(
"In SolutionUserObject::discontinuousPointValue map and map2 have different size");
1009 for (
auto & k : map)
1011 if (map2.find(k.first) == map2.end())
1013 "In SolutionUserObject::discontinuousPointValue map and map2 have differing keys");
1026 Real val = (*_serialized_solution)(dof_index);
1029 Real val2 = (*_serialized_solution2)(dof_index);
1037 const unsigned int local_var_index,
1038 unsigned int func_num,
1039 const std::set<subdomain_id_type> * subdomain_ids)
const 1048 (*_mesh_function)(p, 0.0, output, subdomain_ids);
1051 else if (func_num == 2)
1052 (*_mesh_function2)(p, 0.0, output, subdomain_ids);
1060 if (output.
size() == 0)
1062 std::ostringstream oss;
1064 mooseError(
"Failed to access the data for variable '",
1070 "' SolutionUserObject");
1072 return output(local_var_index);
1075 std::map<const Elem *, Real>
1078 const unsigned int local_var_index,
1079 unsigned int func_num,
1080 const std::set<subdomain_id_type> * subdomain_ids)
const 1083 std::map<const Elem *, DenseVector<Number>> temporary_output;
1089 _mesh_function->discontinuous_value(p, 0.0, temporary_output, subdomain_ids);
1092 else if (func_num == 2)
1093 _mesh_function2->discontinuous_value(p, 0.0, temporary_output, subdomain_ids);
1101 if (temporary_output.size() == 0)
1103 std::ostringstream oss;
1105 mooseError(
"Failed to access the data for variable '",
1111 "' SolutionUserObject");
1115 std::map<const Elem *, Real> output;
1116 for (
auto & k : temporary_output)
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);
1130 const unsigned int local_var_index,
1131 unsigned int func_num,
1132 const std::set<subdomain_id_type> * subdomain_ids)
const 1135 std::vector<Gradient> output;
1144 else if (func_num == 2)
1153 if (output.size() == 0)
1155 std::ostringstream oss;
1157 mooseError(
"Failed to access the data for variable '",
1163 "' SolutionUserObject");
1165 return output[local_var_index];
1168 std::map<const Elem *, RealGradient>
1171 const unsigned int local_var_index,
1172 unsigned int func_num,
1173 const std::set<subdomain_id_type> * subdomain_ids)
const 1176 std::map<const Elem *, std::vector<Gradient>> temporary_output;
1182 _mesh_function->discontinuous_gradient(p, 0.0, temporary_output, subdomain_ids);
1185 else if (func_num == 2)
1186 _mesh_function2->discontinuous_gradient(p, 0.0, temporary_output, subdomain_ids);
1194 if (temporary_output.size() == 0)
1196 std::ostringstream oss;
1198 mooseError(
"Failed to access the data for variable '",
1204 "' SolutionUserObject");
1208 std::map<const Elem *, RealGradient> output;
1209 for (
auto & k : temporary_output)
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];
1220 const std::vector<std::string> &
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);
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;
1252 for (
const auto & it : id_to_block)
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.
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...
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.
virtual dof_id_type maxElemId() const
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.
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.
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.
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...
void readExodusII()
Method for reading an ExodusII file, which is called when 'file_type = exodusII is set in the input f...
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)
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...
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)
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...
bool hasExtension(const std::string &filename, std::string ext, bool strip_exodus_ext=false)
Function tests if the supplied filename as the desired extension.
std::unique_ptr< MeshFunction > _mesh_function2
Pointer to second libMesh::MeshFuntion, used for interpolation.
void paramError(const std::string ¶m, 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
Real scalarValue(Real t, const std::string &var_name) const
Returns a value of a global variable.
SolutionUserObject(const InputParameters ¶meters)
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 ...
virtual ~SolutionUserObject()
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.
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.
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
std::vector< std::string > _nodal_variables
Stores names of nodal variables.
std::vector< std::string > _elemental_variables
Stores names of elemental variables.
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.
std::unique_ptr< MeshBase > _mesh
Pointer the libmesh::mesh object.
Real _rotation1_angle
angle (in degrees) which to rotate through about vector _rotation1_vector
virtual void initialize() override
Called before execute() is ever called so that data can be cleared.
std::vector< Real > _translation
Translation.