Merge branch 'feature/updateBiominExampleCO2Tables' into 'master'

[example][biomineralization][material][co2] use new type of CO2 tables,...

See merge request !3062

examples/biomineralization/.doc_config

@@ -21,7 +21,7 @@
 
     "doc/fluidmaterial.md" : [
         "doc/fluidmaterial_intro.md",
-        "material/co2tableslaboratory.hh",
+        "material/co2tables.hh",
         "material/components/suspendedbiomass.hh",
         "material/fluidsystems/biominsimplechemistry.hh",
         "material/fluidsystems/icpcomplexsalinitybrine.hh"

examples/biomineralization/README.md

@@ -122,7 +122,7 @@ We will discuss the different parts of the code in detail subsequently.
         │   └── icpcomplexsalinitybrine.hh
         ├── solidsystems/
         │   └── biominsolids.hh
-        └── co2tableslaboratory.hh
+        └── co2tables.hh
 ```
 
 In order to define a simulation setup in DuMu<sup>x</sup>, you need to implement compile-time settings,
@@ -149,7 +149,7 @@ Similar strategies might be useful when simulating experimental setups with boun
 
 __Part 4__ discusses the code concerned with the fluid (files in the folder `material/`),
 especially the multi-component fluidsystems.
-The CO<sub>2</sub> properties are stored in the CO<sub>2</sub> tables in the subfolder `material` (`co2tableslaboratory.hh`,`co2valueslaboratory.inc`).
+The CO<sub>2</sub> properties are stored in the CO<sub>2</sub> tables in the subfolder `material` (`co2tables.hh`,`co2values.inc`).
 
 __Part 5__ discusses the code concerned with the solid properties (files in the folder `material/`),
 especially the multi-component solidsystems and variable solid volume fractions.

examples/biomineralization/doc/_intro.md

@@ -120,7 +120,7 @@ We will discuss the different parts of the code in detail subsequently.
         │   └── icpcomplexsalinitybrine.hh
         ├── solidsystems/
         │   └── biominsolids.hh
-        └── co2tableslaboratory.hh
+        └── co2tables.hh
 ```
 
 In order to define a simulation setup in DuMu<sup>x</sup>, you need to implement compile-time settings,
@@ -147,7 +147,7 @@ Similar strategies might be useful when simulating experimental setups with boun
 
 __Part 4__ discusses the code concerned with the fluid (files in the folder `material/`),
 especially the multi-component fluidsystems.
-The CO<sub>2</sub> properties are stored in the CO<sub>2</sub> tables in the subfolder `material` (`co2tableslaboratory.hh`,`co2valueslaboratory.inc`).
+The CO<sub>2</sub> properties are stored in the CO<sub>2</sub> tables in the subfolder `material` (`co2tables.hh`,`co2values.inc`).
 
 __Part 5__ discusses the code concerned with the solid properties (files in the folder `material/`),
 especially the multi-component solidsystems and variable solid volume fractions.

examples/biomineralization/doc/fluidmaterial.md

@@ -30,7 +30,7 @@ For further specialization, the overview over the material subfolder is split in
 ## Fluids in the folder `material`
 
 As this example is about biomineralization involving many components with complex inteactions, some specific fluid material files are necessary.
-A CO_2-Table file provides tabulated CO_2 properties according to @Span1996 in `material/co2tableslaboratory.hh`
+A CO_2-Table file provides tabulated CO_2 properties according to @Span1996 in `material/co2tables.hh`
 In the component subfolder, `material/components/suspendedbiomass.hh` defines the component suspended biomass, which is the mobile form of biomass being transported suspended in the aqueous fluid phase.
 In the fluidsystem subfolder, the biomineralization fluidsystem `material/fluidsystems/biominsimplechemistry.hh` as well as
 the complex salinity brine adapter `material/fluidsystems/icpcomplexsalinitybrine.hh` can be found.
@@ -48,24 +48,174 @@ The subsequent documentation is structured as follows:
 [@Span1996]: https://aip.scitation.org/doi/abs/10.1063/1.555991 "A new equation of state for carbon dioxide covering the fluid region from the triple-point temperature to 1100 K at pressures up to 800 MPa"
 
 
-## The CO2 tables (`co2tableslaboratory.hh`)
 
-This file contains the __co2table class__ which forwards to tabulated properties of CO2 according to Span and Wagner 1996.
-The real work (creating the tables) is done by some external program by Span and Wagner 1996 which provides the ready-to-use tables.
+```cpp
+/*!
+ * \file
+ * \ingroup Components
+ * \brief A reader and tables for CO$_2$ tabulated material laws that depend
+ *        on pressure and temperature.
+ */
+
+
+#include <dune/common/float_cmp.hh>
+```
+
+## The CO2 tables (`co2tables.hh`)
+
+This file contains the __co2table class__ which forwards to tabulated properties of CO2.
+The tables are generated using the NIST (National Institute of Standards
+and Technology) Standard Reference Database Number 69
+(https://doi.org/10.18434/T4D303).
+
+Copyright for NIST Standard Reference Data is governed by the Standard
+Reference Data Act (https://www.nist.gov/srd/public-law).
+
+######################################################################
+In case you are using this the data generated with this script
+please cite the following publications:
+
+P.J. Linstrom and W.G. Mallard, Eds.,
+NIST Chemistry WebBook, NIST Standard Reference Database Number 69,
+National Institute of Standards and Technology, Gaithersburg MD, 20899,
+https://doi.org/10.18434/T4D303, (retrieved [insert date]).
+
+Span, Roland, and Wolfgang Wagner.
+"A new equation of state for carbon dioxide covering
+the fluid region from the triple‐point temperature
+to 1100 K at pressures up to 800 MPa."
+Journal of physical and chemical reference data 25.6 (1996): 1509-1596.
+https://doi.org/10.1063/1.555991
+
+######################################################################
+
+The density and the enthalpy are calculated using the equation of Span and
+Wagner (2009 "A New Equation of State for Carbon Dioxide Covering the Fluid
+Region from the Triple-Point Temperature to 1100 K at Pressures up to 800 MPa").
+Therefore, the maximum pressure limit is the lowest of the following values:
+* 800.0000 MPa
+* The pressure at which a density of 1178.5 kg/m3 is reached.
 
 
 <details open>
-<summary><b>Click to hide/show the file documentation</b> (or inspect the [source code](../material/co2tableslaboratory.hh))</summary>
+<summary><b>Click to hide/show the file documentation</b> (or inspect the [source code](../material/co2tables.hh))</summary>
 
 
 
 ```cpp
 
-#include <assert.h>
-#include <test/porousmediumflow/co2/co2tablereader.hh>
-namespace Dumux::ICP {
-#include "co2valueslaboratory.inc"
-}// end namespace Dumux::ICP
+namespace Dumux::BiomineralizationCO2Tables {
+
+/*!
+ * \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"
+
+
+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
 ```
 
 

examples/biomineralization/doc/fluidmaterial_intro.md

@@ -24,7 +24,7 @@ For further specialization, the overview over the material subfolder is split in
 ## Fluids in the folder `material`
 
 As this example is about biomineralization involving many components with complex inteactions, some specific fluid material files are necessary.
-A CO_2-Table file provides tabulated CO_2 properties according to @Span1996 in `material/co2tableslaboratory.hh`
+A CO_2-Table file provides tabulated CO_2 properties according to @Span1996 in `material/co2tables.hh`
 In the component subfolder, `material/components/suspendedbiomass.hh` defines the component suspended biomass, which is the mobile form of biomass being transported suspended in the aqueous fluid phase.
 In the fluidsystem subfolder, the biomineralization fluidsystem `material/fluidsystems/biominsimplechemistry.hh` as well as
 the complex salinity brine adapter `material/fluidsystems/icpcomplexsalinitybrine.hh` can be found.

examples/biomineralization/doc/setup.md

@@ -835,7 +835,7 @@ We include the necessary material files
 ```cpp
 #include <examples/biomineralization/material/fluidsystems/biominsimplechemistry.hh>
 #include <examples/biomineralization/material/solidsystems/biominsolids.hh>
-#include <examples/biomineralization/material/co2tableslaboratory.hh>
+#include <examples/biomineralization/material/co2tables.hh>
 ```
 
 We include the problem and spatial parameters headers used for this simulation.
@@ -886,7 +886,7 @@ template<class TypeTag>
 struct FluidSystem<TypeTag, TTag::MICPColumnSimpleChemistry>
 {
     using Scalar = GetPropType<TypeTag, Properties::Scalar>;
-    using CO2Tables = Dumux::ICP::CO2Tables;
+    using CO2Tables = BiomineralizationCO2Tables::CO2Tables;
     using H2OTabulated = Components::TabulatedComponent<Components::H2O<Scalar>>;
     using type = Dumux::FluidSystems::BioMinSimpleChemistryFluid<Scalar, CO2Tables, H2OTabulated>;
 };

examples/biomineralization/material/co2tables.hh

+// -*- 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 reader and tables for CO$_2$ tabulated material laws that depend
+ *        on pressure and temperature.
+ */
+#ifndef DUMUX_EXAMPLE_BIOMINERALIZATION_CO2TABLES_HH
+#define DUMUX_EXAMPLE_BIOMINERALIZATION_CO2TABLES_HH
+
+#include <dune/common/float_cmp.hh>
+
+// ## The CO2 tables (`co2tables.hh`)
+//
+// This file contains the __co2table class__ which forwards to tabulated properties of CO2.
+// The tables are generated using the NIST (National Institute of Standards
+// and Technology) Standard Reference Database Number 69
+// (https://doi.org/10.18434/T4D303).
+//
+// Copyright for NIST Standard Reference Data is governed by the Standard
+// Reference Data Act (https://www.nist.gov/srd/public-law).
+//
+// ######################################################################
+// In case you are using this the data generated with this script
+// please cite the following publications:
+//
+// P.J. Linstrom and W.G. Mallard, Eds.,
+// NIST Chemistry WebBook, NIST Standard Reference Database Number 69,
+// National Institute of Standards and Technology, Gaithersburg MD, 20899,
+// https://doi.org/10.18434/T4D303, (retrieved [insert date]).
+//
+// Span, Roland, and Wolfgang Wagner.
+// "A new equation of state for carbon dioxide covering
+// the fluid region from the triple‐point temperature
+// to 1100 K at pressures up to 800 MPa."
+// Journal of physical and chemical reference data 25.6 (1996): 1509-1596.
+// https://doi.org/10.1063/1.555991
+//
+// ######################################################################
+//
+// The density and the enthalpy are calculated using the equation of Span and
+// Wagner (2009 "A New Equation of State for Carbon Dioxide Covering the Fluid
+// Region from the Triple-Point Temperature to 1100 K at Pressures up to 800 MPa").
+// Therefore, the maximum pressure limit is the lowest of the following values:
+// * 800.0000 MPa
+// * The pressure at which a density of 1178.5 kg/m3 is reached.
+//
+// [[content]]
+//
+// [[codeblock]]
+
+namespace Dumux::BiomineralizationCO2Tables {
+
+/*!
+ * \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
+// [[/codeblock]]
+// [[/content]]
+
+#endif

examples/biomineralization/material/co2tableslaboratory.hh

-// -*- 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 2 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/>.   *
- *****************************************************************************/
-
-#ifndef DUMUX_ICP_CO2TABLES_LABORATORY_HH
-#define DUMUX_ICP_CO2TABLES_LABORATORY_HH
-
-// ## The CO2 tables (`co2tableslaboratory.hh`)
-//
-// This file contains the __co2table class__ which forwards to tabulated properties of CO2 according to Span and Wagner 1996.
-// The real work (creating the tables) is done by some external program by Span and Wagner 1996 which provides the ready-to-use tables.
-//
-// [[content]]
-//
-// [[codeblock]]
-
-#include <assert.h>
-#include <test/porousmediumflow/co2/co2tablereader.hh>
-namespace Dumux::ICP {
-#include "co2valueslaboratory.inc"
-}// end namespace Dumux::ICP
-// [[/codeblock]]
-// [[/content]]
-#endif

examples/biomineralization/properties.hh

@@ -47,7 +47,7 @@
 // We include the necessary material files
 #include <examples/biomineralization/material/fluidsystems/biominsimplechemistry.hh>
 #include <examples/biomineralization/material/solidsystems/biominsolids.hh>
-#include <examples/biomineralization/material/co2tableslaboratory.hh>
+#include <examples/biomineralization/material/co2tables.hh>
 
 // We include the problem and spatial parameters headers used for this simulation.
 #include "problem.hh"
@@ -93,7 +93,7 @@ template<class TypeTag>
 struct FluidSystem<TypeTag, TTag::MICPColumnSimpleChemistry>
 {
     using Scalar = GetPropType<TypeTag, Properties::Scalar>;
-    using CO2Tables = Dumux::ICP::CO2Tables;
+    using CO2Tables = BiomineralizationCO2Tables::CO2Tables;
     using H2OTabulated = Components::TabulatedComponent<Components::H2O<Scalar>>;
     using type = Dumux::FluidSystems::BioMinSimpleChemistryFluid<Scalar, CO2Tables, H2OTabulated>;
 };