www.mooseframework.org
SolutionUserObject.h
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 #pragma once
11 
12 // MOOSE includes
13 #include "GeneralUserObject.h"
14 
15 // Forward declarations
16 namespace libMesh
17 {
18 class ExodusII_IO;
19 class EquationSystems;
20 class System;
21 class MeshFunction;
22 template <class T>
23 class NumericVector;
24 }
25 
31 {
32 public:
34 
36  virtual ~SolutionUserObject(); // empty dtor required for unique_ptr with forward declarations
37 
41  virtual void timestepSetup() override;
42 
48  unsigned int getLocalVarIndex(const std::string & var_name) const;
49 
62  const Point & p,
63  const std::string & var_name,
64  const MooseEnum & weighting_type = weightingType(),
65  const std::set<subdomain_id_type> * subdomain_ids = nullptr) const;
66 
77  const Point & p,
78  const unsigned int local_var_index,
79  const std::set<subdomain_id_type> * subdomain_ids = nullptr) const;
80 
91  const Point & p,
92  const std::string & var_name,
93  const std::set<subdomain_id_type> * subdomain_ids = nullptr) const;
94 
107  std::map<const Elem *, Real>
109  Point pt,
110  const unsigned int local_var_index,
111  const std::set<subdomain_id_type> * subdomain_ids = nullptr) const;
112 
125  std::map<const Elem *, Real>
127  const Point & p,
128  const std::string & var_name,
129  const std::set<subdomain_id_type> * subdomain_ids = nullptr) const;
130 
144  const Point & p,
145  const std::string & var_name,
146  const MooseEnum & weighting_type = weightingType(),
147  const std::set<subdomain_id_type> * subdomain_ids = nullptr) const;
148 
160  const Point & p,
161  const std::string & var_name,
162  const std::set<subdomain_id_type> * subdomain_ids = nullptr) const;
163 
175  Point pt,
176  const unsigned int local_var_index,
177  const std::set<subdomain_id_type> * subdomain_ids = nullptr) const;
178 
191  std::map<const Elem *, RealGradient> discontinuousPointValueGradient(
192  Real t,
193  const Point & p,
194  const std::string & var_name,
195  const std::set<subdomain_id_type> * subdomain_ids = nullptr) const;
196 
209  std::map<const Elem *, RealGradient> discontinuousPointValueGradient(
210  Real t,
211  Point pt,
212  const unsigned int local_var_index,
213  const std::set<subdomain_id_type> * subdomain_ids = nullptr) const;
214 
221  Real directValue(const Node * node, const std::string & var_name) const;
222 
229  Real directValue(const Elem * elem, const std::string & var_name) const;
230 
238  Real scalarValue(Real t, const std::string & var_name) const;
239 
240  // Required pure virtual function (not used)
241  virtual void initialize() override;
242 
243  // Required pure virtual function (not used)
244  virtual void finalize() override;
245 
246  // Required pure virtual function (not used)
247  virtual void execute() override;
248 
250  virtual void initialSetup() override;
251 
252  const std::vector<std::string> & variableNames() const;
253 
254  bool isVariableNodal(const std::string & var_name) const;
255 
257  {
258  return MooseEnum("found_first=1 average=2 smallest_element_id=4 largest_element_id=8",
259  "found_first");
260  }
261 
266  unsigned int getMeshFileDimension() const { return _mesh->spatial_dimension(); }
267 
271  const std::string getMeshFileName() const { return _mesh_file; }
272 
278  const std::map<SubdomainName, SubdomainID> & getBlockNamesToIds() const
279  {
280  return _block_name_to_id;
281  }
282 
288  const std::map<SubdomainID, SubdomainName> & getBlockIdsToNames() const
289  {
290  return _block_id_to_name;
291  }
292 
298 
299 protected:
305  void readXda();
306 
311  void readExodusII();
312 
320  virtual Real directValue(dof_id_type dof_index) const;
321 
326  void updateExodusTimeInterpolation(Real time);
327 
332  bool updateExodusBracketingTimeIndices(Real time);
333 
341  Real evalMeshFunction(const Point & p,
342  const unsigned int local_var_index,
343  unsigned int func_num,
344  const std::set<subdomain_id_type> * subdomain_ids = nullptr) const;
345 
354  std::map<const Elem *, Real>
355  evalMultiValuedMeshFunction(const Point & p,
356  const unsigned int local_var_index,
357  unsigned int func_num,
358  const std::set<subdomain_id_type> * subdomain_ids = nullptr) const;
359 
368  evalMeshFunctionGradient(const Point & p,
369  const unsigned int local_var_index,
370  unsigned int func_num,
371  const std::set<subdomain_id_type> * subdomain_ids = nullptr) const;
372 
382  std::map<const Elem *, RealGradient> evalMultiValuedMeshFunctionGradient(
383  const Point & p,
384  const unsigned int local_var_index,
385  unsigned int func_num,
386  const std::set<subdomain_id_type> * subdomain_ids = nullptr) const;
387 
392 
395 
397  std::string _mesh_file;
398 
400  std::string _es_file;
401 
403  std::string _system_name;
404 
406  std::vector<std::string> _system_variables;
407 
409  std::map<std::string, unsigned int> _local_variable_index;
410 
412  std::vector<std::string> _nodal_variables;
413 
415  std::vector<std::string> _elemental_variables;
416 
418  std::vector<std::string> _scalar_variables;
419 
422 
425 
427  std::unique_ptr<MeshBase> _mesh;
428 
430  std::unique_ptr<EquationSystems> _es;
431 
433  System * _system;
434 
436  std::unique_ptr<MeshFunction> _mesh_function;
437 
439  std::unique_ptr<ExodusII_IO> _exodusII_io;
440 
442  std::unique_ptr<NumericVector<Number>> _serialized_solution;
443 
445  std::unique_ptr<EquationSystems> _es2;
446 
448  System * _system2;
449 
451  std::unique_ptr<MeshFunction> _mesh_function2;
452 
454  std::unique_ptr<NumericVector<Number>> _serialized_solution2;
455 
458 
461 
463  const std::vector<Real> * _exodus_times;
464 
467 
470 
472  std::vector<Real> _scale;
473 
475  std::vector<Real> _scale_multiplier;
476 
478  std::vector<Real> _translation;
479 
482 
485 
488 
491 
494 
497 
500 
503 
505  std::map<SubdomainName, SubdomainID> _block_name_to_id;
506 
508  std::map<SubdomainID, SubdomainName> _block_id_to_name;
509 
510 private:
511  static Threads::spin_mutex _solution_user_object_mutex;
512 };
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
std::map< SubdomainID, SubdomainName > _block_id_to_name
Map from block names to block IDs. Read from the ExodusII file.
static InputParameters validParams()
const std::map< SubdomainID, SubdomainName > & getBlockIdsToNames() const
Get the map from block id to block name.
std::vector< Real > _scale_multiplier
scale_multiplier parameter
std::unique_ptr< NumericVector< Number > > _serialized_solution
Pointer to the serial solution vector.
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...
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.
The main MOOSE class responsible for handling user-defined parameters in almost every MOOSE system...
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 ...
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. ...
The following methods are specializations for using the libMesh::Parallel::packed_range_* routines fo...
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.
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
int _exodus_index1
Time index 1, used for interpolation.
Real directValue(const Node *node, const std::string &var_name) const
Return a value directly from a Node.
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 &#39;file_type = exodusII is set in the input f...
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.
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
std::unique_ptr< MeshFunction > _mesh_function2
Pointer to second libMesh::MeshFuntion, used for interpolation.
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.
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.
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 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.
const std::map< SubdomainName, SubdomainID > & getBlockNamesToIds() const
Get the map from block name to block ID.
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
bool _initialized
True if initial_setup has executed.
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.
unsigned int getMeshFileDimension() const
Return the spatial dimension of the mesh file.
std::vector< Real > _scale
Scale parameter.
bool _interpolate_times
Flag for triggering interpolation of ExodusII data.
std::vector< std::string > _nodal_variables
Stores names of nodal variables.
const InputParameters & parameters() const
Get the parameters of the object.
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.
const std::string getMeshFileName() const
Return the name of the mesh file this object read the solution from.
std::string _mesh_file
The XDA or ExodusII file that is being read.
System * _system
Pointer libMesh::System class storing the read solution.
static MooseEnum weightingType()
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
virtual void initialize() override
Called before execute() is ever called so that data can be cleared.
std::vector< Real > _translation
Translation.