From e29dd18d7746cc47bcea6e8868cdb3fc7927d3ae Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Atgeirr=20Fl=C3=B8=20Rasmussen?= <atgeirr@sintef.no>
Date: Mon, 19 Mar 2012 13:55:26 +0100
Subject: [PATCH] Added RockCompressibility class.

---
 Makefile.am                            |  2 +
 opm/core/fluid/RockCompressibility.cpp | 73 ++++++++++++++++++++++++++
 opm/core/fluid/RockCompressibility.hpp | 52 ++++++++++++++++++
 3 files changed, 127 insertions(+)
 create mode 100644 opm/core/fluid/RockCompressibility.cpp
 create mode 100644 opm/core/fluid/RockCompressibility.hpp

diff --git a/Makefile.am b/Makefile.am
index add0ae09..5ec45d06 100644
--- a/Makefile.am
+++ b/Makefile.am
@@ -29,6 +29,7 @@ opm/core/fluid/IncompPropertiesFromDeck.cpp             \
 opm/core/fluid/PvtPropertiesBasic.cpp                   \
 opm/core/fluid/PvtPropertiesIncompFromDeck.cpp          \
 opm/core/fluid/RockBasic.cpp                            \
+opm/core/fluid/RockCompressibility.cpp                  \
 opm/core/fluid/RockFromDeck.cpp                         \
 opm/core/fluid/SaturationPropsBasic.cpp                 \
 opm/core/fluid/SaturationPropsFromDeck.cpp              \
@@ -107,6 +108,7 @@ opm/core/fluid/IncompPropertiesFromDeck.hpp             \
 opm/core/fluid/PvtPropertiesBasic.hpp                   \
 opm/core/fluid/PvtPropertiesIncompFromDeck.hpp          \
 opm/core/fluid/RockBasic.hpp                            \
+opm/core/fluid/RockCompressibility.hpp                  \
 opm/core/fluid/RockFromDeck.hpp                         \
 opm/core/fluid/SaturationPropsBasic.hpp                 \
 opm/core/fluid/SaturationPropsFromDeck.hpp              \
diff --git a/opm/core/fluid/RockCompressibility.cpp b/opm/core/fluid/RockCompressibility.cpp
new file mode 100644
index 00000000..ca4e3ec6
--- /dev/null
+++ b/opm/core/fluid/RockCompressibility.cpp
@@ -0,0 +1,73 @@
+/*
+  Copyright 2012 SINTEF ICT, Applied Mathematics.
+
+  This file is part of the Open Porous Media project (OPM).
+
+  OPM 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.
+
+  OPM 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 OPM.  If not, see <http://www.gnu.org/licenses/>.
+*/
+
+#include <opm/core/fluid/RockCompressibility.hpp>
+#include <opm/core/eclipse/EclipseGridParser.hpp>
+#include <opm/core/utility/ErrorMacros.hpp>
+#include <opm/core/utility/linearInterpolation.hpp>
+
+namespace Opm
+{
+
+    RockCompressibility::RockCompressibility(const EclipseGridParser& deck)
+        : pref_(0.0),
+          rock_comp_(0.0)
+    {
+        if (deck.hasField("ROCKTAB")) {
+            const table_t& rt = deck.getROCKTAB().rocktab_;
+            int n = rt[0].size();
+            p_.resize(n);
+            poromult_.resize(n);
+            for (int i = 0; i < n; ++i) {
+                p_[i] = rt[0][i][0];
+                poromult_[i] = rt[0][i][1];
+            }
+        } else if (deck.hasField("ROCK")) {
+            const ROCK& r = deck.getROCK();
+            pref_ = r.rock_compressibilities_[0][0];
+            rock_comp_ = r.rock_compressibilities_[0][1];
+        } else {
+            std::cout << "**** warning: no rock compressibility data found in deck (ROCK or ROCKTAB)." << std::endl;
+        }
+    }
+
+    double RockCompressibility::poroMult(double pressure)
+    {
+        if (p_.empty()) {
+            // Approximating with a quadratic curve.
+            const double cpnorm = rock_comp_*(pressure - pref_);
+            return (1.0 + cpnorm + 0.5*cpnorm*cpnorm);
+        } else {
+            return Opm::linearInterpolation(p_, poromult_, pressure);
+        }
+    }
+
+    double RockCompressibility::rockComp(double pressure)
+    {
+        if (p_.empty()) {
+            return rock_comp_;
+        } else {
+            const double poromult = Opm::linearInterpolation(p_, poromult_, pressure);
+            const double dporomultdp = Opm::linearInterpolationDerivative(p_, poromult_, pressure);
+            return dporomultdp/poromult;
+        }
+    }
+
+} // namespace Opm
+
diff --git a/opm/core/fluid/RockCompressibility.hpp b/opm/core/fluid/RockCompressibility.hpp
new file mode 100644
index 00000000..746ffeac
--- /dev/null
+++ b/opm/core/fluid/RockCompressibility.hpp
@@ -0,0 +1,52 @@
+/*
+  Copyright 2012 SINTEF ICT, Applied Mathematics.
+
+  This file is part of the Open Porous Media project (OPM).
+
+  OPM 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.
+
+  OPM 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 OPM.  If not, see <http://www.gnu.org/licenses/>.
+*/
+
+#ifndef OPM_ROCKCOMPRESSIBILITY_HEADER_INCLUDED
+#define OPM_ROCKCOMPRESSIBILITY_HEADER_INCLUDED
+
+#include <vector>
+
+namespace Opm
+{
+
+    class EclipseGridParser;
+
+    class RockCompressibility
+    {
+    public:
+        /// Construct from input deck.
+        RockCompressibility(const EclipseGridParser& deck);
+
+        /// Porosity multiplier.
+        double poroMult(double pressure);
+
+        /// Rock compressibility = (d poro / d p)*(1 / poro).
+        double rockComp(double pressure);
+
+    private:
+        std::vector<double> p_;
+        std::vector<double> poromult_;
+        double pref_;
+        double rock_comp_;
+    };
+
+} // namespace Opm
+
+
+#endif // OPM_ROCKCOMPRESSIBILITY_HEADER_INCLUDED