From 9ef601496459f753b575da822db1c07dce14aed5 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Halvor=20M=C3=B8ll=20Nilsen?= Date: Mon, 27 Aug 2012 12:20:03 +0200 Subject: [PATCH 01/23] Addes support for initialisation of three phases using SWOF and SGOF --- opm/core/utility/initState_impl.hpp | 28 ++++++++++++++++++++-------- 1 file changed, 20 insertions(+), 8 deletions(-) diff --git a/opm/core/utility/initState_impl.hpp b/opm/core/utility/initState_impl.hpp index c50db8a91..f66704348 100644 --- a/opm/core/utility/initState_impl.hpp +++ b/opm/core/utility/initState_impl.hpp @@ -512,11 +512,12 @@ namespace Opm State& state) { const int num_phases = props.numPhases(); - if (num_phases != 2) { - THROW("initStateFromDeck(): currently handling only two-phase scenarios."); - } + state.init(grid, num_phases); if (deck.hasField("EQUIL")) { + if (num_phases != 2) { + THROW("initStateFromDeck(): currently handling only two-phase scenarios."); + } // Set saturations depending on oil-water contact. const EQUIL& equil= deck.getEQUIL(); if (equil.equil.size() != 1) { @@ -535,11 +536,22 @@ namespace Opm const std::vector& sw_deck = deck.getFloatingPointValue("SWAT"); const std::vector& p_deck = deck.getFloatingPointValue("PRESSURE"); const int num_cells = grid.number_of_cells; - for (int c = 0; c < num_cells; ++c) { - int c_deck = (grid.global_cell == NULL) ? c : grid.global_cell[c]; - s[2*c] = sw_deck[c_deck]; - s[2*c + 1] = 1.0 - s[2*c]; - p[c] = p_deck[c_deck]; + if(num_phases == 2){ + for (int c = 0; c < num_cells; ++c) { + int c_deck = (grid.global_cell == NULL) ? c : grid.global_cell[c]; + s[2*c] = sw_deck[c_deck]; + s[2*c + 1] = 1.0 - s[2*c]; + p[c] = p_deck[c_deck]; + } + }else{ + const std::vector& sg_deck = deck.getFloatingPointValue("SGAS"); + for (int c = 0; c < num_cells; ++c) { + int c_deck = (grid.global_cell == NULL) ? c : grid.global_cell[c]; + s[2*c] = sw_deck[c_deck]; + s[2*c + 1] = 1.0 - (sw_deck[c_deck] + sg_deck[c_deck]); + s[2*c + 2] = sg_deck[c_deck]; + p[c] = p_deck[c_deck]; + } } } else { THROW("initStateFromDeck(): we must either have EQUIL, or both SWAT and PRESSURE."); From cf9b9cdf20c0150059a944ae3f7a218b2d9ca568 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Halvor=20M=C3=B8ll=20Nilsen?= Date: Mon, 27 Aug 2012 12:22:32 +0200 Subject: [PATCH 02/23] Started work on supporting 3 phases for wellreport. --- opm/core/utility/miscUtilities.cpp | 45 +++++++++++++++++++----------- 1 file changed, 28 insertions(+), 17 deletions(-) diff --git a/opm/core/utility/miscUtilities.cpp b/opm/core/utility/miscUtilities.cpp index 16c59b17e..cc7ff1276 100644 --- a/opm/core/utility/miscUtilities.cpp +++ b/opm/core/utility/miscUtilities.cpp @@ -574,9 +574,10 @@ namespace Opm { int nw = well_bhp.size(); ASSERT(nw == wells.number_of_wells); - if (props.numPhases() != 2) { - THROW("WellReport for now assumes two phase flow."); - } + int np = props.numPhases(); + //if (props.numPhases() != 2) { + // THROW("WellReport for now assumes two phase flow."); + //} const double* visc = props.viscosity(); std::vector data_now; data_now.reserve(1 + 3*nw); @@ -586,7 +587,8 @@ namespace Opm double well_rate_total = 0.0; double well_rate_water = 0.0; for (int perf = wells.well_connpos[w]; perf < wells.well_connpos[w + 1]; ++perf) { - const double perf_rate = well_perfrates[perf]*(unit::day/unit::second); + const double perf_rate = unit::convert::to(well_perfrates[perf], + unit::cubic(unit::meter)/unit::day); well_rate_total += perf_rate; if (perf_rate > 0.0) { // Injection. @@ -594,11 +596,14 @@ namespace Opm } else { // Production. const int cell = wells.well_cells[perf]; - double mob[2]; + double mob[np]; props.relperm(1, &saturation[2*cell], &cell, mob, 0); - mob[0] /= visc[0]; - mob[1] /= visc[1]; - const double fracflow = mob[0]/(mob[0] + mob[1]); + double tmob=0; + for(int i=0; i < np; ++i){ + mob[i] /= visc[i]; + tmob += mob[i]; + } + const double fracflow = mob[0]/tmob; well_rate_water += perf_rate*fracflow; } } @@ -627,9 +632,10 @@ namespace Opm // TODO: refactor, since this is almost identical to the other push(). int nw = well_bhp.size(); ASSERT(nw == wells.number_of_wells); - if (props.numPhases() != 2) { - THROW("WellReport for now assumes two phase flow."); - } + int np = props.numPhases(); + //if (props.numPhases() != 2) { + // THROW("WellReport for now assumes two phase flow."); + //} std::vector data_now; data_now.reserve(1 + 3*nw); data_now.push_back(time/unit::day); @@ -640,20 +646,25 @@ namespace Opm for (int perf = wells.well_connpos[w]; perf < wells.well_connpos[w + 1]; ++perf) { const double perf_rate = well_perfrates[perf]*(unit::day/unit::second); well_rate_total += perf_rate; - if (perf_rate > 0.0) { + if (perf_rate > 0.0) { // Injection. well_rate_water += perf_rate*wells.comp_frac[0]; } else { // Production. const int cell = wells.well_cells[perf]; - double mob[2]; + double mob[np]; props.relperm(1, &s[2*cell], &cell, mob, 0); - double visc[2]; + double visc[np]; props.viscosity(1, &p[cell], &z[2*cell], &cell, visc, 0); - mob[0] /= visc[0]; - mob[1] /= visc[1]; - const double fracflow = mob[0]/(mob[0] + mob[1]); + double tmob=0; + for(int i=0; i < np; ++i){ + mob[i] /= visc[i]; + tmob += mob[i]; + } + const double fracflow = mob[0]/(tmob); well_rate_water += perf_rate*fracflow; + //const double fracflow = mob[0]/(tmob); + //well_rate_water += perf_rate*fracflow; } } data_now.push_back(well_rate_total); From af8a6fa7da151c9a8930610f9c0489c39a070674 Mon Sep 17 00:00:00 2001 From: Xavier Raynaud Date: Mon, 27 Aug 2012 13:19:04 +0200 Subject: [PATCH 03/23] Changed private->protected and added virtual to make class derivation. --- opm/core/pressure/CompressibleTpfa.hpp | 15 ++++++++------- 1 file changed, 8 insertions(+), 7 deletions(-) diff --git a/opm/core/pressure/CompressibleTpfa.hpp b/opm/core/pressure/CompressibleTpfa.hpp index 9fc5c43fd..affada84e 100644 --- a/opm/core/pressure/CompressibleTpfa.hpp +++ b/opm/core/pressure/CompressibleTpfa.hpp @@ -90,16 +90,15 @@ namespace Opm bool singularPressure() const; private: - void computePerSolveDynamicData(const double dt, - const BlackoilState& state, - const WellState& well_state); - void computeWellPotentials(const BlackoilState& state); + virtual void computePerSolveDynamicData(const double dt, + const BlackoilState& state, + const WellState& well_state); void computePerIterationDynamicData(const double dt, const BlackoilState& state, const WellState& well_state); - void computeCellDynamicData(const double dt, - const BlackoilState& state, - const WellState& well_state); + virtual void computeCellDynamicData(const double dt, + const BlackoilState& state, + const WellState& well_state); void computeFaceDynamicData(const double dt, const BlackoilState& state, const WellState& well_state); @@ -114,6 +113,8 @@ namespace Opm double incrementNorm() const; void computeResults(BlackoilState& state, WellState& well_state) const; + protected: + void computeWellPotentials(const BlackoilState& state); // ------ Data that will remain unmodified after construction. ------ const UnstructuredGrid& grid_; From d9cff689b49f382a20ed61c28ef64180a2cea5f6 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Atgeirr=20Fl=C3=B8=20Rasmussen?= Date: Mon, 27 Aug 2012 16:48:21 +0200 Subject: [PATCH 04/23] Added class SinglePvtDead, add parameter 'props_use_spline' to simulators. Recall that the class that used to be called SinglePvtDead has been renamed to SinglePvtDeadSpline. If 'props_use_spline' is true, that class is used (this is the default), which makes a monotone spline that is uniformly, densely sampled. The new class simply uses linear interpolation in the input tables. --- opm/core/fluid/BlackoilPropertiesFromDeck.cpp | 5 +++-- opm/core/fluid/BlackoilPropertiesFromDeck.hpp | 3 ++- opm/core/fluid/blackoil/BlackoilPvtProperties.hpp | 2 +- 3 files changed, 6 insertions(+), 4 deletions(-) diff --git a/opm/core/fluid/BlackoilPropertiesFromDeck.cpp b/opm/core/fluid/BlackoilPropertiesFromDeck.cpp index 004e6f88f..e7244d87c 100644 --- a/opm/core/fluid/BlackoilPropertiesFromDeck.cpp +++ b/opm/core/fluid/BlackoilPropertiesFromDeck.cpp @@ -23,10 +23,11 @@ namespace Opm { BlackoilPropertiesFromDeck::BlackoilPropertiesFromDeck(const EclipseGridParser& deck, - const UnstructuredGrid& grid) + const UnstructuredGrid& grid, + const bool use_spline) { rock_.init(deck, grid); - pvt_.init(deck); + pvt_.init(deck, use_spline); satprops_.init(deck, grid); if (pvt_.numPhases() != satprops_.numPhases()) { THROW("BlackoilPropertiesBasic::BlackoilPropertiesBasic() - Inconsistent number of phases in pvt data (" diff --git a/opm/core/fluid/BlackoilPropertiesFromDeck.hpp b/opm/core/fluid/BlackoilPropertiesFromDeck.hpp index a2a6c3de9..c434b7f4b 100644 --- a/opm/core/fluid/BlackoilPropertiesFromDeck.hpp +++ b/opm/core/fluid/BlackoilPropertiesFromDeck.hpp @@ -43,7 +43,8 @@ namespace Opm /// mapping from cell indices (typically from a processed grid) /// to logical cartesian indices consistent with the deck. BlackoilPropertiesFromDeck(const EclipseGridParser& deck, - const UnstructuredGrid& grid); + const UnstructuredGrid& grid, + const bool use_spline); /// Destructor. virtual ~BlackoilPropertiesFromDeck(); diff --git a/opm/core/fluid/blackoil/BlackoilPvtProperties.hpp b/opm/core/fluid/blackoil/BlackoilPvtProperties.hpp index 627f266b7..ca193a4d4 100644 --- a/opm/core/fluid/blackoil/BlackoilPvtProperties.hpp +++ b/opm/core/fluid/blackoil/BlackoilPvtProperties.hpp @@ -47,7 +47,7 @@ namespace Opm BlackoilPvtProperties(); /// Initialize from deck. - void init(const EclipseGridParser& deck); + void init(const EclipseGridParser& deck, const bool use_spline); /// Number of active phases. int numPhases() const; From 7a79bd1872ff5e3207c96d7dfcb98dae0194521e Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Atgeirr=20Fl=C3=B8=20Rasmussen?= Date: Mon, 27 Aug 2012 17:56:01 +0200 Subject: [PATCH 05/23] Enable choice of spline-smoothed saturation props (or not). --- opm/core/fluid/BlackoilPropertiesFromDeck.cpp | 22 +- opm/core/fluid/BlackoilPropertiesFromDeck.hpp | 3 +- opm/core/fluid/IncompPropertiesFromDeck.hpp | 2 +- opm/core/fluid/SaturationPropsFromDeck.hpp | 45 ++- .../fluid/SaturationPropsFromDeck_impl.hpp | 271 ++++++++++++++++++ opm/core/fluid/SaturationPropsInterface.hpp | 85 ++++++ 6 files changed, 397 insertions(+), 31 deletions(-) create mode 100644 opm/core/fluid/SaturationPropsFromDeck_impl.hpp create mode 100644 opm/core/fluid/SaturationPropsInterface.hpp diff --git a/opm/core/fluid/BlackoilPropertiesFromDeck.cpp b/opm/core/fluid/BlackoilPropertiesFromDeck.cpp index e7244d87c..fdc9519cf 100644 --- a/opm/core/fluid/BlackoilPropertiesFromDeck.cpp +++ b/opm/core/fluid/BlackoilPropertiesFromDeck.cpp @@ -28,10 +28,20 @@ namespace Opm { rock_.init(deck, grid); pvt_.init(deck, use_spline); - satprops_.init(deck, grid); - if (pvt_.numPhases() != satprops_.numPhases()) { + // Unfortunate lack of pointer smartness here... + if (use_spline) { + SaturationPropsFromDeck<>* ptr = new SaturationPropsFromDeck<>(); + ptr->init(deck, grid); + satprops_.reset(ptr); + } else { + SaturationPropsFromDeck* ptr + = new SaturationPropsFromDeck(); + ptr->init(deck, grid); + satprops_.reset(ptr); + } + if (pvt_.numPhases() != satprops_->numPhases()) { THROW("BlackoilPropertiesBasic::BlackoilPropertiesBasic() - Inconsistent number of phases in pvt data (" - << pvt_.numPhases() << ") and saturation-dependent function data (" << satprops_.numPhases() << ")."); + << pvt_.numPhases() << ") and saturation-dependent function data (" << satprops_->numPhases() << ")."); } } @@ -236,7 +246,7 @@ namespace Opm double* kr, double* dkrds) const { - satprops_.relperm(n, s, cells, kr, dkrds); + satprops_->relperm(n, s, cells, kr, dkrds); } @@ -255,7 +265,7 @@ namespace Opm double* pc, double* dpcds) const { - satprops_.capPress(n, s, cells, pc, dpcds); + satprops_->capPress(n, s, cells, pc, dpcds); } @@ -271,7 +281,7 @@ namespace Opm double* smin, double* smax) const { - satprops_.satRange(n, cells, smin, smax); + satprops_->satRange(n, cells, smin, smax); } diff --git a/opm/core/fluid/BlackoilPropertiesFromDeck.hpp b/opm/core/fluid/BlackoilPropertiesFromDeck.hpp index c434b7f4b..60dc1e754 100644 --- a/opm/core/fluid/BlackoilPropertiesFromDeck.hpp +++ b/opm/core/fluid/BlackoilPropertiesFromDeck.hpp @@ -26,6 +26,7 @@ #include #include #include +#include struct UnstructuredGrid; @@ -166,7 +167,7 @@ namespace Opm private: RockFromDeck rock_; BlackoilPvtProperties pvt_; - SaturationPropsFromDeck satprops_; + boost::scoped_ptr satprops_; mutable std::vector B_; mutable std::vector dB_; mutable std::vector R_; diff --git a/opm/core/fluid/IncompPropertiesFromDeck.hpp b/opm/core/fluid/IncompPropertiesFromDeck.hpp index 68623e5cc..5820bb668 100644 --- a/opm/core/fluid/IncompPropertiesFromDeck.hpp +++ b/opm/core/fluid/IncompPropertiesFromDeck.hpp @@ -135,7 +135,7 @@ namespace Opm private: RockFromDeck rock_; PvtPropertiesIncompFromDeck pvt_; - SaturationPropsFromDeck satprops_; + SaturationPropsFromDeck<> satprops_; }; diff --git a/opm/core/fluid/SaturationPropsFromDeck.hpp b/opm/core/fluid/SaturationPropsFromDeck.hpp index b83a12711..a978bbac0 100644 --- a/opm/core/fluid/SaturationPropsFromDeck.hpp +++ b/opm/core/fluid/SaturationPropsFromDeck.hpp @@ -20,8 +20,8 @@ #ifndef OPM_SATURATIONPROPSFROMDECK_HEADER_INCLUDED #define OPM_SATURATIONPROPSFROMDECK_HEADER_INCLUDED +#include #include -#include #include #include @@ -30,7 +30,23 @@ struct UnstructuredGrid; namespace Opm { - class SaturationPropsFromDeck : public BlackoilPhases + + /// Class storing saturation functions in a uniform table, + /// densely sampled from a monotone spline, + /// using linear interpolation. + class SatFuncSetUniform; + + /// Class storing saturation functions in a nonuniform table + /// using linear interpolation. + class SatFuncSetNonuniform; + + + /// Interface to saturation functions from deck. + /// Possible values for template argument: + /// SatFuncSetNonuniform, + /// SatFuncSetUniform (default). + template + class SaturationPropsFromDeck : public SaturationPropsInterface { public: /// Default constructor. @@ -88,30 +104,12 @@ namespace Opm private: PhaseUsage phase_usage_; - class SatFuncSet - { - public: - void init(const EclipseGridParser& deck, const int table_num, PhaseUsage phase_usg); - void evalKr(const double* s, double* kr) const; - void evalKrDeriv(const double* s, double* kr, double* dkrds) const; - void evalPc(const double* s, double* pc) const; - void evalPcDeriv(const double* s, double* pc, double* dpcds) const; - double smin_[PhaseUsage::MaxNumPhases]; - double smax_[PhaseUsage::MaxNumPhases]; - private: - PhaseUsage phase_usage; // A copy of the outer class' phase_usage_. - UniformTableLinear krw_; - UniformTableLinear krow_; - UniformTableLinear pcow_; - UniformTableLinear krg_; - UniformTableLinear krog_; - UniformTableLinear pcog_; - double krocw_; // = krow_(s_wc) - }; std::vector satfuncset_; std::vector cell_to_func_; // = SATNUM - 1 - const SatFuncSet& funcForCell(const int cell) const; + typedef SatFuncSet Funcs; + + const Funcs& funcForCell(const int cell) const; }; @@ -119,6 +117,7 @@ namespace Opm } // namespace Opm +#include #endif // OPM_SATURATIONPROPSFROMDECK_HEADER_INCLUDED diff --git a/opm/core/fluid/SaturationPropsFromDeck_impl.hpp b/opm/core/fluid/SaturationPropsFromDeck_impl.hpp new file mode 100644 index 000000000..f64f12bcf --- /dev/null +++ b/opm/core/fluid/SaturationPropsFromDeck_impl.hpp @@ -0,0 +1,271 @@ +/* + Copyright 2012 SINTEF ICT, Applied Mathematics. + + This file is part of the Open Porous Media project (OPM). + + OPM is free software: you can redistribute it and/or modify + it under the terms of the GNU General Public License as published by + the Free Software Foundation, either version 3 of the License, or + (at your option) any later version. + + OPM is distributed in the hope that it will be useful, + but WITHOUT ANY WARRANTY; without even the implied warranty of + MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + GNU General Public License for more details. + + You should have received a copy of the GNU General Public License + along with OPM. If not, see . +*/ + +#ifndef OPM_SATURATIONPROPSFROMDECK_IMPL_HEADER_INCLUDED +#define OPM_SATURATIONPROPSFROMDECK_IMPL_HEADER_INCLUDED + + +#include +#include +#include +#include + +namespace Opm +{ + + /// Class storing saturation functions in a uniform table, + /// densely sampled from a monotone spline, + /// using linear interpolation. + class SatFuncSetUniform : public BlackoilPhases + { + public: + void init(const EclipseGridParser& deck, const int table_num, PhaseUsage phase_usg); + void evalKr(const double* s, double* kr) const; + void evalKrDeriv(const double* s, double* kr, double* dkrds) const; + void evalPc(const double* s, double* pc) const; + void evalPcDeriv(const double* s, double* pc, double* dpcds) const; + double smin_[PhaseUsage::MaxNumPhases]; + double smax_[PhaseUsage::MaxNumPhases]; + private: + PhaseUsage phase_usage; // A copy of the outer class' phase_usage_. + UniformTableLinear krw_; + UniformTableLinear krow_; + UniformTableLinear pcow_; + UniformTableLinear krg_; + UniformTableLinear krog_; + UniformTableLinear pcog_; + double krocw_; // = krow_(s_wc) + }; + + + + + /// Class storing saturation functions in a nonuniform table + /// using linear interpolation. + class SatFuncSetNonuniform : public BlackoilPhases + { + public: + void init(const EclipseGridParser& deck, const int table_num, PhaseUsage phase_usg); + void evalKr(const double* s, double* kr) const; + void evalKrDeriv(const double* s, double* kr, double* dkrds) const; + void evalPc(const double* s, double* pc) const; + void evalPcDeriv(const double* s, double* pc, double* dpcds) const; + double smin_[PhaseUsage::MaxNumPhases]; + double smax_[PhaseUsage::MaxNumPhases]; + private: + PhaseUsage phase_usage; // A copy of the outer class' phase_usage_. + NonuniformTableLinear krw_; + NonuniformTableLinear krow_; + NonuniformTableLinear pcow_; + NonuniformTableLinear krg_; + NonuniformTableLinear krog_; + NonuniformTableLinear pcog_; + double krocw_; // = krow_(s_wc) + }; + + + + // ----------- Methods of SaturationPropsFromDeck --------- + + + /// Default constructor. + template + SaturationPropsFromDeck::SaturationPropsFromDeck() + { + } + + /// Initialize from deck. + template + void SaturationPropsFromDeck::init(const EclipseGridParser& deck, + const UnstructuredGrid& grid) + { + phase_usage_ = phaseUsageFromDeck(deck); + + // Extract input data. + // Oil phase should be active. + if (!phase_usage_.phase_used[Liquid]) { + THROW("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()); + const int num_cells = grid.number_of_cells; + cell_to_func_.resize(num_cells); + const int* gc = grid.global_cell; + for (int cell = 0; cell < num_cells; ++cell) { + const int deck_pos = (gc == NULL) ? cell : gc[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) { + THROW("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) { + THROW("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) { + THROW("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_); + } + } + + + + + /// \return P, the number of phases. + template + int SaturationPropsFromDeck::numPhases() const + { + return phase_usage_.num_phases; + } + + + + + /// Relative permeability. + /// \param[in] n Number of data points. + /// \param[in] s Array of nP saturation values. + /// \param[in] cells Array of n cell indices to be associated with the s values. + /// \param[out] kr Array of nP relperm values, array must be valid before calling. + /// \param[out] dkrds If non-null: array of nP^2 relperm derivative values, + /// array must be valid before calling. + /// The P^2 derivative matrix is + /// m_{ij} = \frac{dkr_i}{ds^j}, + /// and is output in Fortran order (m_00 m_10 m_20 m01 ...) + template + void SaturationPropsFromDeck::relperm(const int n, + const double* s, + const int* cells, + double* kr, + double* dkrds) const + { + ASSERT (cells != 0); + + const int np = phase_usage_.num_phases; + if (dkrds) { +// #pragma omp parallel for + for (int i = 0; i < n; ++i) { + funcForCell(cells[i]).evalKrDeriv(s + np*i, kr + np*i, dkrds + np*np*i); + } + } else { +// #pragma omp parallel for + for (int i = 0; i < n; ++i) { + funcForCell(cells[i]).evalKr(s + np*i, kr + np*i); + } + } + } + + + + + /// Capillary pressure. + /// \param[in] n Number of data points. + /// \param[in] s Array of nP saturation values. + /// \param[in] cells Array of n cell indices to be associated with the s values. + /// \param[out] pc Array of nP capillary pressure values, array must be valid before calling. + /// \param[out] dpcds If non-null: array of nP^2 derivative values, + /// array must be valid before calling. + /// The P^2 derivative matrix is + /// m_{ij} = \frac{dpc_i}{ds^j}, + /// and is output in Fortran order (m_00 m_10 m_20 m01 ...) + template + void SaturationPropsFromDeck::capPress(const int n, + const double* s, + const int* cells, + double* pc, + double* dpcds) const + { + ASSERT (cells != 0); + + const int np = phase_usage_.num_phases; + if (dpcds) { +// #pragma omp parallel for + for (int i = 0; i < n; ++i) { + funcForCell(cells[i]).evalPcDeriv(s + np*i, pc + np*i, dpcds + np*np*i); + } + } else { +// #pragma omp parallel for + for (int i = 0; i < n; ++i) { + funcForCell(cells[i]).evalPc(s + np*i, pc + np*i); + } + } + } + + + + + /// Obtain the range of allowable saturation values. + /// \param[in] n Number of data points. + /// \param[in] cells Array of n cell indices. + /// \param[out] smin Array of nP minimum s values, array must be valid before calling. + /// \param[out] smax Array of nP maximum s values, array must be valid before calling. + template + void SaturationPropsFromDeck::satRange(const int n, + const int* cells, + double* smin, + double* smax) const + { + ASSERT (cells != 0); + + const int np = phase_usage_.num_phases; + for (int i = 0; i < n; ++i) { + for (int p = 0; p < np; ++p) { + smin[np*i + p] = funcForCell(cells[i]).smin_[p]; + smax[np*i + p] = funcForCell(cells[i]).smax_[p]; + } + } + } + + + // Map the cell number to the correct function set. + template + const typename SaturationPropsFromDeck::Funcs& + SaturationPropsFromDeck::funcForCell(const int cell) const + { + return cell_to_func_.empty() ? satfuncset_[0] : satfuncset_[cell_to_func_[cell]]; + } + + + +} // namespace Opm + +#endif // OPM_SATURATIONPROPSFROMDECK_IMPL_HEADER_INCLUDED diff --git a/opm/core/fluid/SaturationPropsInterface.hpp b/opm/core/fluid/SaturationPropsInterface.hpp new file mode 100644 index 000000000..e55995677 --- /dev/null +++ b/opm/core/fluid/SaturationPropsInterface.hpp @@ -0,0 +1,85 @@ +/* + Copyright 2012 SINTEF ICT, Applied Mathematics. + + This file is part of the Open Porous Media project (OPM). + + OPM is free software: you can redistribute it and/or modify + it under the terms of the GNU General Public License as published by + the Free Software Foundation, either version 3 of the License, or + (at your option) any later version. + + OPM is distributed in the hope that it will be useful, + but WITHOUT ANY WARRANTY; without even the implied warranty of + MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + GNU General Public License for more details. + + You should have received a copy of the GNU General Public License + along with OPM. If not, see . +*/ + +#ifndef OPM_SATURATIONPROPSINTERFACE_HEADER_INCLUDED +#define OPM_SATURATIONPROPSINTERFACE_HEADER_INCLUDED + +#include + + +namespace Opm +{ + + class SaturationPropsInterface : public BlackoilPhases + { + public: + /// Virtual destructor. + virtual ~SaturationPropsInterface() {}; + + /// \return P, the number of phases. + virtual int numPhases() const = 0; + + /// Relative permeability. + /// \param[in] n Number of data points. + /// \param[in] s Array of nP saturation values. + /// \param[out] kr Array of nP relperm values, array must be valid before calling. + /// \param[out] dkrds If non-null: array of nP^2 relperm derivative values, + /// array must be valid before calling. + /// The P^2 derivative matrix is + /// m_{ij} = \frac{dkr_i}{ds^j}, + /// and is output in Fortran order (m_00 m_10 m_20 m01 ...) + virtual void relperm(const int n, + const double* s, + const int* cells, + double* kr, + double* dkrds) const = 0; + + /// Capillary pressure. + /// \param[in] n Number of data points. + /// \param[in] s Array of nP saturation values. + /// \param[out] pc Array of nP capillary pressure values, array must be valid before calling. + /// \param[out] dpcds If non-null: array of nP^2 derivative values, + /// array must be valid before calling. + /// The P^2 derivative matrix is + /// m_{ij} = \frac{dpc_i}{ds^j}, + /// and is output in Fortran order (m_00 m_10 m_20 m01 ...) + virtual void capPress(const int n, + const double* s, + const int* cells, + double* pc, + double* dpcds) const = 0; + + /// Obtain the range of allowable saturation values. + /// \param[in] n Number of data points. + /// \param[out] smin Array of nP minimum s values, array must be valid before calling. + /// \param[out] smax Array of nP maximum s values, array must be valid before calling. + virtual void satRange(const int n, + const int* cells, + double* smin, + double* smax) const = 0; + + }; + + + +} // namespace Opm + + + +#endif // OPM_SATURATIONPROPSINTERFACE_HEADER_INCLUDED From 9aa349b1e50040c20d9a6373e369888026441d62 Mon Sep 17 00:00:00 2001 From: Xavier Raynaud Date: Tue, 28 Aug 2012 13:12:54 +0200 Subject: [PATCH 06/23] Corrected comment. --- .../transport/reorder/TransportModelCompressibleTwophase.cpp | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/opm/core/transport/reorder/TransportModelCompressibleTwophase.cpp b/opm/core/transport/reorder/TransportModelCompressibleTwophase.cpp index 697f94a57..01882053b 100644 --- a/opm/core/transport/reorder/TransportModelCompressibleTwophase.cpp +++ b/opm/core/transport/reorder/TransportModelCompressibleTwophase.cpp @@ -124,7 +124,7 @@ namespace Opm // // [[ incompressible was: r(s) = s - s0 + dt/pv*( influx + outflux*f(s) ) ]] // - // r(s) = s - B*z0 + dt/pv*( influx + outflux*f(s) ) + // r(s) = s - B*z0 + s*(poro - poro0)/poro0 + dt/pv*( influx + outflux*f(s) ) // // @@@ What about the source term // From 94c04f343d5187701bc1461fa90825af69ca9849 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Halvor=20M=C3=B8ll=20Nilsen?= Date: Tue, 28 Aug 2012 14:27:14 +0200 Subject: [PATCH 07/23] Added param initializer on fluids with param to addjust table length. Moved internal class SatFunc to SatFuncStone2. --- opm/core/fluid/BlackoilPropertiesFromDeck.cpp | 16 ++++++++- opm/core/fluid/BlackoilPropertiesFromDeck.hpp | 7 ++-- opm/core/fluid/SaturationPropsFromDeck.hpp | 33 +++++-------------- 3 files changed, 29 insertions(+), 27 deletions(-) diff --git a/opm/core/fluid/BlackoilPropertiesFromDeck.cpp b/opm/core/fluid/BlackoilPropertiesFromDeck.cpp index 004e6f88f..ae6a61306 100644 --- a/opm/core/fluid/BlackoilPropertiesFromDeck.cpp +++ b/opm/core/fluid/BlackoilPropertiesFromDeck.cpp @@ -18,7 +18,7 @@ */ #include - +#include namespace Opm { @@ -34,6 +34,20 @@ namespace Opm } } + BlackoilPropertiesFromDeck::BlackoilPropertiesFromDeck(const EclipseGridParser& deck, + const UnstructuredGrid& grid, + const parameter::ParameterGroup& param) + { + rock_.init(deck, grid); + pvt_.init(deck); + satprops_.init(deck, grid, param); + + if (pvt_.numPhases() != satprops_.numPhases()) { + THROW("BlackoilPropertiesBasic::BlackoilPropertiesBasic() - Inconsistent number of phases in pvt data (" + << pvt_.numPhases() << ") and saturation-dependent function data (" << satprops_.numPhases() << ")."); + } + } + BlackoilPropertiesFromDeck::~BlackoilPropertiesFromDeck() { } diff --git a/opm/core/fluid/BlackoilPropertiesFromDeck.hpp b/opm/core/fluid/BlackoilPropertiesFromDeck.hpp index a2a6c3de9..1d0a0d380 100644 --- a/opm/core/fluid/BlackoilPropertiesFromDeck.hpp +++ b/opm/core/fluid/BlackoilPropertiesFromDeck.hpp @@ -26,6 +26,7 @@ #include #include #include +#include struct UnstructuredGrid; @@ -43,8 +44,10 @@ namespace Opm /// mapping from cell indices (typically from a processed grid) /// to logical cartesian indices consistent with the deck. BlackoilPropertiesFromDeck(const EclipseGridParser& deck, - const UnstructuredGrid& grid); - + const UnstructuredGrid& grid); + BlackoilPropertiesFromDeck(const EclipseGridParser& deck, + const UnstructuredGrid& grid, + const parameter::ParameterGroup& param); /// Destructor. virtual ~BlackoilPropertiesFromDeck(); diff --git a/opm/core/fluid/SaturationPropsFromDeck.hpp b/opm/core/fluid/SaturationPropsFromDeck.hpp index b83a12711..2c6a89a06 100644 --- a/opm/core/fluid/SaturationPropsFromDeck.hpp +++ b/opm/core/fluid/SaturationPropsFromDeck.hpp @@ -19,10 +19,11 @@ #ifndef OPM_SATURATIONPROPSFROMDECK_HEADER_INCLUDED #define OPM_SATURATIONPROPSFROMDECK_HEADER_INCLUDED - +#include #include #include #include +#include #include struct UnstructuredGrid; @@ -44,6 +45,10 @@ namespace Opm void init(const EclipseGridParser& deck, const UnstructuredGrid& grid); + void init(const EclipseGridParser& deck, + const UnstructuredGrid& grid, + const parameter::ParameterGroup& param); + /// \return P, the number of phases. int numPhases() const; @@ -87,31 +92,11 @@ namespace Opm double* smax) const; private: - PhaseUsage phase_usage_; - class SatFuncSet - { - public: - void init(const EclipseGridParser& deck, const int table_num, PhaseUsage phase_usg); - void evalKr(const double* s, double* kr) const; - void evalKrDeriv(const double* s, double* kr, double* dkrds) const; - void evalPc(const double* s, double* pc) const; - void evalPcDeriv(const double* s, double* pc, double* dpcds) const; - double smin_[PhaseUsage::MaxNumPhases]; - double smax_[PhaseUsage::MaxNumPhases]; - private: - PhaseUsage phase_usage; // A copy of the outer class' phase_usage_. - UniformTableLinear krw_; - UniformTableLinear krow_; - UniformTableLinear pcow_; - UniformTableLinear krg_; - UniformTableLinear krog_; - UniformTableLinear pcog_; - double krocw_; // = krow_(s_wc) - }; - std::vector satfuncset_; + PhaseUsage phase_usage_; + std::vector satfuncset_; std::vector cell_to_func_; // = SATNUM - 1 - const SatFuncSet& funcForCell(const int cell) const; + const SatFuncStone2& funcForCell(const int cell) const; }; From 6852be422c458f6deb43ef6fc22eba7b561a02e5 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Halvor=20M=C3=B8ll=20Nilsen?= Date: Tue, 28 Aug 2012 16:41:06 +0200 Subject: [PATCH 08/23] Added new SatFuncSimple fluid. Introduced a simple fluid which has no problem with strange black oil behavior. Intended for testing, but for now it is used in SaturationPropsFromDeck. --- opm/core/fluid/SaturationPropsFromDeck.hpp | 8 +++++--- 1 file changed, 5 insertions(+), 3 deletions(-) diff --git a/opm/core/fluid/SaturationPropsFromDeck.hpp b/opm/core/fluid/SaturationPropsFromDeck.hpp index 2c6a89a06..bd70d8029 100644 --- a/opm/core/fluid/SaturationPropsFromDeck.hpp +++ b/opm/core/fluid/SaturationPropsFromDeck.hpp @@ -24,6 +24,7 @@ #include #include #include +#include #include struct UnstructuredGrid; @@ -92,11 +93,12 @@ namespace Opm double* smax) const; private: - PhaseUsage phase_usage_; - std::vector satfuncset_; + PhaseUsage phase_usage_; + typedef SatFuncSimple satfunc_t; + std::vector satfuncset_; std::vector cell_to_func_; // = SATNUM - 1 - const SatFuncStone2& funcForCell(const int cell) const; + const satfunc_t& funcForCell(const int cell) const; }; From b4a4e297a884632a8992a162239c85f62c72a7f8 Mon Sep 17 00:00:00 2001 From: Xavier Raynaud Date: Wed, 29 Aug 2012 10:54:21 +0200 Subject: [PATCH 09/23] Added functionality to compute porosity for compressible fluid. --- opm/core/utility/miscUtilities.cpp | 19 +++++++++++++++++++ opm/core/utility/miscUtilities.hpp | 11 +++++++++++ 2 files changed, 30 insertions(+) diff --git a/opm/core/utility/miscUtilities.cpp b/opm/core/utility/miscUtilities.cpp index 16c59b17e..ac9552dc3 100644 --- a/opm/core/utility/miscUtilities.cpp +++ b/opm/core/utility/miscUtilities.cpp @@ -70,6 +70,25 @@ namespace Opm } } + /// @brief Computes porosity of all cells in a grid, with rock compressibility effects. + /// @param[in] grid a grid + /// @param[in] porosity_standard array of grid.number_of_cells porosity values (at standard conditions) + /// @param[in] rock_comp rock compressibility properties + /// @param[in] pressure pressure by cell + /// @param[out] porosity porosity (at reservoir condition) + void computePorosity(const UnstructuredGrid& grid, + const double* porosity_standard, + const RockCompressibility& rock_comp, + const std::vector& pressure, + std::vector& porosity) + { + int num_cells = grid.number_of_cells; + porosity.resize(num_cells); + for (int i = 0; i < num_cells; ++i) { + porosity[i] = porosity_standard[i]*rock_comp.poroMult(pressure[i]); + } + } + /// @brief Computes total saturated volumes over all grid cells. /// @param[in] pv the pore volume by cell. diff --git a/opm/core/utility/miscUtilities.hpp b/opm/core/utility/miscUtilities.hpp index 37b252bc2..de452fbac 100644 --- a/opm/core/utility/miscUtilities.hpp +++ b/opm/core/utility/miscUtilities.hpp @@ -54,6 +54,17 @@ namespace Opm const std::vector& pressure, std::vector& porevol); + /// @brief Computes porosity of all cells in a grid, with rock compressibility effects. + /// @param[in] grid a grid + /// @param[in] porosity_standard array of grid.number_of_cells porosity values (at standard conditions) + /// @param[in] rock_comp rock compressibility properties + /// @param[in] pressure pressure by cell + /// @param[out] porosity porosity (at reservoir condition) + void computePorosity(const UnstructuredGrid& grid, + const double* porosity_standard, + const RockCompressibility& rock_comp, + const std::vector& pressure, + std::vector& porosity); /// @brief Computes total saturated volumes over all grid cells. /// @param[in] pv the pore volume by cell. From 2dede29f20824b5c2f089b91950c7fe6d474c7b3 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Halvor=20M=C3=B8ll=20Nilsen?= Date: Fri, 31 Aug 2012 17:01:07 +0200 Subject: [PATCH 10/23] Introduced posibility to change number of sample points for pvt. Did change the PVTW calculation so derivatives are exact. Extended the test functions for pvt and relperm --- opm/core/fluid/BlackoilPropertiesFromDeck.cpp | 3 ++- opm/core/fluid/blackoil/BlackoilPvtProperties.hpp | 2 +- 2 files changed, 3 insertions(+), 2 deletions(-) diff --git a/opm/core/fluid/BlackoilPropertiesFromDeck.cpp b/opm/core/fluid/BlackoilPropertiesFromDeck.cpp index ae6a61306..67943fc76 100644 --- a/opm/core/fluid/BlackoilPropertiesFromDeck.cpp +++ b/opm/core/fluid/BlackoilPropertiesFromDeck.cpp @@ -39,7 +39,8 @@ namespace Opm const parameter::ParameterGroup& param) { rock_.init(deck, grid); - pvt_.init(deck); + int samples = param.getDefault("dead_tab_size", 1025); + pvt_.init(deck, samples); satprops_.init(deck, grid, param); if (pvt_.numPhases() != satprops_.numPhases()) { diff --git a/opm/core/fluid/blackoil/BlackoilPvtProperties.hpp b/opm/core/fluid/blackoil/BlackoilPvtProperties.hpp index 627f266b7..99fd5b8f4 100644 --- a/opm/core/fluid/blackoil/BlackoilPvtProperties.hpp +++ b/opm/core/fluid/blackoil/BlackoilPvtProperties.hpp @@ -47,7 +47,7 @@ namespace Opm BlackoilPvtProperties(); /// Initialize from deck. - void init(const EclipseGridParser& deck); + void init(const EclipseGridParser& deck,const int samples = 16); /// Number of active phases. int numPhases() const; From 0701daa4dc2f8f2bb6a748ab2565e04591840c4b Mon Sep 17 00:00:00 2001 From: Xavier Raynaud Date: Mon, 3 Sep 2012 11:29:48 +0200 Subject: [PATCH 11/23] Improve documentation. --- opm/core/utility/miscUtilities.hpp | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/opm/core/utility/miscUtilities.hpp b/opm/core/utility/miscUtilities.hpp index de452fbac..74cdc6748 100644 --- a/opm/core/utility/miscUtilities.hpp +++ b/opm/core/utility/miscUtilities.hpp @@ -44,7 +44,7 @@ namespace Opm /// @brief Computes pore volume of all cells in a grid, with rock compressibility effects. /// @param[in] grid a grid - /// @param[in] porosity array of grid.number_of_cells porosity values + /// @param[in] porosity array of grid.number_of_cells porosity values (at reference pressure) /// @param[in] rock_comp rock compressibility properties /// @param[in] pressure pressure by cell /// @param[out] porevol the pore volume by cell. From a8cb519d4eb194189d9a9bf9b7fc1cbf25468143 Mon Sep 17 00:00:00 2001 From: Xavier Raynaud Date: Mon, 3 Sep 2012 12:26:22 +0200 Subject: [PATCH 12/23] Fixed comment. --- opm/core/utility/miscUtilities.hpp | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/opm/core/utility/miscUtilities.hpp b/opm/core/utility/miscUtilities.hpp index 74cdc6748..2def319a9 100644 --- a/opm/core/utility/miscUtilities.hpp +++ b/opm/core/utility/miscUtilities.hpp @@ -56,7 +56,7 @@ namespace Opm /// @brief Computes porosity of all cells in a grid, with rock compressibility effects. /// @param[in] grid a grid - /// @param[in] porosity_standard array of grid.number_of_cells porosity values (at standard conditions) + /// @param[in] porosity_standard array of grid.number_of_cells porosity values (at reference presure) /// @param[in] rock_comp rock compressibility properties /// @param[in] pressure pressure by cell /// @param[out] porosity porosity (at reservoir condition) From 913054c47369279b1d074d71f08757fa48ee9536 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Atgeirr=20Fl=C3=B8=20Rasmussen?= Date: Mon, 3 Sep 2012 13:54:50 +0200 Subject: [PATCH 13/23] Added more checks in 3-phase init code. --- opm/core/utility/initState_impl.hpp | 14 +++++++++----- 1 file changed, 9 insertions(+), 5 deletions(-) diff --git a/opm/core/utility/initState_impl.hpp b/opm/core/utility/initState_impl.hpp index f66704348..9a635e93f 100644 --- a/opm/core/utility/initState_impl.hpp +++ b/opm/core/utility/initState_impl.hpp @@ -512,11 +512,10 @@ namespace Opm State& state) { const int num_phases = props.numPhases(); - state.init(grid, num_phases); if (deck.hasField("EQUIL")) { if (num_phases != 2) { - THROW("initStateFromDeck(): currently handling only two-phase scenarios."); + THROW("initStateFromDeck(): EQUIL-based init currently handling only two-phase scenarios."); } // Set saturations depending on oil-water contact. const EQUIL& equil= deck.getEQUIL(); @@ -536,22 +535,27 @@ namespace Opm const std::vector& sw_deck = deck.getFloatingPointValue("SWAT"); const std::vector& p_deck = deck.getFloatingPointValue("PRESSURE"); const int num_cells = grid.number_of_cells; - if(num_phases == 2){ + if (num_phases == 2) { for (int c = 0; c < num_cells; ++c) { int c_deck = (grid.global_cell == NULL) ? c : grid.global_cell[c]; s[2*c] = sw_deck[c_deck]; s[2*c + 1] = 1.0 - s[2*c]; p[c] = p_deck[c_deck]; } - }else{ + } else if (num_phases == 3) { + if (!deck.hasField("SGAS")) { + THROW("initStateFromDeck(): missing SGAS keyword in 3-phase init (only SWAT and PRESSURE found)."); + } const std::vector& sg_deck = deck.getFloatingPointValue("SGAS"); - for (int c = 0; c < num_cells; ++c) { + for (int c = 0; c < num_cells; ++c) { int c_deck = (grid.global_cell == NULL) ? c : grid.global_cell[c]; s[2*c] = sw_deck[c_deck]; s[2*c + 1] = 1.0 - (sw_deck[c_deck] + sg_deck[c_deck]); s[2*c + 2] = sg_deck[c_deck]; p[c] = p_deck[c_deck]; } + } else { + THROW("initStateFromDeck(): init with SWAT etc. only available with 2 or 3 phases."); } } else { THROW("initStateFromDeck(): we must either have EQUIL, or both SWAT and PRESSURE."); From c2d41a6639c0f24663ccdf0c3a6ed38c92853d06 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Atgeirr=20Fl=C3=B8=20Rasmussen?= Date: Mon, 3 Sep 2012 14:04:52 +0200 Subject: [PATCH 14/23] Whitespace cleanup and adding a check for #phases <= 3. --- opm/core/utility/miscUtilities.cpp | 340 ++++++++++++++--------------- 1 file changed, 170 insertions(+), 170 deletions(-) diff --git a/opm/core/utility/miscUtilities.cpp b/opm/core/utility/miscUtilities.cpp index 0b8d20103..e1658b83a 100644 --- a/opm/core/utility/miscUtilities.cpp +++ b/opm/core/utility/miscUtilities.cpp @@ -40,14 +40,14 @@ namespace Opm /// @param[out] porevol the pore volume by cell. void computePorevolume(const UnstructuredGrid& grid, const double* porosity, - std::vector& porevol) + std::vector& porevol) { - int num_cells = grid.number_of_cells; - porevol.resize(num_cells); - std::transform(porosity, porosity + num_cells, - grid.cell_volumes, - porevol.begin(), - std::multiplies()); + int num_cells = grid.number_of_cells; + porevol.resize(num_cells); + std::transform(porosity, porosity + num_cells, + grid.cell_volumes, + porevol.begin(), + std::multiplies()); } @@ -98,20 +98,20 @@ namespace Opm /// For each phase p, we compute /// sat_vol_p = sum_i s_p_i pv_i void computeSaturatedVol(const std::vector& pv, - const std::vector& s, - double* sat_vol) + const std::vector& s, + double* sat_vol) { - const int num_cells = pv.size(); - const int np = s.size()/pv.size(); - if (int(s.size()) != num_cells*np) { - THROW("Sizes of s and pv vectors do not match."); - } - std::fill(sat_vol, sat_vol + np, 0.0); - for (int c = 0; c < num_cells; ++c) { - for (int p = 0; p < np; ++p) { - sat_vol[p] += pv[c]*s[np*c + p]; - } - } + const int num_cells = pv.size(); + const int np = s.size()/pv.size(); + if (int(s.size()) != num_cells*np) { + THROW("Sizes of s and pv vectors do not match."); + } + std::fill(sat_vol, sat_vol + np, 0.0); + for (int c = 0; c < num_cells; ++c) { + for (int p = 0; p < np; ++p) { + sat_vol[p] += pv[c]*s[np*c + p]; + } + } } @@ -123,28 +123,28 @@ namespace Opm /// For each phase p, we compute /// aver_sat_p = (sum_i s_p_i pv_i) / (sum_i pv_i). void computeAverageSat(const std::vector& pv, - const std::vector& s, - double* aver_sat) + const std::vector& s, + double* aver_sat) { - const int num_cells = pv.size(); - const int np = s.size()/pv.size(); - if (int(s.size()) != num_cells*np) { - THROW("Sizes of s and pv vectors do not match."); - } - double tot_pv = 0.0; - // Note that we abuse the output array to accumulate the - // saturated pore volumes. - std::fill(aver_sat, aver_sat + np, 0.0); - for (int c = 0; c < num_cells; ++c) { - tot_pv += pv[c]; - for (int p = 0; p < np; ++p) { - aver_sat[p] += pv[c]*s[np*c + p]; - } - } - // Must divide by pore volumes to get saturations. - for (int p = 0; p < np; ++p) { - aver_sat[p] /= tot_pv; - } + const int num_cells = pv.size(); + const int np = s.size()/pv.size(); + if (int(s.size()) != num_cells*np) { + THROW("Sizes of s and pv vectors do not match."); + } + double tot_pv = 0.0; + // Note that we abuse the output array to accumulate the + // saturated pore volumes. + std::fill(aver_sat, aver_sat + np, 0.0); + for (int c = 0; c < num_cells; ++c) { + tot_pv += pv[c]; + for (int p = 0; p < np; ++p) { + aver_sat[p] += pv[c]*s[np*c + p]; + } + } + // Must divide by pore volumes to get saturations. + for (int p = 0; p < np; ++p) { + aver_sat[p] /= tot_pv; + } } @@ -161,38 +161,38 @@ namespace Opm /// where P = s.size()/src.size(). /// @param[out] produced must also point to a valid array with P elements. void computeInjectedProduced(const IncompPropertiesInterface& props, - const std::vector& s, - const std::vector& src, - const double dt, - double* injected, - double* produced) + const std::vector& s, + const std::vector& src, + const double dt, + double* injected, + double* produced) { - const int num_cells = src.size(); - const int np = s.size()/src.size(); - if (int(s.size()) != num_cells*np) { - THROW("Sizes of s and src vectors do not match."); - } - std::fill(injected, injected + np, 0.0); - std::fill(produced, produced + np, 0.0); - const double* visc = props.viscosity(); - std::vector mob(np); - for (int c = 0; c < num_cells; ++c) { - if (src[c] > 0.0) { - injected[0] += src[c]*dt; - } else if (src[c] < 0.0) { - const double flux = -src[c]*dt; - const double* sat = &s[np*c]; - props.relperm(1, sat, &c, &mob[0], 0); - double totmob = 0.0; - for (int p = 0; p < np; ++p) { - mob[p] /= visc[p]; - totmob += mob[p]; - } - for (int p = 0; p < np; ++p) { - produced[p] += (mob[p]/totmob)*flux; - } - } - } + const int num_cells = src.size(); + const int np = s.size()/src.size(); + if (int(s.size()) != num_cells*np) { + THROW("Sizes of s and src vectors do not match."); + } + std::fill(injected, injected + np, 0.0); + std::fill(produced, produced + np, 0.0); + const double* visc = props.viscosity(); + std::vector mob(np); + for (int c = 0; c < num_cells; ++c) { + if (src[c] > 0.0) { + injected[0] += src[c]*dt; + } else if (src[c] < 0.0) { + const double flux = -src[c]*dt; + const double* sat = &s[np*c]; + props.relperm(1, sat, &c, &mob[0], 0); + double totmob = 0.0; + for (int p = 0; p < np; ++p) { + mob[p] /= visc[p]; + totmob += mob[p]; + } + for (int p = 0; p < np; ++p) { + produced[p] += (mob[p]/totmob)*flux; + } + } + } } @@ -203,9 +203,9 @@ namespace Opm /// @param[in] s saturation values (for all phases) /// @param[out] totmob total mobilities. void computeTotalMobility(const Opm::IncompPropertiesInterface& props, - const std::vector& cells, - const std::vector& s, - std::vector& totmob) + const std::vector& cells, + const std::vector& s, + std::vector& totmob) { std::vector pmobc; @@ -231,10 +231,10 @@ namespace Opm /// @param[out] totmob total mobility /// @param[out] omega fractional-flow weighted fluid densities. void computeTotalMobilityOmega(const Opm::IncompPropertiesInterface& props, - const std::vector& cells, - const std::vector& s, - std::vector& totmob, - std::vector& omega) + const std::vector& cells, + const std::vector& s, + std::vector& totmob, + std::vector& omega) { std::vector pmobc; @@ -331,33 +331,33 @@ namespace Opm /// (+) positive inflow of first phase (water) /// (-) negative total outflow of both phases void computeTransportSource(const UnstructuredGrid& grid, - const std::vector& src, - const std::vector& faceflux, - const double inflow_frac, + const std::vector& src, + const std::vector& faceflux, + const double inflow_frac, const Wells* wells, const std::vector& well_perfrates, - std::vector& transport_src) + std::vector& transport_src) { - int nc = grid.number_of_cells; - transport_src.resize(nc); + int nc = grid.number_of_cells; + transport_src.resize(nc); // Source term and boundary contributions. - for (int c = 0; c < nc; ++c) { - transport_src[c] = 0.0; - transport_src[c] += src[c] > 0.0 ? inflow_frac*src[c] : src[c]; - for (int hf = grid.cell_facepos[c]; hf < grid.cell_facepos[c + 1]; ++hf) { - int f = grid.cell_faces[hf]; - const int* f2c = &grid.face_cells[2*f]; - double bdy_influx = 0.0; - if (f2c[0] == c && f2c[1] == -1) { - bdy_influx = -faceflux[f]; - } else if (f2c[0] == -1 && f2c[1] == c) { - bdy_influx = faceflux[f]; - } - if (bdy_influx != 0.0) { - transport_src[c] += bdy_influx > 0.0 ? inflow_frac*bdy_influx : bdy_influx; - } - } - } + for (int c = 0; c < nc; ++c) { + transport_src[c] = 0.0; + transport_src[c] += src[c] > 0.0 ? inflow_frac*src[c] : src[c]; + for (int hf = grid.cell_facepos[c]; hf < grid.cell_facepos[c + 1]; ++hf) { + int f = grid.cell_faces[hf]; + const int* f2c = &grid.face_cells[2*f]; + double bdy_influx = 0.0; + if (f2c[0] == c && f2c[1] == -1) { + bdy_influx = -faceflux[f]; + } else if (f2c[0] == -1 && f2c[1] == c) { + bdy_influx = faceflux[f]; + } + if (bdy_influx != 0.0) { + transport_src[c] += bdy_influx > 0.0 ? inflow_frac*bdy_influx : bdy_influx; + } + } + } // Well contributions. if (wells) { @@ -392,52 +392,52 @@ namespace Opm /// @param[in] face_flux signed per-face fluxes /// @param[out] cell_velocity the estimated velocities. void estimateCellVelocity(const UnstructuredGrid& grid, - const std::vector& face_flux, - std::vector& cell_velocity) + const std::vector& face_flux, + std::vector& cell_velocity) { - const int dim = grid.dimensions; - cell_velocity.clear(); - cell_velocity.resize(grid.number_of_cells*dim, 0.0); - for (int face = 0; face < grid.number_of_faces; ++face) { - int c[2] = { grid.face_cells[2*face], grid.face_cells[2*face + 1] }; - const double* fc = &grid.face_centroids[face*dim]; - double flux = face_flux[face]; - for (int i = 0; i < 2; ++i) { - if (c[i] >= 0) { - const double* cc = &grid.cell_centroids[c[i]*dim]; - for (int d = 0; d < dim; ++d) { - double v_contrib = fc[d] - cc[d]; - v_contrib *= flux/grid.cell_volumes[c[i]]; - cell_velocity[c[i]*dim + d] += (i == 0) ? v_contrib : -v_contrib; - } - } - } - } + const int dim = grid.dimensions; + cell_velocity.clear(); + cell_velocity.resize(grid.number_of_cells*dim, 0.0); + for (int face = 0; face < grid.number_of_faces; ++face) { + int c[2] = { grid.face_cells[2*face], grid.face_cells[2*face + 1] }; + const double* fc = &grid.face_centroids[face*dim]; + double flux = face_flux[face]; + for (int i = 0; i < 2; ++i) { + if (c[i] >= 0) { + const double* cc = &grid.cell_centroids[c[i]*dim]; + for (int d = 0; d < dim; ++d) { + double v_contrib = fc[d] - cc[d]; + v_contrib *= flux/grid.cell_volumes[c[i]]; + cell_velocity[c[i]*dim + d] += (i == 0) ? v_contrib : -v_contrib; + } + } + } + } } /// Extract a vector of water saturations from a vector of /// interleaved water and oil saturations. void toWaterSat(const std::vector& sboth, - std::vector& sw) + std::vector& sw) { - int num = sboth.size()/2; - sw.resize(num); - for (int i = 0; i < num; ++i) { - sw[i] = sboth[2*i]; - } + int num = sboth.size()/2; + sw.resize(num); + for (int i = 0; i < num; ++i) { + sw[i] = sboth[2*i]; + } } /// Make a a vector of interleaved water and oil saturations from /// a vector of water saturations. void toBothSat(const std::vector& sw, - std::vector& sboth) + std::vector& sboth) { - int num = sw.size(); - sboth.resize(2*num); - for (int i = 0; i < num; ++i) { - sboth[2*i] = sw[i]; - sboth[2*i + 1] = 1.0 - sw[i]; - } + int num = sw.size(); + sboth.resize(2*num); + for (int i = 0; i < num; ++i) { + sboth[2*i] = sw[i]; + sboth[2*i + 1] = 1.0 - sw[i]; + } } @@ -450,30 +450,30 @@ namespace Opm if (np != 2) { THROW("wellsToSrc() requires a 2 phase case."); } - src.resize(num_cells); - for (int w = 0; w < wells.number_of_wells; ++w) { + src.resize(num_cells); + for (int w = 0; w < wells.number_of_wells; ++w) { const int cur = wells.ctrls[w]->current; - if (wells.ctrls[w]->num != 1) { - MESSAGE("In wellsToSrc(): well has more than one control, all but current control will be ignored."); - } - if (wells.ctrls[w]->type[cur] != RESERVOIR_RATE) { - THROW("In wellsToSrc(): well is something other than RESERVOIR_RATE."); - } - if (wells.well_connpos[w+1] - wells.well_connpos[w] != 1) { - THROW("In wellsToSrc(): well has multiple perforations."); - } + if (wells.ctrls[w]->num != 1) { + MESSAGE("In wellsToSrc(): well has more than one control, all but current control will be ignored."); + } + if (wells.ctrls[w]->type[cur] != RESERVOIR_RATE) { + THROW("In wellsToSrc(): well is something other than RESERVOIR_RATE."); + } + if (wells.well_connpos[w+1] - wells.well_connpos[w] != 1) { + THROW("In wellsToSrc(): well has multiple perforations."); + } for (int p = 0; p < np; ++p) { if (wells.ctrls[w]->distr[np*cur + p] != 1.0) { THROW("In wellsToSrc(): well not controlled on total rate."); } } - double flow = wells.ctrls[w]->target[cur]; + double flow = wells.ctrls[w]->target[cur]; if (wells.type[w] == INJECTOR) { flow *= wells.comp_frac[np*w + 0]; // Obtaining water rate for inflow source. } - const double cell = wells.well_cells[wells.well_connpos[w]]; - src[cell] = flow; - } + const double cell = wells.well_cells[wells.well_connpos[w]]; + src[cell] = flow; + } } @@ -594,9 +594,10 @@ namespace Opm int nw = well_bhp.size(); ASSERT(nw == wells.number_of_wells); int np = props.numPhases(); - //if (props.numPhases() != 2) { - // THROW("WellReport for now assumes two phase flow."); - //} + const int max_np = 3; + if (np > max_np) { + THROW("WellReport for now assumes #phases <= " << max_np); + } const double* visc = props.viscosity(); std::vector data_now; data_now.reserve(1 + 3*nw); @@ -606,8 +607,8 @@ namespace Opm double well_rate_total = 0.0; double well_rate_water = 0.0; for (int perf = wells.well_connpos[w]; perf < wells.well_connpos[w + 1]; ++perf) { - const double perf_rate = unit::convert::to(well_perfrates[perf], - unit::cubic(unit::meter)/unit::day); + const double perf_rate = unit::convert::to(well_perfrates[perf], + unit::cubic(unit::meter)/unit::day); well_rate_total += perf_rate; if (perf_rate > 0.0) { // Injection. @@ -615,13 +616,13 @@ namespace Opm } else { // Production. const int cell = wells.well_cells[perf]; - double mob[np]; + double mob[max_np]; props.relperm(1, &saturation[2*cell], &cell, mob, 0); - double tmob=0; - for(int i=0; i < np; ++i){ + double tmob = 0; + for(int i = 0; i < np; ++i) { mob[i] /= visc[i]; tmob += mob[i]; - } + } const double fracflow = mob[0]/tmob; well_rate_water += perf_rate*fracflow; } @@ -652,9 +653,10 @@ namespace Opm int nw = well_bhp.size(); ASSERT(nw == wells.number_of_wells); int np = props.numPhases(); - //if (props.numPhases() != 2) { - // THROW("WellReport for now assumes two phase flow."); - //} + const int max_np = 3; + if (np > max_np) { + THROW("WellReport for now assumes #phases <= " << max_np); + } std::vector data_now; data_now.reserve(1 + 3*nw); data_now.push_back(time/unit::day); @@ -665,25 +667,23 @@ namespace Opm for (int perf = wells.well_connpos[w]; perf < wells.well_connpos[w + 1]; ++perf) { const double perf_rate = well_perfrates[perf]*(unit::day/unit::second); well_rate_total += perf_rate; - if (perf_rate > 0.0) { + if (perf_rate > 0.0) { // Injection. well_rate_water += perf_rate*wells.comp_frac[0]; } else { // Production. const int cell = wells.well_cells[perf]; - double mob[np]; + double mob[max_np]; props.relperm(1, &s[2*cell], &cell, mob, 0); - double visc[np]; + double visc[max_np]; props.viscosity(1, &p[cell], &z[2*cell], &cell, visc, 0); - double tmob=0; - for(int i=0; i < np; ++i){ + double tmob = 0; + for(int i = 0; i < np; ++i) { mob[i] /= visc[i]; tmob += mob[i]; - } + } const double fracflow = mob[0]/(tmob); well_rate_water += perf_rate*fracflow; - //const double fracflow = mob[0]/(tmob); - //well_rate_water += perf_rate*fracflow; } } data_now.push_back(well_rate_total); From 489501b49b3669548a3e6657e51de156692f87a5 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Atgeirr=20Fl=C3=B8=20Rasmussen?= Date: Mon, 3 Sep 2012 14:16:54 +0200 Subject: [PATCH 15/23] Documented new constructor. Fixed layout. --- opm/core/fluid/BlackoilPropertiesFromDeck.hpp | 9 +++++++++ 1 file changed, 9 insertions(+) diff --git a/opm/core/fluid/BlackoilPropertiesFromDeck.hpp b/opm/core/fluid/BlackoilPropertiesFromDeck.hpp index 1d0a0d380..e6b9c5b41 100644 --- a/opm/core/fluid/BlackoilPropertiesFromDeck.hpp +++ b/opm/core/fluid/BlackoilPropertiesFromDeck.hpp @@ -45,9 +45,18 @@ namespace Opm /// to logical cartesian indices consistent with the deck. BlackoilPropertiesFromDeck(const EclipseGridParser& deck, const UnstructuredGrid& grid); + + /// Initialize from deck, grid and parameters. + /// \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. + /// \param param Parameters. Accepted parameters include: + /// tab_size (200) number of uniform sample points for saturation tables. BlackoilPropertiesFromDeck(const EclipseGridParser& deck, const UnstructuredGrid& grid, const parameter::ParameterGroup& param); + /// Destructor. virtual ~BlackoilPropertiesFromDeck(); From 4e1647bb6272fa98ba57936d6761b8b39578f90b Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Atgeirr=20Fl=C3=B8=20Rasmussen?= Date: Mon, 3 Sep 2012 15:07:03 +0200 Subject: [PATCH 16/23] Formatting fixes. --- opm/core/fluid/blackoil/BlackoilPvtProperties.hpp | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/opm/core/fluid/blackoil/BlackoilPvtProperties.hpp b/opm/core/fluid/blackoil/BlackoilPvtProperties.hpp index 99fd5b8f4..191c6fa89 100644 --- a/opm/core/fluid/blackoil/BlackoilPvtProperties.hpp +++ b/opm/core/fluid/blackoil/BlackoilPvtProperties.hpp @@ -47,7 +47,7 @@ namespace Opm BlackoilPvtProperties(); /// Initialize from deck. - void init(const EclipseGridParser& deck,const int samples = 16); + void init(const EclipseGridParser& deck, const int samples = 16); /// Number of active phases. int numPhases() const; From e9c4c2499cc362facef490fae23178c2b72fe6b0 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Atgeirr=20Fl=C3=B8=20Rasmussen?= Date: Mon, 3 Sep 2012 15:09:55 +0200 Subject: [PATCH 17/23] Documented parameters. --- opm/core/fluid/BlackoilPropertiesFromDeck.hpp | 3 ++- 1 file changed, 2 insertions(+), 1 deletion(-) diff --git a/opm/core/fluid/BlackoilPropertiesFromDeck.hpp b/opm/core/fluid/BlackoilPropertiesFromDeck.hpp index e6b9c5b41..b63c0af88 100644 --- a/opm/core/fluid/BlackoilPropertiesFromDeck.hpp +++ b/opm/core/fluid/BlackoilPropertiesFromDeck.hpp @@ -52,7 +52,8 @@ namespace Opm /// mapping from cell indices (typically from a processed grid) /// to logical cartesian indices consistent with the deck. /// \param param Parameters. Accepted parameters include: - /// tab_size (200) number of uniform sample points for saturation tables. + /// dead_tab_size (1025) number of uniform sample points for dead-oil pvt tables. + /// tab_size_kr (200) number of uniform sample points for saturation tables. BlackoilPropertiesFromDeck(const EclipseGridParser& deck, const UnstructuredGrid& grid, const parameter::ParameterGroup& param); From 0101e1f575b0ce50222a5b1d61bcb84d875d4196 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Atgeirr=20Fl=C3=B8=20Rasmussen?= Date: Tue, 4 Sep 2012 11:49:05 +0200 Subject: [PATCH 18/23] Updated doc to match new parameter names. New parameters are (default): pvt_tab_size (200) sat_tab_size (200) threephase_model ("simple") [also accepts "stone2"]. --- opm/core/fluid/BlackoilPropertiesFromDeck.hpp | 5 +++-- 1 file changed, 3 insertions(+), 2 deletions(-) diff --git a/opm/core/fluid/BlackoilPropertiesFromDeck.hpp b/opm/core/fluid/BlackoilPropertiesFromDeck.hpp index 97a29d0af..80a8d8192 100644 --- a/opm/core/fluid/BlackoilPropertiesFromDeck.hpp +++ b/opm/core/fluid/BlackoilPropertiesFromDeck.hpp @@ -53,8 +53,9 @@ namespace Opm /// mapping from cell indices (typically from a processed grid) /// to logical cartesian indices consistent with the deck. /// \param[in] param Parameters. Accepted parameters include: - /// dead_tab_size (1025) number of uniform sample points for dead-oil pvt tables. - /// tab_size_kr (200) number of uniform sample points for saturation tables. + /// 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 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, From 7c7609698f617a4ce640720f6cb09421365474ad Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Atgeirr=20Fl=C3=B8=20Rasmussen?= Date: Tue, 4 Sep 2012 12:52:17 +0200 Subject: [PATCH 19/23] Added proper Doxygen tags. --- opm/core/fluid/BlackoilPropertiesFromDeck.hpp | 10 +++++----- 1 file changed, 5 insertions(+), 5 deletions(-) diff --git a/opm/core/fluid/BlackoilPropertiesFromDeck.hpp b/opm/core/fluid/BlackoilPropertiesFromDeck.hpp index b63c0af88..13f5799d8 100644 --- a/opm/core/fluid/BlackoilPropertiesFromDeck.hpp +++ b/opm/core/fluid/BlackoilPropertiesFromDeck.hpp @@ -39,19 +39,19 @@ namespace Opm { public: /// Initialize from deck and grid. - /// \param deck Deck input parser - /// \param grid Grid to which property object applies, needed for the + /// \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); /// Initialize from deck, grid and parameters. - /// \param deck Deck input parser - /// \param grid Grid to which property object applies, needed for the + /// \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 param Parameters. Accepted parameters include: + /// \param[in] param Parameters. Accepted parameters include: /// dead_tab_size (1025) number of uniform sample points for dead-oil pvt tables. /// tab_size_kr (200) number of uniform sample points for saturation tables. BlackoilPropertiesFromDeck(const EclipseGridParser& deck, From 14aab769f2f9d90875d2568596a2924b0eb421bc Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Atgeirr=20Fl=C3=B8=20Rasmussen?= Date: Tue, 4 Sep 2012 13:06:13 +0200 Subject: [PATCH 20/23] Fix stride bug in SGAS/SWAT init. --- opm/core/utility/initState_impl.hpp | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/opm/core/utility/initState_impl.hpp b/opm/core/utility/initState_impl.hpp index 9a635e93f..7d01ade0e 100644 --- a/opm/core/utility/initState_impl.hpp +++ b/opm/core/utility/initState_impl.hpp @@ -549,9 +549,9 @@ namespace Opm const std::vector& sg_deck = deck.getFloatingPointValue("SGAS"); for (int c = 0; c < num_cells; ++c) { int c_deck = (grid.global_cell == NULL) ? c : grid.global_cell[c]; - s[2*c] = sw_deck[c_deck]; - s[2*c + 1] = 1.0 - (sw_deck[c_deck] + sg_deck[c_deck]); - s[2*c + 2] = sg_deck[c_deck]; + s[3*c] = sw_deck[c_deck]; + s[3*c + 1] = 1.0 - (sw_deck[c_deck] + sg_deck[c_deck]); + s[3*c + 2] = sg_deck[c_deck]; p[c] = p_deck[c_deck]; } } else { From 0f9d572f9a4a148c4feb1cb38731faf27a7c0b87 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Atgeirr=20Fl=C3=B8=20Rasmussen?= Date: Tue, 4 Sep 2012 13:12:06 +0200 Subject: [PATCH 21/23] Fix stride bug in WellReport::push(). Use convert::to(). --- opm/core/utility/miscUtilities.cpp | 7 ++++--- 1 file changed, 4 insertions(+), 3 deletions(-) diff --git a/opm/core/utility/miscUtilities.cpp b/opm/core/utility/miscUtilities.cpp index e1658b83a..aeafbfca7 100644 --- a/opm/core/utility/miscUtilities.cpp +++ b/opm/core/utility/miscUtilities.cpp @@ -665,7 +665,8 @@ namespace Opm double well_rate_total = 0.0; double well_rate_water = 0.0; for (int perf = wells.well_connpos[w]; perf < wells.well_connpos[w + 1]; ++perf) { - const double perf_rate = well_perfrates[perf]*(unit::day/unit::second); + const double perf_rate = unit::convert::to(well_perfrates[perf], + unit::cubic(unit::meter)/unit::day); well_rate_total += perf_rate; if (perf_rate > 0.0) { // Injection. @@ -674,9 +675,9 @@ namespace Opm // Production. const int cell = wells.well_cells[perf]; double mob[max_np]; - props.relperm(1, &s[2*cell], &cell, mob, 0); + props.relperm(1, &s[np*cell], &cell, mob, 0); double visc[max_np]; - props.viscosity(1, &p[cell], &z[2*cell], &cell, visc, 0); + props.viscosity(1, &p[cell], &z[np*cell], &cell, visc, 0); double tmob = 0; for(int i = 0; i < np; ++i) { mob[i] /= visc[i]; From 23ae62f8dc6339eff646cdb1e85f12ac9c76be72 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Atgeirr=20Fl=C3=B8=20Rasmussen?= Date: Tue, 4 Sep 2012 13:15:26 +0200 Subject: [PATCH 22/23] Constified some local variables. --- opm/core/fluid/BlackoilPropertiesFromDeck.cpp | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/opm/core/fluid/BlackoilPropertiesFromDeck.cpp b/opm/core/fluid/BlackoilPropertiesFromDeck.cpp index 67943fc76..5c834fe8c 100644 --- a/opm/core/fluid/BlackoilPropertiesFromDeck.cpp +++ b/opm/core/fluid/BlackoilPropertiesFromDeck.cpp @@ -39,7 +39,7 @@ namespace Opm const parameter::ParameterGroup& param) { rock_.init(deck, grid); - int samples = param.getDefault("dead_tab_size", 1025); + const int samples = param.getDefault("dead_tab_size", 1025); pvt_.init(deck, samples); satprops_.init(deck, grid, param); From 0f91bc6a349a83fdcd3966e8a03eb62950f54085 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Atgeirr=20Fl=C3=B8=20Rasmussen?= Date: Wed, 5 Sep 2012 11:28:54 +0200 Subject: [PATCH 23/23] Massive whitespace cleanup: entire fluid subdir. --- opm/core/fluid/BlackoilPropertiesBasic.cpp | 72 +++---- opm/core/fluid/BlackoilPropertiesBasic.hpp | 20 +- opm/core/fluid/BlackoilPropertiesFromDeck.cpp | 10 +- opm/core/fluid/BlackoilPropertiesFromDeck.hpp | 12 +- .../fluid/BlackoilPropertiesInterface.hpp | 6 +- opm/core/fluid/IncompPropertiesBasic.cpp | 70 +++--- opm/core/fluid/IncompPropertiesBasic.hpp | 34 +-- opm/core/fluid/IncompPropertiesFromDeck.cpp | 42 ++-- opm/core/fluid/IncompPropertiesFromDeck.hpp | 14 +- opm/core/fluid/IncompPropertiesInterface.hpp | 6 +- opm/core/fluid/PvtPropertiesBasic.cpp | 112 +++++----- opm/core/fluid/PvtPropertiesBasic.hpp | 18 +- .../fluid/PvtPropertiesIncompFromDeck.cpp | 70 +++--- .../fluid/PvtPropertiesIncompFromDeck.hpp | 12 +- opm/core/fluid/RockBasic.hpp | 8 +- opm/core/fluid/RockCompressibility.cpp | 6 +- opm/core/fluid/RockFromDeck.cpp | 2 +- opm/core/fluid/RockFromDeck.hpp | 2 +- opm/core/fluid/SaturationPropsBasic.cpp | 204 +++++++++--------- opm/core/fluid/SaturationPropsBasic.hpp | 22 +- opm/core/fluid/SaturationPropsFromDeck.hpp | 10 +- .../fluid/SaturationPropsFromDeck_impl.hpp | 20 +- opm/core/fluid/SaturationPropsInterface.hpp | 4 +- .../fluid/blackoil/BlackoilPvtProperties.hpp | 6 +- 24 files changed, 391 insertions(+), 391 deletions(-) diff --git a/opm/core/fluid/BlackoilPropertiesBasic.cpp b/opm/core/fluid/BlackoilPropertiesBasic.cpp index be8211c15..2496d98ef 100644 --- a/opm/core/fluid/BlackoilPropertiesBasic.cpp +++ b/opm/core/fluid/BlackoilPropertiesBasic.cpp @@ -26,20 +26,20 @@ namespace Opm { BlackoilPropertiesBasic::BlackoilPropertiesBasic(const parameter::ParameterGroup& param, - const int dim, - const int num_cells) + const int dim, + const int num_cells) { - double poro = param.getDefault("porosity", 1.0); - using namespace Opm::unit; - using namespace Opm::prefix; - double perm = param.getDefault("permeability", 100.0)*milli*darcy; + double poro = param.getDefault("porosity", 1.0); + using namespace Opm::unit; + using namespace Opm::prefix; + double perm = param.getDefault("permeability", 100.0)*milli*darcy; rock_.init(dim, num_cells, poro, perm); - pvt_.init(param); + pvt_.init(param); satprops_.init(param); - if (pvt_.numPhases() != satprops_.numPhases()) { - THROW("BlackoilPropertiesBasic::BlackoilPropertiesBasic() - Inconsistent number of phases in pvt data (" - << pvt_.numPhases() << ") and saturation-dependent function data (" << satprops_.numPhases() << ")."); - } + if (pvt_.numPhases() != satprops_.numPhases()) { + THROW("BlackoilPropertiesBasic::BlackoilPropertiesBasic() - Inconsistent number of phases in pvt data (" + << pvt_.numPhases() << ") and saturation-dependent function data (" << satprops_.numPhases() << ")."); + } } BlackoilPropertiesBasic::~BlackoilPropertiesBasic() @@ -90,11 +90,11 @@ namespace Opm /// \param[out] dmudp If non-null: array of nP viscosity derivative values, /// array must be valid before calling. void BlackoilPropertiesBasic::viscosity(const int n, - const double* p, - const double* z, - const int* /*cells*/, - double* mu, - double* dmudp) const + const double* p, + const double* z, + const int* /*cells*/, + double* mu, + double* dmudp) const { if (dmudp) { THROW("BlackoilPropertiesBasic::viscosity() -- derivatives of viscosity not yet implemented."); @@ -114,16 +114,16 @@ namespace Opm /// array must be valid before calling. The matrices are output /// in Fortran order. void BlackoilPropertiesBasic::matrix(const int n, - const double* /*p*/, - const double* /*z*/, - const int* /*cells*/, - double* A, - double* dAdp) const + const double* /*p*/, + const double* /*z*/, + const int* /*cells*/, + double* A, + double* dAdp) const { - const int np = numPhases(); - ASSERT(np <= 2); - double B[2]; // Must be enough since component classes do not handle more than 2. - pvt_.B(1, 0, 0, B); + const int np = numPhases(); + ASSERT(np <= 2); + double B[2]; // Must be enough since component classes do not handle more than 2. + pvt_.B(1, 0, 0, B); // Compute A matrix // #pragma omp parallel for for (int i = 0; i < n; ++i) { @@ -152,8 +152,8 @@ namespace Opm /// of a call to the method matrix(). /// \param[out] rho Array of nP density values, array must be valid before calling. void BlackoilPropertiesBasic::density(const int n, - const double* A, - double* rho) const + const double* A, + double* rho) const { const int np = numPhases(); const double* sdens = pvt_.surfaceDensities(); @@ -186,10 +186,10 @@ namespace Opm /// m_{ij} = \frac{dkr_i}{ds^j}, /// and is output in Fortran order (m_00 m_10 m_20 m01 ...) void BlackoilPropertiesBasic::relperm(const int n, - const double* s, - const int* /*cells*/, - double* kr, - double* dkrds) const + const double* s, + const int* /*cells*/, + double* kr, + double* dkrds) const { satprops_.relperm(n, s, kr, dkrds); } @@ -205,10 +205,10 @@ namespace Opm /// m_{ij} = \frac{dpc_i}{ds^j}, /// and is output in Fortran order (m_00 m_10 m_20 m01 ...) void BlackoilPropertiesBasic::capPress(const int n, - const double* s, - const int* /*cells*/, - double* pc, - double* dpcds) const + const double* s, + const int* /*cells*/, + double* pc, + double* dpcds) const { satprops_.capPress(n, s, pc, dpcds); } @@ -226,7 +226,7 @@ namespace Opm double* smin, double* smax) const { - satprops_.satRange(n, smin, smax); + satprops_.satRange(n, smin, smax); } diff --git a/opm/core/fluid/BlackoilPropertiesBasic.hpp b/opm/core/fluid/BlackoilPropertiesBasic.hpp index d879036ef..44a76013b 100644 --- a/opm/core/fluid/BlackoilPropertiesBasic.hpp +++ b/opm/core/fluid/BlackoilPropertiesBasic.hpp @@ -35,16 +35,16 @@ namespace Opm { public: /// Construct from parameters. - /// The following parameters are accepted (defaults): - /// num_phases (2) Must be 1 or 2. - /// relperm_func ("Linear") Must be "Constant", "Linear" or "Quadratic". - /// rho1 [rho2, rho3] (1.0e3) Density in kg/m^3 - /// mu1 [mu2, mu3] (1.0) Viscosity in cP - /// porosity (1.0) Porosity - /// permeability (100.0) Permeability in mD + /// The following parameters are accepted (defaults): + /// num_phases (2) Must be 1 or 2. + /// relperm_func ("Linear") Must be "Constant", "Linear" or "Quadratic". + /// rho1 [rho2, rho3] (1.0e3) Density in kg/m^3 + /// mu1 [mu2, mu3] (1.0) Viscosity in cP + /// porosity (1.0) Porosity + /// permeability (100.0) Permeability in mD BlackoilPropertiesBasic(const parameter::ParameterGroup& param, - const int dim, - const int num_cells); + const int dim, + const int num_cells); /// Destructor. virtual ~BlackoilPropertiesBasic(); @@ -151,7 +151,7 @@ namespace Opm double* dpcds) const; - /// Obtain the range of allowable saturation values. +/// Obtain the range of allowable saturation values. /// In cell cells[i], saturation of phase p is allowed to be /// in the interval [smin[i*P + p], smax[i*P + p]]. /// \param[in] n Number of data points. diff --git a/opm/core/fluid/BlackoilPropertiesFromDeck.cpp b/opm/core/fluid/BlackoilPropertiesFromDeck.cpp index 7cdae4ae2..e2457f39b 100644 --- a/opm/core/fluid/BlackoilPropertiesFromDeck.cpp +++ b/opm/core/fluid/BlackoilPropertiesFromDeck.cpp @@ -33,10 +33,10 @@ namespace Opm satprops_.reset(ptr); ptr->init(deck, grid, 200); - if (pvt_.numPhases() != satprops_->numPhases()) { - THROW("BlackoilPropertiesFromDeck::BlackoilPropertiesFromDeck() - Inconsistent number of phases in pvt data (" - << pvt_.numPhases() << ") and saturation-dependent function data (" << satprops_->numPhases() << ")."); - } + if (pvt_.numPhases() != satprops_->numPhases()) { + THROW("BlackoilPropertiesFromDeck::BlackoilPropertiesFromDeck() - Inconsistent number of phases in pvt data (" + << pvt_.numPhases() << ") and saturation-dependent function data (" << satprops_->numPhases() << ")."); + } } BlackoilPropertiesFromDeck::BlackoilPropertiesFromDeck(const EclipseGridParser& deck, @@ -319,7 +319,7 @@ namespace Opm double* smin, double* smax) const { - satprops_->satRange(n, cells, smin, smax); + satprops_->satRange(n, cells, smin, smax); } diff --git a/opm/core/fluid/BlackoilPropertiesFromDeck.hpp b/opm/core/fluid/BlackoilPropertiesFromDeck.hpp index aceb534b8..eb8b947b4 100644 --- a/opm/core/fluid/BlackoilPropertiesFromDeck.hpp +++ b/opm/core/fluid/BlackoilPropertiesFromDeck.hpp @@ -41,15 +41,15 @@ namespace Opm public: /// Initialize from deck and grid. /// \param[in] deck Deck input parser - /// \param[in] grid Grid to which property object applies, needed for the + /// \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); + const UnstructuredGrid& grid); /// Initialize from deck, grid and parameters. /// \param[in] deck Deck input parser - /// \param[in] grid Grid to which property object applies, needed for the + /// \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: @@ -167,9 +167,9 @@ namespace Opm double* dpcds) const; - /// Obtain the range of allowable saturation values. - /// In cell cells[i], saturation of phase p is allowed to be - /// in the interval [smin[i*P + p], smax[i*P + p]]. + /// Obtain the range of allowable saturation values. + /// In cell cells[i], saturation of phase p is allowed to be + /// in the interval [smin[i*P + p], smax[i*P + p]]. /// \param[in] n Number of data points. /// \param[in] cells Array of n cell indices. /// \param[out] smin Array of nP minimum s values, array must be valid before calling. diff --git a/opm/core/fluid/BlackoilPropertiesInterface.hpp b/opm/core/fluid/BlackoilPropertiesInterface.hpp index dd8a857fb..501c8a406 100644 --- a/opm/core/fluid/BlackoilPropertiesInterface.hpp +++ b/opm/core/fluid/BlackoilPropertiesInterface.hpp @@ -138,9 +138,9 @@ namespace Opm double* dpcds) const = 0; - /// Obtain the range of allowable saturation values. - /// In cell cells[i], saturation of phase p is allowed to be - /// in the interval [smin[i*P + p], smax[i*P + p]]. + /// Obtain the range of allowable saturation values. + /// In cell cells[i], saturation of phase p is allowed to be + /// in the interval [smin[i*P + p], smax[i*P + p]]. /// \param[in] n Number of data points. /// \param[in] cells Array of n cell indices. /// \param[out] smin Array of nP minimum s values, array must be valid before calling. diff --git a/opm/core/fluid/IncompPropertiesBasic.cpp b/opm/core/fluid/IncompPropertiesBasic.cpp index 0ce3c88aa..ab68017c1 100644 --- a/opm/core/fluid/IncompPropertiesBasic.cpp +++ b/opm/core/fluid/IncompPropertiesBasic.cpp @@ -28,22 +28,22 @@ namespace Opm { IncompPropertiesBasic::IncompPropertiesBasic(const parameter::ParameterGroup& param, - const int dim, - const int num_cells) + const int dim, + const int num_cells) { - double poro = param.getDefault("porosity", 1.0); - using namespace Opm::unit; - using namespace Opm::prefix; - double perm = param.getDefault("permeability", 100.0)*milli*darcy; + double poro = param.getDefault("porosity", 1.0); + using namespace Opm::unit; + using namespace Opm::prefix; + double perm = param.getDefault("permeability", 100.0)*milli*darcy; rock_.init(dim, num_cells, poro, perm); - pvt_.init(param); + pvt_.init(param); satprops_.init(param); - if (pvt_.numPhases() != satprops_.numPhases()) { - THROW("IncompPropertiesBasic::IncompPropertiesBasic() - Inconsistent number of phases in pvt data (" - << pvt_.numPhases() << ") and saturation-dependent function data (" << satprops_.numPhases() << ")."); - } - viscosity_.resize(pvt_.numPhases()); - pvt_.mu(1, 0, 0, &viscosity_[0]); + if (pvt_.numPhases() != satprops_.numPhases()) { + THROW("IncompPropertiesBasic::IncompPropertiesBasic() - Inconsistent number of phases in pvt data (" + << pvt_.numPhases() << ") and saturation-dependent function data (" << satprops_.numPhases() << ")."); + } + viscosity_.resize(pvt_.numPhases()); + pvt_.mu(1, 0, 0, &viscosity_[0]); } IncompPropertiesBasic::IncompPropertiesBasic(const int num_phases, @@ -56,14 +56,14 @@ namespace Opm const int num_cells) { rock_.init(dim, num_cells, por, perm); - pvt_.init(num_phases, rho, mu); + pvt_.init(num_phases, rho, mu); satprops_.init(num_phases, relpermfunc); - if (pvt_.numPhases() != satprops_.numPhases()) { - THROW("IncompPropertiesBasic::IncompPropertiesBasic() - Inconsistent number of phases in pvt data (" - << pvt_.numPhases() << ") and saturation-dependent function data (" << satprops_.numPhases() << ")."); - } - viscosity_.resize(pvt_.numPhases()); - pvt_.mu(1, 0, 0, &viscosity_[0]); + if (pvt_.numPhases() != satprops_.numPhases()) { + THROW("IncompPropertiesBasic::IncompPropertiesBasic() - Inconsistent number of phases in pvt data (" + << pvt_.numPhases() << ") and saturation-dependent function data (" << satprops_.numPhases() << ")."); + } + viscosity_.resize(pvt_.numPhases()); + pvt_.mu(1, 0, 0, &viscosity_[0]); } IncompPropertiesBasic::~IncompPropertiesBasic() @@ -109,7 +109,7 @@ namespace Opm /// \return Array of P viscosity values. const double* IncompPropertiesBasic::viscosity() const { - return &viscosity_[0]; + return &viscosity_[0]; } /// \return Array of P density values. @@ -117,7 +117,7 @@ namespace Opm { // No difference between reservoir and surface densities // modelled by this class. - return pvt_.surfaceDensities(); + return pvt_.surfaceDensities(); } /// \return Array of P density values. @@ -125,7 +125,7 @@ namespace Opm { // No difference between reservoir and surface densities // modelled by this class. - return pvt_.surfaceDensities(); + return pvt_.surfaceDensities(); } /// \param[in] n Number of data points. @@ -138,10 +138,10 @@ namespace Opm /// m_{ij} = \frac{dkr_i}{ds^j}, /// and is output in Fortran order (m_00 m_10 m_20 m_01 ...) void IncompPropertiesBasic::relperm(const int n, - const double* s, - const int* /*cells*/, - double* kr, - double* dkrds) const + const double* s, + const int* /*cells*/, + double* kr, + double* dkrds) const { satprops_.relperm(n, s, kr, dkrds); } @@ -157,10 +157,10 @@ namespace Opm /// m_{ij} = \frac{dpc_i}{ds^j}, /// and is output in Fortran order (m_00 m_10 m_20 m_01 ...) void IncompPropertiesBasic::capPress(const int n, - const double* s, - const int* /*cells*/, - double* pc, - double* dpcds) const + const double* s, + const int* /*cells*/, + double* pc, + double* dpcds) const { satprops_.capPress(n, s, pc, dpcds); } @@ -174,11 +174,11 @@ namespace Opm /// \param[out] smin Array of nP minimum s values, array must be valid before calling. /// \param[out] smax Array of nP maximum s values, array must be valid before calling. void IncompPropertiesBasic::satRange(const int n, - const int* /*cells*/, - double* smin, - double* smax) const + const int* /*cells*/, + double* smin, + double* smax) const { - satprops_.satRange(n, smin, smax); + satprops_.satRange(n, smin, smax); } } // namespace Opm diff --git a/opm/core/fluid/IncompPropertiesBasic.hpp b/opm/core/fluid/IncompPropertiesBasic.hpp index b90b0876f..72c39b1b0 100644 --- a/opm/core/fluid/IncompPropertiesBasic.hpp +++ b/opm/core/fluid/IncompPropertiesBasic.hpp @@ -42,29 +42,29 @@ namespace Opm { public: /// Construct from parameters. - /// The following parameters are accepted (defaults): - /// num_phases (2) Must be 1 or 2. - /// relperm_func ("Linear") Must be "Constant", "Linear" or "Quadratic". - /// rho1 [rho2, rho3] (1.0e3) Density in kg/m^3 - /// mu1 [mu2, mu3] (1.0) Viscosity in cP - /// porosity (1.0) Porosity - /// permeability (100.0) Permeability in mD + /// The following parameters are accepted (defaults): + /// num_phases (2) Must be 1 or 2. + /// relperm_func ("Linear") Must be "Constant", "Linear" or "Quadratic". + /// rho1 [rho2, rho3] (1.0e3) Density in kg/m^3 + /// mu1 [mu2, mu3] (1.0) Viscosity in cP + /// porosity (1.0) Porosity + /// permeability (100.0) Permeability in mD IncompPropertiesBasic(const parameter::ParameterGroup& param, - const int dim, - const int num_cells); + const int dim, + const int num_cells); /// Construct from arguments a basic two phase fluid. IncompPropertiesBasic(const int num_phases, const SaturationPropsBasic::RelPermFunc& relpermfunc, const std::vector& rho, - const std::vector& mu, + const std::vector& mu, const double porosity, const double permeability, const int dim, - const int num_cells); + const int num_cells); - /// Destructor. + /// Destructor. virtual ~IncompPropertiesBasic(); // ---- Rock interface ---- @@ -132,9 +132,9 @@ namespace Opm double* dpcds) const; - /// Obtain the range of allowable saturation values. - /// In cell cells[i], saturation of phase p is allowed to be - /// in the interval [smin[i*P + p], smax[i*P + p]]. + /// Obtain the range of allowable saturation values. + /// In cell cells[i], saturation of phase p is allowed to be + /// in the interval [smin[i*P + p], smax[i*P + p]]. /// \param[in] n Number of data points. /// \param[in] cells Array of n cell indices. /// \param[out] smin Array of nP minimum s values, array must be valid before calling. @@ -145,9 +145,9 @@ namespace Opm double* smax) const; private: RockBasic rock_; - PvtPropertiesBasic pvt_; + PvtPropertiesBasic pvt_; SaturationPropsBasic satprops_; - std::vector viscosity_; + std::vector viscosity_; }; diff --git a/opm/core/fluid/IncompPropertiesFromDeck.cpp b/opm/core/fluid/IncompPropertiesFromDeck.cpp index 45e49791f..9aa690980 100644 --- a/opm/core/fluid/IncompPropertiesFromDeck.cpp +++ b/opm/core/fluid/IncompPropertiesFromDeck.cpp @@ -27,15 +27,15 @@ namespace Opm { IncompPropertiesFromDeck::IncompPropertiesFromDeck(const EclipseGridParser& deck, - const UnstructuredGrid& grid) + const UnstructuredGrid& grid) { rock_.init(deck, grid); - pvt_.init(deck); + pvt_.init(deck); satprops_.init(deck, grid, 200); - if (pvt_.numPhases() != satprops_.numPhases()) { - THROW("IncompPropertiesFromDeck::IncompPropertiesFromDeck() - Inconsistent number of phases in pvt data (" - << pvt_.numPhases() << ") and saturation-dependent function data (" << satprops_.numPhases() << ")."); - } + if (pvt_.numPhases() != satprops_.numPhases()) { + THROW("IncompPropertiesFromDeck::IncompPropertiesFromDeck() - Inconsistent number of phases in pvt data (" + << pvt_.numPhases() << ") and saturation-dependent function data (" << satprops_.numPhases() << ")."); + } } IncompPropertiesFromDeck::~IncompPropertiesFromDeck() @@ -81,19 +81,19 @@ namespace Opm /// \return Array of P viscosity values. const double* IncompPropertiesFromDeck::viscosity() const { - return pvt_.viscosity(); + return pvt_.viscosity(); } /// \return Array of P density values. const double* IncompPropertiesFromDeck::density() const { - return pvt_.reservoirDensities(); + return pvt_.reservoirDensities(); } /// \return Array of P density values. const double* IncompPropertiesFromDeck::surfaceDensity() const { - return pvt_.surfaceDensities(); + return pvt_.surfaceDensities(); } /// \param[in] n Number of data points. @@ -106,10 +106,10 @@ namespace Opm /// m_{ij} = \frac{dkr_i}{ds^j}, /// and is output in Fortran order (m_00 m_10 m_20 m_01 ...) void IncompPropertiesFromDeck::relperm(const int n, - const double* s, - const int* cells, - double* kr, - double* dkrds) const + const double* s, + const int* cells, + double* kr, + double* dkrds) const { satprops_.relperm(n, s, cells, kr, dkrds); } @@ -125,10 +125,10 @@ namespace Opm /// m_{ij} = \frac{dpc_i}{ds^j}, /// and is output in Fortran order (m_00 m_10 m_20 m_01 ...) void IncompPropertiesFromDeck::capPress(const int n, - const double* s, - const int* cells, - double* pc, - double* dpcds) const + const double* s, + const int* cells, + double* pc, + double* dpcds) const { satprops_.capPress(n, s, cells, pc, dpcds); } @@ -142,11 +142,11 @@ namespace Opm /// \param[out] smin Array of nP minimum s values, array must be valid before calling. /// \param[out] smax Array of nP maximum s values, array must be valid before calling. void IncompPropertiesFromDeck::satRange(const int n, - const int* cells, - double* smin, - double* smax) const + const int* cells, + double* smin, + double* smax) const { - satprops_.satRange(n, cells, smin, smax); + satprops_.satRange(n, cells, smin, smax); } } // namespace Opm diff --git a/opm/core/fluid/IncompPropertiesFromDeck.hpp b/opm/core/fluid/IncompPropertiesFromDeck.hpp index 8f917d3b9..6680eaa49 100644 --- a/opm/core/fluid/IncompPropertiesFromDeck.hpp +++ b/opm/core/fluid/IncompPropertiesFromDeck.hpp @@ -47,13 +47,13 @@ namespace Opm public: /// Initialize from deck and grid. /// \param deck Deck input parser - /// \param grid Grid to which property object applies, needed for the + /// \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); + const UnstructuredGrid& grid); - /// Destructor. + /// Destructor. virtual ~IncompPropertiesFromDeck(); // ---- Rock interface ---- @@ -121,9 +121,9 @@ namespace Opm double* dpcds) const; - /// Obtain the range of allowable saturation values. - /// In cell cells[i], saturation of phase p is allowed to be - /// in the interval [smin[i*P + p], smax[i*P + p]]. + /// Obtain the range of allowable saturation values. + /// In cell cells[i], saturation of phase p is allowed to be + /// in the interval [smin[i*P + p], smax[i*P + p]]. /// \param[in] n Number of data points. /// \param[in] cells Array of n cell indices. /// \param[out] smin Array of nP minimum s values, array must be valid before calling. @@ -134,7 +134,7 @@ namespace Opm double* smax) const; private: RockFromDeck rock_; - PvtPropertiesIncompFromDeck pvt_; + PvtPropertiesIncompFromDeck pvt_; SaturationPropsFromDeck satprops_; }; diff --git a/opm/core/fluid/IncompPropertiesInterface.hpp b/opm/core/fluid/IncompPropertiesInterface.hpp index a5025e8a1..4f8bfb516 100644 --- a/opm/core/fluid/IncompPropertiesInterface.hpp +++ b/opm/core/fluid/IncompPropertiesInterface.hpp @@ -109,9 +109,9 @@ namespace Opm double* pc, double* dpcds) const = 0; - /// Obtain the range of allowable saturation values. - /// In cell cells[i], saturation of phase p is allowed to be - /// in the interval [smin[i*P + p], smax[i*P + p]]. + /// Obtain the range of allowable saturation values. + /// In cell cells[i], saturation of phase p is allowed to be + /// in the interval [smin[i*P + p], smax[i*P + p]]. /// \param[in] n Number of data points. /// \param[in] cells Array of n cell indices. /// \param[out] smin Array of nP minimum s values, array must be valid before calling. diff --git a/opm/core/fluid/PvtPropertiesBasic.cpp b/opm/core/fluid/PvtPropertiesBasic.cpp index 3d7d7bb2d..5220502ad 100644 --- a/opm/core/fluid/PvtPropertiesBasic.cpp +++ b/opm/core/fluid/PvtPropertiesBasic.cpp @@ -34,41 +34,41 @@ namespace Opm void PvtPropertiesBasic::init(const parameter::ParameterGroup& param) { - int num_phases = param.getDefault("num_phases", 2); - if (num_phases > 3 || num_phases < 1) { - THROW("PvtPropertiesBasic::init() illegal num_phases: " << num_phases); - } - density_.resize(num_phases); - viscosity_.resize(num_phases); - // We currently do not allow the user to set B. - formation_volume_factor_.clear(); - formation_volume_factor_.resize(num_phases, 1.0); + int num_phases = param.getDefault("num_phases", 2); + if (num_phases > 3 || num_phases < 1) { + THROW("PvtPropertiesBasic::init() illegal num_phases: " << num_phases); + } + density_.resize(num_phases); + viscosity_.resize(num_phases); + // We currently do not allow the user to set B. + formation_volume_factor_.clear(); + formation_volume_factor_.resize(num_phases, 1.0); - // Setting mu and rho from parameters - using namespace Opm::prefix; - using namespace Opm::unit; - const double kgpm3 = kilogram/cubic(meter); - const double cP = centi*Poise; - std::string rname[3] = { "rho1", "rho2", "rho3" }; - double rdefault[3] = { 1.0e3, 1.0e3, 1.0e3 }; - std::string vname[3] = { "mu1", "mu2", "mu3" }; - double vdefault[3] = { 1.0, 1.0, 1.0 }; - for (int phase = 0; phase < num_phases; ++phase) { - density_[phase] = kgpm3*param.getDefault(rname[phase], rdefault[phase]); - viscosity_[phase] = cP*param.getDefault(vname[phase], vdefault[phase]); - } + // Setting mu and rho from parameters + using namespace Opm::prefix; + using namespace Opm::unit; + const double kgpm3 = kilogram/cubic(meter); + const double cP = centi*Poise; + std::string rname[3] = { "rho1", "rho2", "rho3" }; + double rdefault[3] = { 1.0e3, 1.0e3, 1.0e3 }; + std::string vname[3] = { "mu1", "mu2", "mu3" }; + double vdefault[3] = { 1.0, 1.0, 1.0 }; + for (int phase = 0; phase < num_phases; ++phase) { + density_[phase] = kgpm3*param.getDefault(rname[phase], rdefault[phase]); + viscosity_[phase] = cP*param.getDefault(vname[phase], vdefault[phase]); + } } void PvtPropertiesBasic::init(const int num_phases, const std::vector& rho, const std::vector& visc) { - if (num_phases > 3 || num_phases < 1) { - THROW("PvtPropertiesBasic::init() illegal num_phases: " << num_phases); - } - // We currently do not allow the user to set B. - formation_volume_factor_.clear(); - formation_volume_factor_.resize(num_phases, 1.0); + if (num_phases > 3 || num_phases < 1) { + THROW("PvtPropertiesBasic::init() illegal num_phases: " << num_phases); + } + // We currently do not allow the user to set B. + formation_volume_factor_.clear(); + formation_volume_factor_.resize(num_phases, 1.0); density_ = rho; viscosity_ = visc; } @@ -87,69 +87,69 @@ namespace Opm void PvtPropertiesBasic::mu(const int n, - const double* /*p*/, - const double* /*z*/, - double* output_mu) const + const double* /*p*/, + const double* /*z*/, + double* output_mu) const { - const int np = numPhases(); + const int np = numPhases(); for (int phase = 0; phase < np; ++phase) { // #pragma omp parallel for for (int i = 0; i < n; ++i) { output_mu[np*i + phase] = viscosity_[phase]; } - } + } } void PvtPropertiesBasic::B(const int n, - const double* /*p*/, - const double* /*z*/, - double* output_B) const + const double* /*p*/, + const double* /*z*/, + double* output_B) const { - const int np = numPhases(); + const int np = numPhases(); for (int phase = 0; phase < np; ++phase) { // #pragma omp parallel for for (int i = 0; i < n; ++i) { output_B[np*i + phase] = formation_volume_factor_[phase]; } - } + } } void PvtPropertiesBasic::dBdp(const int n, - const double* /*p*/, - const double* /*z*/, - double* output_B, - double* output_dBdp) const + const double* /*p*/, + const double* /*z*/, + double* output_B, + double* output_dBdp) const { - const int np = numPhases(); + const int np = numPhases(); for (int phase = 0; phase < np; ++phase) { // #pragma omp parallel for for (int i = 0; i < n; ++i) { output_B[np*i + phase] = formation_volume_factor_[phase]; output_dBdp[np*i + phase] = 0.0; } - } + } } void PvtPropertiesBasic::R(const int n, - const double* /*p*/, - const double* /*z*/, - double* output_R) const + const double* /*p*/, + const double* /*z*/, + double* output_R) const { - const int np = numPhases(); - std::fill(output_R, output_R + n*np, 0.0); + const int np = numPhases(); + std::fill(output_R, output_R + n*np, 0.0); } void PvtPropertiesBasic::dRdp(const int n, - const double* /*p*/, - const double* /*z*/, - double* output_R, - double* output_dRdp) const + const double* /*p*/, + const double* /*z*/, + double* output_R, + double* output_dRdp) const { - const int np = numPhases(); - std::fill(output_R, output_R + n*np, 0.0); - std::fill(output_dRdp, output_dRdp + n*np, 0.0); + const int np = numPhases(); + std::fill(output_R, output_R + n*np, 0.0); + std::fill(output_dRdp, output_dRdp + n*np, 0.0); } } // namespace Opm diff --git a/opm/core/fluid/PvtPropertiesBasic.hpp b/opm/core/fluid/PvtPropertiesBasic.hpp index d1a3e9f3d..3e30e807c 100644 --- a/opm/core/fluid/PvtPropertiesBasic.hpp +++ b/opm/core/fluid/PvtPropertiesBasic.hpp @@ -38,11 +38,11 @@ namespace Opm PvtPropertiesBasic(); /// Initialize from parameters. - /// The following parameters are accepted (defaults): - /// num_phases (2) Must be 1, 2 or 3. - /// rho1 [rho2, rho3] (1.0e3) Density in kg/m^3 - /// mu1 [mu2, mu3] (1.0) Viscosity in cP - void init(const parameter::ParameterGroup& param); + /// The following parameters are accepted (defaults): + /// num_phases (2) Must be 1, 2 or 3. + /// rho1 [rho2, rho3] (1.0e3) Density in kg/m^3 + /// mu1 [mu2, mu3] (1.0) Viscosity in cP + void init(const parameter::ParameterGroup& param); /// Initialize from arguments. /// Basic multi phase fluid pvt properties. @@ -55,7 +55,7 @@ namespace Opm /// Densities of stock components at surface conditions. /// \return Array of size numPhases(). - const double* surfaceDensities() const; + const double* surfaceDensities() const; /// Viscosity as a function of p and z. void mu(const int n, @@ -90,9 +90,9 @@ namespace Opm double* output_dRdp) const; private: - std::vector density_; - std::vector viscosity_; - std::vector formation_volume_factor_; + std::vector density_; + std::vector viscosity_; + std::vector formation_volume_factor_; }; } diff --git a/opm/core/fluid/PvtPropertiesIncompFromDeck.cpp b/opm/core/fluid/PvtPropertiesIncompFromDeck.cpp index 62cc1b231..34a4ba002 100644 --- a/opm/core/fluid/PvtPropertiesIncompFromDeck.cpp +++ b/opm/core/fluid/PvtPropertiesIncompFromDeck.cpp @@ -38,54 +38,54 @@ namespace Opm { typedef std::vector > > table_t; // If we need multiple regions, this class and the SinglePvt* classes must change. - int region_number = 0; + 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]) { - THROW("PvtPropertiesIncompFromDeck::init() -- must have gas and oil phases (only) in deck input.\n"); - } + if (phase_usage.phase_used[PhaseUsage::Vapour] || + !phase_usage.phase_used[PhaseUsage::Aqua] || + !phase_usage.phase_used[PhaseUsage::Liquid]) { + THROW("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 { - THROW("Input is missing DENSITY\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 { + THROW("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) { - MESSAGE("Compressibility effects in PVTW are ignored."); - } + if (deck.hasField("PVTW")) { + const std::vector& pvtw = deck.getPVTW().pvtw_[region_number]; + if (pvtw[2] != 0.0 || pvtw[4] != 0.0) { + 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; - THROW("Input is missing PVTW\n"); - } + 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; + THROW("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) { - MESSAGE("Compressibility effects in PVCDO are ignored."); - } + if (deck.hasField("PVCDO")) { + const std::vector& pvcdo = deck.getPVCDO().pvcdo_[region_number]; + if (pvcdo[2] != 0.0 || pvcdo[4] != 0.0) { + 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 { - THROW("Input is missing PVCDO\n"); - } + viscosity_[phase_usage.phase_pos[PhaseUsage::Liquid]] = pvcdo[3]; + } else { + THROW("Input is missing PVCDO\n"); + } } const double* PvtPropertiesIncompFromDeck::surfaceDensities() const diff --git a/opm/core/fluid/PvtPropertiesIncompFromDeck.hpp b/opm/core/fluid/PvtPropertiesIncompFromDeck.hpp index acbabaa49..01bbeed1d 100644 --- a/opm/core/fluid/PvtPropertiesIncompFromDeck.hpp +++ b/opm/core/fluid/PvtPropertiesIncompFromDeck.hpp @@ -39,14 +39,14 @@ namespace Opm PvtPropertiesIncompFromDeck(); /// Initialize from deck. - void init(const EclipseGridParser& deck); + void init(const EclipseGridParser& deck); /// Number of active phases. int numPhases() const; /// Densities of stock components at surface conditions. /// \return Array of size numPhases(). - const double* surfaceDensities() const; + const double* surfaceDensities() const; /// Densities of stock components at reservoir conditions. /// Note: a reasonable question to ask is why there can be @@ -58,15 +58,15 @@ namespace Opm /// reporting and using data given in terms of surface values, /// we need to handle this difference. /// \return Array of size numPhases(). - const double* reservoirDensities() const; + const double* reservoirDensities() const; /// Viscosities. const double* viscosity() const; private: - std::tr1::array surface_density_; - std::tr1::array reservoir_density_; - std::tr1::array viscosity_; + std::tr1::array surface_density_; + std::tr1::array reservoir_density_; + std::tr1::array viscosity_; }; } diff --git a/opm/core/fluid/RockBasic.hpp b/opm/core/fluid/RockBasic.hpp index 965637807..2adc3ffd6 100644 --- a/opm/core/fluid/RockBasic.hpp +++ b/opm/core/fluid/RockBasic.hpp @@ -35,9 +35,9 @@ namespace Opm /// Initialize with homogenous porosity and permeability. void init(const int dimensions, - const int num_cells, - const double poro, - const double perm); + const int num_cells, + const double poro, + const double perm); /// \return D, the number of spatial dimensions. int numDimensions() const @@ -66,7 +66,7 @@ namespace Opm } private: - int dimensions_; + int dimensions_; std::vector porosity_; std::vector permeability_; }; diff --git a/opm/core/fluid/RockCompressibility.cpp b/opm/core/fluid/RockCompressibility.cpp index 104861013..af65d852b 100644 --- a/opm/core/fluid/RockCompressibility.cpp +++ b/opm/core/fluid/RockCompressibility.cpp @@ -69,8 +69,8 @@ namespace Opm const double cpnorm = rock_comp_*(pressure - pref_); return (1.0 + cpnorm + 0.5*cpnorm*cpnorm); } else { - // return Opm::linearInterpolation(p_, poromult_, pressure); - return Opm::linearInterpolationExtrap(p_, poromult_, pressure); + // return Opm::linearInterpolation(p_, poromult_, pressure); + return Opm::linearInterpolationExtrap(p_, poromult_, pressure); } } @@ -81,7 +81,7 @@ namespace Opm } else { //const double poromult = Opm::linearInterpolation(p_, poromult_, pressure); //const double dporomultdp = Opm::linearInterpolationDerivative(p_, poromult_, pressure); - const double poromult = Opm::linearInterpolationExtrap(p_, poromult_, pressure); + const double poromult = Opm::linearInterpolationExtrap(p_, poromult_, pressure); const double dporomultdp = Opm::linearInterpolationDerivativeExtrap(p_, poromult_, pressure); return dporomultdp/poromult; diff --git a/opm/core/fluid/RockFromDeck.cpp b/opm/core/fluid/RockFromDeck.cpp index 7ba8903fd..ec7db9a6a 100644 --- a/opm/core/fluid/RockFromDeck.cpp +++ b/opm/core/fluid/RockFromDeck.cpp @@ -51,7 +51,7 @@ namespace Opm /// Initialize from deck and cell mapping. /// \param deck Deck input parser - /// \param grid grid to which property object applies, needed for the + /// \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, diff --git a/opm/core/fluid/RockFromDeck.hpp b/opm/core/fluid/RockFromDeck.hpp index 65027e425..63e71a6e4 100644 --- a/opm/core/fluid/RockFromDeck.hpp +++ b/opm/core/fluid/RockFromDeck.hpp @@ -37,7 +37,7 @@ namespace Opm /// Initialize from deck and grid. /// \param deck Deck input parser - /// \param grid Grid to which property object applies, needed for the + /// \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, diff --git a/opm/core/fluid/SaturationPropsBasic.cpp b/opm/core/fluid/SaturationPropsBasic.cpp index 760f934d4..10ddf23dd 100644 --- a/opm/core/fluid/SaturationPropsBasic.cpp +++ b/opm/core/fluid/SaturationPropsBasic.cpp @@ -29,64 +29,64 @@ namespace Opm namespace { - struct KrFunConstant - { - double kr(double) - { - return 1.0; - } - double dkrds(double) - { - return 0.0; - } - }; + struct KrFunConstant + { + double kr(double) + { + return 1.0; + } + double dkrds(double) + { + return 0.0; + } + }; - struct KrFunLinear - { - double kr(double s) - { - return s; - } - double dkrds(double) - { - return 1.0; - } - }; + struct KrFunLinear + { + double kr(double s) + { + return s; + } + double dkrds(double) + { + return 1.0; + } + }; - struct KrFunQuadratic - { - double kr(double s) - { - return s*s; - } - double dkrds(double s) - { - return 2.0*s; - } - }; + struct KrFunQuadratic + { + double kr(double s) + { + return s*s; + } + double dkrds(double s) + { + return 2.0*s; + } + }; - template - static inline void evalAllKrDeriv(const int n, const int np, - const double* s, double* kr, double* dkrds, Fun fun) - { - if (dkrds == 0) { + template + static inline void evalAllKrDeriv(const int n, const int np, + const double* s, double* kr, double* dkrds, Fun fun) + { + if (dkrds == 0) { // #pragma omp parallel for - for (int i = 0; i < n*np; ++i) { - kr[i] = fun.kr(s[i]); - } - return; - } + for (int i = 0; i < n*np; ++i) { + kr[i] = fun.kr(s[i]); + } + return; + } // #pragma omp parallel for - for (int i = 0; i < n; ++i) { - std::fill(dkrds + i*np*np, dkrds + (i+1)*np*np, 0.0); - for (int phase = 0; phase < np; ++phase) { - kr[i*np + phase] = fun.kr(s[i*np + phase]); - // Only diagonal elements in derivative. - dkrds[i*np*np + phase*np + phase] = fun.dkrds(s[i*np + phase]); - } - } - } + for (int i = 0; i < n; ++i) { + std::fill(dkrds + i*np*np, dkrds + (i+1)*np*np, 0.0); + for (int phase = 0; phase < np; ++phase) { + kr[i*np + phase] = fun.kr(s[i*np + phase]); + // Only diagonal elements in derivative. + dkrds[i*np*np + phase*np + phase] = fun.dkrds(s[i*np + phase]); + } + } + } } // anon namespace @@ -109,25 +109,25 @@ namespace Opm /// Initialize from parameters. void SaturationPropsBasic::init(const parameter::ParameterGroup& param) { - int num_phases = param.getDefault("num_phases", 2); - if (num_phases > 2 || num_phases < 1) { - THROW("SaturationPropsBasic::init() illegal num_phases: " << num_phases); - } + int num_phases = param.getDefault("num_phases", 2); + if (num_phases > 2 || num_phases < 1) { + THROW("SaturationPropsBasic::init() illegal num_phases: " << num_phases); + } num_phases_ = num_phases; - //std::string rpf = param.getDefault("relperm_func", std::string("Unset")); - std::string rpf = param.getDefault("relperm_func", std::string("Linear")); - if (rpf == "Constant") { - relperm_func_ = Constant; - if(num_phases!=1){ - THROW("Constant relperm with more than one phase???"); - } - } else if (rpf == "Linear") { - relperm_func_ = Linear; - } else if (rpf == "Quadratic") { - relperm_func_ = Quadratic; - } else { - THROW("SaturationPropsBasic::init() illegal relperm_func: " << rpf); - } + //std::string rpf = param.getDefault("relperm_func", std::string("Unset")); + std::string rpf = param.getDefault("relperm_func", std::string("Linear")); + if (rpf == "Constant") { + relperm_func_ = Constant; + if(num_phases!=1){ + THROW("Constant relperm with more than one phase???"); + } + } else if (rpf == "Linear") { + relperm_func_ = Linear; + } else if (rpf == "Quadratic") { + relperm_func_ = Quadratic; + } else { + THROW("SaturationPropsBasic::init() illegal relperm_func: " << rpf); + } } @@ -136,7 +136,7 @@ namespace Opm /// \return P, the number of phases. int SaturationPropsBasic::numPhases() const { - return num_phases_; + return num_phases_; } @@ -152,29 +152,29 @@ namespace Opm /// m_{ij} = \frac{dkr_i}{ds^j}, /// and is output in Fortran order (m_00 m_10 m_20 m01 ...) void SaturationPropsBasic::relperm(const int n, - const double* s, - double* kr, - double* dkrds) const + const double* s, + double* kr, + double* dkrds) const { - switch (relperm_func_) { - case Constant: - { - evalAllKrDeriv(n, num_phases_, s, kr, dkrds, KrFunConstant()); - break; - } - case Linear: - { - evalAllKrDeriv(n, num_phases_, s, kr, dkrds, KrFunLinear()); - break; - } - case Quadratic: - { - evalAllKrDeriv(n, num_phases_, s, kr, dkrds, KrFunQuadratic()); - break; - } - default: - THROW("SaturationPropsBasic::relperm() unhandled relperm func type: " << relperm_func_); - } + switch (relperm_func_) { + case Constant: + { + evalAllKrDeriv(n, num_phases_, s, kr, dkrds, KrFunConstant()); + break; + } + case Linear: + { + evalAllKrDeriv(n, num_phases_, s, kr, dkrds, KrFunLinear()); + break; + } + case Quadratic: + { + evalAllKrDeriv(n, num_phases_, s, kr, dkrds, KrFunQuadratic()); + break; + } + default: + THROW("SaturationPropsBasic::relperm() unhandled relperm func type: " << relperm_func_); + } } @@ -190,13 +190,13 @@ namespace Opm /// m_{ij} = \frac{dpc_i}{ds^j}, /// and is output in Fortran order (m_00 m_10 m_20 m01 ...) void SaturationPropsBasic::capPress(const int n, - const double* /*s*/, - double* pc, - double* dpcds) const + const double* /*s*/, + double* pc, + double* dpcds) const { - std::fill(pc, pc + num_phases_*n, 0.0); + std::fill(pc, pc + num_phases_*n, 0.0); if (dpcds) { - std::fill(dpcds, dpcds + num_phases_*num_phases_*n, 0.0); + std::fill(dpcds, dpcds + num_phases_*num_phases_*n, 0.0); } } @@ -207,11 +207,11 @@ namespace Opm /// \param[out] smin Array of nP minimum s values, array must be valid before calling. /// \param[out] smax Array of nP maximum s values, array must be valid before calling. void SaturationPropsBasic::satRange(const int n, - double* smin, - double* smax) const + double* smin, + double* smax) const { - std::fill(smin, smin + num_phases_*n, 0.0); - std::fill(smax, smax + num_phases_*n, 1.0); + std::fill(smin, smin + num_phases_*n, 0.0); + std::fill(smax, smax + num_phases_*n, 1.0); } diff --git a/opm/core/fluid/SaturationPropsBasic.hpp b/opm/core/fluid/SaturationPropsBasic.hpp index 789e0a351..2f90672d8 100644 --- a/opm/core/fluid/SaturationPropsBasic.hpp +++ b/opm/core/fluid/SaturationPropsBasic.hpp @@ -40,16 +40,16 @@ namespace Opm SaturationPropsBasic(); /// Initialize from parameters. - /// The following parameters are accepted (defaults): - /// num_phases (2) Must be 1 or 2. - /// relperm_func ("Linear") Must be "Constant", "Linear" or "Quadratic". + /// The following parameters are accepted (defaults): + /// num_phases (2) Must be 1 or 2. + /// relperm_func ("Linear") Must be "Constant", "Linear" or "Quadratic". void init(const parameter::ParameterGroup& param); - enum RelPermFunc { Constant, Linear, Quadratic }; + enum RelPermFunc { Constant, Linear, Quadratic }; /// Initialize from arguments a basic Saturation property. void init(const int num_phases, - const RelPermFunc& relperm_func) + const RelPermFunc& relperm_func) { num_phases_ = num_phases; relperm_func_ = relperm_func; @@ -86,18 +86,18 @@ namespace Opm double* pc, double* dpcds) const; - /// Obtain the range of allowable saturation values. + /// Obtain the range of allowable saturation values. /// \param[in] n Number of data points. /// \param[out] smin Array of nP minimum s values, array must be valid before calling. /// \param[out] smax Array of nP maximum s values, array must be valid before calling. - void satRange(const int n, - double* smin, - double* smax) const; + void satRange(const int n, + double* smin, + double* smax) const; private: - int num_phases_; - RelPermFunc relperm_func_; + int num_phases_; + RelPermFunc relperm_func_; }; diff --git a/opm/core/fluid/SaturationPropsFromDeck.hpp b/opm/core/fluid/SaturationPropsFromDeck.hpp index 42461985c..0a0599ec8 100644 --- a/opm/core/fluid/SaturationPropsFromDeck.hpp +++ b/opm/core/fluid/SaturationPropsFromDeck.hpp @@ -50,7 +50,7 @@ namespace Opm /// Initialize from deck and grid. /// \param[in] deck Deck input parser - /// \param[in] grid Grid to which property object applies, needed for the + /// \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. @@ -92,14 +92,14 @@ namespace Opm double* pc, double* dpcds) const; - /// Obtain the range of allowable saturation values. + /// Obtain the range of allowable saturation values. /// \param[in] n Number of data points. /// \param[out] smin Array of nP minimum s values, array must be valid before calling. /// \param[out] smax Array of nP maximum s values, array must be valid before calling. - void satRange(const int n, + void satRange(const int n, const int* cells, - double* smin, - double* smax) const; + double* smin, + double* smax) const; private: PhaseUsage phase_usage_; diff --git a/opm/core/fluid/SaturationPropsFromDeck_impl.hpp b/opm/core/fluid/SaturationPropsFromDeck_impl.hpp index a4df7f425..f7a80b453 100644 --- a/opm/core/fluid/SaturationPropsFromDeck_impl.hpp +++ b/opm/core/fluid/SaturationPropsFromDeck_impl.hpp @@ -105,7 +105,7 @@ namespace Opm template int SaturationPropsFromDeck::numPhases() const { - return phase_usage_.num_phases; + return phase_usage_.num_phases; } @@ -191,18 +191,18 @@ namespace Opm template void SaturationPropsFromDeck::satRange(const int n, const int* cells, - double* smin, - double* smax) const + double* smin, + double* smax) const { ASSERT (cells != 0); - const int np = phase_usage_.num_phases; - for (int i = 0; i < n; ++i) { - for (int p = 0; p < np; ++p) { - smin[np*i + p] = funcForCell(cells[i]).smin_[p]; - smax[np*i + p] = funcForCell(cells[i]).smax_[p]; - } - } + const int np = phase_usage_.num_phases; + for (int i = 0; i < n; ++i) { + for (int p = 0; p < np; ++p) { + smin[np*i + p] = funcForCell(cells[i]).smin_[p]; + smax[np*i + p] = funcForCell(cells[i]).smax_[p]; + } + } } diff --git a/opm/core/fluid/SaturationPropsInterface.hpp b/opm/core/fluid/SaturationPropsInterface.hpp index e55995677..30f277cc1 100644 --- a/opm/core/fluid/SaturationPropsInterface.hpp +++ b/opm/core/fluid/SaturationPropsInterface.hpp @@ -65,11 +65,11 @@ namespace Opm double* pc, double* dpcds) const = 0; - /// Obtain the range of allowable saturation values. + /// Obtain the range of allowable saturation values. /// \param[in] n Number of data points. /// \param[out] smin Array of nP minimum s values, array must be valid before calling. /// \param[out] smax Array of nP maximum s values, array must be valid before calling. - virtual void satRange(const int n, + virtual void satRange(const int n, const int* cells, double* smin, double* smax) const = 0; diff --git a/opm/core/fluid/blackoil/BlackoilPvtProperties.hpp b/opm/core/fluid/blackoil/BlackoilPvtProperties.hpp index d3a42e321..8633fd9ee 100644 --- a/opm/core/fluid/blackoil/BlackoilPvtProperties.hpp +++ b/opm/core/fluid/blackoil/BlackoilPvtProperties.hpp @@ -70,7 +70,7 @@ namespace Opm /// Densities of stock components at surface conditions. /// \return Array of size numPhases(). - const double* surfaceDensities() const; + const double* surfaceDensities() const; /// Viscosity as a function of p and z. void mu(const int n, @@ -111,11 +111,11 @@ namespace Opm PhaseUsage phase_usage_; - int region_number_; + int region_number_; std::vector > props_; - double densities_[MaxNumPhases]; + double densities_[MaxNumPhases]; mutable std::vector data1_; mutable std::vector data2_; };