// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*- // vi: set et ts=4 sw=4 sts=4: /***************************************************************************** * See the file COPYING for full copying permissions. * * * * This program is free software: you can redistribute it and/or modify * * it under the terms of the GNU General Public License as published by * * the Free Software Foundation, either version 3 of the License, or * * (at your option) any later version. * * * * This program is distributed in the hope that it will be useful, * * but WITHOUT ANY WARRANTY; without even the implied warranty of * * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the * * GNU General Public License for more details. * * * * You should have received a copy of the GNU General Public License * * along with this program. If not, see <http://www.gnu.org/licenses/>. * *****************************************************************************/ /*! * \file * \ingroup Components * \brief A generic template for tabulated material laws that depend * on two parameters. */ #ifndef DUMUX_EXERCISE_CO2TABLES_HH #define DUMUX_EXERCISE_CO2TABLES_HH #include <dune/common/float_cmp.hh> namespace Dumux::BioMinCO2Tables { /*! * \ingroup Components * \brief A generic template for tabulated material laws that depend * on two parameters. */ template <class Traits> class TabulatedProperties { using Scalar = typename Traits::Scalar; static constexpr auto numTempSteps = Traits::numTempSteps; static constexpr auto numPressSteps = Traits::numPressSteps; public: TabulatedProperties() = default; constexpr Scalar minTemp() const { return Traits::minTemp; } constexpr Scalar maxTemp() const { return Traits::maxTemp; } constexpr Scalar minPress() const { return Traits::minPress; } constexpr Scalar maxPress() const { return Traits::maxPress; } constexpr bool applies(Scalar temperature, Scalar pressure) const { return minTemp() <= temperature && temperature <= maxTemp() && minPress() <= pressure && pressure <= maxPress(); } constexpr Scalar at(Scalar temperature, Scalar pressure) const { if (!applies(temperature, pressure)) { if (temperature<minTemp()) temperature = minTemp(); else if (temperature>maxTemp()) temperature = maxTemp(); if (pressure<minPress()) pressure = minPress(); else if (pressure>maxPress()) pressure = maxPress(); } const int i = findTempIdx_(temperature); const int j = findPressIdx_(pressure); const Scalar tempAtI = temperatureAt_(i); const Scalar tempAtI1 = temperatureAt_(i + 1); const Scalar pressAtI = pressureAt_(j); const Scalar pressAtI1 = pressureAt_(j + 1); const Scalar alpha = (temperature - tempAtI)/(tempAtI1 - tempAtI); const Scalar beta = (pressure - pressAtI)/(pressAtI1 - pressAtI); // bi-linear interpolation const Scalar lowresValue = (1-alpha)*(1-beta)*val(i, j) + (1-alpha)*( beta)*val(i, j + 1) + ( alpha)*(1-beta)*val(i + 1, j) + ( alpha)*( beta)*val(i + 1, j + 1); // return the weighted sum of the low- and high-resolution values return lowresValue; } constexpr Scalar val(int i, int j) const { return Traits::vals[i][j]; } private: constexpr int findTempIdx_(Scalar temperature) const { if (Dune::FloatCmp::eq<Scalar>(temperature, maxTemp())) return numTempSteps - 2; const int result = static_cast<int>((temperature - minTemp())/(maxTemp() - minTemp())*(numTempSteps - 1)); using std::clamp; return clamp(result, 0, numTempSteps - 2); } constexpr int findPressIdx_(Scalar pressure) const { if (Dune::FloatCmp::eq<Scalar>(pressure, maxPress())) return numPressSteps - 2; const int result = static_cast<int>((pressure - minPress())/(maxPress() - minPress())*(numPressSteps - 1)); using std::clamp; return clamp(result, 0, numPressSteps - 2); } constexpr Scalar temperatureAt_(int i) const { return i*(maxTemp() - minTemp())/(numTempSteps - 1) + minTemp(); } constexpr Scalar pressureAt_(int j) const { return j*(maxPress() - minPress())/(numPressSteps - 1) + minPress(); } }; #ifndef DOXYGEN // hide from doxygen // the real work is done by some external program which provides // ready-to-use tables. #include "co2values.inc" #endif using TabulatedDensity = TabulatedProperties<TabulatedDensityTraits>; using TabulatedEnthalpy = TabulatedProperties<TabulatedEnthalpyTraits>; // this class collects all the tabulated quantities in one convenient place struct CO2Tables { static constexpr inline TabulatedEnthalpy tabulatedEnthalpy = {}; static constexpr inline TabulatedDensity tabulatedDensity = {}; }; } // end namespace Dumux::GeneratedCO2Tables #endif