www.mooseframework.org
PiecewiseLinearInterpolationMaterial.C
Go to the documentation of this file.
1 /****************************************************************/
2 /* DO NOT MODIFY THIS HEADER */
3 /* MOOSE - Multiphysics Object Oriented Simulation Environment */
4 /* */
5 /* (c) 2010 Battelle Energy Alliance, LLC */
6 /* ALL RIGHTS RESERVED */
7 /* */
8 /* Prepared by Battelle Energy Alliance, LLC */
9 /* Under Contract No. DE-AC07-05ID14517 */
10 /* With the U. S. Department of Energy */
11 /* */
12 /* See COPYRIGHT for full restrictions */
13 /****************************************************************/
14 
16 
17 // MOOSE includes
18 #include "MooseVariable.h"
19 
20 template <>
23 {
25  params.addClassDescription("Compute a property using a piecewise linear interpolation to define "
26  "its dependence on a variable");
27  params.addRequiredParam<std::string>("property",
28  "The name of the property this material will compute");
29  params.addRequiredCoupledVar(
30  "variable",
31  "The name of the variable whose value is used as the abscissa in the interpolation");
32  params.addParam<std::vector<Real>>("x", "The abscissa values");
33  params.addParam<std::vector<Real>>("y", "The ordinate values");
34  params.addParam<std::vector<Real>>("xy_data",
35  "All function data, supplied in abscissa, ordinate pairs");
36  params.addParam<Real>("scale_factor", 1.0, "Scale factor to be applied to the ordinate values");
37  return params;
38 }
39 
41  const InputParameters & parameters)
43  _prop_name(getParam<std::string>("property")),
44  _coupled_var(coupledValue("variable")),
45  _scale_factor(getParam<Real>("scale_factor"))
46 {
47  std::vector<Real> x;
48  std::vector<Real> y;
49 
50  if ((parameters.isParamValid("x")) || (parameters.isParamValid("y")))
51  {
52  if (!((parameters.isParamValid("x")) && (parameters.isParamValid("y"))))
53  mooseError("In PiecewiseLinearInterpolationMaterial ",
54  _name,
55  ": Both 'x' and 'y' must be specified if either one is specified.");
56 
57  if (parameters.isParamValid("xy_data"))
58  mooseError("In PiecewiseLinearInterpolationMaterial ",
59  _name,
60  ": Cannot specify 'x', 'y', and 'xy_data' together.");
61 
62  x = getParam<std::vector<Real>>("x");
63  y = getParam<std::vector<Real>>("y");
64  }
65  else if (parameters.isParamValid("xy_data"))
66  {
67  std::vector<Real> xy = getParam<std::vector<Real>>("xy_data");
68  unsigned int xy_size = xy.size();
69  if (xy_size % 2 != 0)
70  mooseError("In PiecewiseLinearInterpolationMaterial ",
71  _name,
72  ": Length of data provided in 'xy_data' must be a multiple of 2.");
73 
74  unsigned int x_size = xy_size / 2;
75  x.reserve(x_size);
76  y.reserve(x_size);
77  for (unsigned int i = 0; i < xy_size / 2; ++i)
78  {
79  x.push_back(xy[i * 2]);
80  y.push_back(xy[i * 2 + 1]);
81  }
82  }
83 
84  try
85  {
86  _linear_interp = libmesh_make_unique<LinearInterpolation>(x, y);
87  }
88  catch (std::domain_error & e)
89  {
90  mooseError("In PiecewiseLinearInterpolationMaterial ", _name, ": ", e.what());
91  }
92 
93  _property = &declareProperty<Real>(_prop_name);
94  const VariableName & vname = getVar("variable", 0)->name();
95  _dproperty = &declarePropertyDerivative<Real>(_prop_name, vname);
96 }
97 
98 void
100 {
101  (*_property)[_qp] = _scale_factor * _linear_interp->sample(_coupled_var[_qp]);
102  (*_dproperty)[_qp] = _scale_factor * _linear_interp->sampleDerivative(_coupled_var[_qp]);
103 }
PiecewiseLinearInterpolationMaterial(const InputParameters &parameters)
MaterialProperty< Real > * _dproperty
First derivative of the material property wrt the coupled variable.
const VariableValue & _coupled_var
Value of the coupled variable to be used as the abscissa in the piecewise linear interpolation.
The main MOOSE class responsible for handling user-defined parameters in almost every MOOSE system...
MooseVariable * getVar(const std::string &var_name, unsigned int comp)
Extract pointer to a coupled variable.
Definition: Coupleable.C:117
InputParameters validParams< PiecewiseLinearInterpolationMaterial >()
static PetscErrorCode Vec x
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...
const std::string & name() const
Get the variable name.
unsigned int _qp
Definition: Material.h:220
bool isParamValid(const std::string &name) const
This method returns parameters that have been initialized in one fashion or another, i.e.
MaterialProperty< Real > * _property
Material property to be calculated.
std::unique_ptr< LinearInterpolation > _linear_interp
LinearInterpolation object.
Materials compute MaterialProperties.
Definition: Material.h:53
void addRequiredCoupledVar(const std::string &name, const std::string &doc_string)
This method adds a coupled variable name pair.
const std::string & _name
The name of this object, reference to value stored in InputParameters.
Definition: MooseObject.h:114
Interface class ("Veneer") to provide generator methods for derivative material property names...
virtual void computeQpProperties() override
Users must override this method.
void addClassDescription(const std::string &doc_string)
This method adds a description of the class that will be displayed in the input file syntax dump...
void addParam(const std::string &name, const S &value, const std::string &doc_string)
These methods add an option parameter and a documentation string to the InputParameters object...
void mooseError(Args &&...args) const
Definition: MooseObject.h:80
std::string _prop_name
Name of the property to be computed.
InputParameters validParams< Material >()
Definition: Material.C:27
const Real _scale_factor
Factor to scale the ordinate values by (default = 1)