www.mooseframework.org
EquilibriumConstantAux.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 "EquilibriumConstantAux.h"
11 
12 registerMooseObject("ChemicalReactionsApp", EquilibriumConstantAux);
13 
16 {
18  params.addCoupledVar(
19  "temperature", 298.15, "The temperature of the aqueous phase (K). Default is 298.15K");
20  params.addRequiredParam<std::vector<Real>>(
21  "temperature_points", "Temperature points where log(Keq) data is evaluated (K)");
22  params.addRequiredParam<std::vector<Real>>(
23  "logk_points", "log(Keq) data evaluated at each value of temperature_points");
24  params.addClassDescription(
25  "Equilibrium constant for a given equilibrium species (in form log10(Keq))");
26  return params;
27 }
28 
30  : AuxKernel(parameters),
31  _temperature(coupledValue("temperature")),
32  _temperature_points(getParam<std::vector<Real>>("temperature_points")),
33  _logk_points(getParam<std::vector<Real>>("logk_points"))
34 {
35  // Check that the number of temperature_points and logk_points are equal
36  if (_temperature_points.size() != _logk_points.size())
37  mooseError("The number of temperature_points and logk_points must be equal in ", _name);
38 
39  if (_temperature_points.size() >= 5)
40  {
41  // If there at least 5 values, then use the Maier-Kelley fit
42  _logk = std::make_unique<EquilibriumConstantFit>(_temperature_points, _logk_points);
43  _logk->generate();
44  }
45  else if ((_temperature_points.size() >= 2) && (_temperature_points.size() <= 4))
46  {
47  // If between 2 and 4 values are provided, use a linear fit
48  _linear_logk = std::make_unique<PolynomialFit>(_temperature_points, _logk_points, 1);
49  _linear_logk->generate();
50  }
51 }
52 
53 Real
55 {
56  if (_temperature_points.size() == 1)
57  return -_logk_points[0];
58 
59  if (_temperature_points.size() > 5)
60  return -_logk->sample(_temperature[_qp]);
61 
62  return -_linear_logk->sample(_temperature[_qp]);
63 }
std::unique_ptr< PolynomialFit > _linear_logk
Linear least-squares fit.
const std::vector< Real > & _logk_points
log(Keq) values at each temperature point
const std::vector< Real > & _temperature_points
Temperature points in data set (in K)
Equilibrium constant (in the form log10(Keq)) calculated using a least-squares fit to the data provid...
std::unique_ptr< EquilibriumConstantFit > _logk
Least-squares fit to data.
void addRequiredParam(const std::string &name, const std::string &doc_string)
virtual Real computeValue() override
const std::string _name
void addCoupledVar(const std::string &name, const std::string &doc_string)
EquilibriumConstantAux(const InputParameters &parameters)
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
static InputParameters validParams()
void mooseError(Args &&... args) const
void addClassDescription(const std::string &doc_string)
static InputParameters validParams()
registerMooseObject("ChemicalReactionsApp", EquilibriumConstantAux)
const VariableValue & _temperature
Temperature (in K)