From 3f683e607e7a6cf018aeef2c722b35f765d49314 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.

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

diff --git a/opm/core/fluid/RockCompressibility.cpp b/opm/core/fluid/RockCompressibility.cpp
new file mode 100644
index 000000000..ca4e3ec6f
--- /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 000000000..746ffeac2
--- /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