diff --git a/opm/core/props/pvt/BlackoilPvtProperties.cpp b/opm/core/props/pvt/BlackoilPvtProperties.cpp index 217107cb..0fd7dda4 100644 --- a/opm/core/props/pvt/BlackoilPvtProperties.cpp +++ b/opm/core/props/pvt/BlackoilPvtProperties.cpp @@ -32,6 +32,9 @@ #include #include +#include +#include +#include namespace Opm { @@ -40,7 +43,6 @@ namespace Opm { } - void BlackoilPvtProperties::init(const EclipseGridParser& deck, const int samples) { // If we need multiple regions, this class and the SinglePvt* classes must change. @@ -113,6 +115,83 @@ namespace Opm } } + void BlackoilPvtProperties::init(Opm::DeckConstPtr newParserDeck, int samples) + { + // If we need multiple regions, this class and the SinglePvt* classes must change. + region_number_ = 0; + + phase_usage_ = phaseUsageFromDeck(newParserDeck); + + // Surface densities. Accounting for different orders in eclipse and our code. + if (newParserDeck->hasKeyword("DENSITY")) { + Opm::DeckKeywordConstPtr densityKeyword = newParserDeck->getKeyword("DENSITY"); + const std::vector& d = densityKeyword->getRecord(region_number_)->getItem(0)->getSIDoubleData(); + enum { ECL_oil = 0, ECL_water = 1, ECL_gas = 2 }; + if (phase_usage_.phase_used[Aqua]) { + densities_[phase_usage_.phase_pos[Aqua]] = d[ECL_water]; + } + if (phase_usage_.phase_used[Vapour]) { + densities_[phase_usage_.phase_pos[Vapour]] = d[ECL_gas]; + } + if (phase_usage_.phase_used[Liquid]) { + densities_[phase_usage_.phase_pos[Liquid]] = d[ECL_oil]; + } + } else { + OPM_THROW(std::runtime_error, "Input is missing DENSITY\n"); + } + + // Set the properties. + props_.resize(phase_usage_.num_phases); + // Water PVT + if (phase_usage_.phase_used[Aqua]) { + if (newParserDeck->hasKeyword("PVTW")) { + Opm::PvtwTable pvtwTable(newParserDeck->getKeyword("PVTW")); + + props_[phase_usage_.phase_pos[Aqua]].reset(new SinglePvtConstCompr(pvtwTable)); + } else { + // Eclipse 100 default. + props_[phase_usage_.phase_pos[Aqua]].reset(new SinglePvtConstCompr(0.5*Opm::prefix::centi*Opm::unit::Poise)); + } + } + // Oil PVT + if (phase_usage_.phase_used[Liquid]) { + if (newParserDeck->hasKeyword("PVDO")) { + Opm::PvdoTable pvdoTable(newParserDeck->getKeyword("PVDO"), region_number_); + + props_[phase_usage_.phase_pos[Liquid]].reset(new SinglePvtDeadSpline(pvdoTable, samples)); + } else if (newParserDeck->hasKeyword("PVTO")) { + Opm::PvtoTable pvtoTable(newParserDeck->getKeyword("PVTO"), /*tableIdx=*/0); + + props_[phase_usage_.phase_pos[Liquid]].reset(new SinglePvtLiveOil(pvtoTable)); + } else if (newParserDeck->hasKeyword("PVCDO")) { + Opm::PvdcoTable pvcdoTable(newParserDeck->getKeyword("PVCDO")); + + props_[phase_usage_.phase_pos[Liquid]].reset(new SinglePvtConstCompr(pvcdoTable)); + } else { + OPM_THROW(std::runtime_error, "Input is missing PVDO or PVTO\n"); + } + } + // Gas PVT + if (phase_usage_.phase_used[Vapour]) { + if (newParserDeck->hasKeyword("PVDG")) { + Opm::PvdgTable pvdgTable(newParserDeck->getKeyword("PVDG"), region_number_); + + props_[phase_usage_.phase_pos[Vapour]].reset(new SinglePvtDeadSpline(pvdgTable, samples)); + } else if (newParserDeck->hasKeyword("PVTG")) { + Opm::PvtgTable pvtgTable(newParserDeck->getKeyword("PVTG"), /*tableIdx=*/0); + + props_[phase_usage_.phase_pos[Vapour]].reset(new SinglePvtLiveGas(pvtgTable)); + } else { + OPM_THROW(std::runtime_error, "Input is missing PVDG or PVTG\n"); + } + } + + // Must inform pvt property objects of phase structure. + for (int i = 0; i < phase_usage_.num_phases; ++i) { + props_[i]->setPhaseConfiguration(phase_usage_.num_phases, phase_usage_.phase_pos); + } + } + const double* BlackoilPvtProperties::surfaceDensities() const { return densities_; diff --git a/opm/core/props/pvt/BlackoilPvtProperties.hpp b/opm/core/props/pvt/BlackoilPvtProperties.hpp index 4ab60921..9b466b2d 100644 --- a/opm/core/props/pvt/BlackoilPvtProperties.hpp +++ b/opm/core/props/pvt/BlackoilPvtProperties.hpp @@ -20,11 +20,12 @@ #ifndef OPM_BLACKOILPVTPROPERTIES_HEADER_INCLUDED #define OPM_BLACKOILPVTPROPERTIES_HEADER_INCLUDED - - #include #include #include + +#include + #include #include @@ -55,6 +56,11 @@ namespace Opm /// data without fitting a spline. void init(const EclipseGridParser& deck, const int samples); + + /// Initialize from deck. + /// \param deck An input deck from the opm-parser module. + void init(Opm::DeckConstPtr deck, int samples); + /// \return Object describing the active phases. PhaseUsage phaseUsage() const; diff --git a/opm/core/props/pvt/SinglePvtConstCompr.hpp b/opm/core/props/pvt/SinglePvtConstCompr.hpp index e97a7242..d6c4d088 100644 --- a/opm/core/props/pvt/SinglePvtConstCompr.hpp +++ b/opm/core/props/pvt/SinglePvtConstCompr.hpp @@ -20,9 +20,12 @@ #ifndef OPM_SINGLEPVTCONSTCOMPR_HEADER_INCLUDED #define OPM_SINGLEPVTCONSTCOMPR_HEADER_INCLUDED - #include #include + +#include +#include + #include #include @@ -55,6 +58,32 @@ namespace Opm visc_comp_ = pvtw[region_number][4]; } + SinglePvtConstCompr(const Opm::PvtwTable &pvtwTable) + { + if (pvtwTable.numRows() != 1) + OPM_THROW(std::runtime_error, + "The table specified by the PVTW keyword is required" + "to exhibit exactly one row."); + ref_press_ = pvtwTable.getPressureColumn()[0]; + ref_B_ = pvtwTable.getFormationFactorColumn()[0]; + comp_ = pvtwTable.getCompressibilityColumn()[0]; + viscosity_ = pvtwTable.getViscosityColumn()[0]; + visc_comp_ = pvtwTable.getViscosibilityColumn()[0]; + } + + SinglePvtConstCompr(const Opm::PvdcoTable &pvdcoTable) + { + if (pvdcoTable.numRows() != 1) + OPM_THROW(std::runtime_error, + "The table specified by the PVDCO keyword is required" + "to exhibit exactly one row."); + ref_press_ = pvdcoTable.getPressureColumn()[0]; + ref_B_ = pvdcoTable.getFormationFactorColumn()[0]; + comp_ = pvdcoTable.getCompressibilityColumn()[0]; + viscosity_ = pvdcoTable.getViscosityColumn()[0]; + visc_comp_ = pvdcoTable.getViscosibilityColumn()[0]; + } + SinglePvtConstCompr(double visc) : ref_press_(0.0), ref_B_(1.0), diff --git a/opm/core/props/pvt/SinglePvtDead.cpp b/opm/core/props/pvt/SinglePvtDead.cpp index 2db23dbe..fb11d2de 100644 --- a/opm/core/props/pvt/SinglePvtDead.cpp +++ b/opm/core/props/pvt/SinglePvtDead.cpp @@ -62,6 +62,30 @@ namespace Opm // << "\n\nvisc\n\n" << viscosity_ << std::endl; } + /// Constructor + SinglePvtDead::SinglePvtDead(const Opm::PvdoTable& pvdoTable) + { + // Copy data + const std::vector& press = pvdoTable.getPressureColumn(); + const std::vector& b = pvdoTable.getFormationFactorColumn(); + const std::vector& visc = pvdoTable.getViscosityColumn(); + + const int sz = b.size(); + std::vector bInv(sz); + for (int i = 0; i < sz; ++i) { + bInv[i] = 1.0 / b[i]; + } + b_ = NonuniformTableLinear(press, bInv); + viscosity_ = NonuniformTableLinear(press, visc); + + // Dumping the created tables. +// static int count = 0; +// std::ofstream os((std::string("dump-") + boost::lexical_cast(count++)).c_str()); +// os.precision(15); +// os << "1/B\n\n" << one_over_B_ +// << "\n\nvisc\n\n" << viscosity_ << std::endl; + } + // Destructor SinglePvtDead::~SinglePvtDead() { diff --git a/opm/core/props/pvt/SinglePvtDead.hpp b/opm/core/props/pvt/SinglePvtDead.hpp index 2c40a951..703da0b8 100644 --- a/opm/core/props/pvt/SinglePvtDead.hpp +++ b/opm/core/props/pvt/SinglePvtDead.hpp @@ -23,6 +23,9 @@ #include #include + +#include + #include namespace Opm @@ -40,6 +43,7 @@ namespace Opm public: typedef std::vector > > table_t; SinglePvtDead(const table_t& pvd_table); + SinglePvtDead(const Opm::PvdoTable &pvdoTable); virtual ~SinglePvtDead(); /// Viscosity as a function of p and z. diff --git a/opm/core/props/pvt/SinglePvtDeadSpline.cpp b/opm/core/props/pvt/SinglePvtDeadSpline.cpp index a0873e0b..0dc3b560 100644 --- a/opm/core/props/pvt/SinglePvtDeadSpline.cpp +++ b/opm/core/props/pvt/SinglePvtDeadSpline.cpp @@ -20,6 +20,9 @@ #include "config.h" #include #include + +#include + #include // Extra includes for debug dumping of tables. @@ -29,7 +32,6 @@ namespace Opm { - //------------------------------------------------------------------------ // Member functions //------------------------------------------------------------------------- @@ -54,13 +56,42 @@ namespace Opm } buildUniformMonotoneTable(press, B_inv, samples, b_); buildUniformMonotoneTable(press, visc, samples, viscosity_); + } - // Dumping the created tables. -// static int count = 0; -// std::ofstream os((std::string("dump-") + boost::lexical_cast(count++)).c_str()); -// os.precision(15); -// os << "1/B\n\n" << one_over_B_ -// << "\n\nvisc\n\n" << viscosity_ << std::endl; + SinglePvtDeadSpline::SinglePvtDeadSpline(const Opm::PvdoTable &pvdoTable, int samples) + { + int numRows = pvdoTable.numRows(); + + // Copy data + const std::vector &press = pvdoTable.getPressureColumn(); + const std::vector &B = pvdoTable.getFormationFactorColumn(); + const std::vector &visc = pvdoTable.getViscosityColumn(); + + std::vector B_inv(numRows); + for (int i = 0; i < numRows; ++i) { + B_inv[i] = 1.0 / B[i]; + } + + buildUniformMonotoneTable(press, B_inv, samples, b_); + buildUniformMonotoneTable(press, visc, samples, viscosity_); + } + + SinglePvtDeadSpline::SinglePvtDeadSpline(const Opm::PvdgTable &pvdgTable, int samples) + { + int numRows = pvdgTable.numRows(); + + // Copy data + const std::vector &press = pvdgTable.getPressureColumn(); + const std::vector &B = pvdgTable.getFormationFactorColumn(); + const std::vector &visc = pvdgTable.getViscosityColumn(); + + std::vector B_inv(numRows); + for (int i = 0; i < numRows; ++i) { + B_inv[i] = 1.0 / B[i]; + } + + buildUniformMonotoneTable(press, B_inv, samples, b_); + buildUniformMonotoneTable(press, visc, samples, viscosity_); } // Destructor diff --git a/opm/core/props/pvt/SinglePvtDeadSpline.hpp b/opm/core/props/pvt/SinglePvtDeadSpline.hpp index 337c8a80..65bd2561 100644 --- a/opm/core/props/pvt/SinglePvtDeadSpline.hpp +++ b/opm/core/props/pvt/SinglePvtDeadSpline.hpp @@ -20,9 +20,12 @@ #ifndef OPM_SINGLEPVTDEADSPLINE_HEADER_INCLUDED #define OPM_SINGLEPVTDEADSPLINE_HEADER_INCLUDED - #include #include + +#include +#include + #include namespace Opm @@ -41,6 +44,8 @@ namespace Opm typedef std::vector > > table_t; SinglePvtDeadSpline(const table_t& pvd_table, const int samples); + SinglePvtDeadSpline(const Opm::PvdoTable &pvdoTable, int samples); + SinglePvtDeadSpline(const Opm::PvdgTable &pvdgTable, int samples); virtual ~SinglePvtDeadSpline(); /// Viscosity as a function of p and z. diff --git a/opm/core/props/pvt/SinglePvtLiveGas.cpp b/opm/core/props/pvt/SinglePvtLiveGas.cpp index 6a43f475..43debe7b 100644 --- a/opm/core/props/pvt/SinglePvtLiveGas.cpp +++ b/opm/core/props/pvt/SinglePvtLiveGas.cpp @@ -32,8 +32,8 @@ #include #include #include -#include +#include namespace Opm { @@ -41,6 +41,7 @@ namespace Opm using Opm::linearInterpolation; using Opm::linearInterpolationDerivative; + //------------------------------------------------------------------------ // Member functions //------------------------------------------------------------------------- @@ -81,6 +82,27 @@ namespace Opm } } + SinglePvtLiveGas::SinglePvtLiveGas(const Opm::PvtgTable& pvtgTable) + { + // GAS, PVTG + saturated_gas_table_.resize(4); + saturated_gas_table_[0] = pvtgTable.getOuterTable()->getPressureColumn(); + saturated_gas_table_[1] = pvtgTable.getOuterTable()->getGasFormationFactorColumn(); + saturated_gas_table_[2] = pvtgTable.getOuterTable()->getGasViscosityColumn(); + saturated_gas_table_[3] = pvtgTable.getOuterTable()->getOilSolubilityColumn(); + + int sz = pvtgTable.getOuterTable()->numRows(); + undersat_gas_tables_.resize(sz); + for (int i=0; igetGasViscosityColumn(); + } + } + // Destructor SinglePvtLiveGas::~SinglePvtLiveGas() { diff --git a/opm/core/props/pvt/SinglePvtLiveGas.hpp b/opm/core/props/pvt/SinglePvtLiveGas.hpp index de595393..c09f8ffc 100644 --- a/opm/core/props/pvt/SinglePvtLiveGas.hpp +++ b/opm/core/props/pvt/SinglePvtLiveGas.hpp @@ -21,6 +21,9 @@ #define OPM_SINGLEPVTLIVEGAS_HEADER_INCLUDED #include + +#include + #include namespace Opm @@ -38,6 +41,7 @@ namespace Opm typedef std::vector > > table_t; SinglePvtLiveGas(const table_t& pvto); + SinglePvtLiveGas(const Opm::PvtgTable& pvtgTable); virtual ~SinglePvtLiveGas(); /// Viscosity as a function of p and z. diff --git a/opm/core/props/pvt/SinglePvtLiveOil.cpp b/opm/core/props/pvt/SinglePvtLiveOil.cpp index d18252ee..66f1036e 100644 --- a/opm/core/props/pvt/SinglePvtLiveOil.cpp +++ b/opm/core/props/pvt/SinglePvtLiveOil.cpp @@ -21,8 +21,8 @@ #include #include #include -#include +#include namespace Opm { @@ -31,7 +31,6 @@ namespace Opm using Opm::linearInterpolationDerivative; using Opm::tableIndex; - //------------------------------------------------------------------------ // Member functions //------------------------------------------------------------------------- @@ -102,7 +101,73 @@ namespace Opm undersat_oil_tables_[i][2].push_back(mu); } } + } + SinglePvtLiveOil::SinglePvtLiveOil(const Opm::PvtoTable &pvtoTable) + { + const auto saturatedPvto = pvtoTable.getOuterTable(); + + // OIL, PVTO + saturated_oil_table_.resize(4); + const int sz = saturatedPvto->numRows(); + for (int k=0; k<4; ++k) { + saturated_oil_table_[k].resize(sz); + } + for (int i=0; igetPressureColumn()[i]; // p + saturated_oil_table_[1][i] = 1.0/saturatedPvto->getOilFormationFactorColumn()[i]; // 1/Bo + saturated_oil_table_[2][i] = saturatedPvto->getOilViscosityColumn()[i]; // mu_o + saturated_oil_table_[3][i] = saturatedPvto->getGasSolubilityColumn()[i]; // Rs + } + + undersat_oil_tables_.resize(sz); + for (int i=0; inumRows(); + undersat_oil_tables_[i][0].resize(tsize); + undersat_oil_tables_[i][1].resize(tsize); + undersat_oil_tables_[i][2].resize(tsize); + for (int j=0; jgetPressureColumn()[j]; // p + undersat_oil_tables_[i][1][j] = 1.0/undersaturatedPvto->getOilFormationFactorColumn()[j]; // 1/Bo + undersat_oil_tables_[i][2][j] = undersaturatedPvto->getOilViscosityColumn()[j]; // mu_o + } + } + + // Complete undersaturated tables by extrapolating from existing data + // as is done in Eclipse and Mrst + int iNext = -1; + for (int i=0; i 1) { + continue; + } + // Look ahead for next record containing undersaturated data + if (iNext < i) { + iNext = i+1; + while (iNext > >::size_type sz_t; + for (sz_t j=1; j + +#include + #include namespace Opm @@ -39,6 +41,7 @@ namespace Opm typedef std::vector > > table_t; SinglePvtLiveOil(const table_t& pvto); + SinglePvtLiveOil(const Opm::PvtoTable &pvtoTable); virtual ~SinglePvtLiveOil(); /// Viscosity as a function of p and z.