diff --git a/doc/doxygen/dumux.bib b/doc/doxygen/dumux.bib index b98d485626718c079c6685e3af7e247baee802d9..ad5f39df176dba77f587472ee2ab69dfb8471530 100644 --- a/doc/doxygen/dumux.bib +++ b/doc/doxygen/dumux.bib @@ -2061,3 +2061,18 @@ author = {F. Fichot and F. Duval and N. Trégourès and C. Béchaud and M. Quint ISSN = {0065-2717}, publisher = {Elsevier} } + +@Article{Shomate1954, + author = {Shomate, C. Howard}, + title = {A Method for Evaluating and Correlating Thermodynamic Data}, + journal = {The Journal of Physical Chemistry}, + year = {1954}, + month = {Apr}, + day = {01}, + publisher = {American Chemical Society}, + volume = {58}, + number = {4}, + pages = {368-372}, + issn = {0022-3654}, + doi = {10.1021/j150514a018} +} diff --git a/dumux/material/components/shomate.hh b/dumux/material/components/shomate.hh new file mode 100644 index 0000000000000000000000000000000000000000..b55c06a1e8861823b553ba0e30771a6b8732a7bb --- /dev/null +++ b/dumux/material/components/shomate.hh @@ -0,0 +1,152 @@ +// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*- +// vi: set et ts=4 sw=4 sts=4: +// +// SPDX-FileCopyrightInfo: Copyright © DuMux Project contributors, see AUTHORS.md in root folder +// SPDX-License-Identifier: GPL-3.0-or-later +// +/*! + * \file + * \ingroup Components + * \brief Apply Shomate equations for enthalpy and heat capacity. + * + * The Shomate equations were inroduced in + * "A method for evaluating and correlating thermodynamic data" by C. Howard Shomate in 1954, \cite Shomate1954. + * They offer analytical equations for computing the heat capacity as well as enthalpy. + * Both equations use a set of component-specific parameters \f$A,B,C,D,E,F,G,H\f$. + * For the heat capacity \f$C_p^0\f$, one obtains + \f[ + C_p^0 = A + Bt + Ct^2 + Dt^3 + E/t^2, + \f] + * while for the enthalpy with a reference state at \f$T=298.15\text{K}\f$, one uses + \f[ + H^0-H_{298.15}^0 = At + \frac{Bt^2}{2} + \frac{Ct^3}{3} + \frac{Dt^4}{4} -\frac{E}{t} + F - H + \f] + *where: + * * \f$ C_p \f$ is the heat capacity in J/(mol*K), + * * \f$ H^0 \f$ represents the standard enthalpy in kJ/mol, + * * \f$ t \f$ is the temperature in K divided by 1000. + * + */ + +#ifndef DUMUX_MATERIAL_COMPONENTS_SHOMATE_HH +#define DUMUX_MATERIAL_COMPONENTS_SHOMATE_HH +#include <vector> +#include <stdexcept> +#include <iostream> + +namespace Dumux { +template<class Scalar> +class ShomateMethod +{ +public: + /*! + * \brief Constructor for Shomate method class + * + * \param tempRange contains the temperature boundaries of the intervals defined by the Shomate method, the number of intervals can vary + * \param coeffs contains the sets of coefficients for each temperature interval bounded by the entries of tempRange + * + */ + ShomateMethod(std::vector<Scalar> tempRange, std::vector<std::array<Scalar,8>> coeffs) + : tempRange_(tempRange) + , coeffs_(coeffs) + { + if (tempRange_.size()-1 != coeffs_.size()) + { + DUNE_THROW(Dune::InvalidStateException,"Temperature range and coefficients must have the same size."); + } + for(size_t i=0; i<coeffs_.size(); i++) + { + if (coeffs_[i].size() != 8) + { + DUNE_THROW(Dune::InvalidStateException,"Shomate coeffs must have 8 elements."); + } + } + } + + std::array<Scalar,8> paramsAtTemperature(const Scalar& T) const + { + // check if T is smaller or higher than allowed min/max T + if (T < tempRange_[0]) + { + if (!warningPrinted_) + { + std::cerr << "Temperature "<< T << " [K] is out of range. Enthalpy values are extrapolated." << std::endl; + warningPrinted_ = true; + } + return coeffs_[0]; + } + + for (size_t i = 0; i < tempRange_.size()-1; i++) + { + if (T >= tempRange_[i] && T < tempRange_[i+1]) + { + return coeffs_[i]; + } + } + + if (!warningPrinted_) + { + std::cerr << "Temperature "<< T << " [K] is out of range. Enthalpy values are extrapolated." << std::endl; + warningPrinted_ = true; + } + return coeffs_[tempRange_.size()-1]; + } + + /** + * @brief return enthalpy in kJ/mol + * + * @param temperature in K + * @param pressure in Pa + * @return Scalar + */ + Scalar enthalpy(const Scalar& temperature, + const Scalar& pressure) const + { + const auto& correctCoeffs = paramsAtTemperature(temperature); + + Scalar A = correctCoeffs[0]; + Scalar B = correctCoeffs[1]; + Scalar C = correctCoeffs[2]; + Scalar D = correctCoeffs[3]; + Scalar E = correctCoeffs[4]; + Scalar F = correctCoeffs[5]; + Scalar H = correctCoeffs[7]; + //convert temperature to correct format + Scalar t = temperature/1000.0; + //calculate standard enthalpy difference according to Shomate + Scalar standardEnthalpyDiff = A*t + B*t*t/2.0 + C*t*t*t/3.0 + D*t*t*t*t/4.0 - E/t + F - H; + //unit is kJ/mol + return standardEnthalpyDiff; + } + + /** + * @brief return heat capacity in J/(mol*K) + * + * @param temperature in K + * @param pressure in Pa + * @return Scalar + */ + Scalar heatCapacity(const Scalar& temperature, + const Scalar& pressure) const + { + auto correctCoeffs = paramsAtTemperature(temperature); + Scalar A = correctCoeffs[0]; + Scalar B = correctCoeffs[1]; + Scalar C = correctCoeffs[2]; + Scalar D = correctCoeffs[3]; + Scalar E = correctCoeffs[4]; + Scalar t = temperature/1000.0; + Scalar heatCapacity = A + B*t + C*t*t + D*t*t*t + E/(t*t); + return heatCapacity; + } + +private: + std::vector<Scalar> tempRange_; + std::vector<std::array<Scalar,8>> coeffs_; + static bool warningPrinted_; +}; + +template<class Scalar> +bool ShomateMethod<Scalar>::warningPrinted_ = false; +} // namespace Dumux +#endif // DUMUX_SHOMATE_HH