ECLPropTable.hpp

Go to the documentation of this file.

/*
  Copyright 2017 Statoil ASA.
 
  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_ECLPROPTABLE_HEADER_INCLUDED
#define OPM_ECLPROPTABLE_HEADER_INCLUDED
 
#include <opm/utility/ECLPiecewiseLinearInterpolant.hpp>
#include <opm/utility/ECLTableInterpolation1D.hpp>
 
#include <functional>
#include <vector>
 
 
namespace Opm {
 
    struct ECLPropTableRawData
    {
        using DataVector = std::vector<double>;
 
        using SizeType = DataVector::size_type;
 
        using ElementIterator = DataVector::const_iterator;
 
        DataVector data{};
 
        DataVector primaryKey{};
 
        SizeType numPrimary{0};
 
        SizeType numRows{0};
 
        SizeType numCols{0};
 
        SizeType numTables{0};
    };
 
    template <class Interpolant>
    struct MakeInterpolants
    {
        template <class Factory>
        static std::vector<Interpolant>
        fromRawData(const ECLPropTableRawData& raw,
                    Factory&&                  construct)
        {
            auto interp = std::vector<Interpolant>{};
 
            const auto numInterp = raw.numTables * raw.numPrimary;
 
            // Table format: numRows*numInterp values of first column
            // (indep. var) followed by numCols-1 dependent variable
            // (function value result) columns of numRows*numInterp values
            // each, one column at a time.
            const auto colStride = raw.numRows * numInterp;
 
            // Position column iterators (independent variable and results
            // respectively) at beginning of each pertinent table column.
            auto xBegin = std::begin(raw.data);
            auto colIt  = std::vector<decltype(xBegin)> {
                xBegin + colStride
            };
 
            for (auto col = 0*raw.numCols + 1;
                      col <   raw.numCols - 1; ++col)
            {
                colIt.push_back(colIt.back() + colStride);
            }
 
            // Construct actual interpolants by invoking the
            // constructor/factory function on each sub-table.
            for (auto i = 0*numInterp;
                      i <   numInterp; ++i, xBegin += raw.numRows)
            {
                auto xEnd = xBegin + raw.numRows;
 
                // Layering violation:
                //    The constructor is expected to advance the result
                //    column iterators across 'numRows' entries.
                interp.push_back(construct(xBegin, xEnd, colIt));
            }
 
            return interp;
        }
    };
 
    class SatFuncInterpolant
    {
    public:
        struct ConvertUnits
        {
            using Converter = std::function<double(const double)>;
 
            Converter indep;
 
            std::vector<Converter> column;
        };
 
        SatFuncInterpolant(const ECLPropTableRawData& raw,
                           const ConvertUnits&        convert);
 
        struct InTable {
            ECLPropTableRawData::SizeType i;
        };
 
        struct ResultColumn {
            ECLPropTableRawData::SizeType i;
        };
 
        std::vector<double>
        interpolate(const InTable&             t,
                    const ResultColumn&        c,
                    const std::vector<double>& x) const;
 
        std::vector<double> connateSat() const;
 
        std::vector<double> criticalSat(const ResultColumn& c) const;
 
        std::vector<double> maximumSat() const;
 
        const std::vector<double>& saturationPoints(const InTable& t) const;
 
        ECLPropTableRawData::SizeType numTables() const
        {
            return this->table_.size();
        }
 
    private:
        class SingleTable
        {
        public:
            using ElmIt = ECLPropTableRawData::ElementIterator;
 
            SingleTable(ElmIt               xBegin,
                        ElmIt               xEnd,
                        const ConvertUnits& convert,
                        std::vector<ElmIt>& colIt);
 
            std::vector<double>
            interpolate(const ResultColumn&        c,
                        const std::vector<double>& x) const;
 
            double connateSat() const;
 
            double criticalSat(const ResultColumn& c) const;
 
            double maximumSat() const;
 
            const std::vector<double>& saturationPoints() const;
 
        private:
            using Extrap = ::Opm::Interp1D::PiecewisePolynomial::
                ExtrapolationPolicy::Constant;
 
            using Backend = ::Opm::Interp1D::
                PiecewisePolynomial::Linear<Extrap>;
 
            Backend interp_;
        };
 
        ECLPropTableRawData::SizeType nResCols_;
 
        std::vector<SingleTable> table_;
    };
 
} // namespace Opm
 
#endif // OPM_ECLPROPTABLE_HEADER_INCLUDED