From f360562aee461116b8227bbf0dae029e0f9f3bf0 Mon Sep 17 00:00:00 2001 From: Andreas Lauser Date: Thu, 17 Apr 2014 12:04:31 +0200 Subject: [PATCH] remove EclipseGridParser compatibility methods from all classes --- opm/core/props/BlackoilPropertiesFromDeck.cpp | 18 - opm/core/props/BlackoilPropertiesFromDeck.hpp | 82 +--- opm/core/props/IncompPropertiesFromDeck.cpp | 12 - opm/core/props/IncompPropertiesFromDeck.hpp | 8 - opm/core/props/phaseUsageFromDeck.hpp | 40 +- .../props/pvt/PvtPropertiesIncompFromDeck.cpp | 55 --- .../props/pvt/PvtPropertiesIncompFromDeck.hpp | 4 - opm/core/props/rock/RockCompressibility.cpp | 28 -- opm/core/props/rock/RockCompressibility.hpp | 5 - opm/core/props/rock/RockFromDeck.cpp | 270 +--------- opm/core/props/rock/RockFromDeck.hpp | 28 -- .../props/satfunc/SaturationPropsFromDeck.hpp | 53 +- .../satfunc/SaturationPropsFromDeck_impl.hpp | 463 +----------------- opm/core/simulator/initState.hpp | 49 +- opm/core/simulator/initStateEquil.hpp | 254 +--------- opm/core/simulator/initStateEquil_impl.hpp | 22 +- opm/core/simulator/initState_impl.hpp | 179 ------- opm/core/wells/WellCollection.cpp | 44 -- opm/core/wells/WellCollection.hpp | 14 - opm/core/wells/WellsGroup.cpp | 95 ---- opm/core/wells/WellsGroup.hpp | 7 - opm/core/wells/WellsManager.cpp | 1 - opm/core/wells/WellsManager.hpp | 2 - tests/test_equil.cpp | 12 +- 24 files changed, 27 insertions(+), 1718 deletions(-) diff --git a/opm/core/props/BlackoilPropertiesFromDeck.cpp b/opm/core/props/BlackoilPropertiesFromDeck.cpp index b672cc30b..1f4a0bb71 100644 --- a/opm/core/props/BlackoilPropertiesFromDeck.cpp +++ b/opm/core/props/BlackoilPropertiesFromDeck.cpp @@ -23,15 +23,6 @@ namespace Opm { - BlackoilPropertiesFromDeck::BlackoilPropertiesFromDeck(const EclipseGridParser& deck, - const UnstructuredGrid& grid, - bool init_rock) - { - init(deck, grid.number_of_cells, grid.global_cell, grid.cartdims, grid.cell_centroids, - grid.dimensions, init_rock); - } - - BlackoilPropertiesFromDeck::BlackoilPropertiesFromDeck(Opm::DeckConstPtr newParserDeck, const UnstructuredGrid& grid, bool init_rock) @@ -40,15 +31,6 @@ namespace Opm grid.cell_centroids, grid.dimensions, init_rock); } - BlackoilPropertiesFromDeck::BlackoilPropertiesFromDeck(const EclipseGridParser& deck, - const UnstructuredGrid& grid, - const parameter::ParameterGroup& param, - bool init_rock) - { - init(deck, grid.number_of_cells, grid.global_cell, grid.cartdims, grid.cell_centroids, - grid.dimensions, param, init_rock); - } - BlackoilPropertiesFromDeck::BlackoilPropertiesFromDeck(Opm::DeckConstPtr newParserDeck, const UnstructuredGrid& grid, const parameter::ParameterGroup& param, diff --git a/opm/core/props/BlackoilPropertiesFromDeck.hpp b/opm/core/props/BlackoilPropertiesFromDeck.hpp index e521e28aa..77d7b6d1b 100644 --- a/opm/core/props/BlackoilPropertiesFromDeck.hpp +++ b/opm/core/props/BlackoilPropertiesFromDeck.hpp @@ -25,7 +25,6 @@ #include #include #include -#include #include #include @@ -42,14 +41,6 @@ namespace Opm class BlackoilPropertiesFromDeck : public BlackoilPropertiesInterface { public: - /// Initialize from deck and grid. - /// \param[in] deck Deck input parser - /// \param[in] grid Grid to which property object applies, needed for the - /// mapping from cell indices (typically from a processed grid) - /// to logical cartesian indices consistent with the deck. - BlackoilPropertiesFromDeck(const EclipseGridParser &deck, - const UnstructuredGrid& grid, bool init_rock=true ); - /// Initialize from deck and grid. /// \param[in] deck Deck input parser /// \param[in] grid Grid to which property object applies, needed for the @@ -58,22 +49,6 @@ namespace Opm BlackoilPropertiesFromDeck(Opm::DeckConstPtr newParserDeck, const UnstructuredGrid& grid, bool init_rock=true ); - /// Initialize from deck, grid and parameters. - /// \param[in] deck Deck input parser - /// \param[in] grid Grid to which property object applies, needed for the - /// mapping from cell indices (typically from a processed grid) - /// to logical cartesian indices consistent with the deck. - /// \param[in] param Parameters. Accepted parameters include: - /// pvt_tab_size (200) number of uniform sample points for dead-oil pvt tables. - /// sat_tab_size (200) number of uniform sample points for saturation tables. - /// threephase_model("simple") three-phase relperm model (accepts "simple" and "stone2"). - /// For both size parameters, a 0 or negative value indicates that no spline fitting is to - /// be done, and the input fluid data used directly for linear interpolation. - BlackoilPropertiesFromDeck(const EclipseGridParser& deck, - const UnstructuredGrid& grid, - const parameter::ParameterGroup& param, - bool init_rock=true); - /// Initialize from deck, grid and parameters. /// \param[in] deck Deck input parser /// \param[in] grid Grid to which property object applies, needed for the @@ -90,41 +65,22 @@ namespace Opm const parameter::ParameterGroup& param, bool init_rock=true); - template - BlackoilPropertiesFromDeck(const EclipseGridParser& deck, - int number_of_cells, - const int* global_cell, - const int* cart_dims, - T begin_cell_centroids, - int dimension, - bool init_rock=true); - - template - BlackoilPropertiesFromDeck(const EclipseGridParser& deck, - int number_of_cells, - const int* global_cell, - const int* cart_dims, - T begin_cell_centroids, - int dimension, - const parameter::ParameterGroup& param, - bool init_rock=true); - - template + template BlackoilPropertiesFromDeck(Opm::DeckConstPtr newParserDeck, int number_of_cells, const int* global_cell, const int* cart_dims, - T begin_cell_centroids, + const CentroidIterator& begin_cell_centroids, int dimension, bool init_rock=true); - template + template BlackoilPropertiesFromDeck(Opm::DeckConstPtr newParserDeck, int number_of_cells, const int* global_cell, const int* cart_dims, - T begin_cell_centroids, + const CentroidIterator& begin_cell_centroids, int dimension, const parameter::ParameterGroup& param, bool init_rock=true); @@ -250,38 +206,20 @@ namespace Opm double* smax) const; private: - template - void init(const EclipseGridParser& deck, + template + void init(Opm::DeckConstPtr deck, int number_of_cells, const int* global_cell, const int* cart_dims, - T begin_cell_centroids, + const CentroidIterator& begin_cell_centroids, int dimension, bool init_rock); - template - void init(const EclipseGridParser& deck, + template + void init(Opm::DeckConstPtr deck, int number_of_cells, const int* global_cell, const int* cart_dims, - T begin_cell_centroids, - int dimension, - const parameter::ParameterGroup& param, - bool init_rock); - - template - void init(Opm::DeckConstPtr newParserDeck, - int number_of_cells, - const int* global_cell, - const int* cart_dims, - T begin_cell_centroids, - int dimension, - bool init_rock); - template - void init(Opm::DeckConstPtr newParserDeck, - int number_of_cells, - const int* global_cell, - const int* cart_dims, - T begin_cell_centroids, + const CentroidIterator& begin_cell_centroids, int dimension, const parameter::ParameterGroup& param, bool init_rock); diff --git a/opm/core/props/IncompPropertiesFromDeck.cpp b/opm/core/props/IncompPropertiesFromDeck.cpp index 04306b285..933e8d527 100644 --- a/opm/core/props/IncompPropertiesFromDeck.cpp +++ b/opm/core/props/IncompPropertiesFromDeck.cpp @@ -26,18 +26,6 @@ namespace Opm { - IncompPropertiesFromDeck::IncompPropertiesFromDeck(const EclipseGridParser& deck, - const UnstructuredGrid& grid) - { - rock_.init(deck, grid); - pvt_.init(deck); - satprops_.init(deck, grid, 200); - if (pvt_.numPhases() != satprops_.numPhases()) { - OPM_THROW(std::runtime_error, "IncompPropertiesFromDeck::IncompPropertiesFromDeck() - Inconsistent number of phases in pvt data (" - << pvt_.numPhases() << ") and saturation-dependent function data (" << satprops_.numPhases() << ")."); - } - } - IncompPropertiesFromDeck::IncompPropertiesFromDeck(Opm::DeckConstPtr newParserDeck, const UnstructuredGrid& grid) { diff --git a/opm/core/props/IncompPropertiesFromDeck.hpp b/opm/core/props/IncompPropertiesFromDeck.hpp index 83718aac9..696990171 100644 --- a/opm/core/props/IncompPropertiesFromDeck.hpp +++ b/opm/core/props/IncompPropertiesFromDeck.hpp @@ -45,14 +45,6 @@ namespace Opm class IncompPropertiesFromDeck : public IncompPropertiesInterface { public: - /// Initialize from deck and grid. - /// \param deck Deck input parser - /// \param grid Grid to which property object applies, needed for the - /// mapping from cell indices (typically from a processed grid) - /// to logical cartesian indices consistent with the deck. - IncompPropertiesFromDeck(const EclipseGridParser& deck, - const UnstructuredGrid& grid); - /// Initialize from deck and grid. /// \param deck Deck input parser /// \param grid Grid to which property object applies, needed for the diff --git a/opm/core/props/phaseUsageFromDeck.hpp b/opm/core/props/phaseUsageFromDeck.hpp index 6ee66075e..f69af0e0a 100644 --- a/opm/core/props/phaseUsageFromDeck.hpp +++ b/opm/core/props/phaseUsageFromDeck.hpp @@ -20,9 +20,8 @@ #ifndef OPM_PHASEUSAGEFROMDECK_HEADER_INCLUDED #define OPM_PHASEUSAGEFROMDECK_HEADER_INCLUDED - -#include #include +#include #include #include @@ -68,43 +67,6 @@ namespace Opm return pu; } - /// Looks at presence of WATER, OIL and GAS keywords in deck - /// to determine active phases. - inline PhaseUsage phaseUsageFromDeck(const EclipseGridParser& deck) - { - PhaseUsage pu; - std::fill(pu.phase_used, pu.phase_used + BlackoilPhases::MaxNumPhases, 0); - - // Discover phase usage. - if (deck.hasField("WATER")) { - pu.phase_used[BlackoilPhases::Aqua] = 1; - } - if (deck.hasField("OIL")) { - pu.phase_used[BlackoilPhases::Liquid] = 1; - } - if (deck.hasField("GAS")) { - pu.phase_used[BlackoilPhases::Vapour] = 1; - } - pu.num_phases = 0; - for (int i = 0; i < BlackoilPhases::MaxNumPhases; ++i) { - pu.phase_pos[i] = pu.num_phases; - pu.num_phases += pu.phase_used[i]; - } - - // Only 2 or 3 phase systems handled. - if (pu.num_phases < 2 || pu.num_phases > 3) { - OPM_THROW(std::runtime_error, "Cannot handle cases with " << pu.num_phases << " phases."); - } - - // We need oil systems, since we do not support the keywords needed for - // water-gas systems. - if (!pu.phase_used[BlackoilPhases::Liquid]) { - OPM_THROW(std::runtime_error, "Cannot handle cases with no OIL, i.e. water-gas systems."); - } - - return pu; - } - /// Looks at presence of WATER, OIL and GAS keywords in deck /// to determine active phases. inline PhaseUsage phaseUsageFromDeck(Opm::DeckConstPtr newParserDeck) diff --git a/opm/core/props/pvt/PvtPropertiesIncompFromDeck.cpp b/opm/core/props/pvt/PvtPropertiesIncompFromDeck.cpp index 330941266..c7b974577 100644 --- a/opm/core/props/pvt/PvtPropertiesIncompFromDeck.cpp +++ b/opm/core/props/pvt/PvtPropertiesIncompFromDeck.cpp @@ -21,7 +21,6 @@ #include "config.h" #include #include -#include #include #include #include @@ -34,60 +33,6 @@ namespace Opm { } - - void PvtPropertiesIncompFromDeck::init(const EclipseGridParser& deck) - { - // If we need multiple regions, this class and the SinglePvt* classes must change. - int region_number = 0; - - PhaseUsage phase_usage = phaseUsageFromDeck(deck); - if (phase_usage.phase_used[PhaseUsage::Vapour] || - !phase_usage.phase_used[PhaseUsage::Aqua] || - !phase_usage.phase_used[PhaseUsage::Liquid]) { - OPM_THROW(std::runtime_error, "PvtPropertiesIncompFromDeck::init() -- must have gas and oil phases (only) in deck input.\n"); - } - - // Surface densities. Accounting for different orders in eclipse and our code. - if (deck.hasField("DENSITY")) { - const std::vector& d = deck.getDENSITY().densities_[region_number]; - enum { ECL_oil = 0, ECL_water = 1, ECL_gas = 2 }; - surface_density_[phase_usage.phase_pos[PhaseUsage::Aqua]] = d[ECL_water]; - surface_density_[phase_usage.phase_pos[PhaseUsage::Liquid]] = d[ECL_oil]; - } else { - OPM_THROW(std::runtime_error, "Input is missing DENSITY\n"); - } - - // Make reservoir densities the same as surface densities initially. - // We will modify them with formation volume factors if found. - reservoir_density_ = surface_density_; - - // Water viscosity. - if (deck.hasField("PVTW")) { - const std::vector& pvtw = deck.getPVTW().pvtw_[region_number]; - if (pvtw[2] != 0.0 || pvtw[4] != 0.0) { - OPM_MESSAGE("Compressibility effects in PVTW are ignored."); - } - reservoir_density_[phase_usage.phase_pos[PhaseUsage::Aqua]] /= pvtw[1]; - viscosity_[phase_usage.phase_pos[PhaseUsage::Aqua]] = pvtw[3]; - } else { - // Eclipse 100 default. - // viscosity_[phase_usage.phase_pos[PhaseUsage::Aqua]] = 0.5*Opm::prefix::centi*Opm::unit::Poise; - OPM_THROW(std::runtime_error, "Input is missing PVTW\n"); - } - - // Oil viscosity. - if (deck.hasField("PVCDO")) { - const std::vector& pvcdo = deck.getPVCDO().pvcdo_[region_number]; - if (pvcdo[2] != 0.0 || pvcdo[4] != 0.0) { - OPM_MESSAGE("Compressibility effects in PVCDO are ignored."); - } - reservoir_density_[phase_usage.phase_pos[PhaseUsage::Liquid]] /= pvcdo[1]; - viscosity_[phase_usage.phase_pos[PhaseUsage::Liquid]] = pvcdo[3]; - } else { - OPM_THROW(std::runtime_error, "Input is missing PVCDO\n"); - } - } - void PvtPropertiesIncompFromDeck::init(Opm::DeckConstPtr newParserDeck ) { // If we need multiple regions, this class and the SinglePvt* classes must change. diff --git a/opm/core/props/pvt/PvtPropertiesIncompFromDeck.hpp b/opm/core/props/pvt/PvtPropertiesIncompFromDeck.hpp index 4ffe9864a..6d8d3527d 100644 --- a/opm/core/props/pvt/PvtPropertiesIncompFromDeck.hpp +++ b/opm/core/props/pvt/PvtPropertiesIncompFromDeck.hpp @@ -20,7 +20,6 @@ #ifndef OPM_PVTPROPERTIESINCOMPFROMDECK_HEADER_INCLUDED #define OPM_PVTPROPERTIESINCOMPFROMDECK_HEADER_INCLUDED -#include #include #include @@ -38,9 +37,6 @@ namespace Opm /// Default constructor. PvtPropertiesIncompFromDeck(); - /// Initialize from deck. - void init(const EclipseGridParser& deck); - /// Initialize from deck. void init(Opm::DeckConstPtr newParserDeck); diff --git a/opm/core/props/rock/RockCompressibility.cpp b/opm/core/props/rock/RockCompressibility.cpp index d1729533c..5ccf49183 100644 --- a/opm/core/props/rock/RockCompressibility.cpp +++ b/opm/core/props/rock/RockCompressibility.cpp @@ -19,7 +19,6 @@ #include "config.h" #include -#include #include #include #include @@ -41,33 +40,6 @@ namespace Opm rock_comp_ = param.getDefault("rock_compressibility", 0.0)/unit::barsa; } - RockCompressibility::RockCompressibility(const EclipseGridParser& deck) - : pref_(0.0), - rock_comp_(0.0) - { - if (deck.hasField("ROCKTAB")) { - const table_t& rt = deck.getROCKTAB().rocktab_; - if (rt.size() != 1) { - OPM_THROW(std::runtime_error, "Can only handle a single region in ROCKTAB."); - } - const int n = rt[0][0].size(); - p_.resize(n); - poromult_.resize(n); - transmult_.resize(n); - for (int i = 0; i < n; ++i) { - p_[i] = rt[0][0][i]; - poromult_[i] = rt[0][1][i]; - transmult_[i] = rt[0][2][i]; - } - } 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; - } - } - RockCompressibility::RockCompressibility(Opm::DeckConstPtr newParserDeck) : pref_(0.0), rock_comp_(0.0) diff --git a/opm/core/props/rock/RockCompressibility.hpp b/opm/core/props/rock/RockCompressibility.hpp index bd375e1ca..071c144e7 100644 --- a/opm/core/props/rock/RockCompressibility.hpp +++ b/opm/core/props/rock/RockCompressibility.hpp @@ -27,16 +27,11 @@ namespace Opm { - class EclipseGridParser; namespace parameter { class ParameterGroup; } class RockCompressibility { public: - /// Construct from input deck. - /// Looks for the keywords ROCK and ROCKTAB. - RockCompressibility(const EclipseGridParser& deck); - /// Construct from input deck. /// Looks for the keywords ROCK and ROCKTAB. RockCompressibility(Opm::DeckConstPtr newParserDeck); diff --git a/opm/core/props/rock/RockFromDeck.cpp b/opm/core/props/rock/RockFromDeck.cpp index e34d784da..3aa9bb077 100644 --- a/opm/core/props/rock/RockFromDeck.cpp +++ b/opm/core/props/rock/RockFromDeck.cpp @@ -21,6 +21,7 @@ #include "config.h" #include #include +#include #include @@ -34,12 +35,8 @@ namespace Opm { enum PermeabilityKind { ScalarPerm, DiagonalPerm, TensorPerm, None, Invalid }; - PermeabilityKind classifyPermeability(const EclipseGridParser& parser); void setScalarPermIfNeeded(std::array& kmap, int i, int j, int k); - PermeabilityKind fillTensor(const EclipseGridParser& parser, - std::vector*>& tensor, - std::array& kmap); PermeabilityKind fillTensor(Opm::DeckConstPtr newParserDeck, std::vector*>& tensor, std::array& kmap); @@ -56,34 +53,6 @@ namespace Opm { } - - /// Initialize from deck and cell mapping. - /// \param deck Deck input parser - /// \param grid grid to which property object applies, needed for the - /// mapping from cell indices (typically from a processed grid) - /// to logical cartesian indices consistent with the deck. - void RockFromDeck::init(const EclipseGridParser& deck, - const UnstructuredGrid& grid) - { - init(deck, grid.number_of_cells, grid.global_cell, grid.cartdims); - } - - /// Initialize from deck and cell mapping. - /// \param deck Deck input parser - /// \param number_of_cells The number of cells in the grid. - /// \param global_cell The mapping fom local to global cell indices. - /// global_cell[i] is the corresponding global index of i. - /// \param cart_dims The size of the underlying cartesian grid. - void RockFromDeck::init(const EclipseGridParser& deck, - int number_of_cells, const int* global_cell, - const int* cart_dims) - { - assignPorosity(deck, number_of_cells, global_cell); - permfield_valid_.assign(number_of_cells, false); - const double perm_threshold = 0.0; // Maybe turn into parameter? - assignPermeability(deck, number_of_cells, global_cell, cart_dims, perm_threshold); - } - void RockFromDeck::init(Opm::DeckConstPtr newParserDeck, int number_of_cells, const int* global_cell, const int* cart_dims) @@ -95,20 +64,6 @@ namespace Opm perm_threshold); } - - void RockFromDeck::assignPorosity(const EclipseGridParser& parser, - int number_of_cells, const int* global_cell) - { - porosity_.assign(number_of_cells, 1.0); - if (parser.hasField("PORO")) { - const std::vector& poro = parser.getFloatingPointValue("PORO"); - for (int c = 0; c < int(porosity_.size()); ++c) { - const int deck_pos = (global_cell == NULL) ? c : global_cell[c]; - porosity_[c] = poro[deck_pos]; - } - } - } - void RockFromDeck::assignPorosity(Opm::DeckConstPtr newParserDeck, int number_of_cells, const int* global_cell) { @@ -124,61 +79,6 @@ namespace Opm } } - - void RockFromDeck::assignPermeability(const EclipseGridParser& parser, - int number_of_cells, - const int* global_cell, - const int* cartdims, - double perm_threshold) - { - const int dim = 3; - const int num_global_cells = cartdims[0]*cartdims[1]*cartdims[2]; - - assert(num_global_cells > 0); - - permeability_.assign(dim * dim * number_of_cells, 0.0); - - std::vector*> tensor; - tensor.reserve(10); - - const std::vector zero(num_global_cells, 0.0); - tensor.push_back(&zero); - - std::array kmap; - PermeabilityKind pkind = fillTensor(parser, tensor, kmap); - if (pkind == Invalid) { - OPM_THROW(std::runtime_error, "Invalid permeability field."); - } - - // Assign permeability values only if such values are - // given in the input deck represented by 'parser'. In - // other words: Don't set any (arbitrary) default values. - // It is infinitely better to experience a reproducible - // crash than subtle errors resulting from a (poorly - // chosen) default value... - // - if (tensor.size() > 1) { - int off = 0; - - for (int c = 0; c < number_of_cells; ++c, off += dim*dim) { - // SharedPermTensor K(dim, dim, &permeability_[off]); - int kix = 0; - const int glob = (global_cell == NULL) ? c : global_cell[c]; - - for (int i = 0; i < dim; ++i) { - for (int j = 0; j < dim; ++j, ++kix) { - // K(i,j) = (*tensor[kmap[kix]])[glob]; - permeability_[off + kix] = (*tensor[kmap[kix]])[glob]; - } - // K(i,i) = std::max(K(i,i), perm_threshold); - permeability_[off + 3*i + i] = std::max(permeability_[off + 3*i + i], perm_threshold); - } - - permfield_valid_[c] = std::vector::value_type(1); - } - } - } - void RockFromDeck::assignPermeability(Opm::DeckConstPtr newParserDeck, int number_of_cells, const int* global_cell, @@ -236,73 +136,6 @@ namespace Opm } namespace { - - /// @brief - /// Classify and verify a given permeability specification - /// from a structural point of view. In particular, we - /// verify that there are no off-diagonal permeability - /// components such as @f$k_{xy}@f$ unless the - /// corresponding diagonal components are known as well. - /// - /// @param parser [in] - /// An Eclipse data parser capable of answering which - /// permeability components are present in a given input - /// deck. - /// - /// @return - /// An enum value with the following possible values: - /// ScalarPerm only one component was given. - /// DiagonalPerm more than one component given. - /// TensorPerm at least one cross-component given. - /// None no components given. - /// Invalid invalid set of components given. - PermeabilityKind classifyPermeability(const EclipseGridParser& parser) - { - const bool xx = parser.hasField("PERMX" ); - const bool xy = parser.hasField("PERMXY"); - const bool xz = parser.hasField("PERMXZ"); - - const bool yx = parser.hasField("PERMYX"); - const bool yy = parser.hasField("PERMY" ); - const bool yz = parser.hasField("PERMYZ"); - - const bool zx = parser.hasField("PERMZX"); - const bool zy = parser.hasField("PERMZY"); - const bool zz = parser.hasField("PERMZ" ); - - int num_cross_comp = xy + xz + yx + yz + zx + zy; - int num_comp = xx + yy + zz + num_cross_comp; - PermeabilityKind retval = None; - if (num_cross_comp > 0) { - retval = TensorPerm; - } else { - if (num_comp == 1) { - retval = ScalarPerm; - } else if (num_comp >= 2) { - retval = DiagonalPerm; - } - } - - bool ok = true; - if (num_comp > 0) { - // At least one tensor component specified on input. - // Verify that any remaining components are OK from a - // structural point of view. In particular, there - // must not be any cross-components (e.g., k_{xy}) - // unless the corresponding diagonal component (e.g., - // k_{xx}) is present as well... - // - ok = xx || !(xy || xz || yx || zx) ; - ok = ok && (yy || !(yx || yz || xy || zy)); - ok = ok && (zz || !(zx || zy || xz || yz)); - } - if (!ok) { - retval = Invalid; - } - - return retval; - } - /// @brief /// Classify and verify a given permeability specification /// from a structural point of view. In particular, we @@ -397,107 +230,6 @@ namespace Opm if (kmap[k] == 0) { kmap[k] = kmap[i]; } } - - /// @brief - /// Extract pointers to appropriate tensor components from - /// input deck. The permeability tensor is, generally, - /// @code - /// [ kxx kxy kxz ] - /// K = [ kyx kyy kyz ] - /// [ kzx kzy kzz ] - /// @endcode - /// We store these values in a linear array using natural - /// ordering with the column index cycling the most rapidly. - /// In particular we use the representation - /// @code - /// [ 0 1 2 3 4 5 6 7 8 ] - /// K = [ kxx, kxy, kxz, kyx, kyy, kyz, kzx, kzy, kzz ] - /// @endcode - /// Moreover, we explicitly enforce symmetric tensors by - /// assigning - /// @code - /// 3 1 6 2 7 5 - /// kyx = kxy, kzx = kxz, kzy = kyz - /// @endcode - /// However, we make no attempt at enforcing positive - /// definite tensors. - /// - /// @param [in] parser - /// An Eclipse data parser capable of answering which - /// permeability components are present in a given input - /// deck as well as retrieving the numerical value of each - /// permeability component in each grid cell. - /// - /// @param [out] tensor - /// @param [out] kmap - PermeabilityKind fillTensor(const EclipseGridParser& parser, - std::vector*>& tensor, - std::array& kmap) - { - PermeabilityKind kind = classifyPermeability(parser); - if (kind == Invalid) { - OPM_THROW(std::runtime_error, "Invalid set of permeability fields given."); - } - assert(tensor.size() == 1); - for (int i = 0; i < 9; ++i) { kmap[i] = 0; } - - enum { xx, xy, xz, // 0, 1, 2 - yx, yy, yz, // 3, 4, 5 - zx, zy, zz }; // 6, 7, 8 - - // ----------------------------------------------------------- - // 1st row: [kxx, kxy, kxz] - if (parser.hasField("PERMX" )) { - kmap[xx] = tensor.size(); - tensor.push_back(&parser.getFloatingPointValue("PERMX" )); - - setScalarPermIfNeeded(kmap, xx, yy, zz); - } - if (parser.hasField("PERMXY")) { - kmap[xy] = kmap[yx] = tensor.size(); // Enforce symmetry. - tensor.push_back(&parser.getFloatingPointValue("PERMXY")); - } - if (parser.hasField("PERMXZ")) { - kmap[xz] = kmap[zx] = tensor.size(); // Enforce symmetry. - tensor.push_back(&parser.getFloatingPointValue("PERMXZ")); - } - - // ----------------------------------------------------------- - // 2nd row: [kyx, kyy, kyz] - if (parser.hasField("PERMYX")) { - kmap[yx] = kmap[xy] = tensor.size(); // Enforce symmetry. - tensor.push_back(&parser.getFloatingPointValue("PERMYX")); - } - if (parser.hasField("PERMY" )) { - kmap[yy] = tensor.size(); - tensor.push_back(&parser.getFloatingPointValue("PERMY" )); - - setScalarPermIfNeeded(kmap, yy, zz, xx); - } - if (parser.hasField("PERMYZ")) { - kmap[yz] = kmap[zy] = tensor.size(); // Enforce symmetry. - tensor.push_back(&parser.getFloatingPointValue("PERMYZ")); - } - - // ----------------------------------------------------------- - // 3rd row: [kzx, kzy, kzz] - if (parser.hasField("PERMZX")) { - kmap[zx] = kmap[xz] = tensor.size(); // Enforce symmetry. - tensor.push_back(&parser.getFloatingPointValue("PERMZX")); - } - if (parser.hasField("PERMZY")) { - kmap[zy] = kmap[yz] = tensor.size(); // Enforce symmetry. - tensor.push_back(&parser.getFloatingPointValue("PERMZY")); - } - if (parser.hasField("PERMZ" )) { - kmap[zz] = tensor.size(); - tensor.push_back(&parser.getFloatingPointValue("PERMZ" )); - - setScalarPermIfNeeded(kmap, zz, xx, yy); - } - return kind; - } - /// @brief /// Extract pointers to appropriate tensor components from /// input deck. The permeability tensor is, generally, diff --git a/opm/core/props/rock/RockFromDeck.hpp b/opm/core/props/rock/RockFromDeck.hpp index 5d8c44aef..34eacb029 100644 --- a/opm/core/props/rock/RockFromDeck.hpp +++ b/opm/core/props/rock/RockFromDeck.hpp @@ -20,9 +20,6 @@ #ifndef OPM_ROCKFROMDECK_HEADER_INCLUDED #define OPM_ROCKFROMDECK_HEADER_INCLUDED - -#include - #include #include @@ -38,24 +35,7 @@ namespace Opm /// Default constructor. RockFromDeck(); - /// Initialize from deck and grid. - /// \param deck Deck input parser - /// \param grid Grid to which property object applies, needed for the - /// mapping from cell indices (typically from a processed grid) - /// to logical cartesian indices consistent with the deck. - void init(const EclipseGridParser& deck, - const UnstructuredGrid& grid); - /// Initialize from deck and cell mapping. - /// \param deck Deck input parser - /// \param number_of_cells The number of cells in the grid. - /// \param global_cell The mapping fom local to global cell indices. - /// global_cell[i] is the corresponding global index of i. - /// \param cart_dims The size of the underlying cartesian grid. - void init(const EclipseGridParser& deck, - int number_of_cells, const int* global_cell, - const int* cart_dims); - /// Initialize from deck and cell mapping. /// \param newParserDeck Deck produced by the opm-parser code /// \param number_of_cells The number of cells in the grid. /// \param global_cell The mapping fom local to global cell indices. @@ -92,17 +72,9 @@ namespace Opm } private: - void assignPorosity(const EclipseGridParser& parser, - int number_of_cells, - const int* global_cell); void assignPorosity(Opm::DeckConstPtr newParserDeck, int number_of_cells, const int* global_cell); - void assignPermeability(const EclipseGridParser& parser, - int number_of_cells, - const int* global_cell, - const int* cart_dims, - const double perm_threshold); void assignPermeability(Opm::DeckConstPtr newParserDeck, int number_of_cells, const int* global_cell, diff --git a/opm/core/props/satfunc/SaturationPropsFromDeck.hpp b/opm/core/props/satfunc/SaturationPropsFromDeck.hpp index ad9c2e702..070d5c2b0 100644 --- a/opm/core/props/satfunc/SaturationPropsFromDeck.hpp +++ b/opm/core/props/satfunc/SaturationPropsFromDeck.hpp @@ -22,7 +22,6 @@ #include #include -#include #include #include #include @@ -52,17 +51,6 @@ namespace Opm /// Default constructor. SaturationPropsFromDeck(); - /// Initialize from deck and grid. - /// \param[in] deck Deck input parser - /// \param[in] grid Grid to which property object applies, needed for the - /// mapping from cell indices (typically from a processed grid) - /// to logical cartesian indices consistent with the deck. - /// \param[in] samples Number of uniform sample points for saturation tables. - /// NOTE: samples will only be used with the SatFuncSetUniform template argument. - void init(const EclipseGridParser& deck, - const UnstructuredGrid& grid, - const int samples); - /// Initialize from deck and grid. /// \param[in] deck Deck input parser /// \param[in] grid Grid to which property object applies, needed for the @@ -74,29 +62,11 @@ namespace Opm const UnstructuredGrid& grid, const int samples); - /// Initialize from deck and grid. - /// \param[in] deck Deck input parser + /// \param[in] newParserDeck Deck input parser /// \param[in] number_of_cells The number of cells of the grid to which property /// object applies, needed for the /// mapping from cell indices (typically from a processed - /// grid) to logical cartesian indices consistent with the - /// deck. - /// \param[in] global_cell The mapping from local cell indices of the grid to - /// global cell indices used in the deck. - /// \param[in] begin_cell_centroids Pointer to the first cell_centroid of the grid. - /// \param[in] dimensions The dimensions of the grid. - /// \param[in] samples Number of uniform sample points for saturation tables. - /// \tparam T The iterator Type for the cell centroids. - /// NOTE: samples will only be used with the SatFuncSetUniform template argument. - template - void init(const EclipseGridParser& deck, - int number_of_cells, - const int* global_cell, - const T& begin_cell_centroids, - int dimensions, - const int samples); - /// grid) to logical cartesian indices consistent with the /// deck. /// \param[in] global_cell The mapping from local cell indices of the grid to @@ -178,27 +148,6 @@ namespace Opm typedef SatFuncSet Funcs; const Funcs& funcForCell(const int cell) const; - - template - void initEPS(const EclipseGridParser& deck, - int number_of_cells, - const int* global_cell, - const T& begin_cell_centroids, - int dimensions); - template - void initEPSHyst(const EclipseGridParser& deck, - int number_of_cells, - const int* global_cell, - const T& begin_cell_centroids, - int dimensions); - template - void initEPSKey(const EclipseGridParser& deck, - int number_of_cells, - const int* global_cell, - const T& begin_cell_centroids, - int dimensions, - const std::string& keyword, - std::vector& scaleparam); template void initEPS(Opm::DeckConstPtr newParserDeck, int number_of_cells, diff --git a/opm/core/props/satfunc/SaturationPropsFromDeck_impl.hpp b/opm/core/props/satfunc/SaturationPropsFromDeck_impl.hpp index 330997ee0..67611aa2f 100644 --- a/opm/core/props/satfunc/SaturationPropsFromDeck_impl.hpp +++ b/opm/core/props/satfunc/SaturationPropsFromDeck_impl.hpp @@ -47,126 +47,15 @@ namespace Opm { } - /// Initialize from deck. - template - void SaturationPropsFromDeck::init(const EclipseGridParser& deck, - const UnstructuredGrid& grid, - const int samples) - { - init(deck, grid.number_of_cells, grid.global_cell, grid.cell_centroids, - grid.dimensions, samples); - } - - /// Initialize from deck. template void SaturationPropsFromDeck::init(Opm::DeckConstPtr newParserDeck, const UnstructuredGrid& grid, const int samples) { - init(newParserDeck, grid.number_of_cells, - grid.global_cell, grid.cell_centroids, - grid.dimensions, samples); - } - - /// Initialize from deck. - template - template< class T> - void SaturationPropsFromDeck::init(const EclipseGridParser& deck, - int number_of_cells, - const int* global_cell, - const T& begin_cell_centroids, - int dimensions, - const int samples) - { - phase_usage_ = phaseUsageFromDeck(deck); - - // Extract input data. - // Oil phase should be active. - if (!phase_usage_.phase_used[Liquid]) { - OPM_THROW(std::runtime_error, "SaturationPropsFromDeck::init() -- oil phase must be active."); - } - - // Obtain SATNUM, if it exists, and create cell_to_func_. - // Otherwise, let the cell_to_func_ mapping be just empty. - int satfuncs_expected = 1; - if (deck.hasField("SATNUM")) { - const std::vector& satnum = deck.getIntegerValue("SATNUM"); - satfuncs_expected = *std::max_element(satnum.begin(), satnum.end()); - cell_to_func_.resize(number_of_cells); - for (int cell = 0; cell < number_of_cells; ++cell) { - const int deck_pos = (global_cell == NULL) ? cell : global_cell[cell]; - cell_to_func_[cell] = satnum[deck_pos] - 1; - } - } - - // Find number of tables, check for consistency. - enum { Uninitialized = -1 }; - int num_tables = Uninitialized; - if (phase_usage_.phase_used[Aqua]) { - const SWOF::table_t& swof_table = deck.getSWOF().swof_; - num_tables = swof_table.size(); - if (num_tables < satfuncs_expected) { - OPM_THROW(std::runtime_error, "Found " << num_tables << " SWOF tables, SATNUM specifies at least " << satfuncs_expected); - } - } - if (phase_usage_.phase_used[Vapour]) { - const SGOF::table_t& sgof_table = deck.getSGOF().sgof_; - int num_sgof_tables = sgof_table.size(); - if (num_sgof_tables < satfuncs_expected) { - OPM_THROW(std::runtime_error, "Found " << num_tables << " SGOF tables, SATNUM specifies at least " << satfuncs_expected); - } - if (num_tables == Uninitialized) { - num_tables = num_sgof_tables; - } else if (num_tables != num_sgof_tables) { - OPM_THROW(std::runtime_error, "Inconsistent number of tables in SWOF and SGOF."); - } - } - - // Initialize tables. - satfuncset_.resize(num_tables); - for (int table = 0; table < num_tables; ++table) { - satfuncset_[table].init(deck, table, phase_usage_, samples); - } - - // Saturation table scaling - do_hyst_ = false; - do_eps_ = false; - do_3pt_ = false; - if (deck.hasField("ENDSCALE")) { - //if (!phase_usage_.phase_used[Aqua] || !phase_usage_.phase_used[Liquid] || phase_usage_.phase_used[Vapour]) { - // OPM_THROW(std::runtime_error, "Currently endpoint-scaling limited to oil-water systems without gas."); - //} - if (deck.getENDSCALE().dir_switch_ != std::string("NODIR")) { - OPM_THROW(std::runtime_error, "SaturationPropsFromDeck::init() -- ENDSCALE: Currently only 'NODIR' accepted."); - } - if (deck.getENDSCALE().revers_switch_ != std::string("REVERS")) { - OPM_THROW(std::runtime_error, "SaturationPropsFromDeck::init() -- ENDSCALE: Currently only 'REVERS' accepted."); - } - if (deck.hasField("SCALECRS")) { - if (deck.getSCALECRS().scalecrs_ == std::string("YES")) { - do_3pt_ = true; - } - } - do_eps_ = true; - initEPS(deck, number_of_cells, global_cell, begin_cell_centroids, dimensions); - - // For now, a primitive detection of hysteresis. TODO: SATOPTS HYSTER/ and EHYSTR - do_hyst_ = deck.hasField("ISWL") || deck.hasField("ISWU") || deck.hasField("ISWCR") || deck.hasField("ISGL") || - deck.hasField("ISGU") || deck.hasField("ISGCR") || deck.hasField("ISOWCR") || deck.hasField("ISOGCR"); - if (do_hyst_) { - if (deck.hasField("KRW") || deck.hasField("KRG") || deck.hasField("KRO") || deck.hasField("KRWR") || - deck.hasField("KRGR") || deck.hasField("KRORW") || deck.hasField("KRORG") || - deck.hasField("IKRW") || deck.hasField("IKRG") || deck.hasField("IKRO") || deck.hasField("IKRWR") || - deck.hasField("IKRGR") || deck.hasField("IKRORW") || deck.hasField("IKRORG") ) { - OPM_THROW(std::runtime_error, "SaturationPropsFromDeck::init() -- ENDSCALE: Currently hysteresis and relperm value scaling can not be combined."); - } - initEPSHyst(deck, number_of_cells, global_cell, begin_cell_centroids, - dimensions); - } - - //OPM_THROW(std::runtime_error, "SaturationPropsFromDeck::init() -- ENDSCALE: Under construction ..."); - } + this->init(newParserDeck, grid.number_of_cells, + grid.global_cell, grid.cell_centroids, + grid.dimensions, samples); } /// Initialize from deck. @@ -532,89 +421,6 @@ namespace Opm return cell_to_func_.empty() ? satfuncset_[0] : satfuncset_[cell_to_func_[cell]]; } - // Initialize saturation scaling parameters - template - template - void SaturationPropsFromDeck::initEPS(const EclipseGridParser& deck, - int number_of_cells, - const int* global_cell, - const T& begin_cell_centroid, - int dimensions) - { - std::vector swl, swcr, swu, sgl, sgcr, sgu, sowcr, sogcr; - std::vector krw, krg, kro, krwr, krgr, krorw, krorg; - // Initialize saturation scaling parameter - initEPSKey(deck, number_of_cells, global_cell, begin_cell_centroid, dimensions, - std::string("SWL"), swl); - initEPSKey(deck, number_of_cells, global_cell, begin_cell_centroid, dimensions, - std::string("SWU"), swu); - initEPSKey(deck, number_of_cells, global_cell, begin_cell_centroid, dimensions, - std::string("SWCR"), swcr); - initEPSKey(deck, number_of_cells, global_cell, begin_cell_centroid, dimensions, - std::string("SGL"), sgl); - initEPSKey(deck, number_of_cells, global_cell, begin_cell_centroid, dimensions, - std::string("SGU"), sgu); - initEPSKey(deck, number_of_cells, global_cell, begin_cell_centroid, dimensions, - std::string("SGCR"), sgcr); - initEPSKey(deck, number_of_cells, global_cell, begin_cell_centroid, dimensions, - std::string("SOWCR"), sowcr); - initEPSKey(deck, number_of_cells, global_cell, begin_cell_centroid, dimensions, - std::string("SOGCR"), sogcr); - initEPSKey(deck, number_of_cells, global_cell, begin_cell_centroid, dimensions, - std::string("KRW"), krw); - initEPSKey(deck, number_of_cells, global_cell, begin_cell_centroid, dimensions, - std::string("KRG"), krg); - initEPSKey(deck, number_of_cells, global_cell, begin_cell_centroid, dimensions, - std::string("KRO"), kro); - initEPSKey(deck, number_of_cells, global_cell, begin_cell_centroid, dimensions, - std::string("KRWR"), krwr); - initEPSKey(deck, number_of_cells, global_cell, begin_cell_centroid, dimensions, - std::string("KRGR"), krgr); - initEPSKey(deck, number_of_cells, global_cell, begin_cell_centroid, dimensions, - std::string("KRORW"), krorw); - initEPSKey(deck, number_of_cells, global_cell, begin_cell_centroid, dimensions, - std::string("KRORG"), krorg); - - eps_transf_.resize(number_of_cells); - - const int wpos = phase_usage_.phase_pos[BlackoilPhases::Aqua]; - const int gpos = phase_usage_.phase_pos[BlackoilPhases::Vapour]; - const bool oilWater = phase_usage_.phase_used[Aqua] && phase_usage_.phase_used[Liquid] && !phase_usage_.phase_used[Vapour]; - const bool oilGas = !phase_usage_.phase_used[Aqua] && phase_usage_.phase_used[Liquid] && phase_usage_.phase_used[Vapour]; - const bool threephase = phase_usage_.phase_used[Aqua] && phase_usage_.phase_used[Liquid] && phase_usage_.phase_used[Vapour]; - - for (int cell = 0; cell < number_of_cells; ++cell) { - if (oilWater) { - // ### krw - initEPSParam(cell, eps_transf_[cell].wat, false, funcForCell(cell).smin_[wpos], funcForCell(cell).swcr_, funcForCell(cell).smax_[wpos], - funcForCell(cell).sowcr_, -1.0, funcForCell(cell).krwr_, funcForCell(cell).krwmax_, swl, swcr, swu, sowcr, sgl, krwr, krw); - // ### krow - initEPSParam(cell, eps_transf_[cell].watoil, true, 0.0, funcForCell(cell).sowcr_, funcForCell(cell).smin_[wpos], - funcForCell(cell).swcr_, -1.0, funcForCell(cell).krorw_, funcForCell(cell).kromax_, swl, sowcr, swl, swcr, sgl, krorw, kro); - } else if (oilGas) { - // ### krg - initEPSParam(cell, eps_transf_[cell].gas, false, funcForCell(cell).smin_[gpos], funcForCell(cell).sgcr_, funcForCell(cell).smax_[gpos], - funcForCell(cell).sogcr_, -1.0, funcForCell(cell).krgr_, funcForCell(cell).krgmax_, sgl, sgcr, sgu, sogcr, swl, krgr, krg); - // ### krog - initEPSParam(cell, eps_transf_[cell].gasoil, true, 0.0, funcForCell(cell).sogcr_, funcForCell(cell).smin_[gpos], - funcForCell(cell).sgcr_, -1.0, funcForCell(cell).krorg_, funcForCell(cell).kromax_, sgl, sogcr, sgl, sgcr, swl, krorg, kro); - } else if (threephase) { - // ### krw - initEPSParam(cell, eps_transf_[cell].wat, false, funcForCell(cell).smin_[wpos], funcForCell(cell).swcr_, funcForCell(cell).smax_[wpos], funcForCell(cell).sowcr_, - funcForCell(cell).smin_[gpos], funcForCell(cell).krwr_, funcForCell(cell).krwmax_, swl, swcr, swu, sowcr, sgl, krwr, krw); - // ### krow - initEPSParam(cell, eps_transf_[cell].watoil, true, 0.0, funcForCell(cell).sowcr_, funcForCell(cell).smin_[wpos], funcForCell(cell).swcr_, - funcForCell(cell).smin_[gpos], funcForCell(cell).krorw_, funcForCell(cell).kromax_, swl, sowcr, swl, swcr, sgl, krorw, kro); - // ### krg - initEPSParam(cell, eps_transf_[cell].gas, false, funcForCell(cell).smin_[gpos], funcForCell(cell).sgcr_, funcForCell(cell).smax_[gpos], funcForCell(cell).sogcr_, - funcForCell(cell).smin_[wpos], funcForCell(cell).krgr_, funcForCell(cell).krgmax_, sgl, sgcr, sgu, sogcr, swl, krgr, krg); - // ### krog - initEPSParam(cell, eps_transf_[cell].gasoil, true, 0.0, funcForCell(cell).sogcr_, funcForCell(cell).smin_[gpos], funcForCell(cell).sgcr_, - funcForCell(cell).smin_[wpos], funcForCell(cell).krorg_, funcForCell(cell).kromax_, sgl, sogcr, sgl, sgcr, swl, krorg, kro); - } - } - } - // Initialize saturation scaling parameters template template @@ -698,91 +504,6 @@ namespace Opm } } - // Initialize hysteresis saturation scaling parameters - template - template - void SaturationPropsFromDeck::initEPSHyst(const EclipseGridParser& deck, - int number_of_cells, - const int* global_cell, - const T& begin_cell_centroid, - int dimensions) - { - std::vector iswl, iswcr, iswu, isgl, isgcr, isgu, isowcr, isogcr; - std::vector ikrw, ikrg, ikro, ikrwr, ikrgr, ikrorw, ikrorg; - // Initialize hysteresis saturation scaling parameters - initEPSKey(deck, number_of_cells, global_cell, begin_cell_centroid, dimensions, - std::string("ISWL"), iswl); - initEPSKey(deck, number_of_cells, global_cell, begin_cell_centroid, dimensions, - std::string("ISWU"), iswu); - initEPSKey(deck, number_of_cells, global_cell, begin_cell_centroid, dimensions, - std::string("ISWCR"), iswcr); - initEPSKey(deck, number_of_cells, global_cell, begin_cell_centroid, dimensions, - std::string("ISGL"), isgl); - initEPSKey(deck, number_of_cells, global_cell, begin_cell_centroid, dimensions, - std::string("ISGU"), isgu); - initEPSKey(deck, number_of_cells, global_cell, begin_cell_centroid, dimensions, - std::string("ISGCR"), isgcr); - initEPSKey(deck, number_of_cells, global_cell, begin_cell_centroid, dimensions, - std::string("ISOWCR"), isowcr); - initEPSKey(deck, number_of_cells, global_cell, begin_cell_centroid, dimensions, - std::string("ISOGCR"), isogcr); - initEPSKey(deck, number_of_cells, global_cell, begin_cell_centroid, dimensions, - std::string("IKRW"), ikrw); - initEPSKey(deck, number_of_cells, global_cell, begin_cell_centroid, dimensions, - std::string("IKRG"), ikrg); - initEPSKey(deck, number_of_cells, global_cell, begin_cell_centroid, dimensions, - std::string("IKRO"), ikro); - initEPSKey(deck, number_of_cells, global_cell, begin_cell_centroid, dimensions, - std::string("IKRWR"), ikrwr); - initEPSKey(deck, number_of_cells, global_cell, begin_cell_centroid, dimensions, - std::string("IKRGR"), ikrgr); - initEPSKey(deck, number_of_cells, global_cell, begin_cell_centroid, dimensions, - std::string("IKRORW"), ikrorw); - initEPSKey(deck, number_of_cells, global_cell, begin_cell_centroid, dimensions, - std::string("IKRORG"), ikrorg); - - - eps_transf_hyst_.resize(number_of_cells); - sat_hyst_.resize(number_of_cells); - - const int wpos = phase_usage_.phase_pos[BlackoilPhases::Aqua]; - const int gpos = phase_usage_.phase_pos[BlackoilPhases::Vapour]; - const bool oilWater = phase_usage_.phase_used[Aqua] && phase_usage_.phase_used[Liquid] && !phase_usage_.phase_used[Vapour]; - const bool oilGas = !phase_usage_.phase_used[Aqua] && phase_usage_.phase_used[Liquid] && phase_usage_.phase_used[Vapour]; - const bool threephase = phase_usage_.phase_used[Aqua] && phase_usage_.phase_used[Liquid] && phase_usage_.phase_used[Vapour]; - - for (int cell = 0; cell < number_of_cells; ++cell) { - if (oilWater) { - // ### krw - initEPSParam(cell, eps_transf_hyst_[cell].wat, false, funcForCell(cell).smin_[wpos], funcForCell(cell).swcr_, funcForCell(cell).smax_[wpos], - funcForCell(cell).sowcr_, -1.0, funcForCell(cell).krwr_, funcForCell(cell).krwmax_, iswl, iswcr, iswu, isowcr, isgl, ikrwr, ikrw); - // ### krow - initEPSParam(cell, eps_transf_hyst_[cell].watoil, true, 0.0, funcForCell(cell).sowcr_, funcForCell(cell).smin_[wpos], - funcForCell(cell).swcr_, -1.0, funcForCell(cell).krorw_, funcForCell(cell).kromax_, iswl, isowcr, iswl, iswcr, isgl, ikrorw, ikro); - } else if (oilGas) { - // ### krg - initEPSParam(cell, eps_transf_hyst_[cell].gas, false, funcForCell(cell).smin_[gpos], funcForCell(cell).sgcr_, funcForCell(cell).smax_[gpos], - funcForCell(cell).sogcr_, -1.0, funcForCell(cell).krgr_, funcForCell(cell).krgmax_, isgl, isgcr, isgu, isogcr, iswl, ikrgr, ikrg); - // ### krog - initEPSParam(cell, eps_transf_hyst_[cell].gasoil, true, 0.0, funcForCell(cell).sogcr_, funcForCell(cell).smin_[gpos], - funcForCell(cell).sgcr_, -1.0, funcForCell(cell).krorg_, funcForCell(cell).kromax_, isgl, isogcr, isgl, isgcr, iswl, ikrorg, ikro); - } else if (threephase) { - // ### krw - initEPSParam(cell, eps_transf_hyst_[cell].wat, false, funcForCell(cell).smin_[wpos], funcForCell(cell).swcr_, funcForCell(cell).smax_[wpos], funcForCell(cell).sowcr_, - funcForCell(cell).smin_[gpos], funcForCell(cell).krwr_, funcForCell(cell).krwmax_, iswl, iswcr, iswu, isowcr, isgl, ikrwr, ikrw); - // ### krow - initEPSParam(cell, eps_transf_hyst_[cell].watoil, true, 0.0, funcForCell(cell).sowcr_, funcForCell(cell).smin_[wpos], funcForCell(cell).swcr_, - funcForCell(cell).smin_[gpos], funcForCell(cell).krorw_, funcForCell(cell).kromax_, iswl, isowcr, iswl, iswcr, isgl, ikrorw, ikro); - // ### krg - initEPSParam(cell, eps_transf_hyst_[cell].gas, false, funcForCell(cell).smin_[gpos], funcForCell(cell).sgcr_, funcForCell(cell).smax_[gpos], funcForCell(cell).sogcr_, - funcForCell(cell).smin_[wpos], funcForCell(cell).krgr_, funcForCell(cell).krgmax_, isgl, isgcr, isgu, isogcr, iswl, ikrgr, ikrg); - // ### krog - initEPSParam(cell, eps_transf_hyst_[cell].gasoil, true, 0.0, funcForCell(cell).sogcr_, funcForCell(cell).smin_[gpos], funcForCell(cell).sgcr_, - funcForCell(cell).smin_[wpos], funcForCell(cell).krorg_, funcForCell(cell).kromax_, isgl, isogcr, isgl, isgcr, iswl, ikrorg, ikro); - } - } - } - // Initialize hysteresis saturation scaling parameters template template @@ -867,184 +588,6 @@ namespace Opm } } - // Initialize saturation scaling parameter - template - template - void SaturationPropsFromDeck::initEPSKey(const EclipseGridParser& deck, - int number_of_cells, - const int* global_cell, - const T& begin_cell_centroid, - int dimensions, - const std::string& keyword, - std::vector& scaleparam) - { - const bool useAqua = phase_usage_.phase_used[Aqua]; - const bool useLiquid = phase_usage_.phase_used[Liquid]; - const bool useVapour = phase_usage_.phase_used[Vapour]; - bool useKeyword = deck.hasField(keyword); - bool hasENPTVD = deck.hasField("ENPTVD"); - bool hasENKRVD = deck.hasField("ENKRVD"); - int itab = 0; - std::vector > > table_dummy; - std::vector > >& table = table_dummy; - - // Active keyword assigned default values for each cell (in case of possible box-wise assignment) - int phase_pos_aqua = phase_usage_.phase_pos[BlackoilPhases::Aqua]; - int phase_pos_vapour = phase_usage_.phase_pos[BlackoilPhases::Vapour]; - if ((keyword[0] == 'S' && (useKeyword || hasENPTVD)) || (keyword[1] == 'S' && useKeyword) ) { - if (keyword == std::string("SWL") || keyword == std::string("ISWL") ) { - if (useAqua && (useKeyword || deck.getENPTVD().mask_[0])) { - itab = 1; - scaleparam.resize(number_of_cells); - for (int i=0; i 0) { - table = deck.getENPTVD().table_; - } - } else if ((keyword[0] == 'K' && (useKeyword || hasENKRVD)) || (keyword[1] == 'K' && useKeyword) ) { - if (keyword == std::string("KRW") || keyword == std::string("IKRW") ) { - if (useAqua && (useKeyword || deck.getENKRVD().mask_[0])) { - itab = 1; - scaleparam.resize(number_of_cells); - for (int i=0; i 0) { - table = deck.getENKRVD().table_; - } - } - - if (scaleparam.empty()) { - return; - } else if (useKeyword) { - // Keyword values from deck - std::cout << "--- Scaling parameter '" << keyword << "' assigned." << std::endl; - const std::vector& val = deck.getFloatingPointValue(keyword); - for (int c = 0; c < int(scaleparam.size()); ++c) { - const int deck_pos = (global_cell == NULL) ? c : global_cell[c]; - scaleparam[c] = val[deck_pos]; - } - } else { - std::cout << "--- Scaling parameter '" << keyword << "' assigned via "; - if (keyword[0] == 'S') - deck.getENPTVD().write(std::cout); - else - deck.getENKRVD().write(std::cout); - for (int cell = 0; cell < number_of_cells; ++cell) { - int jtab = cell_to_func_.empty() ? 0 : cell_to_func_[cell]; - if (table[itab][jtab][0] != -1.0) { - std::vector& depth = table[0][jtab]; - std::vector& val = table[itab][jtab]; - double zc = UgGridHelpers - ::getCoordinate(UgGridHelpers::increment(begin_cell_centroid, cell, - dimensions), - dimensions-1); - if (zc >= depth.front() && zc <= depth.back()) { //don't want extrap outside depth interval - scaleparam[cell] = linearInterpolation(depth, val, zc); - } - } - } - } - } - // Initialize saturation scaling parameter template template diff --git a/opm/core/simulator/initState.hpp b/opm/core/simulator/initState.hpp index 4c3a8cc6a..7f9ed9809 100644 --- a/opm/core/simulator/initState.hpp +++ b/opm/core/simulator/initState.hpp @@ -20,7 +20,7 @@ #ifndef OPM_INITSTATE_HEADER_INCLUDED #define OPM_INITSTATE_HEADER_INCLUDED -#include // DeckConstPtr +#include struct UnstructuredGrid; @@ -28,7 +28,6 @@ namespace Opm { namespace parameter { class ParameterGroup; } - class EclipseGridParser; class IncompPropertiesInterface; class BlackoilPropertiesInterface; @@ -168,27 +167,7 @@ namespace Opm template void initStateFromDeck(const UnstructuredGrid& grid, const Props& props, - const EclipseGridParser& deck, - const double gravity, - State& state); - - /// Initialize a two-phase state from input deck. - /// If EQUIL is present: - /// - saturation is set according to the water-oil contact, - /// - pressure is set to hydrostatic equilibrium. - /// Otherwise: - /// - saturation is set according to SWAT, - /// - pressure is set according to PRESSURE. - template - void initStateFromDeck(int number_of_cells, - const int* global_cell, - int number_of_faces, - FaceCells face_cells, - FCI begin_face_centroids, - CCI begin_cell_centroids, - int dimensions, - const Props& props, - const EclipseGridParser& deck, + Opm::DeckConstPtr newParserDeck, const double gravity, State& state); @@ -203,31 +182,9 @@ namespace Opm template void initBlackoilStateFromDeck(const UnstructuredGrid& grid, const Props& props, - const EclipseGridParser& deck, + Opm::DeckConstPtr newParserDeck, const double gravity, State& state); - - /// Initialize a two-phase water-oil blackoil state from input deck. - /// If EQUIL is present: - /// - saturation is set according to the water-oil contact, - /// - pressure is set to hydrostatic equilibrium. - /// Otherwise: - /// - saturation is set according to SWAT, - /// - pressure is set according to PRESSURE. - /// In addition, this function sets surfacevol. - template - void initBlackoilStateFromDeck(int number_of_cells, - const int* global_cell, - int number_of_faces, - FaceCells face_cells, - FCI begin_face_centroids, - CCI begin_cell_centroids, - int dimensions, - const Props& props, - const EclipseGridParser& deck, - const double gravity, - State& state); - /// Initialize a blackoil state from input deck. template void initBlackoilStateFromDeck(int number_of_cells, diff --git a/opm/core/simulator/initStateEquil.hpp b/opm/core/simulator/initStateEquil.hpp index 98ebadbb3..fd8c35a13 100644 --- a/opm/core/simulator/initStateEquil.hpp +++ b/opm/core/simulator/initStateEquil.hpp @@ -22,7 +22,6 @@ #include #include -#include #include #include #include @@ -60,13 +59,6 @@ namespace Opm * \param[in] deck Simulation deck, used to obtain EQUIL and related data. * \param[in] gravity Acceleration of gravity, assumed to be in Z direction. */ - void initStateEquil(const UnstructuredGrid& grid, - const BlackoilPropertiesInterface& props, - const EclipseGridParser& deck, - const double gravity, - BlackoilState& state); - - void initStateEquil(const UnstructuredGrid& grid, const BlackoilPropertiesInterface& props, const Opm::DeckConstPtr newParserDeck, @@ -188,54 +180,6 @@ namespace Opm const std::vector gas_saturation); namespace DeckDependent { - - inline - std::vector - getEquil(const EclipseGridParser& deck) - { - if (deck.hasField("EQUIL")) { - const EQUIL& eql = deck.getEQUIL(); - - typedef std::vector::size_type sz_t; - const sz_t nrec = eql.equil.size(); - - std::vector ret; - ret.reserve(nrec); - for (sz_t r = 0; r < nrec; ++r) { - const EquilLine& rec = eql.equil[r]; - - EquilRecord record = - { - { rec.datum_depth_ , - rec.datum_depth_pressure_ } - , - { rec.water_oil_contact_depth_ , - rec.oil_water_cap_pressure_ } - , - { rec.gas_oil_contact_depth_ , - rec.gas_oil_cap_pressure_ } - , - rec.live_oil_table_index_ - , - rec.wet_gas_table_index_ - , - rec.N_ - }; - if (record.N != 0) { - OPM_THROW(std::domain_error, - "kw EQUIL, item 9: Only N=0 supported."); - } - ret.push_back(record); - } - - return ret; - } - else { - OPM_THROW(std::domain_error, - "Deck does not provide equilibration data."); - } - } - inline std::vector getEquil(const Opm::DeckConstPtr newParserDeck) @@ -282,29 +226,6 @@ namespace Opm } } - - inline - std::vector - equilnum(const EclipseGridParser& deck, - const UnstructuredGrid& G ) - { - std::vector eqlnum; - if (deck.hasField("EQLNUM")) { - const std::vector& e = deck.getIntegerValue("EQLNUM"); - eqlnum.reserve(e.size()); - std::transform(e.begin(), e.end(), std::back_inserter(eqlnum), - std::bind2nd(std::minus(), 1)); - } - else { - // No explicit equilibration region. - // All cells in region zero. - eqlnum.assign(G.number_of_cells, 0); - } - - return eqlnum; - } - - inline std::vector equilnum(const Opm::DeckConstPtr newParserDeck, @@ -331,180 +252,7 @@ namespace Opm } - template - class InitialStateComputer; - - template <> - class InitialStateComputer { - public: - InitialStateComputer(const BlackoilPropertiesInterface& props, - const EclipseGridParser& deck , - const UnstructuredGrid& G , - const double grav = unit::gravity) - : pp_(props.numPhases(), - std::vector(G.number_of_cells)), - sat_(props.numPhases(), - std::vector(G.number_of_cells)), - rs_(G.number_of_cells), - rv_(G.number_of_cells) - { - // Get the equilibration records. - const std::vector rec = getEquil(deck); - - // Create (inverse) region mapping. - const RegionMapping<> eqlmap(equilnum(deck, G)); - - // Create Rs functions. - rs_func_.reserve(rec.size()); - if (deck.hasField("DISGAS")) { - for (size_t i = 0; i < rec.size(); ++i) { - const int cell = *(eqlmap.cells(i).begin()); - if (rec[i].live_oil_table_index > 0) { - if (deck.hasField("RSVD")) { - // TODO When this kw is actually parsed, also check for proper number of available tables - // For now, just use dummy ... - std::vector depth; depth.push_back(0.0); depth.push_back(100.0); - std::vector rs; rs.push_back(0.0); rs.push_back(100.0); - rs_func_.push_back(std::make_shared(props, cell, depth, rs)); - } else { - OPM_THROW(std::runtime_error, "Cannot initialise: RSVD table " << (rec[i].live_oil_table_index) << " not available."); - } - } else { - if (rec[i].goc.depth != rec[i].main.depth) { - OPM_THROW(std::runtime_error, - "Cannot initialise: when no explicit RSVD table is given, \n" - "datum depth must be at the gas-oil-contact. " - "In EQUIL region " << (i + 1) << " (counting from 1), this does not hold."); - } - const double p_contact = rec[i].main.press; - rs_func_.push_back(std::make_shared(props, cell, p_contact)); - } - } - } else { - for (size_t i = 0; i < rec.size(); ++i) { - rs_func_.push_back(std::make_shared()); - } - } - - rv_func_.reserve(rec.size()); - if (deck.hasField("VAPOIL")) { - for (size_t i = 0; i < rec.size(); ++i) { - const int cell = *(eqlmap.cells(i).begin()); - if (rec[i].wet_gas_table_index > 0) { - if (deck.hasField("RVVD")) { - // TODO When this kw is actually parsed, also check for proper number of available tables - // For now, just use dummy ... - std::vector depth; depth.push_back(0.0); depth.push_back(100.0); - std::vector rv; rv.push_back(0.0); rv.push_back(0.0001); - rv_func_.push_back(std::make_shared(props, cell, depth, rv)); - } else { - OPM_THROW(std::runtime_error, "Cannot initialise: RVVD table " << (rec[i].wet_gas_table_index) << " not available."); - } - } else { - if (rec[i].goc.depth != rec[i].main.depth) { - OPM_THROW(std::runtime_error, - "Cannot initialise: when no explicit RVVD table is given, \n" - "datum depth must be at the gas-oil-contact. " - "In EQUIL region " << (i + 1) << " (counting from 1), this does not hold."); - } - const double p_contact = rec[i].main.press + rec[i].goc.press; - rv_func_.push_back(std::make_shared(props, cell, p_contact)); - } - } - } else { - for (size_t i = 0; i < rec.size(); ++i) { - rv_func_.push_back(std::make_shared()); - } - } - - // Compute pressures, saturations, rs and rv factors. - calcPressSatRsRv(eqlmap, rec, props, G, grav); - - // Modify oil pressure in no-oil regions so that the pressures of present phases can - // be recovered from the oil pressure and capillary relations. - } - - typedef std::vector Vec; - typedef std::vector PVec; // One per phase. - - const PVec& press() const { return pp_; } - const PVec& saturation() const { return sat_; } - const Vec& rs() const { return rs_; } - const Vec& rv() const { return rv_; } - - private: - typedef DensityCalculator RhoCalc; - typedef EquilReg EqReg; - - std::vector< std::shared_ptr > rs_func_; - std::vector< std::shared_ptr > rv_func_; - - PVec pp_; - PVec sat_; - Vec rs_; - Vec rv_; - - template - void - calcPressSatRsRv(const RMap& reg , - const std::vector< EquilRecord >& rec , - const Opm::BlackoilPropertiesInterface& props, - const UnstructuredGrid& G , - const double grav) - { - typedef Miscibility::NoMixing NoMix; - - for (typename RMap::RegionId - r = 0, nr = reg.numRegions(); - r < nr; ++r) - { - const typename RMap::CellRange cells = reg.cells(r); - - const int repcell = *cells.begin(); - const RhoCalc calc(props, repcell); - const EqReg eqreg(rec[r], calc, - rs_func_[r], rv_func_[r], - props.phaseUsage()); - - PVec press = phasePressures(G, eqreg, cells, grav); - - const PVec sat = phaseSaturations(eqreg, cells, props, press); - - const int np = props.numPhases(); - for (int p = 0; p < np; ++p) { - copyFromRegion(press[p], cells, pp_[p]); - copyFromRegion(sat[p], cells, sat_[p]); - } - if (props.phaseUsage().phase_used[BlackoilPhases::Liquid] - && props.phaseUsage().phase_used[BlackoilPhases::Vapour]) { - const int oilpos = props.phaseUsage().phase_pos[BlackoilPhases::Liquid]; - const int gaspos = props.phaseUsage().phase_pos[BlackoilPhases::Vapour]; - const Vec rs = computeRs(G, cells, press[oilpos], *(rs_func_[r]), sat[gaspos]); - const Vec rv = computeRs(G, cells, press[gaspos], *(rv_func_[r]), sat[oilpos]); - copyFromRegion(rs, cells, rs_); - copyFromRegion(rv, cells, rv_); - } - } - } - - template - void copyFromRegion(const Vec& source, - const CellRangeType& cells, - Vec& destination) - { - auto s = source.begin(); - auto c = cells.begin(); - const auto e = cells.end(); - for (; c != e; ++c, ++s) { - destination[*c] = *s; - } - } - - }; - - - template <> - class InitialStateComputer { + class InitialStateComputer { public: InitialStateComputer(const BlackoilPropertiesInterface& props, const Opm::DeckConstPtr newParserDeck, diff --git a/opm/core/simulator/initStateEquil_impl.hpp b/opm/core/simulator/initStateEquil_impl.hpp index d227d54d1..cd08a98ec 100644 --- a/opm/core/simulator/initStateEquil_impl.hpp +++ b/opm/core/simulator/initStateEquil_impl.hpp @@ -735,33 +735,13 @@ namespace Opm * \param[in] deck Simulation deck, used to obtain EQUIL and related data. * \param[in] gravity Acceleration of gravity, assumed to be in Z direction. */ - void initStateEquil(const UnstructuredGrid& grid, - const BlackoilPropertiesInterface& props, - const EclipseGridParser& deck, - const double gravity, - BlackoilState& state) - { - typedef Equil::DeckDependent::InitialStateComputer ISC; - ISC isc(props, deck, grid, gravity); - const auto pu = props.phaseUsage(); - const int ref_phase = pu.phase_used[BlackoilPhases::Liquid] - ? pu.phase_pos[BlackoilPhases::Liquid] - : pu.phase_pos[BlackoilPhases::Aqua]; - state.pressure() = isc.press()[ref_phase]; - state.saturation() = convertSats(isc.saturation()); - state.gasoilratio() = isc.rs(); - state.rv() = isc.rv(); - // TODO: state.surfacevol() must be computed from s, rs, rv. - } - - void initStateEquil(const UnstructuredGrid& grid, const BlackoilPropertiesInterface& props, const Opm::DeckConstPtr newParserDeck, const double gravity, BlackoilState& state) { - typedef Equil::DeckDependent::InitialStateComputer ISC; + typedef Equil::DeckDependent::InitialStateComputer ISC; ISC isc(props, newParserDeck, grid, gravity); const auto pu = props.phaseUsage(); const int ref_phase = pu.phase_used[BlackoilPhases::Liquid] diff --git a/opm/core/simulator/initState_impl.hpp b/opm/core/simulator/initState_impl.hpp index 7f9435ff2..cf8cd984e 100644 --- a/opm/core/simulator/initState_impl.hpp +++ b/opm/core/simulator/initState_impl.hpp @@ -24,7 +24,6 @@ #include #include #include -#include #include #include #include @@ -673,129 +672,6 @@ namespace Opm begin_cell_centroids, state); } - /// Initialize a state from input deck. - template - void initStateFromDeck(const UnstructuredGrid& grid, - const Props& props, - const EclipseGridParser& deck, - const double gravity, - State& state) - { - initStateFromDeck(grid.number_of_cells, grid.global_cell, - grid.number_of_faces, UgGridHelpers::faceCells(grid), - grid.face_centroids, grid.cell_centroids, - grid.dimensions, props, deck, - gravity, state); - } - - template - void initStateFromDeck(int number_of_cells, - const int* global_cell, - int number_of_faces, - FaceCells face_cells, - FCI begin_face_centroids, - CCI begin_cell_centroids, - int dimensions, - const Props& props, - const EclipseGridParser& deck, - const double gravity, - State& state) - { - const int num_phases = props.numPhases(); - const PhaseUsage pu = phaseUsageFromDeck(deck); - if (num_phases != pu.num_phases) { - OPM_THROW(std::runtime_error, "initStateFromDeck(): user specified property object with " << num_phases << " phases, " - "found " << pu.num_phases << " phases in deck."); - } - if (deck.hasField("EQUIL") && deck.hasField("PRESSURE")) { - OPM_THROW(std::runtime_error, "initStateFromDeck(): The deck must either specify the initial " - "condition using the PRESSURE _or_ the EQUIL keyword (currently it has both)"); - } - state.init(number_of_cells, number_of_faces, num_phases); - if (deck.hasField("EQUIL")) { - if (num_phases != 2) { - OPM_THROW(std::runtime_error, "initStateFromDeck(): EQUIL-based init currently handling only two-phase scenarios."); - } - if (pu.phase_used[BlackoilPhases::Vapour]) { - OPM_THROW(std::runtime_error, "initStateFromDeck(): EQUIL-based init currently handling only oil-water scenario (no gas)."); - } - // Set saturations depending on oil-water contact. - const EQUIL& equil= deck.getEQUIL(); - if (equil.equil.size() != 1) { - OPM_THROW(std::runtime_error, "initStateFromDeck(): No region support yet."); - } - const double woc = equil.equil[0].water_oil_contact_depth_; - initWaterOilContact(number_of_cells, begin_cell_centroids, dimensions, - props, woc, WaterBelow, state); - // Set pressure depending on densities and depths. - const double datum_z = equil.equil[0].datum_depth_; - const double datum_p = equil.equil[0].datum_depth_pressure_; - initHydrostaticPressure(number_of_cells, begin_cell_centroids, dimensions, - props, woc, gravity, datum_z, datum_p, state); - } else if (deck.hasField("PRESSURE")) { - // Set saturations from SWAT/SGAS, pressure from PRESSURE. - std::vector& s = state.saturation(); - std::vector& p = state.pressure(); - const std::vector& p_deck = deck.getFloatingPointValue("PRESSURE"); - if (num_phases == 2) { - if (!pu.phase_used[BlackoilPhases::Aqua]) { - // oil-gas: we require SGAS - if (!deck.hasField("SGAS")) { - OPM_THROW(std::runtime_error, "initStateFromDeck(): missing SGAS keyword in 2-phase init"); - } - const std::vector& sg_deck = deck.getFloatingPointValue("SGAS"); - const int gpos = pu.phase_pos[BlackoilPhases::Vapour]; - const int opos = pu.phase_pos[BlackoilPhases::Liquid]; - for (int c = 0; c < number_of_cells; ++c) { - int c_deck = (global_cell == NULL) ? c : global_cell[c]; - s[2*c + gpos] = sg_deck[c_deck]; - s[2*c + opos] = 1.0 - sg_deck[c_deck]; - p[c] = p_deck[c_deck]; - } - } else { - // water-oil or water-gas: we require SWAT - if (!deck.hasField("SWAT")) { - OPM_THROW(std::runtime_error, "initStateFromDeck(): missing SWAT keyword in 2-phase init"); - } - const std::vector& sw_deck = deck.getFloatingPointValue("SWAT"); - const int wpos = pu.phase_pos[BlackoilPhases::Aqua]; - const int nwpos = (wpos + 1) % 2; - for (int c = 0; c < number_of_cells; ++c) { - int c_deck = (global_cell == NULL) ? c : global_cell[c]; - s[2*c + wpos] = sw_deck[c_deck]; - s[2*c + nwpos] = 1.0 - sw_deck[c_deck]; - p[c] = p_deck[c_deck]; - } - } - } else if (num_phases == 3) { - const bool has_swat_sgas = deck.hasField("SWAT") && deck.hasField("SGAS"); - if (!has_swat_sgas) { - OPM_THROW(std::runtime_error, "initStateFromDeck(): missing SGAS or SWAT keyword in 3-phase init."); - } - const int wpos = pu.phase_pos[BlackoilPhases::Aqua]; - const int gpos = pu.phase_pos[BlackoilPhases::Vapour]; - const int opos = pu.phase_pos[BlackoilPhases::Liquid]; - const std::vector& sw_deck = deck.getFloatingPointValue("SWAT"); - const std::vector& sg_deck = deck.getFloatingPointValue("SGAS"); - for (int c = 0; c < number_of_cells; ++c) { - int c_deck = (global_cell == NULL) ? c : global_cell[c]; - s[3*c + wpos] = sw_deck[c_deck]; - s[3*c + opos] = 1.0 - (sw_deck[c_deck] + sg_deck[c_deck]); - s[3*c + gpos] = sg_deck[c_deck]; - p[c] = p_deck[c_deck]; - } - } else { - OPM_THROW(std::runtime_error, "initStateFromDeck(): init with SWAT etc. only available with 2 or 3 phases."); - } - } else { - OPM_THROW(std::runtime_error, "initStateFromDeck(): we must either have EQUIL, or PRESSURE and SWAT/SOIL/SGAS."); - } - - // Finally, init face pressures. - initFacePressure(dimensions, number_of_faces, face_cells, begin_face_centroids, - begin_cell_centroids, state); - } - /// Initialize surface volume from pressure and saturation by z = As. /// Here saturation is used as an intial guess for z in the /// computation of A. @@ -975,60 +851,6 @@ namespace Opm } } - - /// Initialize a blackoil state from input deck. - template - void initBlackoilStateFromDeck(const UnstructuredGrid& grid, - const Props& props, - const EclipseGridParser& deck, - const double gravity, - State& state) - { - initBlackoilStateFromDeck(grid.number_of_cells, grid.global_cell, - grid.number_of_faces, UgGridHelpers::faceCells(grid), - grid.face_centroids, grid.cell_centroids,grid.dimensions, - props, deck, gravity, state); - } - - template - void initBlackoilStateFromDeck(int number_of_cells, - const int* global_cell, - int number_of_faces, - FaceCells face_cells, - FCI begin_face_centroids, - CCI begin_cell_centroids, - int dimensions, - const Props& props, - const EclipseGridParser& deck, - const double gravity, - State& state) - { - initStateFromDeck(number_of_cells, global_cell, number_of_faces, - face_cells, begin_face_centroids, begin_cell_centroids, - dimensions, props, deck, gravity, state); - if (deck.hasField("RS")) { - const std::vector& rs_deck = deck.getFloatingPointValue("RS"); - for (int c = 0; c < number_of_cells; ++c) { - int c_deck = (global_cell == NULL) ? c : global_cell[c]; - state.gasoilratio()[c] = rs_deck[c_deck]; - } - initBlackoilSurfvolUsingRSorRV(number_of_cells, props, state); - computeSaturation(props,state); - } else if (deck.hasField("RV")){ - const std::vector& rv_deck = deck.getFloatingPointValue("RV"); - for (int c = 0; c < number_of_cells; ++c) { - int c_deck = (global_cell == NULL) ? c : global_cell[c]; - state.rv()[c] = rv_deck[c_deck]; - } - initBlackoilSurfvolUsingRSorRV(number_of_cells, props, state); - computeSaturation(props,state); - } - - else { - OPM_THROW(std::runtime_error, "Temporarily, we require the RS or the RV field."); - } - } - /// Initialize a blackoil state from input deck. template void initBlackoilStateFromDeck(const UnstructuredGrid& grid, @@ -1084,7 +906,6 @@ namespace Opm } } - } // namespace Opm diff --git a/opm/core/wells/WellCollection.cpp b/opm/core/wells/WellCollection.cpp index 9218db06e..e7688c974 100644 --- a/opm/core/wells/WellCollection.cpp +++ b/opm/core/wells/WellCollection.cpp @@ -79,50 +79,6 @@ namespace Opm child->setParent(parent); } - - void WellCollection::addChild(const std::string& child_name, - const std::string& parent_name, - const EclipseGridParser& deck) - { - WellsGroupInterface* parent = findNode(parent_name); - if (!parent) { - roots_.push_back(createWellsGroup(parent_name, deck)); - parent = roots_[roots_.size() - 1].get(); - } - - std::shared_ptr child; - - for (size_t i = 0; i < roots_.size(); ++i) { - if (roots_[i]->name() == child_name) { - child = roots_[i]; - // We've found a new parent to the previously thought root, need to remove it - for(size_t j = i; j < roots_.size() - 1; ++j) { - roots_[j] = roots_[j+1]; - } - - roots_.resize(roots_.size()-1); - break; - } - } - - if (!child.get()) { - child = createWellsGroup(child_name, deck); - } - - WellsGroup* parent_as_group = static_cast (parent); - if (!parent_as_group) { - OPM_THROW(std::runtime_error, "Trying to add child to group named " << parent_name << ", but it's not a group."); - } - parent_as_group->addChild(child); - - if(child->isLeafNode()) { - leaf_nodes_.push_back(static_cast(child.get())); - } - - child->setParent(parent); - } - - const std::vector& WellCollection::getLeafNodes() const { return leaf_nodes_; } diff --git a/opm/core/wells/WellCollection.hpp b/opm/core/wells/WellCollection.hpp index f084a04d8..b597c4b64 100644 --- a/opm/core/wells/WellCollection.hpp +++ b/opm/core/wells/WellCollection.hpp @@ -16,9 +16,6 @@ You should have received a copy of the GNU General Public License along with OPM. If not, see . */ - - - #ifndef OPM_WELLCOLLECTION_HPP #define OPM_WELLCOLLECTION_HPP @@ -27,7 +24,6 @@ #include #include -#include #include #include @@ -47,16 +43,6 @@ namespace Opm void addGroup(GroupConstPtr groupChild, std::string parent_name, size_t timeStep, const PhaseUsage& phaseUsage); - /// Adds and creates if necessary the child to the collection - /// and appends it to parent's children. Also adds and creates the parent - /// if necessary. - /// \param[in] child name of child node - /// \param[in] parent name of parent node - /// \param[in] deck deck from which we will extract group control data - void addChild(const std::string& child, - const std::string& parent, - const EclipseGridParser& deck); - /// Adds the child to the collection /// and appends it to parent's children. /// \param[in] child the child node diff --git a/opm/core/wells/WellsGroup.cpp b/opm/core/wells/WellsGroup.cpp index 050e18712..c4e1a1243 100644 --- a/opm/core/wells/WellsGroup.cpp +++ b/opm/core/wells/WellsGroup.cpp @@ -1045,101 +1045,6 @@ namespace Opm } } // anonymous namespace - std::shared_ptr createWellsGroup(const std::string& name, - const EclipseGridParser& deck) - { - PhaseUsage phase_usage = phaseUsageFromDeck(deck); - - std::shared_ptr return_value; - // First we need to determine whether it's a group or just a well: - bool isWell = false; - if (deck.hasField("WELSPECS")) { - WELSPECS wspecs = deck.getWELSPECS(); - for (size_t i = 0; i < wspecs.welspecs.size(); i++) { - if (wspecs.welspecs[i].name_ == name) { - isWell = true; - break; - } - } - } - // For now, assume that if it isn't a well, it's a group - - if (isWell) { - ProductionSpecification production_specification; - InjectionSpecification injection_specification; - if (deck.hasField("WCONINJE")) { - WCONINJE wconinje = deck.getWCONINJE(); - for (size_t i = 0; i < wconinje.wconinje.size(); i++) { - if (wconinje.wconinje[i].well_ == name) { - WconinjeLine line = wconinje.wconinje[i]; - injection_specification.BHP_limit_ = line.BHP_limit_; - injection_specification.injector_type_ = toInjectorType(line.injector_type_); - injection_specification.control_mode_ = toInjectionControlMode(line.control_mode_); - injection_specification.surface_flow_max_rate_ = line.surface_flow_max_rate_; - injection_specification.reservoir_flow_max_rate_ = line.reservoir_flow_max_rate_; - production_specification.guide_rate_ = 0.0; // We know we're not a producer - } - } - } - - if (deck.hasField("WCONPROD")) { - WCONPROD wconprod = deck.getWCONPROD(); - for (size_t i = 0; i < wconprod.wconprod.size(); i++) { - if (wconprod.wconprod[i].well_ == name) { - WconprodLine line = wconprod.wconprod[i]; - production_specification.BHP_limit_ = line.BHP_limit_; - production_specification.reservoir_flow_max_rate_ = line.reservoir_flow_max_rate_; - production_specification.oil_max_rate_ = line.oil_max_rate_; - production_specification.control_mode_ = toProductionControlMode(line.control_mode_); - production_specification.water_max_rate_ = line.water_max_rate_; - injection_specification.guide_rate_ = 0.0; // we know we're not an injector - } - } - } - return_value.reset(new WellNode(name, production_specification, injection_specification, phase_usage)); - } else { - InjectionSpecification injection_specification; - if (deck.hasField("GCONINJE")) { - GCONINJE gconinje = deck.getGCONINJE(); - for (size_t i = 0; i < gconinje.gconinje.size(); i++) { - if (gconinje.gconinje[i].group_ == name) { - GconinjeLine line = gconinje.gconinje[i]; - injection_specification.injector_type_ = toInjectorType(line.injector_type_); - injection_specification.control_mode_ = toInjectionControlMode(line.control_mode_); - injection_specification.surface_flow_max_rate_ = line.surface_flow_max_rate_; - injection_specification.reservoir_flow_max_rate_ = line.resv_flow_max_rate_; - injection_specification.reinjection_fraction_target_ = line.reinjection_fraction_target_; - injection_specification.voidage_replacment_fraction_ = line.voidage_replacement_fraction_; - } - } - } - - ProductionSpecification production_specification; - if (deck.hasField("GCONPROD")) { - std::cout << "Searching in gconprod " << std::endl; - std::cout << "name= " << name << std::endl; - GCONPROD gconprod = deck.getGCONPROD(); - for (size_t i = 0; i < gconprod.gconprod.size(); i++) { - if (gconprod.gconprod[i].group_ == name) { - GconprodLine line = gconprod.gconprod[i]; - production_specification.oil_max_rate_ = line.oil_max_rate_; - std::cout << "control_mode = " << line.control_mode_ << std::endl; - production_specification.control_mode_ = toProductionControlMode(line.control_mode_); - production_specification.water_max_rate_ = line.water_max_rate_; - production_specification.gas_max_rate_ = line.gas_max_rate_; - production_specification.liquid_max_rate_ = line.liquid_max_rate_; - production_specification.procedure_ = toProductionProcedure(line.procedure_); - production_specification.reservoir_flow_max_rate_ = line.resv_max_rate_; - } - } - } - - return_value.reset(new WellsGroup(name, production_specification, injection_specification, phase_usage)); - } - - return return_value; - } - std::shared_ptr createGroupWellsGroup(GroupConstPtr group, size_t timeStep, const PhaseUsage& phase_usage ) { InjectionSpecification injection_specification; diff --git a/opm/core/wells/WellsGroup.hpp b/opm/core/wells/WellsGroup.hpp index 03a61d5e1..697a23cbc 100644 --- a/opm/core/wells/WellsGroup.hpp +++ b/opm/core/wells/WellsGroup.hpp @@ -22,7 +22,6 @@ #include #include -#include #include #include @@ -403,12 +402,6 @@ namespace Opm bool shut_well_; }; - /// Creates the WellsGroupInterface for the given name - /// \param[in] name the name of the wells group. - /// \param[in] deck the deck from which to fetch information. - std::shared_ptr createWellsGroup(const std::string& name, - const EclipseGridParser& deck); - /// Creates the WellsGroupInterface for the given well /// \param[in] well the Well to construct object for /// \param[in] timeStep the time step in question diff --git a/opm/core/wells/WellsManager.cpp b/opm/core/wells/WellsManager.cpp index 8166801ea..2688c4fc8 100644 --- a/opm/core/wells/WellsManager.cpp +++ b/opm/core/wells/WellsManager.cpp @@ -21,7 +21,6 @@ #include -#include #include #include #include diff --git a/opm/core/wells/WellsManager.hpp b/opm/core/wells/WellsManager.hpp index 44112e72e..8e7220ffe 100644 --- a/opm/core/wells/WellsManager.hpp +++ b/opm/core/wells/WellsManager.hpp @@ -33,8 +33,6 @@ struct UnstructuredGrid; namespace Opm { - - class EclipseGridParser; struct WellData { WellType type; diff --git a/tests/test_equil.cpp b/tests/test_equil.cpp index 05f8df8d7..e94f5d267 100644 --- a/tests/test_equil.cpp +++ b/tests/test_equil.cpp @@ -337,7 +337,7 @@ BOOST_AUTO_TEST_CASE (DeckAllDead) Opm::ParserPtr parser(new Opm::Parser() ); Opm::DeckConstPtr deck = parser->parseFile("deadfluids.DATA"); Opm::BlackoilPropertiesFromDeck props(deck, *grid, false); - Opm::Equil::DeckDependent::InitialStateComputer comp(props, deck, *grid, 10.0); + Opm::Equil::DeckDependent::InitialStateComputer comp(props, deck, *grid, 10.0); const auto& pressures = comp.press(); BOOST_REQUIRE(pressures.size() == 3); BOOST_REQUIRE(int(pressures[0].size()) == grid->number_of_cells); @@ -416,7 +416,7 @@ BOOST_AUTO_TEST_CASE (DeckWithCapillary) Opm::DeckConstPtr deck = parser->parseFile("capillary.DATA"); Opm::BlackoilPropertiesFromDeck props(deck, grid, false); - Opm::Equil::DeckDependent::InitialStateComputer comp(props, deck, grid, 10.0); + Opm::Equil::DeckDependent::InitialStateComputer comp(props, deck, grid, 10.0); const auto& pressures = comp.press(); BOOST_REQUIRE(pressures.size() == 3); BOOST_REQUIRE(int(pressures[0].size()) == grid.number_of_cells); @@ -455,7 +455,7 @@ BOOST_AUTO_TEST_CASE (DeckWithCapillaryOverlap) Opm::DeckConstPtr deck = parser->parseFile("capillary_overlap.DATA"); Opm::BlackoilPropertiesFromDeck props(deck, grid, false); - Opm::Equil::DeckDependent::InitialStateComputer comp(props, deck, grid, 9.80665); + Opm::Equil::DeckDependent::InitialStateComputer comp(props, deck, grid, 9.80665); const auto& pressures = comp.press(); BOOST_REQUIRE(pressures.size() == 3); BOOST_REQUIRE(int(pressures[0].size()) == grid.number_of_cells); @@ -516,7 +516,7 @@ BOOST_AUTO_TEST_CASE (DeckWithLiveOil) Opm::DeckConstPtr deck = parser->parseFile("equil_liveoil.DATA"); Opm::BlackoilPropertiesFromDeck props(deck, grid, false); - Opm::Equil::DeckDependent::InitialStateComputer comp(props, deck, grid, 9.80665); + Opm::Equil::DeckDependent::InitialStateComputer comp(props, deck, grid, 9.80665); const auto& pressures = comp.press(); BOOST_REQUIRE(pressures.size() == 3); BOOST_REQUIRE(int(pressures[0].size()) == grid.number_of_cells); @@ -594,7 +594,7 @@ BOOST_AUTO_TEST_CASE (DeckWithLiveGas) Opm::DeckConstPtr deck = parser->parseFile("equil_livegas.DATA"); Opm::BlackoilPropertiesFromDeck props(deck, grid, false); - Opm::Equil::DeckDependent::InitialStateComputer comp(props, deck, grid, 9.80665); + Opm::Equil::DeckDependent::InitialStateComputer comp(props, deck, grid, 9.80665); const auto& pressures = comp.press(); BOOST_REQUIRE(pressures.size() == 3); BOOST_REQUIRE(int(pressures[0].size()) == grid.number_of_cells); @@ -675,7 +675,7 @@ BOOST_AUTO_TEST_CASE (DeckWithRSVDAndRVVD) Opm::DeckConstPtr deck = parser->parseFile("equil_rsvd_and_rvvd.DATA"); Opm::BlackoilPropertiesFromDeck props(deck, grid, false); - Opm::Equil::DeckDependent::InitialStateComputer comp(props, deck, grid, 9.80665); + Opm::Equil::DeckDependent::InitialStateComputer comp(props, deck, grid, 9.80665); const auto& pressures = comp.press(); BOOST_REQUIRE(pressures.size() == 3); BOOST_REQUIRE(int(pressures[0].size()) == grid.number_of_cells);