OilfieldToolsNet/opm-simulators
4
0
Fork 0
mirror of https://github.com/OPM/opm-simulators.git synced 2025-01-18 21:52:57 -06:00
Code Issues Packages Projects Releases Wiki Activity

Merge branch 'master' of github.com:OPM/opm-core

Browse Source
This commit is contained in:
babrodtk 2016-08-30 14:53:11 +02:00
commit 273c99b36f
48 changed files with 1088 additions and 780 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

6
examples/compute_initial_state.cpp
View File

@ -36,6 +36,8 @@
#include <boost/filesystem.hpp>
#include <fstream>
namespace
{
void warnIfUnusedParams(const Opm::parameter::ParameterGroup& param)
@ -88,9 +90,9 @@ try
Opm::ParseContext parseContext;
Opm::ParserPtr parser(new Opm::Parser() );
Opm::DeckConstPtr deck = parser->parseFile(deck_filename , parseContext);
Opm::EclipseStateConstPtr eclipseState(new Opm::EclipseState(deck, parseContext));
Opm::EclipseStateConstPtr eclipseState(new Opm::EclipseState(*deck, parseContext));
const double grav = param.getDefault("gravity", unit::gravity);
GridManager gm(deck);
GridManager gm(eclipseState->getInputGrid());
const UnstructuredGrid& grid = *gm.c_grid();
BlackoilPropertiesFromDeck props(deck, eclipseState, grid, param);
warnIfUnusedParams(param);

2
examples/compute_tof_from_files.cpp
View File

@ -54,7 +54,7 @@
#include <vector>
#include <numeric>
#include <iterator>
#include <fstream>
namespace
{

14
examples/diagnose_relperm.cpp
View File

@ -76,9 +76,9 @@ try
{ ParseContext::INTERNAL_ERROR_UNINITIALIZED_THPRES, InputError::IGNORE}
});
Opm::DeckConstPtr deck(parser->parseFile(eclipseFilename, parseContext));
eclState.reset(new EclipseState(deck, parseContext));
eclState.reset(new EclipseState(*deck, parseContext));
GridManager gm(deck);
GridManager gm(eclState->getInputGrid());
const UnstructuredGrid& grid = *gm.c_grid();
using boost::filesystem::path;
path fpath(eclipseFilename);
@ -92,13 +92,13 @@ try
std::string logFile = baseName + ".SATFUNCLOG";
std::shared_ptr<EclipsePRTLog> prtLog = std::make_shared<EclipsePRTLog>(logFile, Log::DefaultMessageTypes);
OpmLog::addBackend( "ECLIPSEPRTLOG" , prtLog );
Opm::time::StopWatch timer;
timer.start();
prtLog->setMessageFormatter(std::make_shared<SimpleMessageFormatter>(true, false));
std::shared_ptr<StreamLog> streamLog = std::make_shared<EclipsePRTLog>(std::cout, Log::DefaultMessageTypes);
OpmLog::addBackend( "STREAMLOG" , streamLog );
streamLog->setMessageLimiter(std::make_shared<MessageLimiter>(10));
streamLog->setMessageFormatter(std::make_shared<SimpleMessageFormatter>(true, true));
RelpermDiagnostics diagnostic;
diagnostic.diagnosis(eclState, deck, grid);
timer.stop();
double tt = timer.secsSinceStart();
std::cout << "relperm diagnostics: " << tt << " seconds." << std::endl;
}
catch (const std::exception &e) {
std::cerr << "Program threw an exception: " << e.what() << "\n";

13
opm/core/linalg/LinearSolverFactory.cpp
View File

@ -61,8 +61,19 @@ namespace Opm
LinearSolverFactory::LinearSolverFactory(const parameter::ParameterGroup& param)
{
#if HAVE_SUITESPARSE_UMFPACK_H
std::string default_solver = "umfpack";
#elif HAVE_DUNE_ISTL
std::string default_solver = "istl";
#elif HAVE_PETSC
std::string default_solver = "petsc";
#else
std::string default_solver = "no_solver_available";
OPM_THROW(std::runtime_error, "No linear solver available, you must have UMFPACK , dune-istl or Petsc installed to use LinearSolverFactory.");
#endif
const std::string ls =
param.getDefault<std::string>("linsolver", "umfpack");
param.getDefault("linsolver", default_solver);
if (ls == "umfpack") {
#if HAVE_SUITESPARSE_UMFPACK_H

51
opm/core/linalg/ParallelIstlInformation.hpp
View File

@ -179,12 +179,24 @@ public:
}
return ownerMask_;
}
/// \brief Get the owner Mask.
///
/// \return A vector with entries 0, and 1. 0 marks an index that we cannot
/// compute correct results for. 1 marks an index that this process
/// is responsible for and computes correct results in parallel.
const std::vector<double>& getOwnerMask() const
{
return ownerMask_;
}
/// \brief Compute one or more global reductions.
///
/// This function can either be used with a container, an operator, and an initial value
/// to compute a reduction. Or with tuples of them to compute multiple reductions with only
/// one global communication.
/// The possible functors needed can be constructed with Opm::Reduction::makeGlobalMaxFunctor(),
/// Opm::Reduction::makeLInfinityNormFunctor(),
/// Opm::Reduction::makeGlobalMinFunctor(), and
/// Opm::Reduction::makeGlobalSumFunctor().
/// \tparam type of the container or the tuple of containers.
@ -574,6 +586,45 @@ private:
(std::pointer_to_binary_function<const T&,const T&,const T&>
((const T&(*)(const T&, const T&))std::max<T>));
}
namespace detail
{
/// \brief Computes the maximum of the absolute values of two values.
template<typename T, typename Enable = void>
struct MaxAbsFunctor
{
using result_type = T;
result_type operator()(const T& t1,
const T& t2)
{
return std::max(std::abs(t1), std::abs(t2));
}
};
// Specialization for unsigned integers. They need their own
// version since abs(x) is ambiguous (as well as somewhat
// meaningless).
template<typename T>
struct MaxAbsFunctor<T, typename std::enable_if<std::is_unsigned<T>::value>::type>
{
using result_type = T;
result_type operator()(const T& t1,
const T& t2)
{
return std::max(t1, t2);
}
};
}
/// \brief Create a functor for computing a global L infinity norm
///
/// To be used with ParallelISTLInformation::computeReduction.
template<class T>
MaskIDOperator<detail::MaxAbsFunctor<T> >
makeLInfinityNormFunctor()
{
return MaskIDOperator<detail::MaxAbsFunctor<T> >();
}
/// \brief Create a functor for computing a global minimum.
///
/// To be used with ParallelISTLInformation::computeReduction.

3
opm/core/props/BlackoilPhases.hpp
View File

@ -60,6 +60,9 @@ namespace Opm
void setFreeOil () { insert(BlackoilPhases::Liquid); }
void setFreeGas () { insert(BlackoilPhases::Vapour); }
bool operator==(const PhasePresence& other) const { return present_ == other.present_; }
bool operator!=(const PhasePresence& other) const { return !this->operator==(other); }
private:
unsigned char present_;

188
opm/core/props/BlackoilPropertiesFromDeck.cpp
View File

@ -244,15 +244,15 @@ namespace Opm
const auto& pu = phaseUsage();
const int np = numPhases();
typedef Opm::LocalAd::Evaluation<double, /*size=*/1> LadEval;
typedef Opm::DenseAd::Evaluation<double, /*size=*/1> Eval;
LadEval pLad = 0.0;
LadEval TLad = 0.0;
LadEval RsLad = 0.0;
LadEval RvLad = 0.0;
LadEval muLad = 0.0;
Eval pEval = 0.0;
Eval TEval = 0.0;
Eval RsEval = 0.0;
Eval RvEval = 0.0;
Eval muEval = 0.0;
pLad.derivatives[0] = 1.0;
pEval.derivatives[0] = 1.0;
R_.resize(n*np);
this->compute_R_(n, p, T, z, cells, &R_[0]);
@ -260,30 +260,30 @@ namespace Opm
for (int i = 0; i < n; ++ i) {
int cellIdx = cells[i];
int pvtRegionIdx = cellPvtRegionIdx_[cellIdx];
pLad.value = p[i];
TLad.value = T[i];
pEval.value = p[i];
TEval.value = T[i];
if (pu.phase_used[BlackoilPhases::Aqua]) {
muLad = waterPvt_.viscosity(pvtRegionIdx, TLad, pLad);
muEval = waterPvt_.viscosity(pvtRegionIdx, TEval, pEval);
int offset = pu.num_phases*cellIdx + pu.phase_pos[BlackoilPhases::Aqua];
mu[offset] = muLad.value;
dmudp[offset] = muLad.derivatives[0];
mu[offset] = muEval.value;
dmudp[offset] = muEval.derivatives[0];
}
if (pu.phase_used[BlackoilPhases::Liquid]) {
RsLad.value = R_[i*np + pu.phase_pos[BlackoilPhases::Liquid]];
muLad = oilPvt_.viscosity(pvtRegionIdx, TLad, pLad, RsLad);
RsEval.value = R_[i*np + pu.phase_pos[BlackoilPhases::Liquid]];
muEval = oilPvt_.viscosity(pvtRegionIdx, TEval, pEval, RsEval);
int offset = pu.num_phases*cellIdx + pu.phase_pos[BlackoilPhases::Liquid];
mu[offset] = muLad.value;
dmudp[offset] = muLad.derivatives[0];
mu[offset] = muEval.value;
dmudp[offset] = muEval.derivatives[0];
}
if (pu.phase_used[BlackoilPhases::Vapour]) {
RvLad.value = R_[i*np + pu.phase_pos[BlackoilPhases::Vapour]];
muLad = gasPvt_.viscosity(pvtRegionIdx, TLad, pLad, RvLad);
RvEval.value = R_[i*np + pu.phase_pos[BlackoilPhases::Vapour]];
muEval = gasPvt_.viscosity(pvtRegionIdx, TEval, pEval, RvEval);
int offset = pu.num_phases*cellIdx + pu.phase_pos[BlackoilPhases::Vapour];
mu[offset] = muLad.value;
dmudp[offset] = muLad.derivatives[0];
mu[offset] = muEval.value;
dmudp[offset] = muEval.derivatives[0];
}
}
}
@ -396,27 +396,27 @@ namespace Opm
{
const auto& pu = phaseUsage();
typedef double LadEval;
typedef double Eval;
LadEval pLad = 0.0;
LadEval TLad = 0.0;
LadEval RsLad = 0.0;
LadEval RvLad = 0.0;
Eval pEval = 0.0;
Eval TEval = 0.0;
Eval RsEval = 0.0;
Eval RvEval = 0.0;
for (int i = 0; i < n; ++ i) {
int cellIdx = cells[i];
int pvtRegionIdx = cellPvtRegionIdx_[cellIdx];
pLad = p[i];
TLad = T[i];
pEval = p[i];
TEval = T[i];
int oilOffset = pu.num_phases*i + pu.phase_pos[BlackoilPhases::Liquid];
int gasOffset = pu.num_phases*i + pu.phase_pos[BlackoilPhases::Vapour];
int waterOffset = pu.num_phases*i + pu.phase_pos[BlackoilPhases::Aqua];
if (pu.phase_used[BlackoilPhases::Aqua]) {
LadEval BLad = 1.0/waterPvt_.inverseFormationVolumeFactor(pvtRegionIdx, TLad, pLad);
Eval BEval = 1.0/waterPvt_.inverseFormationVolumeFactor(pvtRegionIdx, TEval, pEval);
B[waterOffset] = BLad;
B[waterOffset] = BEval;
}
if (pu.phase_used[BlackoilPhases::Liquid]) {
@ -424,18 +424,18 @@ namespace Opm
double maxRs = 0.0;
if (pu.phase_used[BlackoilPhases::Vapour]) {
currentRs = (z[oilOffset] == 0.0) ? 0.0 : z[gasOffset]/z[oilOffset];
maxRs = oilPvt_.saturatedGasDissolutionFactor(pvtRegionIdx, TLad, pLad);
maxRs = oilPvt_.saturatedGasDissolutionFactor(pvtRegionIdx, TEval, pEval);
}
LadEval BLad;
Eval BEval;
if (currentRs >= maxRs) {
BLad = 1.0/oilPvt_.saturatedInverseFormationVolumeFactor(pvtRegionIdx, TLad, pLad);
BEval = 1.0/oilPvt_.saturatedInverseFormationVolumeFactor(pvtRegionIdx, TEval, pEval);
}
else {
RsLad = currentRs;
BLad = 1.0/oilPvt_.inverseFormationVolumeFactor(pvtRegionIdx, TLad, pLad, RsLad);
RsEval = currentRs;
BEval = 1.0/oilPvt_.inverseFormationVolumeFactor(pvtRegionIdx, TEval, pEval, RsEval);
}
B[oilOffset] = BLad;
B[oilOffset] = BEval;
}
if (pu.phase_used[BlackoilPhases::Vapour]) {
@ -443,18 +443,18 @@ namespace Opm
double maxRv = 0.0;
if (pu.phase_used[BlackoilPhases::Liquid]) {
currentRv = (z[gasOffset] == 0.0) ? 0.0 : z[oilOffset]/z[gasOffset];
maxRv = gasPvt_.saturatedOilVaporizationFactor(pvtRegionIdx, TLad, pLad);
maxRv = gasPvt_.saturatedOilVaporizationFactor(pvtRegionIdx, TEval, pEval);
}
LadEval BLad;
Eval BEval;
if (currentRv >= maxRv) {
BLad = 1.0/gasPvt_.saturatedInverseFormationVolumeFactor(pvtRegionIdx, TLad, pLad);
BEval = 1.0/gasPvt_.saturatedInverseFormationVolumeFactor(pvtRegionIdx, TEval, pEval);
}
else {
RvLad = currentRv;
BLad = 1.0/gasPvt_.inverseFormationVolumeFactor(pvtRegionIdx, TLad, pLad, RvLad);
RvEval = currentRv;
BEval = 1.0/gasPvt_.inverseFormationVolumeFactor(pvtRegionIdx, TEval, pEval, RvEval);
}
B[gasOffset] = BLad;
B[gasOffset] = BEval;
}
}
}
@ -469,30 +469,30 @@ namespace Opm
{
const auto& pu = phaseUsage();
typedef Opm::LocalAd::Evaluation<double, /*size=*/1> LadEval;
typedef Opm::DenseAd::Evaluation<double, /*size=*/1> Eval;
LadEval pLad = 0.0;
LadEval TLad = 0.0;
LadEval RsLad = 0.0;
LadEval RvLad = 0.0;
Eval pEval = 0.0;
Eval TEval = 0.0;
Eval RsEval = 0.0;
Eval RvEval = 0.0;
pLad.derivatives[0] = 1.0;
pEval.derivatives[0] = 1.0;
for (int i = 0; i < n; ++ i) {
int cellIdx = cells[i];
int pvtRegionIdx = cellPvtRegionIdx_[cellIdx];
pLad.value = p[i];
TLad.value = T[i];
pEval.value = p[i];
TEval.value = T[i];
int oilOffset = pu.num_phases*i + pu.phase_pos[BlackoilPhases::Liquid];
int gasOffset = pu.num_phases*i + pu.phase_pos[BlackoilPhases::Vapour];
int waterOffset = pu.num_phases*i + pu.phase_pos[BlackoilPhases::Aqua];
if (pu.phase_used[BlackoilPhases::Aqua]) {
LadEval BLad = 1.0/waterPvt_.inverseFormationVolumeFactor(pvtRegionIdx, TLad, pLad);
Eval BEval = 1.0/waterPvt_.inverseFormationVolumeFactor(pvtRegionIdx, TEval, pEval);
B[waterOffset] = BLad.value;
dBdp[waterOffset] = BLad.derivatives[0];
B[waterOffset] = BEval.value;
dBdp[waterOffset] = BEval.derivatives[0];
}
if (pu.phase_used[BlackoilPhases::Liquid]) {
@ -500,19 +500,19 @@ namespace Opm
double maxRs = 0.0;
if (pu.phase_used[BlackoilPhases::Vapour]) {
currentRs = (z[oilOffset] == 0.0) ? 0.0 : z[gasOffset]/z[oilOffset];
maxRs = oilPvt_.saturatedGasDissolutionFactor(pvtRegionIdx, TLad.value, pLad.value);
maxRs = oilPvt_.saturatedGasDissolutionFactor(pvtRegionIdx, TEval.value, pEval.value);
}
LadEval BLad;
Eval BEval;
if (currentRs >= maxRs) {
BLad = 1.0/oilPvt_.saturatedInverseFormationVolumeFactor(pvtRegionIdx, TLad, pLad);
BEval = 1.0/oilPvt_.saturatedInverseFormationVolumeFactor(pvtRegionIdx, TEval, pEval);
}
else {
RsLad.value = currentRs;
BLad = 1.0/oilPvt_.inverseFormationVolumeFactor(pvtRegionIdx, TLad, pLad, RsLad);
RsEval.value = currentRs;
BEval = 1.0/oilPvt_.inverseFormationVolumeFactor(pvtRegionIdx, TEval, pEval, RsEval);
}
B[oilOffset] = BLad.value;
dBdp[oilOffset] = BLad.derivatives[0];
B[oilOffset] = BEval.value;
dBdp[oilOffset] = BEval.derivatives[0];
}
if (pu.phase_used[BlackoilPhases::Vapour]) {
@ -520,19 +520,19 @@ namespace Opm
double maxRv = 0.0;
if (pu.phase_used[BlackoilPhases::Liquid]) {
currentRv = (z[gasOffset] == 0.0) ? 0.0 : z[oilOffset]/z[gasOffset];
maxRv = gasPvt_.saturatedOilVaporizationFactor(pvtRegionIdx, TLad.value, pLad.value);
maxRv = gasPvt_.saturatedOilVaporizationFactor(pvtRegionIdx, TEval.value, pEval.value);
}
LadEval BLad;
Eval BEval;
if (currentRv >= maxRv) {
BLad = 1.0/gasPvt_.saturatedInverseFormationVolumeFactor(pvtRegionIdx, TLad, pLad);
BEval = 1.0/gasPvt_.saturatedInverseFormationVolumeFactor(pvtRegionIdx, TEval, pEval);
}
else {
RvLad.value = currentRv;
BLad = 1.0/gasPvt_.inverseFormationVolumeFactor(pvtRegionIdx, TLad, pLad, RvLad);
RvEval.value = currentRv;
BEval = 1.0/gasPvt_.inverseFormationVolumeFactor(pvtRegionIdx, TEval, pEval, RvEval);
}
B[gasOffset] = BLad.value;
dBdp[gasOffset] = BLad.derivatives[0];
B[gasOffset] = BEval.value;
dBdp[gasOffset] = BEval.derivatives[0];
}
}
}
@ -546,16 +546,16 @@ namespace Opm
{
const auto& pu = phaseUsage();
typedef double LadEval;
typedef double Eval;
LadEval pLad = 0.0;
LadEval TLad = 0.0;
Eval pEval = 0.0;
Eval TEval = 0.0;
for (int i = 0; i < n; ++ i) {
int cellIdx = cells[i];
int pvtRegionIdx = cellPvtRegionIdx_[cellIdx];
pLad = p[i];
TLad = T[i];
pEval = p[i];
TEval = T[i];
int oilOffset = pu.num_phases*i + pu.phase_pos[BlackoilPhases::Liquid];
int gasOffset = pu.num_phases*i + pu.phase_pos[BlackoilPhases::Vapour];
@ -566,29 +566,29 @@ namespace Opm
}
if (pu.phase_used[BlackoilPhases::Liquid]) {
LadEval RsSatLad = oilPvt_.saturatedGasDissolutionFactor(pvtRegionIdx, TLad, pLad);
Eval RsSatEval = oilPvt_.saturatedGasDissolutionFactor(pvtRegionIdx, TEval, pEval);
double currentRs = 0.0;
if (pu.phase_used[BlackoilPhases::Vapour]) {
currentRs = (z[oilOffset] == 0.0) ? 0.0 : z[gasOffset]/z[oilOffset];
}
RsSatLad = std::min(RsSatLad, currentRs);
RsSatEval = std::min(RsSatEval, currentRs);
R[oilOffset] = RsSatLad;
R[oilOffset] = RsSatEval;
}
if (pu.phase_used[BlackoilPhases::Vapour]) {
LadEval RvSatLad = gasPvt_.saturatedOilVaporizationFactor(pvtRegionIdx, TLad, pLad);
Eval RvSatEval = gasPvt_.saturatedOilVaporizationFactor(pvtRegionIdx, TEval, pEval);
double currentRv = 0.0;
if (pu.phase_used[BlackoilPhases::Liquid]) {
currentRv = (z[gasOffset] == 0.0) ? 0.0 : z[oilOffset]/z[gasOffset];
}
RvSatLad = std::min(RvSatLad, currentRv);
RvSatEval = std::min(RvSatEval, currentRv);
R[gasOffset] = RvSatLad;
R[gasOffset] = RvSatEval;
}
}
}
@ -603,19 +603,19 @@ namespace Opm
{
const auto& pu = phaseUsage();
typedef Opm::LocalAd::Evaluation<double, /*size=*/1> LadEval;
typedef Opm::MathToolbox<LadEval> Toolbox;
typedef Opm::DenseAd::Evaluation<double, /*size=*/1> Eval;
typedef Opm::MathToolbox<Eval> Toolbox;
LadEval pLad = 0.0;
LadEval TLad = 0.0;
Eval pEval = 0.0;
Eval TEval = 0.0;
pLad.derivatives[0] = 1.0;
pEval.derivatives[0] = 1.0;
for (int i = 0; i < n; ++ i) {
int cellIdx = cells[i];
int pvtRegionIdx = cellPvtRegionIdx_[cellIdx];
pLad.value = p[i];
TLad.value = T[i];
pEval.value = p[i];
TEval.value = T[i];
int oilOffset = pu.num_phases*i + pu.phase_pos[BlackoilPhases::Liquid];
int gasOffset = pu.num_phases*i + pu.phase_pos[BlackoilPhases::Vapour];
@ -626,31 +626,31 @@ namespace Opm
}
if (pu.phase_used[BlackoilPhases::Liquid]) {
LadEval RsSatLad = oilPvt_.saturatedGasDissolutionFactor(pvtRegionIdx, TLad, pLad);
Eval RsSatEval = oilPvt_.saturatedGasDissolutionFactor(pvtRegionIdx, TEval, pEval);
LadEval currentRs = 0.0;
Eval currentRs = 0.0;
if (pu.phase_used[BlackoilPhases::Vapour]) {
currentRs = (z[oilOffset] == 0.0) ? 0.0 : z[gasOffset]/z[oilOffset];
}
RsSatLad = Toolbox::min(RsSatLad, currentRs);
RsSatEval = Toolbox::min(RsSatEval, currentRs);
R[oilOffset] = RsSatLad.value;
dRdp[oilOffset] = RsSatLad.derivatives[0];
R[oilOffset] = RsSatEval.value;
dRdp[oilOffset] = RsSatEval.derivatives[0];
}
if (pu.phase_used[BlackoilPhases::Vapour]) {
LadEval RvSatLad = gasPvt_.saturatedOilVaporizationFactor(pvtRegionIdx, TLad, pLad);
Eval RvSatEval = gasPvt_.saturatedOilVaporizationFactor(pvtRegionIdx, TEval, pEval);
LadEval currentRv = 0.0;
Eval currentRv = 0.0;
if (pu.phase_used[BlackoilPhases::Liquid]) {
currentRv = (z[gasOffset] == 0.0) ? 0.0 : z[oilOffset]/z[gasOffset];
}
RvSatLad = Toolbox::min(RvSatLad, currentRv);
RvSatEval = Toolbox::min(RvSatEval, currentRv);
R[gasOffset] = RvSatLad.value;
dRdp[gasOffset] = RvSatLad.derivatives[0];
R[gasOffset] = RvSatEval.value;
dRdp[gasOffset] = RvSatEval.derivatives[0];
}
}
}

10
opm/core/props/phaseUsageFromDeck.hpp
View File

@ -88,8 +88,14 @@ namespace Opm
}
pu.num_phases = 0;
for (int i = 0; i < BlackoilPhases::MaxNumPhases; ++i) {
pu.phase_pos[i] = pu.num_phases;
pu.num_phases += pu.phase_used[i];
if (pu.phase_used[i]) {
pu.phase_pos[i] = pu.num_phases;
pu.num_phases += 1;
}
else {
//Set to ridiculous value on purpose: should never be used
pu.phase_pos[i] = 2000000000;
}
}
// Only 2 or 3 phase systems handled.

14
opm/core/props/rock/RockCompressibility.cpp
View File

@ -22,6 +22,7 @@
#include <opm/core/utility/parameters/ParameterGroup.hpp>
#include <opm/core/utility/Units.hpp>
#include <opm/common/ErrorMacros.hpp>
#include <opm/common/OpmLog/OpmLog.hpp>
#include <opm/core/utility/linearInterpolation.hpp>
#include <opm/parser/eclipse/EclipseState/Tables/RocktabTable.hpp>
@ -64,14 +65,18 @@ namespace Opm
if (rockKeyword.size() != 1) {
// here it would be better not to use std::cout directly but to add the
// warning to some "warning list"...
std::cout << "Can only handle a single region in ROCK ("<<rockKeyword.size()<<" regions specified)."
<< " Ignoring all except for the first.\n";
OpmLog::warning("Can only handle a single region in ROCK ("
+ std::to_string(rockKeyword.size())
+ " regions specified)."
+ " Ignoring all except for the first.\n"
+ "In file " + rockKeyword.getFileName()
+ ", line " + std::to_string(rockKeyword.getLineNumber()) + "\n");
}
pref_ = rockKeyword.getRecord(0).getItem("PREF").getSIDouble(0);
rock_comp_ = rockKeyword.getRecord(0).getItem("COMPRESSIBILITY").getSIDouble(0);
} else {
std::cout << "**** warning: no rock compressibility data found in deck (ROCK or ROCKTAB)." << std::endl;
OpmLog::warning("No rock compressibility data found in deck (ROCK or ROCKTAB).");
}
}
@ -96,7 +101,8 @@ namespace Opm
if (p_.empty()) {
// Approximating poro multiplier with a quadratic curve,
// we must use its derivative.
return rock_comp_ + 2 * rock_comp_ * rock_comp_ * (pressure - pref_);
const double cpnorm = rock_comp_*(pressure - pref_);
return rock_comp_ + cpnorm*rock_comp_;
} else {
return Opm::linearInterpolationDerivative(p_, poromult_, pressure);
}

243
opm/core/props/satfunc/RelpermDiagnostics.cpp
View File

@ -24,29 +24,6 @@
namespace Opm{
RelpermDiagnostics::Counter::Counter()
:error(0)
,warning(0)
,problem(0)
,bug(0)
{
}
std::vector<std::string>
RelpermDiagnostics::getMessages() const
{
return messages_;
}
void RelpermDiagnostics::phaseCheck_(DeckConstPtr deck)
{
bool hasWater = deck->hasKeyword("WATER");
@ -56,31 +33,26 @@ namespace Opm{
if (hasWater && hasGas && !hasOil && !hasSolvent) {
const std::string msg = "System: Water-Gas system.";
std::cout << msg << std::endl;
OpmLog::info(msg);
fluidSystem_ = FluidSystem::WaterGas;
}
if (hasWater && hasOil && !hasGas && !hasSolvent) {
const std::string msg = "System: Oil-Water system.";
std::cout << msg << std::endl;
OpmLog::info(msg);
fluidSystem_ = FluidSystem::OilWater;
}
if (hasOil && hasGas && !hasWater && !hasSolvent) {
const std::string msg = "System: Oil-Gas system.";
std::cout << msg << std::endl;
OpmLog::info(msg);
fluidSystem_ = FluidSystem::OilGas;
}
if (hasOil && hasWater && hasGas && !hasSolvent) {
const std::string msg = "System: Black-oil system.";
std::cout << msg << std::endl;
OpmLog::info(msg);
fluidSystem_ = FluidSystem::BlackOil;
}
if (hasSolvent) {
const std::string msg = "System: Solvent model.";
std::cout << msg << std::endl;
OpmLog::info(msg);
fluidSystem_ = FluidSystem::Solvent;
}
@ -107,31 +79,25 @@ namespace Opm{
bool family2 = ((!swfnTables.empty() && !sgfnTables.empty()) || !sgwfnTables.empty()) && (!sof3Tables.empty() || !sof2Tables.empty());
if (family1 && family2) {
const std::string msg = "-- Error: Saturation families should not be mixed.\n Use either SGOF and SWOF or SGFN, SWFN and SOF3.";
messages_.push_back(msg);
const std::string msg = "Saturation families should not be mixed.\n Use either SGOF and SWOF or SGFN, SWFN and SOF3.";
OpmLog::error(msg);
counter_.error += 1;
}
if (!family1 && !family2) {
const std::string msg = "-- Error, Saturations function must be specified using either \n \
const std::string msg = "Saturations function must be specified using either \n \
family 1 or family 2 keywords \n \
Use either SGOF and SWOF or SGFN, SWFN and SOF3.";
messages_.push_back(msg);
OpmLog::error(msg);
counter_.error += 1;
}
if (family1 && !family2) {
satFamily_ = SaturationFunctionFamily::FamilyI;
const std::string msg = "Relative permeability input format: Saturation Family I.";
std::cout << msg << std::endl;
OpmLog::info(msg);
}
if (!family1 && family2) {
satFamily_ = SaturationFunctionFamily::FamilyII;
const std::string msg = "Relative permeability input format: Saturation Family II.";
std::cout << msg << std::endl;
OpmLog::info(msg);
}
}
@ -144,7 +110,6 @@ namespace Opm{
{
const int numSatRegions = deck->getKeyword("TABDIMS").getRecord(0).getItem("NTSFUN").get< int >(0);
const std::string msg = "Number of saturation regions: " + std::to_string(numSatRegions) + "\n";
std::cout << msg << std::endl;
OpmLog::info(msg);
const auto& tableManager = eclState->getTableManager();
const TableContainer& swofTables = tableManager.getSwofTables();
@ -217,32 +182,24 @@ namespace Opm{
const std::string regionIdx = std::to_string(satnumIdx);
//Check sw column.
if (sw.front() < 0.0 || sw.back() > 1.0) {
const std::string msg = "-- Error: In SWOF table SATNUM = "+ regionIdx + ", saturation should be in range [0,1].";
messages_.push_back(msg);
const std::string msg = "In SWOF table SATNUM = "+ regionIdx + ", saturation should be in range [0,1].";
OpmLog::error(msg);
counter_.error += 1;
}
//TODO check endpoint sw.back() == 1. - Sor.
//Check krw column.
if (krw.front() != 0.0) {
const std::string msg = "-- Error: In SWOF table SATNUM = " + regionIdx + ", first value of krw should be 0.";
messages_.push_back(msg);
const std::string msg = "In SWOF table SATNUM = " + regionIdx + ", first value of krw should be 0.";
OpmLog::error(msg);
counter_.error += 1;
}
if (krw.front() < 0.0 || krw.back() > 1.0) {
const std::string msg = "-- Error: In SWOF table SATNUM = " + regionIdx + ", krw should be in range [0,1].";
messages_.push_back(msg);
const std::string msg = "In SWOF table SATNUM = " + regionIdx + ", krw should be in range [0,1].";
OpmLog::error(msg);
counter_.error += 1;
}
///Check krow column.
if (krow.front() > 1.0 || krow.back() < 0.0) {
const std::string msg = "-- Error: In SWOF table SATNUM = "+ regionIdx + ", krow should be in range [0, 1].";
messages_.push_back(msg);
const std::string msg = "In SWOF table SATNUM = "+ regionIdx + ", krow should be in range [0, 1].";
OpmLog::error(msg);
counter_.error += 1;
}
///TODO check if run with gas.
}
@ -260,38 +217,28 @@ namespace Opm{
const std::string regionIdx = std::to_string(satnumIdx);
//Check sw column.
if (sg.front() < 0.0 || sg.back() > 1.0) {
const std::string msg = "-- Error: In SGOF table SATNUM = " + regionIdx + ", saturation should be in range [0,1].";
messages_.push_back(msg);
const std::string msg = "In SGOF table SATNUM = " + regionIdx + ", saturation should be in range [0,1].";
OpmLog::error(msg);
counter_.error += 1;
}
if (sg.front() != 0.0) {
const std::string msg = "-- Error: In SGOF table SATNUM = " + regionIdx + ", first value of sg should be 0.";
messages_.push_back(msg);
const std::string msg = "In SGOF table SATNUM = " + regionIdx + ", first value of sg should be 0.";
OpmLog::error(msg);
counter_.error += 1;
}
//TODO check endpoint sw.back() == 1. - Sor.
//Check krw column.
if (krg.front() != 0.0) {
const std::string msg = "-- Error: In SGOF table SATNUM = " + regionIdx + ", first value of krg should be 0.";
messages_.push_back(msg);
const std::string msg = "In SGOF table SATNUM = " + regionIdx + ", first value of krg should be 0.";
OpmLog::error(msg);
counter_.error += 1;
}
if (krg.front() < 0.0 || krg.back() > 1.0) {
const std::string msg = "-- Error: In SGOF table SATNUM = " + regionIdx + ", krg should be in range [0,1].";
messages_.push_back(msg);
const std::string msg = "In SGOF table SATNUM = " + regionIdx + ", krg should be in range [0,1].";
OpmLog::error(msg);
counter_.error += 1;
}
//Check krow column.
if (krog.front() > 1.0 || krog.back() < 0.0) {
const std::string msg = "-- Error: In SGOF table SATNUM = " + regionIdx + ", krog should be in range [0, 1].";
messages_.push_back(msg);
const std::string msg = "In SGOF table SATNUM = " + regionIdx + ", krog should be in range [0, 1].";
OpmLog::error(msg);
counter_.error += 1;
}
//TODO check if run with water.
}
@ -306,36 +253,26 @@ namespace Opm{
//Check sl column.
//TODO first value means sl = swco + sor
if (sl.front() < 0.0 || sl.back() > 1.0) {
const std::string msg = "-- Error: In SLGOF table SATNUM = " + regionIdx + ", saturation should be in range [0,1].";
messages_.push_back(msg);
const std::string msg = "In SLGOF table SATNUM = " + regionIdx + ", saturation should be in range [0,1].";
OpmLog::error(msg);
counter_.error += 1;
}
if (sl.back() != 1.0) {
const std::string msg = "-- Error: In SLGOF table SATNUM = " + regionIdx + ", last value of sl should be 1.";
messages_.push_back(msg);
const std::string msg = "In SLGOF table SATNUM = " + regionIdx + ", last value of sl should be 1.";
OpmLog::error(msg);
counter_.error += 1;
}
if (krg.front() > 1.0 || krg.back() < 0) {
const std::string msg = "-- Error: In SLGOF table SATNUM = " + regionIdx + ", krg shoule be in range [0, 1].";
messages_.push_back(msg);
const std::string msg = "In SLGOF table SATNUM = " + regionIdx + ", krg shoule be in range [0, 1].";
OpmLog::error(msg);
counter_.error += 1;
}
if (krg.back() != 0.0) {
const std::string msg = "-- Error: In SLGOF table SATNUM = " + regionIdx + ", last value of krg hould be 0.";
messages_.push_back(msg);
const std::string msg = "In SLGOF table SATNUM = " + regionIdx + ", last value of krg hould be 0.";
OpmLog::error(msg);
counter_.error += 1;
}
if (krog.front() < 0.0 || krog.back() > 1.0) {
const std::string msg = "-- Error: In SLGOF table SATNUM = " + regionIdx + ", krog shoule be in range [0, 1].";
messages_.push_back(msg);
const std::string msg = "In SLGOF table SATNUM = " + regionIdx + ", krog shoule be in range [0, 1].";
OpmLog::error(msg);
counter_.error += 1;
}
}
@ -351,25 +288,19 @@ namespace Opm{
const std::string regionIdx = std::to_string(satnumIdx);
//Check sw column.
if (sw.front() < 0.0 || sw.back() > 1.0) {
const std::string msg = "-- Error: In SWFN table SATNUM = " + regionIdx + ", saturation should be in range [0,1].";
messages_.push_back(msg);
const std::string msg = "In SWFN table SATNUM = " + regionIdx + ", saturation should be in range [0,1].";
OpmLog::error(msg);
counter_.error += 1;
}
//Check krw column.
if (krw.front() < 0.0 || krw.back() > 1.0) {
const std::string msg = "-- Error: In SWFN table SATNUM = " + regionIdx + ", krw should be in range [0,1].";
messages_.push_back(msg);
const std::string msg = "In SWFN table SATNUM = " + regionIdx + ", krw should be in range [0,1].";
OpmLog::error(msg);
counter_.error += 1;
}
if (krw.front() != 0.0) {
const std::string msg = "-- Error: In SWFN table SATNUM = " + regionIdx + ", first value of krw should be 0.";
messages_.push_back(msg);
const std::string msg = "In SWFN table SATNUM = " + regionIdx + ", first value of krw should be 0.";
OpmLog::error(msg);
counter_.error += 1;
}
}
@ -385,24 +316,18 @@ namespace Opm{
const std::string regionIdx = std::to_string(satnumIdx);
//Check sg column.
if (sg.front() < 0.0 || sg.back() > 1.0) {
const std::string msg = "-- Error: In SGFN table SATNUM = " + regionIdx + ", saturation should be in range [0,1].";
messages_.push_back(msg);
const std::string msg = "In SGFN table SATNUM = " + regionIdx + ", saturation should be in range [0,1].";
OpmLog::error(msg);
counter_.error += 1;
}
//Check krg column.
if (krg.front() < 0.0 || krg.back() > 1.0) {
const std::string msg = "-- Error: In SGFN table SATNUM = " + regionIdx + ", krg should be in range [0,1].";
messages_.push_back(msg);
const std::string msg = "In SGFN table SATNUM = " + regionIdx + ", krg should be in range [0,1].";
OpmLog::error(msg);
counter_.error += 1;
}
if (krg.front() != 0.0) {
const std::string msg = "-- Error: In SGFN table SATNUM = " + regionIdx + ", first value of krg should be 0.";
messages_.push_back(msg);
const std::string msg = "In SGFN table SATNUM = " + regionIdx + ", first value of krg should be 0.";
OpmLog::error(msg);
counter_.error += 1;
}
}
@ -420,46 +345,34 @@ namespace Opm{
//Check so column.
//TODO: The max so = 1 - Swco
if (so.front() < 0.0 || so.back() > 1.0) {
const std::string msg = "-- Error: In SOF3 table SATNUM = " + regionIdx + ", saturation should be in range [0,1].";
messages_.push_back(msg);
const std::string msg = "In SOF3 table SATNUM = " + regionIdx + ", saturation should be in range [0,1].";
OpmLog::error(msg);
counter_.error += 1;
}
//Check krow column.
if (krow.front() < 0.0 || krow.back() > 1.0) {
const std::string msg = "-- Error: In SOF3 table SATNUM = " + regionIdx + ", krow should be in range [0,1].";
messages_.push_back(msg);
const std::string msg = "In SOF3 table SATNUM = " + regionIdx + ", krow should be in range [0,1].";
OpmLog::error(msg);
counter_.error += 1;
}
if (krow.front() != 0.0) {
const std::string msg = "-- Error: In SOF3 table SATNUM = " + regionIdx + ", first value of krow should be 0.";
messages_.push_back(msg);
const std::string msg = "In SOF3 table SATNUM = " + regionIdx + ", first value of krow should be 0.";
OpmLog::error(msg);
counter_.error += 1;
}
//Check krog column.
if (krog.front() < 0.0 || krog.back() > 1.0) {
const std::string msg = "-- Error: In SOF3 table SATNUM = " + regionIdx + ", krog should be in range [0,1].";
messages_.push_back(msg);
const std::string msg = "In SOF3 table SATNUM = " + regionIdx + ", krog should be in range [0,1].";
OpmLog::error(msg);
counter_.error += 1;
}
if (krog.front() != 0.0) {
const std::string msg = "-- Error: In SOF3 table SATNUM = " + regionIdx + ", first value of krog should be 0.";
messages_.push_back(msg);
const std::string msg = "In SOF3 table SATNUM = " + regionIdx + ", first value of krog should be 0.";
OpmLog::error(msg);
counter_.error += 1;
}
if (krog.back() != krow.back()) {
const std::string msg = "-- Error: In SOF3 table SATNUM = " + regionIdx + ", max value of krog and krow should be the same.";
messages_.push_back(msg);
const std::string msg = "In SOF3 table SATNUM = " + regionIdx + ", max value of krog and krow should be the same.";
OpmLog::error(msg);
counter_.error += 1;
}
}
@ -476,24 +389,18 @@ namespace Opm{
//Check so column.
//TODO: The max so = 1 - Swco
if (so.front() < 0.0 || so.back() > 1.0) {
const std::string msg = "-- Error: In SOF2 table SATNUM = " + regionIdx + ", saturation should be in range [0,1].";
messages_.push_back(msg);
const std::string msg = "In SOF2 table SATNUM = " + regionIdx + ", saturation should be in range [0,1].";
OpmLog::error(msg);
counter_.error += 1;
}
//Check krow column.
if (kro.front() < 0.0 || kro.back() > 1.0) {
const std::string msg = "-- Error: In SOF2 table SATNUM = " + regionIdx + ", krow should be in range [0,1].";
messages_.push_back(msg);
const std::string msg = "In SOF2 table SATNUM = " + regionIdx + ", krow should be in range [0,1].";
OpmLog::error(msg);
counter_.error += 1;
}
if (kro.front() != 0.0) {
const std::string msg = "-- Error: In SOF2 table SATNUM = " + regionIdx + ", first value of krow should be 0.";
messages_.push_back(msg);
const std::string msg = "In SOF2 table SATNUM = " + regionIdx + ", first value of krow should be 0.";
OpmLog::error(msg);
counter_.error += 1;
}
}
@ -510,39 +417,29 @@ namespace Opm{
const std::string regionIdx = std::to_string(satnumIdx);
//Check sg column.
if (sg.front() < 0.0 || sg.back() > 1.0) {
const std::string msg = "-- Error: In SGWFN table SATNUM = " + regionIdx + ", saturation should be in range [0,1].";
messages_.push_back(msg);
const std::string msg = "In SGWFN table SATNUM = " + regionIdx + ", saturation should be in range [0,1].";
OpmLog::error(msg);
counter_.error += 1;
}
//Check krg column.
if (krg.front() < 0.0 || krg.back() > 1.0) {
const std::string msg = "-- Error: In SGWFN table SATNUM = " + regionIdx + ", krg should be in range [0,1].";
messages_.push_back(msg);
const std::string msg = "In SGWFN table SATNUM = " + regionIdx + ", krg should be in range [0,1].";
OpmLog::error(msg);
counter_.error += 1;
}
if (krg.front() != 0.0) {
const std::string msg = "-- Error: In SGWFN table SATNUM = " + regionIdx + ", first value of krg should be 0.";
messages_.push_back(msg);
const std::string msg = "In SGWFN table SATNUM = " + regionIdx + ", first value of krg should be 0.";
OpmLog::error(msg);
counter_.error += 1;
}
//Check krgw column.
//TODO check saturation sw = 1. - sg
if (krgw.front() > 1.0 || krgw.back() < 0.0) {
const std::string msg = "-- Error: In SGWFN table SATNUM = " + regionIdx + ", krgw should be in range [0,1].";
messages_.push_back(msg);
const std::string msg = "In SGWFN table SATNUM = " + regionIdx + ", krgw should be in range [0,1].";
OpmLog::error(msg);
counter_.error += 1;
}
if (krgw.back() != 0.0) {
const std::string msg = "-- Error: In SGWFN table SATNUM = " + regionIdx + ", last value of krgw should be 0.";
messages_.push_back(msg);
const std::string msg = "In SGWFN table SATNUM = " + regionIdx + ", last value of krgw should be 0.";
OpmLog::error(msg);
counter_.error += 1;
}
}
@ -556,18 +453,14 @@ namespace Opm{
const std::string regionIdx = std::to_string(satnumIdx);
//Check sw column.
if (sw.front() < 0.0 || sw.back() > 1.0) {
const std::string msg = "-- Error: In SGCWMIS table SATNUM = " + regionIdx + ", saturation should be in range [0,1].";
messages_.push_back(msg);
const std::string msg = "In SGCWMIS table SATNUM = " + regionIdx + ", saturation should be in range [0,1].";
OpmLog::error(msg);
counter_.error += 1;
}
//Check critical gas column.
if (sgc.front() < 0.0 || sgc.back() > 1.0) {
const std::string msg = "-- Error: In SGCWMIS table SATNUM = " + regionIdx + ", critical gas saturation should be in range [0,1].";
messages_.push_back(msg);
const std::string msg = "In SGCWMIS table SATNUM = " + regionIdx + ", critical gas saturation should be in range [0,1].";
OpmLog::error(msg);
counter_.error += 1;
}
}
@ -583,18 +476,14 @@ namespace Opm{
const std::string regionIdx = std::to_string(satnumIdx);
//Check sw column.
if (sw.front() < 0.0 || sw.back() > 1.0) {
const std::string msg = "-- Error: In SORWMIS table SATNUM = " + regionIdx + ", saturation should be in range [0,1].";
messages_.push_back(msg);
const std::string msg = "In SORWMIS table SATNUM = " + regionIdx + ", saturation should be in range [0,1].";
OpmLog::error(msg);
counter_.error += 1;
}
//Check critical oil column.
if (sor.front() < 0.0 || sor.back() > 1.0) {
const std::string msg = "-- Error: In SORWMIS table SATNUM = " + regionIdx + ", critical oil saturation should be in range [0,1].";
messages_.push_back(msg);
const std::string msg = "In SORWMIS table SATNUM = " + regionIdx + ", critical oil saturation should be in range [0,1].";
OpmLog::error(msg);
counter_.error += 1;
}
}
@ -610,26 +499,20 @@ namespace Opm{
const std::string regionIdx = std::to_string(satnumIdx);
//Check phase fraction column.
if (frac.front() < 0.0 || frac.back() > 1.0) {
const std::string msg = "-- Error: In SSFN table SATNUM = " + regionIdx + ", phase fraction should be in range [0,1].";
messages_.push_back(msg);
const std::string msg = "In SSFN table SATNUM = " + regionIdx + ", phase fraction should be in range [0,1].";
OpmLog::error(msg);
counter_.error += 1;
}
//Check gas relperm multiplier column.
if (krgm.front() < 0.0 || krgm.back() > 1.0) {
const std::string msg = "-- Error: In SSFN table SATNUM = " + regionIdx + ", gas relative permeability multiplier should be in range [0,1].";
messages_.push_back(msg);
const std::string msg = "In SSFN table SATNUM = " + regionIdx + ", gas relative permeability multiplier should be in range [0,1].";
OpmLog::error(msg);
counter_.error += 1;
}
//Check solvent relperm multiplier column.
if (krsm.front() < 0.0 || krsm.back() > 1.0) {
const std::string msg = "-- Error: In SSFN table SATNUM = " + regionIdx + ", solvent relative permeability multiplier should be in range [0,1].";
messages_.push_back(msg);
const std::string msg = "In SSFN table SATNUM = " + regionIdx + ", solvent relative permeability multiplier should be in range [0,1].";
OpmLog::error(msg);
counter_.error += 1;
}
}
@ -647,18 +530,14 @@ namespace Opm{
const std::string regionIdx = std::to_string(satnumIdx);
//Check phase fraction column.
if (frac.front() < 0.0 || frac.back() > 1.0) {
const std::string msg = "-- Error: In MISC table SATNUM = " + regionIdx + ", phase fraction should be in range [0,1].";
messages_.push_back(msg);
const std::string msg = "In MISC table SATNUM = " + regionIdx + ", phase fraction should be in range [0,1].";
OpmLog::error(msg);
counter_.error += 1;
}
//Check miscibility column.
if (misc.front() < 0.0 || misc.back() > 1.0) {
const std::string msg = "-- Error: In MISC table SATNUM = " + regionIdx + ", miscibility should be in range [0,1].";
messages_.push_back(msg);
const std::string msg = "In MISC table SATNUM = " + regionIdx + ", miscibility should be in range [0,1].";
OpmLog::error(msg);
counter_.error += 1;
}
}
@ -676,26 +555,20 @@ namespace Opm{
const std::string regionIdx = std::to_string(satnumIdx);
//Check phase fraction column.
if (frac.front() < 0.0 || frac.back() > 1.0) {
const std::string msg = "-- Error: In MSFN table SATNUM = " + regionIdx + ", total gas fraction should be in range [0,1].";
messages_.push_back(msg);
const std::string msg = "In MSFN table SATNUM = " + regionIdx + ", total gas fraction should be in range [0,1].";
OpmLog::error(msg);
counter_.error += 1;
}
//Check gas_solvent relperm multiplier column.
if (krgsm.front() < 0.0 || krgsm.back() > 1.0) {
const std::string msg = "-- Error: In MSFN table SATNUM = " + regionIdx + ", gas+solvent relative permeability multiplier should be in range [0,1].";
messages_.push_back(msg);
const std::string msg = "In MSFN table SATNUM = " + regionIdx + ", gas+solvent relative permeability multiplier should be in range [0,1].";
OpmLog::error(msg);
counter_.error += 1;
}
//Check oil relperm multiplier column.
if (krom.front() > 1.0 || krom.back() < 0.0) {
const std::string msg = "-- Error: In MSFN table SATNUM = " + regionIdx + ", oil relative permeability multiplier should be in range [0,1].";
messages_.push_back(msg);
const std::string msg = "In MSFN table SATNUM = " + regionIdx + ", oil relative permeability multiplier should be in range [0,1].";
OpmLog::error(msg);
counter_.error += 1;
}
}
@ -722,16 +595,12 @@ namespace Opm{
const std::string regionIdx = std::to_string(satnumIdx + 1);
///Consistency check.
if (unscaledEpsInfo_[satnumIdx].Sgu > (1. - unscaledEpsInfo_[satnumIdx].Swl)) {
const std::string msg = "-- Warning: In saturation table SATNUM = " + regionIdx + ", Sgmax should not exceed 1-Swco.";
messages_.push_back(msg);
const std::string msg = "In saturation table SATNUM = " + regionIdx + ", Sgmax should not exceed 1-Swco.";
OpmLog::warning(msg);
counter_.warning += 1;
}
if (unscaledEpsInfo_[satnumIdx].Sgl > (1. - unscaledEpsInfo_[satnumIdx].Swu)) {
const std::string msg = "-- Warning: In saturation table SATNUM = " + regionIdx + ", Sgco should not exceed 1-Swmax.";
messages_.push_back(msg);
const std::string msg = "In saturation table SATNUM = " + regionIdx + ", Sgco should not exceed 1-Swmax.";
OpmLog::warning(msg);
counter_.warning += 1;
}
//Krow(Sou) == Krog(Sou) for three-phase
@ -762,25 +631,19 @@ namespace Opm{
krog_value = table.evaluate("KROG" , Sou);
}
if (krow_value != krog_value) {
const std::string msg = "-- Warning: In saturation table SATNUM = " + regionIdx + ", Krow(Somax) should be equal to Krog(Somax).";
messages_.push_back(msg);
const std::string msg = "In saturation table SATNUM = " + regionIdx + ", Krow(Somax) should be equal to Krog(Somax).";
OpmLog::warning(msg);
counter_.warning += 1;
}
}
//Krw(Sw=0)=Krg(Sg=0)=Krow(So=0)=Krog(So=0)=0.
//Mobile fluid requirements
if (((unscaledEpsInfo_[satnumIdx].Sowcr + unscaledEpsInfo_[satnumIdx].Swcr)-1) >= 0) {
const std::string msg = "-- Warning: In saturation table SATNUM = " + regionIdx + ", Sowcr + Swcr should be less than 1.";
messages_.push_back(msg);
const std::string msg = "In saturation table SATNUM = " + regionIdx + ", Sowcr + Swcr should be less than 1.";
OpmLog::warning(msg);
counter_.warning += 1;
}
if (((unscaledEpsInfo_[satnumIdx].Sogcr + unscaledEpsInfo_[satnumIdx].Sgcr + unscaledEpsInfo_[satnumIdx].Swl) - 1 ) > 0) {
const std::string msg = "-- Warning: In saturation table SATNUM = " + regionIdx + ", Sogcr + Sgcr + Swco should be less than 1.";
messages_.push_back(msg);
const std::string msg = "In saturation table SATNUM = " + regionIdx + ", Sogcr + Sgcr + Swco should be less than 1.";
OpmLog::warning(msg);
counter_.warning += 1;
}
}
}

15
opm/core/props/satfunc/RelpermDiagnostics.hpp
View File

@ -62,8 +62,6 @@ namespace Opm {
DeckConstPtr deck,
const GridT& grid);
std::vector<std::string> getMessages() const;
private:
enum FluidSystem {
OilWater,
@ -83,22 +81,9 @@ namespace Opm {
SaturationFunctionFamily satFamily_;
struct Counter {
Counter();
int error;
int warning;
int problem;
int bug;
};
Counter counter_;
std::vector<Opm::EclEpsScalingPointsInfo<double> > unscaledEpsInfo_;
std::vector<Opm::EclEpsScalingPointsInfo<double> > scaledEpsInfo_;
std::vector<std::string> messages_;
///Store scaled information.
std::vector<std::string> scaled_messages_;
///Check the phase that used.
void phaseCheck_(DeckConstPtr deck);

94
opm/core/props/satfunc/RelpermDiagnostics_impl.hpp
View File

@ -33,41 +33,12 @@ namespace Opm {
Opm::DeckConstPtr deck,
const GridT& grid)
{
std::cout << "\n\n***************Saturation Functions Diagnostics***************\n\n";
OpmLog::info("\n***************Saturation Functions Diagnostics***************");
phaseCheck_(deck);
satFamilyCheck_(eclState);
tableCheck_(eclState, deck);
unscaledEndPointsCheck_(deck, eclState);
scaledEndPointsCheck_(deck, eclState, grid);
if (!messages_.empty()) {
std::sort(messages_.begin(), messages_.end());
auto it = std::unique(messages_.begin(), messages_.end());
messages_.erase(it, messages_.end());
std::cout << std::endl;
for (const auto& x : messages_) {
std::cout << x << std::endl;
}
}
int limits = 0;
if (!scaled_messages_.empty()) {
std::cout << std::endl;
for (const auto& x : scaled_messages_) {
if (limits < 10) {
std::cout << x << std::endl;
limits++;
} else {
std::cout << "\nMore inconsistencies exist. Check saturation function input and LOGFILE!" << std::endl;
break;
}
}
}
const std::string summary_msg = "\n\nError summary:" +
std::string("\nWarnings " + std::to_string(counter_.warning)) +
std::string("\nProblems " + std::to_string(counter_.problem)) +
std::string("\nErrors " + std::to_string(counter_.error)) +
std::string("\nBugs " + std::to_string(counter_.bug))+ "\n";
std::cout << summary_msg << std::endl;
}
template <class GridT>
@ -83,7 +54,8 @@ namespace Opm {
EclEpsGridProperties epsGridProperties;
epsGridProperties.initFromDeck(deck, eclState, /*imbibition=*/false);
const auto& satnum = eclState->get3DProperties().getIntGridProperty("SATNUM");
const std::string tag = "Scaled endpoints";
for (int c = 0; c < nc; ++c) {
const int cartIdx = compressedToCartesianIdx[c];
const std::string satnumIdx = std::to_string(satnum.iget(cartIdx));
@ -98,82 +70,62 @@ namespace Opm {
// SGU <= 1.0 - SWL
if (scaledEpsInfo_[c].Sgu > (1.0 - scaledEpsInfo_[c].Swl)) {
const std::string msg = "-- Warning: For scaled endpoints input, cell" + cellIdx + " SATNUM = " + satnumIdx + ", SGU exceed 1.0 - SWL";
scaled_messages_.push_back(msg);
OpmLog::warning(msg);
counter_.warning += 1;
const std::string msg = "For scaled endpoints input, cell" + cellIdx + " SATNUM = " + satnumIdx + ", SGU exceed 1.0 - SWL";
OpmLog::warning(tag, msg);
}
// SGL <= 1.0 - SWU
if (scaledEpsInfo_[c].Sgl > (1.0 - scaledEpsInfo_[c].Swu)) {
const std::string msg = "-- Warning: For scaled endpoints input, cell" + cellIdx + " SATNUM = " + satnumIdx + ", SGL exceed 1.0 - SWU";
scaled_messages_.push_back(msg);
OpmLog::warning(msg);
counter_.warning += 1;
const std::string msg = "For scaled endpoints input, cell" + cellIdx + " SATNUM = " + satnumIdx + ", SGL exceed 1.0 - SWU";
OpmLog::warning(tag, msg);
}
if (deck->hasKeyword("SCALECRS") && fluidSystem_ == FluidSystem::BlackOil) {
// Mobilility check.
if ((scaledEpsInfo_[c].Sowcr + scaledEpsInfo_[c].Swcr) >= 1.0) {
const std::string msg = "-- Warning: For scaled endpoints input, cell" + cellIdx + " SATNUM = " + satnumIdx + ", SOWCR + SWCR exceed 1.0";
scaled_messages_.push_back(msg);
OpmLog::warning(msg);
counter_.warning += 1;
const std::string msg = "For scaled endpoints input, cell" + cellIdx + " SATNUM = " + satnumIdx + ", SOWCR + SWCR exceed 1.0";
OpmLog::warning(tag, msg);
}
if ((scaledEpsInfo_[c].Sogcr + scaledEpsInfo_[c].Sgcr + scaledEpsInfo_[c].Swl) >= 1.0) {
const std::string msg = "-- Warning: For scaled endpoints input, cell" + cellIdx + " SATNUM = " + satnumIdx + ", SOGCR + SGCR + SWL exceed 1.0";
scaled_messages_.push_back(msg);
OpmLog::warning(msg);
counter_.warning += 1;
const std::string msg = "For scaled endpoints input, cell" + cellIdx + " SATNUM = " + satnumIdx + ", SOGCR + SGCR + SWL exceed 1.0";
OpmLog::warning(tag, msg);
}
}
///Following rules come from NEXUS.
if (fluidSystem_ != FluidSystem::WaterGas) {
if (scaledEpsInfo_[c].Swl > scaledEpsInfo_[c].Swcr) {
const std::string msg = "-- Warning: For scaled endpoints input, cell" + cellIdx + " SATNUM = " + satnumIdx + ", SWL > SWCR";
scaled_messages_.push_back(msg);
OpmLog::warning(msg);
counter_.warning += 1;
const std::string msg = "For scaled endpoints input, cell" + cellIdx + " SATNUM = " + satnumIdx + ", SWL > SWCR";
OpmLog::warning(tag, msg);
}
if (scaledEpsInfo_[c].Swcr > scaledEpsInfo_[c].Sowcr) {
const std::string msg = "-- Warning: For scaled endpoints input, cell" + cellIdx + " SATNUM = " + satnumIdx + ", SWCR > SOWCR";
scaled_messages_.push_back(msg);
OpmLog::warning(msg);
counter_.warning += 1;
const std::string msg = "For scaled endpoints input, cell" + cellIdx + " SATNUM = " + satnumIdx + ", SWCR > SOWCR";
OpmLog::warning(tag, msg);
}
if (scaledEpsInfo_[c].Sowcr > scaledEpsInfo_[c].Swu) {
const std::string msg = "-- Warning: For scaled endpoints input, cell" + cellIdx + " SATNUM = " + satnumIdx + ", SOWCR > SWU";
scaled_messages_.push_back(msg);
OpmLog::warning(msg);
counter_.warning += 1;
const std::string msg = "For scaled endpoints input, cell" + cellIdx + " SATNUM = " + satnumIdx + ", SOWCR > SWU";
OpmLog::warning(tag, msg);
}
}
if (fluidSystem_ != FluidSystem::OilWater) {
if (scaledEpsInfo_[c].Sgl > scaledEpsInfo_[c].Sgcr) {
const std::string msg = "-- Warning: For scaled endpoints input, cell" + cellIdx + " SATNUM = " + satnumIdx + ", SGL > SGCR";
scaled_messages_.push_back(msg);
OpmLog::warning(msg);
counter_.warning += 1;
const std::string msg = "For scaled endpoints input, cell" + cellIdx + " SATNUM = " + satnumIdx + ", SGL > SGCR";
OpmLog::warning(tag, msg);
}
}
if (fluidSystem_ != FluidSystem::BlackOil) {
if (scaledEpsInfo_[c].Sgcr > scaledEpsInfo_[c].Sogcr) {
const std::string msg = "-- Warning: For scaled endpoints input, cell" + cellIdx + " SATNUM = " + satnumIdx + ", SGCR > SOGCR";
scaled_messages_.push_back(msg);
OpmLog::warning(msg);
counter_.warning += 1;
const std::string msg = "For scaled endpoints input, cell" + cellIdx + " SATNUM = " + satnumIdx + ", SGCR > SOGCR";
OpmLog::warning(tag, msg);
}
if (scaledEpsInfo_[c].Sogcr > scaledEpsInfo_[c].Sgu) {
const std::string msg = "-- Warning: For scaled endpoints input, cell" + cellIdx + " SATNUM = " + satnumIdx + ", SOGCR > SGU";
scaled_messages_.push_back(msg);
OpmLog::warning(msg);
counter_.warning += 1;
const std::string msg = "For scaled endpoints input, cell" + cellIdx + " SATNUM = " + satnumIdx + ", SOGCR > SGU";
OpmLog::warning(tag, msg);
}
}
}

43
opm/core/props/satfunc/SaturationPropsFromDeck.cpp
View File

@ -136,9 +136,11 @@ namespace Opm
double* pc,
double* dpcds) const
{
assert(cells != 0);
assert(cells != 0);
assert(phaseUsage_.phase_used[BlackoilPhases::Liquid]);
const int np = numPhases();
if (dpcds) {
ExplicitArraysSatDerivativesFluidState fluidState(phaseUsage_);
typedef ExplicitArraysSatDerivativesFluidState::Evaluation Evaluation;
@ -151,12 +153,26 @@ namespace Opm
MaterialLaw::capillaryPressures(capillaryPressures, params, fluidState);
// copy the values calculated using opm-material to the target arrays
for (int pcPhaseIdx = 0; pcPhaseIdx < np; ++pcPhaseIdx) {
double sign = (pcPhaseIdx == BlackoilPhases::Aqua)? -1.0 : 1.0;
pc[np*i + pcPhaseIdx] = sign*capillaryPressures[pcPhaseIdx].value;
for (int canonicalPhaseIdx = 0; canonicalPhaseIdx < BlackoilPhases::MaxNumPhases; ++canonicalPhaseIdx) {
// skip unused phases
if ( ! phaseUsage_.phase_used[canonicalPhaseIdx]) {
continue;
}
const int pcPhaseIdx = phaseUsage_.phase_pos[canonicalPhaseIdx];
for (int satPhaseIdx = 0; satPhaseIdx < np; ++satPhaseIdx)
dpcds[np*np*i + satPhaseIdx*np + pcPhaseIdx] = sign*capillaryPressures[pcPhaseIdx].derivatives[satPhaseIdx];
const double sign = (canonicalPhaseIdx == BlackoilPhases::Aqua)? -1.0 : 1.0;
// in opm-material the wetting phase is the reference phase
// for two-phase problems i.e water for oil-water system,
// but for flow it is always oil. Add oil (liquid) capillary pressure value
// to shift the reference phase to oil
pc[np*i + pcPhaseIdx] = capillaryPressures[BlackoilPhases::Liquid].value + sign * capillaryPressures[canonicalPhaseIdx].value;
for (int canonicalSatPhaseIdx = 0; canonicalSatPhaseIdx < BlackoilPhases::MaxNumPhases; ++canonicalSatPhaseIdx) {
if ( ! phaseUsage_.phase_used[canonicalSatPhaseIdx])
continue;
const int satPhaseIdx = phaseUsage_.phase_pos[canonicalSatPhaseIdx];
dpcds[np*np*i + satPhaseIdx*np + pcPhaseIdx] = capillaryPressures[BlackoilPhases::Liquid].derivatives[canonicalSatPhaseIdx] + sign * capillaryPressures[canonicalPhaseIdx].derivatives[canonicalSatPhaseIdx];
}
}
}
} else {
@ -170,9 +186,18 @@ namespace Opm
MaterialLaw::capillaryPressures(capillaryPressures, params, fluidState);
// copy the values calculated using opm-material to the target arrays
for (int pcPhaseIdx = 0; pcPhaseIdx < np; ++pcPhaseIdx) {
double sign = (pcPhaseIdx == BlackoilPhases::Aqua)? -1.0 : 1.0;
pc[np*i + pcPhaseIdx] = sign*capillaryPressures[pcPhaseIdx];
for (int canonicalPhaseIdx = 0; canonicalPhaseIdx < BlackoilPhases::MaxNumPhases; ++canonicalPhaseIdx) {
// skip unused phases
if ( ! phaseUsage_.phase_used[canonicalPhaseIdx])
continue;
const int pcPhaseIdx = phaseUsage_.phase_pos[canonicalPhaseIdx];
double sign = (canonicalPhaseIdx == BlackoilPhases::Aqua)? -1.0 : 1.0;
// in opm-material the wetting phase is the reference phase
// for two-phase problems i.e water for oil-water system,
// but for flow it is always oil. Add oil (liquid) capillary pressure value
// to shift the reference phase to oil
pc[np*i + pcPhaseIdx] = capillaryPressures[BlackoilPhases::Liquid] + sign * capillaryPressures[canonicalPhaseIdx];
}
}
}

7
opm/core/simulator/BlackoilState.cpp
View File

@ -1,6 +1,8 @@
#include "BlackoilState.hpp"
#include <opm/common/util/numeric/cmp.hpp>
#include <opm/core/props/BlackoilPropertiesInterface.hpp>
#include <opm/core/simulator/WellState.hpp>
#include <opm/output/Wells.hpp>
using namespace Opm;
@ -21,15 +23,18 @@ BlackoilState::BlackoilState( size_t num_cells , size_t num_faces , size_t num_p
}
BlackoilState::BlackoilState( const BlackoilState& other )
: SimulationDataContainer(other)
: SimulationDataContainer(other),
hydrocarbonstate_(other.hydroCarbonState())
{
setBlackoilStateReferencePointers();
}
BlackoilState& BlackoilState::operator=( const BlackoilState& other )
{
SimulationDataContainer::operator=(other);
setBlackoilStateReferencePointers();
hydrocarbonstate_ = other.hydroCarbonState();
return *this;
}

11
opm/core/simulator/BlackoilState.hpp
View File

@ -30,6 +30,12 @@
namespace Opm
{
enum HydroCarbonState {
GasOnly = 0,
GasAndOil = 1,
OilOnly = 2
};
/// Simulator state for a blackoil simulator.
class BlackoilState : public SimulationDataContainer
{
@ -62,10 +68,12 @@ namespace Opm
std::vector<double>& surfacevol () { return *surfacevol_ref_; }
std::vector<double>& gasoilratio () { return *gasoilratio_ref_; }
std::vector<double>& rv () { return *rv_ref_; }
std::vector<HydroCarbonState>& hydroCarbonState() { return hydrocarbonstate_; }
const std::vector<double>& surfacevol () const { return *surfacevol_ref_; }
const std::vector<double>& gasoilratio () const { return *gasoilratio_ref_; }
const std::vector<double>& rv () const { return *rv_ref_; }
const std::vector<HydroCarbonState>& hydroCarbonState() const { return hydrocarbonstate_; }
private:
void setBlackoilStateReferencePointers();
@ -73,6 +81,9 @@ namespace Opm
std::vector<double>* gasoilratio_ref_;
std::vector<double>* rv_ref_;
// A vector storing the hydro carbon state.
std::vector<HydroCarbonState> hydrocarbonstate_;
};
} // namespace Opm

6
opm/core/simulator/ExplicitArraysSatDerivativesFluidState.hpp
View File

@ -20,8 +20,8 @@
#ifndef OPM_EXPLICIT_ARRAYS_SAT_DERIVATIVES_FLUID_STATE_HEADER_INCLUDED
#define OPM_EXPLICIT_ARRAYS_SAT_DERIVATIVES_FLUID_STATE_HEADER_INCLUDED
#include <opm/material/localad/Evaluation.hpp>
#include <opm/material/localad/Math.hpp>
#include <opm/material/densead/Evaluation.hpp>
#include <opm/material/densead/Math.hpp>
#include <opm/core/props/BlackoilPhases.hpp>
@ -44,7 +44,7 @@ public:
enum { numPhases = BlackoilPhases::MaxNumPhases };
enum { numComponents = 3 };
typedef Opm::LocalAd::Evaluation<double, numPhases> Evaluation;
typedef Opm::DenseAd::Evaluation<double, numPhases> Evaluation;
typedef Evaluation Scalar;
ExplicitArraysSatDerivativesFluidState(const PhaseUsage& phaseUsage)

2
opm/core/simulator/SimulatorReport.cpp
View File

@ -27,6 +27,7 @@ namespace Opm
: pressure_time(0.0),
transport_time(0.0),
total_time(0.0),
total_linearizations( 0 ),
total_newton_iterations( 0 ),
total_linear_iterations( 0 ),
verbose_(verbose)
@ -61,6 +62,7 @@ namespace Opm
{
os << "Total time taken (seconds): " << total_time
<< "\nSolver time (seconds): " << pressure_time
<< "\nOverall Linearizations: " << total_linearizations
<< "\nOverall Newton Iterations: " << total_newton_iterations
<< "\nOverall Linear Iterations: " << total_linear_iterations
<< std::endl;

1
opm/core/simulator/SimulatorReport.hpp
View File

@ -32,6 +32,7 @@ namespace Opm
double transport_time;
double total_time;
unsigned int total_linearizations;
unsigned int total_newton_iterations;
unsigned int total_linear_iterations;

47
opm/core/simulator/WellState.hpp
View File

@ -22,8 +22,11 @@
#include <opm/core/wells.h>
#include <opm/core/well_controls.h>
#include <opm/output/Wells.hpp>
#include <array>
#include <map>
#include <memory>
#include <string>
#include <vector>
#include <cassert>
@ -49,6 +52,7 @@ namespace Opm
{
// clear old name mapping
wellMap_.clear();
wells_.reset( clone_wells( wells ) );
if (wells) {
const int nw = wells->number_of_wells;
@ -208,6 +212,43 @@ namespace Opm
return wellRates().size() / numWells();
}
virtual data::Wells report() const
{
return { { /* WellState offers no completion data, so that has to be added later */ },
this->bhp(),
this->temperature(),
this->wellRates(),
this->perfPress(),
this->perfRates() };
}
virtual ~WellState() {}
WellState() = default;
WellState( const WellState& rhs ) :
bhp_( rhs.bhp_ ),
thp_( rhs.thp_ ),
temperature_( rhs.temperature_ ),
wellrates_( rhs.wellrates_ ),
perfrates_( rhs.perfrates_ ),
perfpress_( rhs.perfpress_ ),
wellMap_( rhs.wellMap_ ),
wells_( clone_wells( rhs.wells_.get() ) )
{}
WellState& operator=( const WellState& rhs ) {
this->bhp_ = rhs.bhp_;
this->thp_ = rhs.thp_;
this->temperature_ = rhs.temperature_;
this->wellrates_ = rhs.wellrates_;
this->perfrates_ = rhs.perfrates_;
this->perfpress_ = rhs.perfpress_;
this->wellMap_ = rhs.wellMap_;
this->wells_.reset( clone_wells( rhs.wells_.get() ) );
return *this;
}
private:
std::vector<double> bhp_;
std::vector<double> thp_;
@ -217,6 +258,12 @@ namespace Opm
std::vector<double> perfpress_;
WellMapType wellMap_;
protected:
struct wdel {
void operator()( Wells* w ) { destroy_wells( w ); }
};
std::unique_ptr< Wells, wdel > wells_;
};
} // namespace Opm

19
opm/core/simulator/initStateEquil.hpp
View File

@ -37,6 +37,7 @@
#include <opm/parser/eclipse/EclipseState/Tables/TableManager.hpp>
#include <opm/parser/eclipse/EclipseState/Tables/RsvdTable.hpp>
#include <opm/parser/eclipse/EclipseState/Tables/RvvdTable.hpp>
#include <opm/common/OpmLog/OpmLog.hpp>
#include <array>
#include <cassert>
@ -199,7 +200,7 @@ namespace Opm
std::vector<EquilRecord>
getEquil(const Opm::EclipseState& state)
{
const auto& init = *state.getInitConfig();
const auto& init = state.getInitConfig();
if( !init.hasEquil() ) {
OPM_THROW(std::domain_error, "Deck does not provide equilibration data.");
@ -264,6 +265,11 @@ namespace Opm
if (deck->hasKeyword("DISGAS")) {
const TableContainer& rsvdTables = tables.getRsvdTables();
for (size_t i = 0; i < rec.size(); ++i) {
if (eqlmap.cells(i).empty())
{
rs_func_.push_back(std::shared_ptr<Miscibility::RsVD>());
continue;
}
const int cell = *(eqlmap.cells(i).begin());
if (!rec[i].liveOilInitConstantRs()) {
if (rsvdTables.size() <= 0 ) {
@ -297,6 +303,11 @@ namespace Opm
if (deck->hasKeyword("VAPOIL")) {
const TableContainer& rvvdTables = tables.getRvvdTables();
for (size_t i = 0; i < rec.size(); ++i) {
if (eqlmap.cells(i).empty())
{
rv_func_.push_back(std::shared_ptr<Miscibility::RvVD>());
continue;
}
const int cell = *(eqlmap.cells(i).begin());
if (!rec[i].wetGasInitConstantRv()) {
if (rvvdTables.size() <= 0) {
@ -381,6 +392,12 @@ namespace Opm
{
for (const auto& r : reg.activeRegions()) {
const auto& cells = reg.cells(r);
if (cells.empty())
{
OpmLog::warning("Equilibration region " + std::to_string(r + 1)
+ " has no active cells");
continue;
}
const int repcell = *cells.begin();
const RhoCalc calc(props, repcell);

50
opm/core/simulator/initStateEquil_impl.hpp
View File

@ -767,28 +767,34 @@ namespace Opm
double sat[BlackoilPhases::MaxNumPhases];
double threshold_sat = 1.0e-6;
sat[waterpos] = smax[waterpos];
sat[gaspos] = smax[gaspos];
sat[oilpos] = 1.0 - sat[waterpos] - sat[gaspos];
if (sw > smax[waterpos]-threshold_sat ) {
sat[waterpos] = smax[waterpos];
props.capPress(1, sat, &cell, pc, 0);
phase_pressures[oilpos][local_index] = phase_pressures[waterpos][local_index] + pc[waterpos];
} else if (sg > smax[gaspos]-threshold_sat) {
sat[gaspos] = smax[gaspos];
props.capPress(1, sat, &cell, pc, 0);
phase_pressures[oilpos][local_index] = phase_pressures[gaspos][local_index] - pc[gaspos];
}
if (sg < smin[gaspos]+threshold_sat) {
sat[gaspos] = smin[gaspos];
props.capPress(1, sat, &cell, pc, 0);
phase_pressures[gaspos][local_index] = phase_pressures[oilpos][local_index] + pc[gaspos];
}
if (sw < smin[waterpos]+threshold_sat) {
sat[waterpos] = smin[waterpos];
props.capPress(1, sat, &cell, pc, 0);
phase_pressures[waterpos][local_index] = phase_pressures[oilpos][local_index] - pc[waterpos];
}
sat[oilpos] = 1.0;
if (water) {
sat[waterpos] = smax[waterpos];
sat[oilpos] -= sat[waterpos];
}
if (gas) {
sat[gaspos] = smax[gaspos];
sat[oilpos] -= sat[gaspos];
}
if (water && sw > smax[waterpos]-threshold_sat ) {
sat[waterpos] = smax[waterpos];
props.capPress(1, sat, &cell, pc, 0);
phase_pressures[oilpos][local_index] = phase_pressures[waterpos][local_index] + pc[waterpos];
} else if (gas && sg > smax[gaspos]-threshold_sat) {
sat[gaspos] = smax[gaspos];
props.capPress(1, sat, &cell, pc, 0);
phase_pressures[oilpos][local_index] = phase_pressures[gaspos][local_index] - pc[gaspos];
}
if (gas && sg < smin[gaspos]+threshold_sat) {
sat[gaspos] = smin[gaspos];
props.capPress(1, sat, &cell, pc, 0);
phase_pressures[gaspos][local_index] = phase_pressures[oilpos][local_index] + pc[gaspos];
}
if (water && sw < smin[waterpos]+threshold_sat) {
sat[waterpos] = smin[waterpos];
props.capPress(1, sat, &cell, pc, 0);
phase_pressures[waterpos][local_index] = phase_pressures[oilpos][local_index] - pc[waterpos];
}
}
return phase_saturations;
}

15
opm/core/simulator/initState_impl.hpp
View File

@ -741,7 +741,7 @@ namespace Opm
}
/// Initialize surface volume from pressure and saturation by z = As.
/// Here saturation is used as an intial guess for z in the
/// Here saturation is used as an initial guess for z in the
/// computation of A.
template <class Props, class State>
void initBlackoilSurfvol(const UnstructuredGrid& grid,
@ -800,10 +800,6 @@ namespace Opm
const Props& props,
State& state)
{
if (props.numPhases() != 3) {
OPM_THROW(std::runtime_error, "initBlackoilSurfvol() is only supported in three-phase simulations.");
}
const std::vector<double>& rs = state.gasoilratio();
const std::vector<double>& rv = state.rv();
@ -828,15 +824,6 @@ namespace Opm
std::vector<double> capPressures(number_of_cells*np);
props.capPress(number_of_cells,&state.saturation()[0],&allcells[0],&capPressures[0],NULL);
std::vector<double> Pw(number_of_cells);
std::vector<double> Pg(number_of_cells);
for (int c = 0; c < number_of_cells; ++c){
Pw[c] = state.pressure()[c] + capPressures[c*np + BlackoilPhases::Aqua];
Pg[c] = state.pressure()[c] + capPressures[c*np + BlackoilPhases::Vapour];
}
double z_tmp;
// Water phase

44
opm/core/utility/initHydroCarbonState.hpp Normal file
View File

@ -0,0 +1,44 @@
#ifndef INITHYDROCARBONSTATE_HPP
#define INITHYDROCARBONSTATE_HPP
#include "opm/core/simulator/BlackoilState.hpp"
namespace Opm
{
void initHydroCarbonState(BlackoilState& state, const PhaseUsage& pu, const int num_cells, const bool has_disgas, const bool has_vapoil) {
enum { Oil = BlackoilPhases::Liquid, Gas = BlackoilPhases::Vapour, Water = BlackoilPhases::Aqua };
// hydrocarbonstate is only used when gas and oil is present
assert(pu.phase_used[Oil]);
std::vector<HydroCarbonState>& hydroCarbonState = state.hydroCarbonState();
hydroCarbonState.resize(num_cells);
if (!pu.phase_used[Gas]) {
// hydroCarbonState should only be used when oil and gas is present. Return OilOnly to avoid potential trouble.
std::fill(hydroCarbonState.begin(), hydroCarbonState.end(), HydroCarbonState::OilOnly);
return;
}
const int np = pu.num_phases;
std::fill(hydroCarbonState.begin(), hydroCarbonState.end(), HydroCarbonState::GasAndOil);
// set hydrocarbon state
const double epsilon = std::sqrt(std::numeric_limits<double>::epsilon());
const std::vector<double>& saturation = state.saturation();
for (int c = 0; c < num_cells; ++c) {
if (pu.phase_used[Water]) {
if ( saturation[c*np + pu.phase_pos[ Water ]] > (1.0 - epsilon)) {
continue; // cases (almost) filled with water is treated as GasAndOil case;
}
}
if ( saturation[c*np + pu.phase_pos[ Gas ]] == 0.0 && has_disgas) {
hydroCarbonState[c] = HydroCarbonState::OilOnly;
continue;
}
if ( saturation[c*np + pu.phase_pos[ Oil ]] == 0.0 && has_vapoil) {
hydroCarbonState[c] = HydroCarbonState::GasOnly;
}
}
}
} // namespace Opm
#endif // INITHYDROCARBONSTATE_HPP

92
opm/core/wells/DynamicListEconLimited.hpp Normal file
View File

@ -0,0 +1,92 @@
/*
Copyright 2016 SINTEF ICT, Applied Mathematics.
This file is part of the Open Porous Media project (OPM).
OPM is free software: you can redistribute it and/or modify
it under the terms of the GNU General Public License as published by
the Free Software Foundation, either version 3 of the License, or
(at your option) any later version.
OPM is distributed in the hope that it will be useful,
but WITHOUT ANY WARRANTY; without even the implied warranty of
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
GNU General Public License for more details.
You should have received a copy of the GNU General Public License
along with OPM. If not, see <http://www.gnu.org/licenses/>.
*/
#ifndef OPM_DYNAMICLISTECONLIMITED_HPP
#define OPM_DYNAMICLISTECONLIMITED_HPP
#include <vector>
#include <string>
#include <map>
#include <cassert>
namespace Opm
{
/// to handle the wells and connections violating economic limits.
class DynamicListEconLimited
{
public:
DynamicListEconLimited() {
}
bool wellShutEconLimited(const std::string& well_name) const {
return std::find(m_shut_wells.begin(), m_shut_wells.end(), well_name) != m_shut_wells.end();
}
void addShutWell(const std::string& well_name) {
assert( !wellShutEconLimited(well_name) );
assert( !wellStoppedEconLimited(well_name) );
m_shut_wells.push_back(well_name);
}
bool wellStoppedEconLimited(const std::string& well_name) const {
return std::find(m_stopped_wells.begin(), m_stopped_wells.end(), well_name) != m_stopped_wells.end();
}
void addStoppedWell(const std::string& well_name) {
assert( !wellShutEconLimited(well_name) );
assert( !wellStoppedEconLimited(well_name) );
m_stopped_wells.push_back(well_name);
}
// TODO: maybe completion better here
bool anyConnectionClosedForWell(const std::string& well_name) const {
return (m_cells_closed_connections.find(well_name) != m_cells_closed_connections.end());
}
const std::vector<int>& getClosedConnectionsForWell(const std::string& well_name) const {
return (m_cells_closed_connections.find(well_name)->second);
}
void addClosedConnectionsForWell(const std::string& well_name,
const int cell_closed_connection) {
if (!anyConnectionClosedForWell(well_name)) {
// first time adding a connection for the well
std::vector<int> vector_cells = {cell_closed_connection};
m_cells_closed_connections[well_name] = vector_cells;
} else {
std::vector<int>& closed_connections = m_cells_closed_connections.find(well_name)->second;
closed_connections.push_back(cell_closed_connection);
}
}
private:
std::vector <std::string> m_shut_wells;
std::vector <std::string> m_stopped_wells;
// using grid cell number to indicate the location of the connections
std::map<std::string, std::vector<int>> m_cells_closed_connections;
};
} // namespace Opm
#endif /* OPM_DYNAMICLISTECONLIMITED_HPP */

61
opm/core/wells/InjectionSpecification.cpp
View File

@ -1,4 +1,26 @@
/*
Copyright 2012 Sintef.
Copyright 2016 Statoil.
This file is part of the Open Porous Media project (OPM).
OPM is free software: you can redistribute it and/or modify it under the terms
of the GNU General Public License as published by the Free Software
Foundation, either version 3 of the License, or (at your option) any later
version.
OPM is distributed in the hope that it will be useful, but WITHOUT ANY
WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR
A PARTICULAR PURPOSE. See the GNU General Public License for more details.
You should have received a copy of the GNU General Public License along with
OPM. If not, see <http://www.gnu.org/licenses/>.
*/
#include <stdexcept>
#include "config.h"
#include <opm/common/ErrorMacros.hpp>
#include <opm/core/wells/InjectionSpecification.hpp>
namespace Opm
@ -18,4 +40,43 @@ namespace Opm
}
std::string
InjectionSpecification::toString(const ControlMode& mode)
{
switch(mode) {
case ControlMode::NONE: return "NONE";
case ControlMode::RATE: return "RATE";
case ControlMode::RESV: return "RESV";
case ControlMode::BHP : return "BHP" ;
case ControlMode::THP : return "THP" ;
case ControlMode::REIN: return "REIN";
case ControlMode::VREP: return "VREP";
case ControlMode::GRUP: return "GRUP";
case ControlMode::FLD : return "FLD" ;
}
OPM_THROW(std::domain_error, "Unknown control mode " << mode << " encountered in injection specification");
}
std::string
InjectionSpecification::toString(const InjectorType& type)
{
switch(type) {
case InjectorType::WATER: return "WATER";
case InjectorType::OIL : return "OIL" ;
case InjectorType::GAS : return "GAS" ;
}
OPM_THROW(std::domain_error, "Unknown injector type " << type << " encountered in injection specification");
}
std::string
InjectionSpecification::toString(const GuideRateType& type)
{
switch(type) {
case GuideRateType::RAT : return "RAT" ;
case GuideRateType::NONE_GRT: return "NONE_GRT";
}
OPM_THROW(std::domain_error, "Unknown guide rate type " << type << " encountered in injection specification");
}
} // namespace Opm

5
opm/core/wells/InjectionSpecification.hpp
View File

@ -2,6 +2,7 @@
#define OPM_INJECTORSPECIFICATION_HPP
#include <opm/core/wells.h>
#include <string>
namespace Opm
{
@ -25,7 +26,9 @@ namespace Opm
};
InjectionSpecification();
static std::string toString(const ControlMode& mode);
static std::string toString(const InjectorType& type);
static std::string toString(const GuideRateType& type);
InjectorType injector_type_;
ControlMode control_mode_;
double surface_flow_max_rate_;

69
opm/core/wells/ProductionSpecification.cpp
View File

@ -1,4 +1,26 @@
/*
Copyright 2012 Sintef.
Copyright 2016 Statoil.
This file is part of the Open Porous Media project (OPM).
OPM is free software: you can redistribute it and/or modify it under the terms
of the GNU General Public License as published by the Free Software
Foundation, either version 3 of the License, or (at your option) any later
version.
OPM is distributed in the hope that it will be useful, but WITHOUT ANY
WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR
A PARTICULAR PURPOSE. See the GNU General Public License for more details.
You should have received a copy of the GNU General Public License along with
OPM. If not, see <http://www.gnu.org/licenses/>.
*/
#include <stdexcept>
#include "config.h"
#include <opm/common/ErrorMacros.hpp>
#include <opm/core/wells/ProductionSpecification.hpp>
namespace Opm
@ -19,4 +41,51 @@ namespace Opm
{
}
std::string
ProductionSpecification::toString(const ControlMode& mode)
{
switch(mode) {
case ControlMode::NONE: return "NONE";
case ControlMode::ORAT: return "ORAT";
case ControlMode::WRAT: return "WRAT";
case ControlMode::GRAT: return "GRAT";
case ControlMode::LRAT: return "LRAT";
case ControlMode::CRAT: return "CRAT";
case ControlMode::RESV: return "RESV";
case ControlMode::PRBL: return "RPBL";
case ControlMode::BHP : return "BHP" ;
case ControlMode::THP : return "THP" ;
case ControlMode::GRUP: return "GRUP";
case ControlMode::FLD : return "FLD" ;
}
OPM_THROW(std::domain_error, "Unknown control mode " << mode << " encountered in production specification");
}
std::string
ProductionSpecification::toString(const Procedure& type)
{
switch(type) {
case Procedure::NONE_P: return "NONE_P";
case Procedure::RATE : return "RATE" ;
case Procedure::WELL : return "WELL" ;
}
OPM_THROW(std::domain_error, "Unknown procedure " << type << " encountered in production specification");
}
std::string
ProductionSpecification::toString(const GuideRateType& type)
{
switch(type) {
case GuideRateType::OIL : return "OIL" ;
case GuideRateType::GAS : return "GAS" ;
case GuideRateType::WATER : return "WATER" ;
case GuideRateType::NONE_GRT: return "NONE_GRT";
}
OPM_THROW(std::domain_error, "Unknown guide rate type " << type << " encountered in production specification");
}
}

4
opm/core/wells/ProductionSpecification.hpp
View File

@ -2,6 +2,7 @@
#define OPM_PRODUCTIONSPECIFICATION_HPP
#include <opm/core/wells.h>
#include <string>
namespace Opm
{
@ -25,6 +26,9 @@ namespace Opm
};
ProductionSpecification();
static std::string toString(const ControlMode& mode);
static std::string toString(const Procedure& type);
static std::string toString(const GuideRateType& type);
ControlMode control_mode_;
Procedure procedure_;

6
opm/core/wells/WellCollection.cpp
View File

@ -29,7 +29,7 @@
namespace Opm
{
void WellCollection::addField(GroupConstPtr fieldGroup, size_t timeStep, const PhaseUsage& phaseUsage) {
void WellCollection::addField(const Group* fieldGroup, size_t timeStep, const PhaseUsage& phaseUsage) {
WellsGroupInterface* fieldNode = findNode(fieldGroup->name());
if (fieldNode) {
OPM_THROW(std::runtime_error, "Trying to add FIELD node, but this already exists. Can only have one FIELD node.");
@ -38,7 +38,7 @@ namespace Opm
roots_.push_back(createGroupWellsGroup(fieldGroup, timeStep, phaseUsage));
}
void WellCollection::addGroup(GroupConstPtr groupChild, std::string parent_name,
void WellCollection::addGroup(const Group* groupChild, std::string parent_name,
size_t timeStep, const PhaseUsage& phaseUsage) {
WellsGroupInterface* parent = findNode(parent_name);
if (!parent) {
@ -60,7 +60,7 @@ namespace Opm
child->setParent(parent);
}
void WellCollection::addWell(WellConstPtr wellChild, size_t timeStep, const PhaseUsage& phaseUsage) {
void WellCollection::addWell(const Well* wellChild, size_t timeStep, const PhaseUsage& phaseUsage) {
if (wellChild->getStatus(timeStep) == WellCommon::SHUT) {
//SHUT wells are not added to the well collection
return;

6
opm/core/wells/WellCollection.hpp
View File

@ -36,11 +36,11 @@ namespace Opm
{
public:
void addField(GroupConstPtr fieldGroup, size_t timeStep, const PhaseUsage& phaseUsage);
void addField(const Group* fieldGroup, size_t timeStep, const PhaseUsage& phaseUsage);
void addWell(WellConstPtr wellChild, size_t timeStep, const PhaseUsage& phaseUsage);
void addWell(const Well* wellChild, size_t timeStep, const PhaseUsage& phaseUsage);
void addGroup(GroupConstPtr groupChild, std::string parent_name,
void addGroup(const Group* groupChild, std::string parent_name,
size_t timeStep, const PhaseUsage& phaseUsage);
/// Adds the child to the collection

36
opm/core/wells/WellsGroup.cpp
View File

@ -345,9 +345,10 @@ namespace Opm
mode);
if (my_rate > target_rate) {
std::cout << "Group " << mode<<" target not met for group " << name() << std::endl;
std::cout << "target = " << target_rate << '\n'
<< "rate = " << my_rate << std::endl;
OpmLog::warning("Group " + InjectionSpecification::toString(mode)
+ " target not met for group " + name() + "\n"
+ "target = " + std::to_string(target_rate) + "\n"
+ "rate = " + std::to_string(my_rate));
applyInjGroupControl(mode, target_rate, true);
injSpec().control_mode_ = mode;
return false;
@ -377,9 +378,10 @@ namespace Opm
child_phases_summed.surf_prod_rates,
mode);
if (std::fabs(my_rate) > target_rate) {
std::cout << "Group" << mode << " target not met for group " << name() << std::endl;
std::cout << "target = " << target_rate << '\n'
<< "rate = " << my_rate << std::endl;
OpmLog::warning("Group" + ProductionSpecification::toString(mode)
+ " target not met for group " + name() + "\n"
+ "target = " + std::to_string(target_rate) + '\n'
+ "rate = " + std::to_string(my_rate));
production_violated = true;
production_mode_violated = mode;
break;
@ -677,9 +679,9 @@ namespace Opm
ctrl_violated = is_producer ? (my_target_bhp > my_well_bhp)
: (my_target_bhp < my_well_bhp);
if (ctrl_violated) {
std::cout << "BHP limit violated for well " << name() << ":\n";
std::cout << "BHP limit = " << my_target_bhp << std::endl;
std::cout << "BHP = " << my_well_bhp << std::endl;
OpmLog::info("BHP limit violated for well " + name() + ":\n"
+ "BHP limit = " + std::to_string(my_target_bhp)
+ "BHP = " + std::to_string(my_well_bhp));
}
break;
}
@ -698,9 +700,9 @@ namespace Opm
const double my_rate_target = well_controls_iget_target(ctrls , ctrl_index);
ctrl_violated = std::fabs(my_rate) - std::fabs(my_rate_target)> std::max(std::abs(my_rate), std::abs(my_rate_target))*1e-6;
if (ctrl_violated) {
std::cout << "RESERVOIR_RATE limit violated for well " << name() << ":\n";
std::cout << "rate limit = " << my_rate_target << std::endl;
std::cout << "rate = " << my_rate << std::endl;
OpmLog::info("RESERVOIR_RATE limit violated for well " + name() + ":\n"
+ "rate limit = " + std::to_string(my_rate_target)
+ "rate = " + std::to_string(my_rate));
}
break;
}
@ -714,9 +716,9 @@ namespace Opm
const double my_rate_target = well_controls_iget_target(ctrls , ctrl_index);
ctrl_violated = std::fabs(my_rate) > std::fabs(my_rate_target);
if (ctrl_violated) {
std::cout << "SURFACE_RATE limit violated for well " << name() << ":\n";
std::cout << "rate limit = " << my_rate_target << std::endl;
std::cout << "rate = " << my_rate << std::endl;
OpmLog::info("SURFACE_RATE limit violated for well " + name() + ":\n"
+ "rate limit = " + std::to_string(my_rate_target)
+ "rate = " + std::to_string(my_rate));
}
break;
}
@ -1062,7 +1064,7 @@ namespace Opm
}
} // anonymous namespace
std::shared_ptr<WellsGroupInterface> createGroupWellsGroup(GroupConstPtr group, size_t timeStep, const PhaseUsage& phase_usage )
std::shared_ptr<WellsGroupInterface> createGroupWellsGroup(const Group* group, size_t timeStep, const PhaseUsage& phase_usage )
{
InjectionSpecification injection_specification;
ProductionSpecification production_specification;
@ -1097,7 +1099,7 @@ namespace Opm
'CMODE_UNDEFINED' - we do not carry that over the specification
objects here.
*/
std::shared_ptr<WellsGroupInterface> createWellWellsGroup(WellConstPtr well, size_t timeStep, const PhaseUsage& phase_usage )
std::shared_ptr<WellsGroupInterface> createWellWellsGroup(const Well* well, size_t timeStep, const PhaseUsage& phase_usage )
{
InjectionSpecification injection_specification;
ProductionSpecification production_specification;

4
opm/core/wells/WellsGroup.hpp
View File

@ -406,14 +406,14 @@ namespace Opm
/// \param[in] well the Well to construct object for
/// \param[in] timeStep the time step in question
/// \param[in] the phase usage
std::shared_ptr<WellsGroupInterface> createWellWellsGroup(WellConstPtr well, size_t timeStep,
std::shared_ptr<WellsGroupInterface> createWellWellsGroup(const Well* well, size_t timeStep,
const PhaseUsage& phase_usage );
/// Creates the WellsGroupInterface for the given Group
/// \param[in] group the Group to construct object for
/// \param[in] timeStep the time step in question
/// \param[in] the phase usage
std::shared_ptr<WellsGroupInterface> createGroupWellsGroup(GroupConstPtr group, size_t timeStep,
std::shared_ptr<WellsGroupInterface> createGroupWellsGroup(const Group* group, size_t timeStep,
const PhaseUsage& phase_usage );
}
#endif /* OPM_WELLSGROUP_HPP */

31
opm/core/wells/WellsManager.cpp
View File

@ -333,11 +333,14 @@ namespace Opm
: w_(0), is_parallel_run_(false)
{
std::vector<double> dummy_well_potentials;
// TODO: not sure about the usage of this WellsManager constructor
// TODO: not sure whether this is the correct thing to do here.
DynamicListEconLimited dummy_list_econ_limited;
init(eclipseState, timeStep, UgGridHelpers::numCells(grid),
UgGridHelpers::globalCell(grid), UgGridHelpers::cartDims(grid),
UgGridHelpers::dimensions(grid),
UgGridHelpers::cell2Faces(grid), UgGridHelpers::beginFaceCentroids(grid),
permeability, dummy_well_potentials);
permeability, dummy_list_econ_limited, dummy_well_potentials);
}
@ -414,9 +417,10 @@ namespace Opm
void WellsManager::setupWellControls(std::vector<WellConstPtr>& wells, size_t timeStep,
void WellsManager::setupWellControls(std::vector< const Well* >& wells, size_t timeStep,
std::vector<std::string>& well_names, const PhaseUsage& phaseUsage,
const std::vector<int>& wells_on_proc) {
const std::vector<int>& wells_on_proc,
const DynamicListEconLimited& list_econ_limited) {
int well_index = 0;
auto well_on_proc = wells_on_proc.begin();
@ -427,18 +431,23 @@ namespace Opm
continue;
}
WellConstPtr well = (*wellIter);
if (well->getStatus(timeStep) == WellCommon::STOP) {
// STOPed wells are kept in the well list but marked as stopped.
well_controls_stop_well(w_->ctrls[well_index]);
}
const auto* well = (*wellIter);
if (well->getStatus(timeStep) == WellCommon::SHUT) {
//SHUT wells are not added to the well list
continue;
}
if (list_econ_limited.wellShutEconLimited(well->name())) {
continue;
}
if (well->getStatus(timeStep) == WellCommon::STOP || list_econ_limited.wellStoppedEconLimited(well->name())) {
// Stopped wells are kept in the well list but marked as stopped.
well_controls_stop_well(w_->ctrls[well_index]);
}
if (well->isInjector(timeStep)) {
const WellInjectionProperties& injectionProperties = well->getInjectionProperties(timeStep);
int ok = 1;
@ -728,12 +737,12 @@ namespace Opm
}
void WellsManager::setupGuideRates(std::vector<WellConstPtr>& wells, const size_t timeStep, std::vector<WellData>& well_data, std::map<std::string, int>& well_names_to_index,
void WellsManager::setupGuideRates(std::vector< const Well* >& wells, const size_t timeStep, std::vector<WellData>& well_data, std::map<std::string, int>& well_names_to_index,
const PhaseUsage& phaseUsage, const std::vector<double>& well_potentials)
{
const int np = phaseUsage.num_phases;
for (auto wellIter = wells.begin(); wellIter != wells.end(); ++wellIter ) {
WellConstPtr well = *wellIter;
const auto* well = *wellIter;
if (well->getStatus(timeStep) == WellCommon::SHUT) {
//SHUT wells does not need guide rates

15
opm/core/wells/WellsManager.hpp
View File

@ -25,6 +25,7 @@
#include <opm/core/wells/WellCollection.hpp>
#include <opm/core/wells/WellsGroup.hpp>
#include <opm/core/wells/DynamicListEconLimited.hpp>
#include <opm/parser/eclipse/EclipseState/Schedule/GroupTree.hpp>
#include <opm/core/utility/CompressedPropertyAccess.hpp>
@ -86,6 +87,7 @@ namespace Opm
const F2C& f2c,
FC begin_face_centroids,
const double* permeability,
const DynamicListEconLimited& list_econ_limited,
bool is_parallel_run=false,
const std::vector<double>& well_potentials={});
@ -155,17 +157,19 @@ namespace Opm
const C2F& cell_to_faces,
FC begin_face_centroids,
const double* permeability,
const DynamicListEconLimited& list_econ_limited,
const std::vector<double>& well_potentials);
// Disable copying and assignment.
WellsManager(const WellsManager& other);
WellsManager& operator=(const WellsManager& other);
static void setupCompressedToCartesian(const int* global_cell, int number_of_cells, std::map<int,int>& cartesian_to_compressed );
void setupWellControls(std::vector<WellConstPtr>& wells, size_t timeStep,
void setupWellControls(std::vector<const Well*>& wells, size_t timeStep,
std::vector<std::string>& well_names, const PhaseUsage& phaseUsage,
const std::vector<int>& wells_on_proc);
const std::vector<int>& wells_on_proc,
const DynamicListEconLimited& list_econ_limited);
template<class C2F, class FC, class NTG>
void createWellsFromSpecs( std::vector<WellConstPtr>& wells, size_t timeStep,
void createWellsFromSpecs( std::vector<const Well*>& wells, size_t timeStep,
const C2F& cell_to_faces,
const int* cart_dims,
FC begin_face_centroids,
@ -178,10 +182,11 @@ namespace Opm
const std::map<int,int>& cartesian_to_compressed,
const double* permeability,
const NTG& ntg,
std::vector<int>& wells_on_proc);
std::vector<int>& wells_on_proc,
const DynamicListEconLimited& list_econ_limited);
void addChildGroups(GroupTreeNodeConstPtr parentNode, std::shared_ptr< const Schedule > schedule, size_t timeStep, const PhaseUsage& phaseUsage);
void setupGuideRates(std::vector<WellConstPtr>& wells, const size_t timeStep, std::vector<WellData>& well_data, std::map<std::string, int>& well_names_to_index,
void setupGuideRates(std::vector<const Well*>& wells, const size_t timeStep, std::vector<WellData>& well_data, std::map<std::string, int>& well_names_to_index,
const PhaseUsage& phaseUsage, const std::vector<double>& well_potentials);
// Data
Wells* w_;

50
opm/core/wells/WellsManager_impl.hpp
View File

@ -2,6 +2,7 @@
#include <opm/core/grid/GridHelpers.hpp>
#include <opm/common/ErrorMacros.hpp>
#include <opm/common/OpmLog/OpmLog.hpp>
#include <opm/core/utility/compressedToCartesian.hpp>
#include <opm/parser/eclipse/EclipseState/Grid/EclipseGrid.hpp>
#include <opm/parser/eclipse/EclipseState/Schedule/CompletionSet.hpp>
@ -106,7 +107,7 @@ getCubeDim(const C2F& c2f,
namespace Opm
{
template<class C2F, class FC, class NTG>
void WellsManager::createWellsFromSpecs(std::vector<WellConstPtr>& wells, size_t timeStep,
void WellsManager::createWellsFromSpecs(std::vector<const Well*>& wells, size_t timeStep,
const C2F& c2f,
const int* cart_dims,
FC begin_face_centroids,
@ -119,7 +120,8 @@ void WellsManager::createWellsFromSpecs(std::vector<WellConstPtr>& wells, size_t
const std::map<int,int>& cartesian_to_compressed,
const double* permeability,
const NTG& ntg,
std::vector<int>& wells_on_proc)
std::vector<int>& wells_on_proc,
const DynamicListEconLimited& list_econ_limited)
{
if (dimensions != 3) {
OPM_THROW(std::domain_error,
@ -137,12 +139,21 @@ void WellsManager::createWellsFromSpecs(std::vector<WellConstPtr>& wells, size_t
// the index of the well according to the eclipse state
int well_index_on_proc = 0;
for (auto wellIter= wells.begin(); wellIter != wells.end(); ++wellIter) {
WellConstPtr well = (*wellIter);
const auto* well = (*wellIter);
if (well->getStatus(timeStep) == WellCommon::SHUT) {
continue;
}
if (list_econ_limited.wellShutEconLimited(well->name())) {
continue;
}
std::vector<int> cells_connection_closed;
if (list_econ_limited.anyConnectionClosedForWell(well->name())) {
cells_connection_closed = list_econ_limited.getClosedConnectionsForWell(well->name());
}
{ // COMPDAT handling
auto completionSet = well->getCompletions(timeStep);
// shut completions and open ones stored in this process will have 1 others 0.
@ -174,6 +185,16 @@ void WellsManager::createWellsFromSpecs(std::vector<WellConstPtr>& wells, size_t
else
{
int cell = cgit->second;
// check if the connection is closed due to economic limits
if (!cells_connection_closed.empty()) {
const bool connection_found = std::find(cells_connection_closed.begin(),
cells_connection_closed.end(), cell)
!= cells_connection_closed.end();
if (connection_found) {
continue;
}
}
PerfData pd;
pd.cell = cell;
{
@ -233,9 +254,9 @@ void WellsManager::createWellsFromSpecs(std::vector<WellConstPtr>& wells, size_t
// Check that the complete well is on this process
if ( sum_completions_on_proc < completionSet->size() )
{
std::cout<< "Well "<< well->name() << " does not seem to be"
<< "completely in the disjoint partition of "
<< "process. Therefore we deactivate it here." << std::endl;
OpmLog::warning("Well " + well->name() + " does not seem to be"
+ "completely in the disjoint partition of "
+ "process. Therefore we deactivate it here.");
// Mark well as not existent on this process
wells_on_proc[wellIter-wells.begin()] = 0;
wellperf_data[well_index_on_proc].clear();
@ -326,13 +347,14 @@ WellsManager(const Opm::EclipseStateConstPtr eclipseState,
const C2F& cell_to_faces,
FC begin_face_centroids,
const double* permeability,
const DynamicListEconLimited& list_econ_limited,
bool is_parallel_run,
const std::vector<double>& well_potentials)
: w_(0), is_parallel_run_(is_parallel_run)
{
init(eclipseState, timeStep, number_of_cells, global_cell,
cart_dims, dimensions,
cell_to_faces, begin_face_centroids, permeability, well_potentials);
cell_to_faces, begin_face_centroids, permeability, list_econ_limited, well_potentials);
}
/// Construct wells from deck.
@ -347,6 +369,7 @@ WellsManager::init(const Opm::EclipseStateConstPtr eclipseState,
const C2F& cell_to_faces,
FC begin_face_centroids,
const double* permeability,
const DynamicListEconLimited& list_econ_limited,
const std::vector<double>& well_potentials)
{
if (dimensions != 3) {
@ -378,8 +401,8 @@ WellsManager::init(const Opm::EclipseStateConstPtr eclipseState,
std::map<std::string, int> well_names_to_index;
auto schedule = eclipseState->getSchedule();
std::vector<WellConstPtr> wells = schedule->getWells(timeStep);
std::vector<int> wells_on_proc;
auto wells = schedule->getWells(timeStep);
std::vector<int> wells_on_proc;
well_names.reserve(wells.size());
well_data.reserve(wells.size());
@ -390,7 +413,7 @@ WellsManager::init(const Opm::EclipseStateConstPtr eclipseState,
DoubleArray ntg_glob(eclipseState, "NTG", 1.0);
NTGArray ntg(ntg_glob, global_cell);
EclipseGridConstPtr eclGrid = eclipseState->getEclipseGrid();
EclipseGridConstPtr eclGrid = eclipseState->getInputGrid();
// use cell thickness (dz) from eclGrid
// dz overwrites values calculated by WellDetails::getCubeDim
@ -409,16 +432,15 @@ WellsManager::init(const Opm::EclipseStateConstPtr eclipseState,
dz,
well_names, well_data, well_names_to_index,
pu, cartesian_to_compressed, permeability, ntg,
wells_on_proc);
wells_on_proc, list_econ_limited);
setupWellControls(wells, timeStep, well_names, pu, wells_on_proc);
setupWellControls(wells, timeStep, well_names, pu, wells_on_proc, list_econ_limited);
{
GroupTreeNodeConstPtr fieldNode =
schedule->getGroupTree(timeStep)->getNode("FIELD");
GroupConstPtr fieldGroup =
schedule->getGroup(fieldNode->name());
const auto* fieldGroup = schedule->getGroup(fieldNode->name());
well_collection_.addField(fieldGroup, timeStep, pu);
addChildGroups(fieldNode, schedule, timeStep, pu);

10
opm/core/wells/wells.c
View File

@ -545,6 +545,16 @@ bool
wells_equal(const struct Wells *W1, const struct Wells *W2 , bool verbose)
/* ---------------------------------------------------------------------- */
{
// Cater the case where W1 and W2 are the same (null) pointers.
if( W1 == W2 )
{
return true;
}
if( W1 == NULL || W2 == NULL)
{
return false;
}
bool are_equal = true;
are_equal = (W1->number_of_wells == W2->number_of_wells);
are_equal = are_equal && (W1->number_of_phases == W2->number_of_phases);

8
tests/norne_pvt.data
View File

@ -15,6 +15,14 @@ DIMENS
1 1 1 /
GRID
DX
1 /
DY
1 /
DZ
1 /
TOPS
1 /
PROPS

64
tests/test_blackoilstate.cpp
View File

@ -28,25 +28,36 @@ using namespace Opm;
using namespace std;
// ACTNUM 1 998*2 3
std::vector<int> get_testBlackoilStateActnum() {
std::vector<int> actnum(10 * 10 * 10, 2);
actnum.front() = 1;
actnum.back() = 3;
return actnum;
}
BOOST_AUTO_TEST_CASE(EqualsDifferentDeckReturnFalse) {
ParseContext parseContext;
const string filename1 = "testBlackoilState1.DATA";
const string filename2 = "testBlackoilState2.DATA";
Opm::ParserPtr parser(new Opm::Parser());
Opm::DeckConstPtr deck1(parser->parseFile(filename1, parseContext));
Opm::DeckConstPtr deck2(parser->parseFile(filename2, parseContext));
GridManager gridManager1(deck1);
const auto es1 = Opm::Parser::parse(filename1);
auto eg1 = es1.getInputGridCopy();
std::vector<int> actnum = get_testBlackoilStateActnum();
eg1->resetACTNUM(actnum.data());
const auto es2 = Opm::Parser::parse(filename2);
auto eg2 = es2.getInputGrid();
GridManager gridManager1(eg1);
const UnstructuredGrid& grid1 = *gridManager1.c_grid();
GridManager gridManager2(deck2);
GridManager gridManager2(eg2);
const UnstructuredGrid& grid2 = *gridManager2.c_grid();
BlackoilState state1( UgGridHelpers::numCells( grid1 ) , UgGridHelpers::numFaces( grid1 ) , 3);
BlackoilState state2( UgGridHelpers::numCells( grid2 ) , UgGridHelpers::numFaces( grid2 ) , 3);
BOOST_CHECK_EQUAL( false , state1.equal(state2) );
BOOST_CHECK( ! state1.equal(state2) );
}
@ -55,44 +66,47 @@ BOOST_AUTO_TEST_CASE(EqualsDifferentDeckReturnFalse) {
BOOST_AUTO_TEST_CASE(EqualsNumericalDifferenceReturnFalse) {
const string filename = "testBlackoilState1.DATA";
Opm::ParseContext parseContext;
Opm::ParserPtr parser(new Opm::Parser());
Opm::DeckConstPtr deck(parser->parseFile(filename , parseContext));
GridManager gridManager(deck);
const auto es = Opm::Parser::parse(filename);
auto eg = es.getInputGridCopy();
std::vector<int> actnum = get_testBlackoilStateActnum();
eg->resetACTNUM(actnum.data());
GridManager gridManager(eg);
const UnstructuredGrid& grid = *gridManager.c_grid();
BlackoilState state1( UgGridHelpers::numCells( grid ) , UgGridHelpers::numFaces( grid ) , 3);
BlackoilState state2( UgGridHelpers::numCells( grid ) , UgGridHelpers::numFaces( grid ) , 3);
BOOST_CHECK_EQUAL( true , state1.equal(state2) );
BOOST_CHECK( state1.equal(state2) );
{
std::vector<double>& p1 = state1.pressure();
std::vector<double>& p2 = state2.pressure();
p1[0] = p1[0] * 2 + 1;
BOOST_CHECK_EQUAL( false , state1.equal(state2) );
BOOST_CHECK( ! state1.equal(state2) );
p1[0] = p2[0];
BOOST_CHECK_EQUAL( true , state1.equal(state2) );
BOOST_CHECK( state1.equal(state2) );
}
{
std::vector<double>& gor1 = state1.gasoilratio();
std::vector<double>& gor2 = state2.gasoilratio();
gor1[0] = gor1[0] * 2 + 1;
BOOST_CHECK_EQUAL( false , state1.equal(state2) );
BOOST_CHECK( ! state1.equal(state2) );
gor1[0] = gor2[0];
BOOST_CHECK_EQUAL( true , state1.equal(state2) );
BOOST_CHECK( state1.equal(state2) );
}
{
std::vector<double>& p1 = state1.facepressure();
std::vector<double>& p2 = state2.facepressure();
p1[0] = p1[0] * 2 + 1;
BOOST_CHECK_EQUAL( false , state1.equal(state2) );
BOOST_CHECK( ! state1.equal(state2) );
p1[0] = p2[0];
BOOST_CHECK_EQUAL( true , state1.equal(state2) );
BOOST_CHECK( state1.equal(state2) );
}
{
@ -101,9 +115,9 @@ BOOST_AUTO_TEST_CASE(EqualsNumericalDifferenceReturnFalse) {
if (f1.size() > 0 ) {
f1[0] = f1[0] * 2 + 1;
BOOST_CHECK_EQUAL( false , state1.equal(state2) );
BOOST_CHECK( ! state1.equal(state2) );
f1[0] = f2[0];
BOOST_CHECK_EQUAL( true , state1.equal(state2) );
BOOST_CHECK( state1.equal(state2) );
}
}
{
@ -112,9 +126,9 @@ BOOST_AUTO_TEST_CASE(EqualsNumericalDifferenceReturnFalse) {
if (sv1.size() > 0) {
sv1[0] = sv1[0] * 2 + 1;
BOOST_CHECK_EQUAL( false , state1.equal(state2) );
BOOST_CHECK( ! state1.equal(state2) );
sv1[0] = sv2[0];
BOOST_CHECK_EQUAL( true , state1.equal(state2) );
BOOST_CHECK( state1.equal(state2) );
}
}
{
@ -122,8 +136,8 @@ BOOST_AUTO_TEST_CASE(EqualsNumericalDifferenceReturnFalse) {
std::vector<double>& sat2 = state2.saturation();
sat1[0] = sat1[0] * 2 + 1;
BOOST_CHECK_EQUAL( false , state1.equal(state2) );
BOOST_CHECK( ! state1.equal(state2) );
sat1[0] = sat2[0];
BOOST_CHECK_EQUAL( true , state1.equal(state2) );
BOOST_CHECK( state1.equal(state2) );
}
}

14
tests/test_equil.cpp
View File

@ -366,7 +366,7 @@ BOOST_AUTO_TEST_CASE (DeckAllDead)
Opm::ParseContext parseContext;
Opm::ParserPtr parser(new Opm::Parser() );
Opm::DeckConstPtr deck = parser->parseFile("deadfluids.DATA" , parseContext);
Opm::EclipseStateConstPtr eclipseState(new Opm::EclipseState(deck, parseContext));
Opm::EclipseStateConstPtr eclipseState(new Opm::EclipseState(*deck, parseContext));
Opm::BlackoilPropertiesFromDeck props(deck, eclipseState, *grid, false);
Opm::EQUIL::DeckDependent::InitialStateComputer comp(props, deck, eclipseState, *grid, 10.0);
const auto& pressures = comp.press();
@ -394,7 +394,7 @@ BOOST_AUTO_TEST_CASE (CapillaryInversion)
Opm::ParserPtr parser(new Opm::Parser() );
Opm::ParseContext parseContext;
Opm::DeckConstPtr deck = parser->parseFile("capillary.DATA" , parseContext);
Opm::EclipseStateConstPtr eclipseState(new Opm::EclipseState(deck , parseContext));
Opm::EclipseStateConstPtr eclipseState(new Opm::EclipseState(*deck , parseContext));
Opm::BlackoilPropertiesFromDeck props(deck, eclipseState, grid, false);
// Test the capillary inversion for oil-water.
@ -448,7 +448,7 @@ BOOST_AUTO_TEST_CASE (DeckWithCapillary)
Opm::ParserPtr parser(new Opm::Parser() );
Opm::ParseContext parseContext;
Opm::DeckConstPtr deck = parser->parseFile("capillary.DATA" , parseContext);
Opm::EclipseStateConstPtr eclipseState(new Opm::EclipseState(deck , parseContext));
Opm::EclipseStateConstPtr eclipseState(new Opm::EclipseState(*deck , parseContext));
Opm::BlackoilPropertiesFromDeck props(deck, eclipseState, grid, false);
Opm::EQUIL::DeckDependent::InitialStateComputer comp(props, deck, eclipseState, grid, 10.0);
@ -489,7 +489,7 @@ BOOST_AUTO_TEST_CASE (DeckWithCapillaryOverlap)
Opm::ParserPtr parser(new Opm::Parser() );
Opm::ParseContext parseContext;
Opm::DeckConstPtr deck = parser->parseFile("capillary_overlap.DATA" , parseContext);
Opm::EclipseStateConstPtr eclipseState(new Opm::EclipseState(deck , parseContext));
Opm::EclipseStateConstPtr eclipseState(new Opm::EclipseState(*deck , parseContext));
Opm::BlackoilPropertiesFromDeck props(deck, eclipseState, grid, false);
Opm::EQUIL::DeckDependent::InitialStateComputer comp(props, deck, eclipseState, grid, 9.80665);
@ -552,7 +552,7 @@ BOOST_AUTO_TEST_CASE (DeckWithLiveOil)
Opm::ParserPtr parser(new Opm::Parser() );
Opm::ParseContext parseContext;
Opm::DeckConstPtr deck = parser->parseFile("equil_liveoil.DATA" , parseContext);
Opm::EclipseStateConstPtr eclipseState(new Opm::EclipseState(deck , parseContext));
Opm::EclipseStateConstPtr eclipseState(new Opm::EclipseState(*deck , parseContext));
Opm::BlackoilPropertiesFromDeck props(deck, eclipseState, grid, false);
Opm::EQUIL::DeckDependent::InitialStateComputer comp(props, deck, eclipseState, grid, 9.80665);
@ -632,7 +632,7 @@ BOOST_AUTO_TEST_CASE (DeckWithLiveGas)
Opm::ParserPtr parser(new Opm::Parser() );
Opm::ParseContext parseContext;
Opm::DeckConstPtr deck = parser->parseFile("equil_livegas.DATA" , parseContext);
Opm::EclipseStateConstPtr eclipseState(new Opm::EclipseState(deck , parseContext));
Opm::EclipseStateConstPtr eclipseState(new Opm::EclipseState(*deck , parseContext));
Opm::BlackoilPropertiesFromDeck props(deck, eclipseState, grid, false);
Opm::EQUIL::DeckDependent::InitialStateComputer comp(props, deck, eclipseState, grid, 9.80665);
@ -715,7 +715,7 @@ BOOST_AUTO_TEST_CASE (DeckWithRSVDAndRVVD)
Opm::ParserPtr parser(new Opm::Parser() );
Opm::ParseContext parseContext;
Opm::DeckConstPtr deck = parser->parseFile("equil_rsvd_and_rvvd.DATA", parseContext);
Opm::EclipseStateConstPtr eclipseState(new Opm::EclipseState(deck , parseContext));
Opm::EclipseStateConstPtr eclipseState(new Opm::EclipseState(*deck , parseContext));
Opm::BlackoilPropertiesFromDeck props(deck, eclipseState, grid, false);
Opm::EQUIL::DeckDependent::InitialStateComputer comp(props, deck, eclipseState, grid, 9.80665);

2
tests/test_norne_pvt.cpp
View File

@ -281,7 +281,7 @@ BOOST_AUTO_TEST_CASE( Test_Norne_PVT) {
std::shared_ptr<const Deck> deck;
deck = parser->parseFile("norne_pvt.data", parseContext);
Opm::EclipseStateConstPtr eclState(new EclipseState(deck, parseContext));
Opm::EclipseStateConstPtr eclState(new EclipseState(*deck, parseContext));
verify_norne_oil_pvt_region1( deck, eclState );
verify_norne_oil_pvt_region2( deck, eclState );

13
tests/test_parallelistlinformation.cpp
View File

@ -33,6 +33,7 @@
#include <functional>
#ifdef HAVE_DUNE_ISTL
template<typename T>
void runSumMaxMinTest(const T offset)
{
@ -45,12 +46,13 @@ void runSumMaxMinTest(const T offset)
assert(comm.indexSet()->size()==x.size());
for(auto it=comm.indexSet()->begin(), itend=comm.indexSet()->end(); it!=itend; ++it)
x[it->local()]=it->global()+offset;
auto containers = std::make_tuple(x, x, x, x);
auto containers = std::make_tuple(x, x, x, x, x);
auto operators = std::make_tuple(Opm::Reduction::makeGlobalSumFunctor<T>(),
Opm::Reduction::makeGlobalMaxFunctor<T>(),
Opm::Reduction::makeGlobalMinFunctor<T>(),
Opm::Reduction::makeInnerProductFunctor<T>());
auto values = std::tuple<T,T,T,T>(0,0,100000, 0);
Opm::Reduction::makeInnerProductFunctor<T>(),
Opm::Reduction::makeLInfinityNormFunctor<T>());
auto values = std::tuple<T,T,T,T,T>(0,0,100000, 0, 0);
auto oldvalues = values;
start = offset;
end = start+N;
@ -59,10 +61,14 @@ void runSumMaxMinTest(const T offset)
BOOST_CHECK(std::get<1>(values)==std::max(N+offset-1, std::get<1>(oldvalues)));
BOOST_CHECK(std::get<2>(values)==std::min(offset, std::get<2>(oldvalues)));
BOOST_CHECK(std::get<3>(values)==((end-1)*end*(2*end-1)-(start-1)*start*(2*start-1))/6+std::get<3>(oldvalues));
// Must avoid std::abs() directly to prevent ambiguity with unsigned integers.
Opm::Reduction::detail::MaxAbsFunctor<T> maxabsfunc;
BOOST_CHECK(std::get<4>(values)==maxabsfunc(offset, N+offset-1));
}
BOOST_AUTO_TEST_CASE(tupleReductionTestInt)
{
runSumMaxMinTest<int>(-200);
runSumMaxMinTest<int>(0);
runSumMaxMinTest<int>(20);
runSumMaxMinTest<int>(-20);
@ -75,6 +81,7 @@ BOOST_AUTO_TEST_CASE(tupleReductionTestUnsignedInt)
}
BOOST_AUTO_TEST_CASE(tupleReductionTestFloat)
{
runSumMaxMinTest<float>(-200);
runSumMaxMinTest<float>(0);
runSumMaxMinTest<float>(20);
runSumMaxMinTest<float>(-20);

8
tests/test_pinchprocessor.cpp
View File

@ -54,9 +54,9 @@ BOOST_AUTO_TEST_CASE(Processing)
Opm::ParserPtr parser(new Opm::Parser());
Opm::ParseContext parseContext({{ ParseContext::PARSE_RANDOM_SLASH , InputError::IGNORE }});
Opm::DeckConstPtr deck = parser->parseFile(filename, parseContext);
std::shared_ptr<EclipseState> eclstate (new Opm::EclipseState(deck, parseContext));
std::shared_ptr<EclipseState> eclstate (new Opm::EclipseState(*deck, parseContext));
const auto& porv = eclstate->get3DProperties().getDoubleGridProperty("PORV").getData();
EclipseGridConstPtr eclgrid = eclstate->getEclipseGrid();
EclipseGridConstPtr eclgrid = eclstate->getInputGrid();
BOOST_CHECK_EQUAL(eclgrid->getMinpvMode(), MinpvMode::EclSTD);
@ -92,10 +92,10 @@ BOOST_AUTO_TEST_CASE(Processing)
std::vector<double> htrans(Opm::UgGridHelpers::numCellFaces(grid));
Grid* ug = const_cast<Grid*>(& grid);
tpfa_htrans_compute(ug, rock.permeability(), htrans.data());
auto transMult = eclstate->getTransMult();
const auto& transMult = eclstate->getTransMult();
std::vector<double> multz(nc, 0.0);
for (int i = 0; i < nc; ++i) {
multz[i] = transMult->getMultiplier(global_cell[i], Opm::FaceDir::ZPlus);
multz[i] = transMult.getMultiplier(global_cell[i], Opm::FaceDir::ZPlus);
}
Opm::NNC nnc(deck, eclgrid);
pinch.process(grid, htrans, actnum, multz, porv, nnc);

16
tests/test_relpermdiagnostics.cpp
View File

@ -34,7 +34,7 @@
#include <boost/test/floating_point_comparison.hpp>
#include <opm/common/utility/platform_dependent/reenable_warnings.h>
#include <opm/common/OpmLog/OpmLog.hpp>
#include <opm/common/OpmLog/EclipsePRTLog.hpp>
#include <opm/common/OpmLog/CounterLog.hpp>
#include <opm/core/grid.h>
#include <opm/core/grid/cart_grid.h>
@ -57,17 +57,13 @@ BOOST_AUTO_TEST_CASE(diagnosis)
{ ParseContext::PARSE_RANDOM_TEXT, InputError::IGNORE}
});
Opm::DeckConstPtr deck(parser->parseFile("../tests/relpermDiagnostics.DATA", parseContext));
eclState.reset(new EclipseState(deck, parseContext));
GridManager gm(deck);
eclState.reset(new EclipseState(*deck, parseContext));
GridManager gm(eclState->getInputGrid());
const UnstructuredGrid& grid = *gm.c_grid();
std::string logFile = "LOGFILE.txt";
std::shared_ptr<EclipsePRTLog> prtLog = std::make_shared<EclipsePRTLog>(logFile, Log::DefaultMessageTypes);
OpmLog::addBackend( "ECLIPSEPRTLOG" , prtLog );
std::shared_ptr<CounterLog> counterLog = std::make_shared<CounterLog>(Log::DefaultMessageTypes);
OpmLog::addBackend( "COUNTERLOG" , counterLog );
RelpermDiagnostics diagnostics;
diagnostics.diagnosis(eclState, deck, grid);
auto msg = diagnostics.getMessages();
BOOST_CHECK(!msg.empty());
BOOST_CHECK_EQUAL(msg.size(), 1);
BOOST_CHECK_EQUAL(1, counterLog->numMessages(Log::MessageType::Warning));
}
BOOST_AUTO_TEST_SUITE_END()

10
tests/test_satfunc.cpp
View File

@ -69,7 +69,7 @@ BOOST_AUTO_TEST_CASE (GwsegStandard)
Opm::ParserPtr parser(new Opm::Parser() );
Opm::ParseContext parseContext;
Opm::DeckConstPtr deck = parser->parseFile("satfuncStandard.DATA", parseContext);
Opm::EclipseStateConstPtr eclipseState(new Opm::EclipseState(deck , parseContext));
Opm::EclipseStateConstPtr eclipseState(new Opm::EclipseState(*deck , parseContext));
Opm::BlackoilPropertiesFromDeck props(deck, eclipseState, grid, param, false);
const int np = 3;
@ -155,7 +155,7 @@ BOOST_AUTO_TEST_CASE (GwsegEPSBase)
Opm::ParserPtr parser(new Opm::Parser() );
Opm::ParseContext parseContext;
Opm::DeckConstPtr deck = parser->parseFile("satfuncEPSBase.DATA" , parseContext);
Opm::EclipseStateConstPtr eclipseState(new Opm::EclipseState(deck , parseContext));
Opm::EclipseStateConstPtr eclipseState(new Opm::EclipseState(*deck , parseContext));
Opm::BlackoilPropertiesFromDeck props(deck, eclipseState, grid, param, false);
const int np = 3;
@ -241,7 +241,7 @@ BOOST_AUTO_TEST_CASE (GwsegEPS_A)
Opm::ParseContext parseContext;
Opm::ParserPtr parser(new Opm::Parser() );
Opm::DeckConstPtr deck = parser->parseFile("satfuncEPS_A.DATA" , parseContext);
Opm::EclipseStateConstPtr eclipseState(new Opm::EclipseState(deck , parseContext));
Opm::EclipseStateConstPtr eclipseState(new Opm::EclipseState(*deck , parseContext));
Opm::BlackoilPropertiesFromDeck props(deck, eclipseState, grid, param, false);
const int np = 3;
@ -493,7 +493,7 @@ BOOST_AUTO_TEST_CASE (GwsegEPS_C)
Opm::ParserPtr parser(new Opm::Parser() );
Opm::ParseContext parseContext;
Opm::DeckConstPtr deck = parser->parseFile("satfuncEPS_C.DATA", parseContext);
Opm::EclipseStateConstPtr eclipseState(new Opm::EclipseState(deck , parseContext));
Opm::EclipseStateConstPtr eclipseState(new Opm::EclipseState(*deck , parseContext));
Opm::BlackoilPropertiesFromDeck props(deck, eclipseState, grid, param, false);
const int np = 3;
@ -596,7 +596,7 @@ BOOST_AUTO_TEST_CASE (GwsegEPS_D)
Opm::ParserPtr parser(new Opm::Parser() );
Opm::ParseContext parseContext;
Opm::DeckConstPtr deck = parser->parseFile("satfuncEPS_D.DATA" , parseContext);
Opm::EclipseStateConstPtr eclipseState(new Opm::EclipseState(deck , parseContext));
Opm::EclipseStateConstPtr eclipseState(new Opm::EclipseState(*deck , parseContext));
Opm::BlackoilPropertiesFromDeck props(deck, eclipseState, grid, param, false);
const int np = 3;

4
tests/test_stoppedwells.cpp
View File

@ -47,8 +47,8 @@ BOOST_AUTO_TEST_CASE(TestStoppedWells)
Opm::ParseContext parseContext;
Opm::DeckConstPtr deck(parser->parseFile(filename , parseContext));
Opm::EclipseStateConstPtr eclipseState(new Opm::EclipseState(deck , parseContext));
Opm::GridManager gridManager(deck);
Opm::EclipseStateConstPtr eclipseState(new Opm::EclipseState(*deck , parseContext));
Opm::GridManager gridManager(eclipseState->getInputGrid());
double target_surfacerate_inj;
double target_surfacerate_prod;

16
tests/test_wellcollection.cpp
View File

@ -44,7 +44,7 @@ BOOST_AUTO_TEST_CASE(AddWellsAndGroupToCollection) {
std::string scheduleFile("wells_group.data");
ParseContext parseContext;
DeckConstPtr deck = parser->parseFile(scheduleFile, parseContext);
EclipseStateConstPtr eclipseState(new EclipseState(deck, parseContext));
EclipseStateConstPtr eclipseState(new EclipseState(*deck, parseContext));
PhaseUsage pu = phaseUsageFromDeck(eclipseState);
GroupTreeNodePtr field=eclipseState->getSchedule()->getGroupTree(2)->getNode("FIELD");
@ -54,24 +54,24 @@ BOOST_AUTO_TEST_CASE(AddWellsAndGroupToCollection) {
WellCollection collection;
// Add groups to WellCollection
GroupConstPtr fieldGroup = eclipseState->getSchedule()->getGroup(field->name());
const auto* fieldGroup = eclipseState->getSchedule()->getGroup(field->name());
collection.addField(fieldGroup, 2, pu);
for (auto iter = field->begin(); iter != field->end(); ++iter) {
GroupConstPtr childGroupNode = eclipseState->getSchedule()->getGroup((*iter).second->name());
const auto* childGroupNode = eclipseState->getSchedule()->getGroup((*iter).second->name());
collection.addGroup(childGroupNode, fieldGroup->name(), 2, pu);
}
GroupConstPtr g1Group = eclipseState->getSchedule()->getGroup(g1->name());
const auto* g1Group = eclipseState->getSchedule()->getGroup(g1->name());
for (auto iter = g1->begin(); iter != g1->end(); ++iter) {
GroupConstPtr childGroupNode = eclipseState->getSchedule()->getGroup((*iter).second->name());
const auto* childGroupNode = eclipseState->getSchedule()->getGroup((*iter).second->name());
collection.addGroup(childGroupNode, g1Group->name(), 2, pu);
}
GroupConstPtr g2Group = eclipseState->getSchedule()->getGroup(g2->name());
const auto* g2Group = eclipseState->getSchedule()->getGroup(g2->name());
for (auto iter = g2->begin(); iter != g2->end(); ++iter) {
GroupConstPtr childGroupNode = eclipseState->getSchedule()->getGroup((*iter).second->name());
const auto* childGroupNode = eclipseState->getSchedule()->getGroup((*iter).second->name());
collection.addGroup(childGroupNode, g2Group->name(), 2, pu);
}
@ -81,7 +81,7 @@ BOOST_AUTO_TEST_CASE(AddWellsAndGroupToCollection) {
// Add wells to WellCollection
WellCollection wellCollection;
std::vector<WellConstPtr> wells = eclipseState->getSchedule()->getWells();
auto wells = eclipseState->getSchedule()->getWells();
for (size_t i=0; i<wells.size(); i++) {
collection.addWell(wells[i], 2, pu);
}

10
tests/test_wellsgroup.cpp
View File

@ -52,13 +52,13 @@ BOOST_AUTO_TEST_CASE(ConstructGroupFromWell) {
std::string scheduleFile("wells_group.data");
ParseContext parseContext;
DeckConstPtr deck = parser->parseFile(scheduleFile, parseContext);
EclipseStateConstPtr eclipseState(new EclipseState(deck , parseContext));
EclipseStateConstPtr eclipseState(new EclipseState(*deck , parseContext));
PhaseUsage pu = phaseUsageFromDeck(eclipseState);
std::vector<WellConstPtr> wells = eclipseState->getSchedule()->getWells();
auto wells = eclipseState->getSchedule()->getWells();
for (size_t i=0; i<wells.size(); i++) {
WellConstPtr well = wells[i];
const auto* well = wells[i];
std::shared_ptr<WellsGroupInterface> wellsGroup = createWellWellsGroup(well, 2, pu);
BOOST_CHECK_EQUAL(well->name(), wellsGroup->name());
if (well->isInjector(2)) {
@ -85,13 +85,13 @@ BOOST_AUTO_TEST_CASE(ConstructGroupFromGroup) {
ParseContext parseContext;
std::string scheduleFile("wells_group.data");
DeckConstPtr deck = parser->parseFile(scheduleFile, parseContext);
EclipseStateConstPtr eclipseState(new EclipseState(deck , parseContext));
EclipseStateConstPtr eclipseState(new EclipseState(*deck , parseContext));
PhaseUsage pu = phaseUsageFromDeck(eclipseState);
std::vector<GroupTreeNodeConstPtr> nodes = eclipseState->getSchedule()->getGroupTree(2)->getNodes();
for (size_t i=0; i<nodes.size(); i++) {
GroupConstPtr group = eclipseState->getSchedule()->getGroup(nodes[i]->name());
const auto* group = eclipseState->getSchedule()->getGroup(nodes[i]->name());
std::shared_ptr<WellsGroupInterface> wellsGroup = createGroupWellsGroup(group, 2, pu);
BOOST_CHECK_EQUAL(group->name(), wellsGroup->name());
if (group->isInjectionGroup(2)) {

416
tests/test_wellsmanager.cpp
View File

@ -29,273 +29,255 @@
#include <boost/test/unit_test.hpp>
#include <opm/parser/eclipse/Parser/Parser.hpp>
#include <opm/parser/eclipse/Parser/ParseContext.hpp>
#include <opm/parser/eclipse/EclipseState/Schedule/Schedule.hpp>
#include <opm/parser/eclipse/Parser/ParseContext.hpp>
#include <opm/parser/eclipse/EclipseState/Schedule/Schedule.hpp>
#include <opm/core/wells/WellsManager.hpp>
#include <opm/core/wells.h>
#include <opm/core/well_controls.h>
#include <opm/core/wells/WellsManager.hpp>
#include <opm/core/wells.h>
#include <opm/core/well_controls.h>
#include <opm/core/grid/GridManager.hpp>
#include <opm/core/grid/GridManager.hpp>
void wells_static_check(const Wells* wells) {
BOOST_CHECK_EQUAL(2, wells->number_of_wells);
BOOST_CHECK_EQUAL(3, wells->number_of_phases);
void wells_static_check(const Wells* wells) {
BOOST_CHECK_EQUAL(2 , wells->number_of_wells);
BOOST_CHECK_EQUAL(3 , wells->number_of_phases);
BOOST_CHECK_EQUAL("INJ1", wells->name[0]);
BOOST_CHECK_EQUAL("PROD1", wells->name[1]);
BOOST_CHECK_EQUAL("INJ1" , wells->name[0]);
BOOST_CHECK_EQUAL("PROD1" , wells->name[1]);
/* The mapping from well number into the wells->WI and wells->well_cells arrays. */
BOOST_CHECK_EQUAL(0, wells->well_connpos[0]);
BOOST_CHECK_EQUAL(1, wells->well_connpos[1]);
BOOST_CHECK_EQUAL(2, wells->well_connpos[2]);
/* The mapping from well number into the wells->WI and wells->well_cells arrays. */
BOOST_CHECK_EQUAL(0 , wells->well_connpos[0]);
BOOST_CHECK_EQUAL(1 , wells->well_connpos[1]);
BOOST_CHECK_EQUAL(2 , wells->well_connpos[2]);
/* Connection factor */
BOOST_CHECK_CLOSE(1.2279166666666664e-12, wells->WI[0], 0.001);
BOOST_CHECK_CLOSE(1.2279166666666664e-12, wells->WI[1], 0.001);
/* Connection factor */
BOOST_CHECK_CLOSE(1.2279166666666664e-12 , wells->WI[0] , 0.001);
BOOST_CHECK_CLOSE(1.2279166666666664e-12 , wells->WI[1] , 0.001);
/* Completed cells */
BOOST_CHECK_EQUAL(0, wells->well_cells[0]);
BOOST_CHECK_EQUAL(9 + 2 * 10 + 2 * 10 * 10, wells->well_cells[1]);
}
/* Completed cells */
BOOST_CHECK_EQUAL(0 , wells->well_cells[0]);
BOOST_CHECK_EQUAL(9 + 2*10 + 2*10*10 , wells->well_cells[1]);
/*
The number of controls is determined by looking at which elements
have been given explicit - non-default - values in the WCONxxxx
keyword. Is that at all interesting?
*/
void check_controls_epoch0(struct WellControls ** ctrls) {
// The injector
{
const struct WellControls * ctrls0 = ctrls[0];
BOOST_CHECK_EQUAL(3, well_controls_get_num(ctrls0)); // The number of controls for the injector == 3??
BOOST_CHECK_EQUAL(SURFACE_RATE, well_controls_iget_type(ctrls0, 0));
BOOST_CHECK_EQUAL(RESERVOIR_RATE, well_controls_iget_type(ctrls0, 1));
BOOST_CHECK_EQUAL(BHP, well_controls_iget_type(ctrls0, 2));
// The different targets
BOOST_CHECK_EQUAL(100.0 / 86400, well_controls_iget_target(ctrls0, 0));
BOOST_CHECK_EQUAL(200.0 / 86400, well_controls_iget_target(ctrls0, 1));
BOOST_CHECK_EQUAL(400 * 100000, well_controls_iget_target(ctrls0, 2));
// Which control is active
BOOST_CHECK_EQUAL(0, well_controls_get_current(ctrls0));
// The phase distribution in the active target
{
const double * distr = well_controls_iget_distr(ctrls0, 0);
BOOST_CHECK_EQUAL(0, distr[0]); // Water
BOOST_CHECK_EQUAL(0, distr[1]); // Oil
BOOST_CHECK_EQUAL(1, distr[2]); // Gas
}
}
// The producer
{
const struct WellControls * ctrls1 = ctrls[1];
BOOST_CHECK_EQUAL(2, well_controls_get_num(ctrls1)); // The number of controls for the producer == 2??
BOOST_CHECK_EQUAL(SURFACE_RATE, well_controls_iget_type(ctrls1, 0));
BOOST_CHECK_EQUAL(BHP, well_controls_iget_type(ctrls1, 1));
/*
The number of controls is determined by looking at which elements
have been given explicit - non-default - values in the WCONxxxx
keyword. Is that at all interesting?
*/
// The different targets
BOOST_CHECK_EQUAL(-20000.0 / 86400, well_controls_iget_target(ctrls1, 0));
BOOST_CHECK_EQUAL(1000 * 100000, well_controls_iget_target(ctrls1, 1));
// Which control is active
BOOST_CHECK_EQUAL(0, well_controls_get_current(ctrls1));
void check_controls_epoch0( struct WellControls ** ctrls) {
// The injector
{
const struct WellControls * ctrls0 = ctrls[0];
BOOST_CHECK_EQUAL( 3 , well_controls_get_num(ctrls0)); // The number of controls for the injector == 3??
BOOST_CHECK_EQUAL( SURFACE_RATE , well_controls_iget_type(ctrls0 , 0) );
BOOST_CHECK_EQUAL( RESERVOIR_RATE , well_controls_iget_type(ctrls0 , 1) );
BOOST_CHECK_EQUAL( BHP , well_controls_iget_type(ctrls0 , 2) );
// The different targets
BOOST_CHECK_EQUAL( 100.0 / 86400 , well_controls_iget_target(ctrls0,0));
BOOST_CHECK_EQUAL( 200.0 / 86400 , well_controls_iget_target(ctrls0,1));
BOOST_CHECK_EQUAL( 400 * 100000 , well_controls_iget_target(ctrls0,2));
// Which control is active
BOOST_CHECK_EQUAL( 0 , well_controls_get_current(ctrls0) );
// The phase distribution in the active target
{
const double * distr = well_controls_iget_distr( ctrls0 , 0 );
BOOST_CHECK_EQUAL( 0 , distr[0] ); // Water
BOOST_CHECK_EQUAL( 0 , distr[1] ); // Oil
BOOST_CHECK_EQUAL( 1 , distr[2] ); // Gas
}
}
// The producer
{
const struct WellControls * ctrls1 = ctrls[1];
BOOST_CHECK_EQUAL( 2 , well_controls_get_num( ctrls1 )); // The number of controls for the producer == 2??
BOOST_CHECK_EQUAL( SURFACE_RATE , well_controls_iget_type(ctrls1 , 0) );
BOOST_CHECK_EQUAL( BHP , well_controls_iget_type(ctrls1 , 1) );
// The different targets
BOOST_CHECK_EQUAL( -20000.0 / 86400 , well_controls_iget_target(ctrls1,0));
BOOST_CHECK_EQUAL( 1000 * 100000 , well_controls_iget_target(ctrls1,1));
// Which control is active
BOOST_CHECK_EQUAL( 0 , well_controls_get_current(ctrls1));
// The phase distribution in the active target
{
const double * distr = well_controls_iget_distr( ctrls1 , 0 );
BOOST_CHECK_EQUAL( 0 , distr[0] ); // Water
BOOST_CHECK_EQUAL( 1 , distr[1] ); // Oil
BOOST_CHECK_EQUAL( 0 , distr[2] ); // Gas
}
}
// The phase distribution in the active target
{
const double * distr = well_controls_iget_distr(ctrls1, 0);
BOOST_CHECK_EQUAL(0, distr[0]); // Water
BOOST_CHECK_EQUAL(1, distr[1]); // Oil
BOOST_CHECK_EQUAL(0, distr[2]); // Gas
}
}
}
void check_controls_epoch1(struct WellControls ** ctrls) {
// The injector
{
const struct WellControls * ctrls0 = ctrls[0];
BOOST_CHECK_EQUAL(3, well_controls_get_num(ctrls0)); // The number of controls for the injector == 3??
BOOST_CHECK_EQUAL(SURFACE_RATE, well_controls_iget_type(ctrls0, 0));
BOOST_CHECK_EQUAL(RESERVOIR_RATE, well_controls_iget_type(ctrls0, 1));
BOOST_CHECK_EQUAL(BHP, well_controls_iget_type(ctrls0, 2));
// The different targets
BOOST_CHECK_CLOSE(10.0 / 86400, well_controls_iget_target(ctrls0, 0), 0.001);
BOOST_CHECK_CLOSE(20.0 / 86400, well_controls_iget_target(ctrls0, 1), 0.001);
BOOST_CHECK_CLOSE(40 * 100000, well_controls_iget_target(ctrls0, 2), 0.001);
void check_controls_epoch1( struct WellControls ** ctrls) {
// The injector
{
const struct WellControls * ctrls0 = ctrls[0];
BOOST_CHECK_EQUAL( 3 , well_controls_get_num(ctrls0)); // The number of controls for the injector == 3??
// Which control is active
BOOST_CHECK_EQUAL(1, well_controls_get_current(ctrls0));
BOOST_CHECK_EQUAL( SURFACE_RATE , well_controls_iget_type(ctrls0 , 0 ));
BOOST_CHECK_EQUAL( RESERVOIR_RATE , well_controls_iget_type(ctrls0 , 1 ));
BOOST_CHECK_EQUAL( BHP , well_controls_iget_type(ctrls0 , 2 ));
// The different targets
BOOST_CHECK_CLOSE( 10.0 / 86400 , well_controls_iget_target(ctrls0 , 0) , 0.001);
BOOST_CHECK_CLOSE( 20.0 / 86400 , well_controls_iget_target(ctrls0 , 1) , 0.001);
BOOST_CHECK_CLOSE( 40 * 100000 , well_controls_iget_target(ctrls0 , 2) , 0.001);
// Which control is active
BOOST_CHECK_EQUAL( 1 , well_controls_get_current(ctrls0));
{
const double * distr = well_controls_iget_distr( ctrls0 , 1 );
BOOST_CHECK_EQUAL( 1 , distr[0] ); // Water
BOOST_CHECK_EQUAL( 0 , distr[1] ); // Oil
BOOST_CHECK_EQUAL( 0 , distr[2] ); // Gas
}
}
// The producer
{
const struct WellControls * ctrls1 = ctrls[1];
BOOST_CHECK_EQUAL( 3 , well_controls_get_num(ctrls1)); // The number of controls for the producer - now 3.
BOOST_CHECK_EQUAL( SURFACE_RATE , well_controls_iget_type(ctrls1 , 0) );
BOOST_CHECK_EQUAL( RESERVOIR_RATE , well_controls_iget_type(ctrls1 , 1) );
BOOST_CHECK_EQUAL( BHP , well_controls_iget_type(ctrls1 , 2) );
// The different targets
BOOST_CHECK_CLOSE( -999.0 / 86400 , well_controls_iget_target(ctrls1 , 0), 0.001);
BOOST_CHECK_CLOSE( -123.0 / 86400 , well_controls_iget_target(ctrls1 , 1), 0.001);
BOOST_CHECK_CLOSE( 100 * 100000 , well_controls_iget_target(ctrls1 , 2), 0.001);
// Which control is active
BOOST_CHECK_EQUAL( 1 , well_controls_get_current(ctrls1) );
{
const double * distr = well_controls_iget_distr( ctrls1 , 1 );
BOOST_CHECK_EQUAL( 1 , distr[0] ); // Water
BOOST_CHECK_EQUAL( 1 , distr[1] ); // Oil
BOOST_CHECK_EQUAL( 1 , distr[2] ); // Gas
}
}
{
const double * distr = well_controls_iget_distr(ctrls0, 1);
BOOST_CHECK_EQUAL(1, distr[0]); // Water
BOOST_CHECK_EQUAL(0, distr[1]); // Oil
BOOST_CHECK_EQUAL(0, distr[2]); // Gas
}
}
// The producer
{
const struct WellControls * ctrls1 = ctrls[1];
BOOST_CHECK_EQUAL(3, well_controls_get_num(ctrls1)); // The number of controls for the producer - now 3.
BOOST_CHECK_EQUAL(SURFACE_RATE, well_controls_iget_type(ctrls1, 0));
BOOST_CHECK_EQUAL(RESERVOIR_RATE, well_controls_iget_type(ctrls1, 1));
BOOST_CHECK_EQUAL(BHP, well_controls_iget_type(ctrls1, 2));
void check_controls_epoch3( struct WellControls ** ctrls) {
// The new producer
const struct WellControls * ctrls1 = ctrls[1];
// Note: controls include default (1 atm) BHP control.
BOOST_CHECK_EQUAL( 6 , well_controls_get_num(ctrls1));
// The different targets
BOOST_CHECK_CLOSE(-999.0 / 86400, well_controls_iget_target(ctrls1, 0), 0.001);
BOOST_CHECK_CLOSE(-123.0 / 86400, well_controls_iget_target(ctrls1, 1), 0.001);
BOOST_CHECK_CLOSE(100 * 100000, well_controls_iget_target(ctrls1, 2), 0.001);
// Which control is active
BOOST_CHECK_EQUAL(1, well_controls_get_current(ctrls1));
{
const double * distr = well_controls_iget_distr(ctrls1, 1);
BOOST_CHECK_EQUAL(1, distr[0]); // Water
BOOST_CHECK_EQUAL(1, distr[1]); // Oil
BOOST_CHECK_EQUAL(1, distr[2]); // Gas
}
}
}
void check_controls_epoch3(struct WellControls ** ctrls) {
// The new producer
const struct WellControls * ctrls1 = ctrls[1];
// Note: controls include default (1 atm) BHP control.
BOOST_CHECK_EQUAL(6, well_controls_get_num(ctrls1));
}
BOOST_AUTO_TEST_CASE(New_Constructor_Works) {
const std::string filename = "wells_manager_data.data";
Opm::ParserPtr parser(new Opm::Parser());
Opm::ParseContext parseContext;
Opm::DeckConstPtr deck(parser->parseFile(filename, parseContext));
BOOST_AUTO_TEST_CASE(New_Constructor_Works) {
Opm::EclipseStateConstPtr eclipseState(new Opm::EclipseState(*deck, parseContext));
Opm::GridManager gridManager(eclipseState->getInputGrid());
const std::string filename = "wells_manager_data.data";
Opm::ParserPtr parser(new Opm::Parser());
Opm::ParseContext parseContext;
Opm::DeckConstPtr deck(parser->parseFile(filename , parseContext));
{
Opm::WellsManager wellsManager(eclipseState, 0, *gridManager.c_grid(), NULL);
wells_static_check(wellsManager.c_wells());
check_controls_epoch0(wellsManager.c_wells()->ctrls);
}
Opm::EclipseStateConstPtr eclipseState(new Opm::EclipseState(deck , parseContext));
Opm::GridManager gridManager(deck);
{
Opm::WellsManager wellsManager(eclipseState, 1, *gridManager.c_grid(), NULL);
wells_static_check(wellsManager.c_wells());
check_controls_epoch1(wellsManager.c_wells()->ctrls);
}
{
Opm::WellsManager wellsManager(eclipseState, 0, *gridManager.c_grid(), NULL);
wells_static_check( wellsManager.c_wells() );
check_controls_epoch0( wellsManager.c_wells()->ctrls );
}
{
Opm::WellsManager wellsManager(eclipseState, 3, *gridManager.c_grid(), NULL);
const Wells* wells = wellsManager.c_wells();
{
Opm::WellsManager wellsManager(eclipseState, 1, *gridManager.c_grid(), NULL);
wells_static_check( wellsManager.c_wells() );
check_controls_epoch1( wellsManager.c_wells()->ctrls );
}
// There is 3 wells in total in the deck at the 3rd schedule step.
// PROD1 is shut and should therefore not be counted.
// The new well is therefore the secound well.
BOOST_CHECK_EQUAL(2, wells->number_of_wells);
BOOST_CHECK_EQUAL(wells->name[0], "INJ1");
BOOST_CHECK_EQUAL(wells->name[1], "NEW");
check_controls_epoch3(wellsManager.c_wells()->ctrls);
}
{
Opm::WellsManager wellsManager(eclipseState, 3, *gridManager.c_grid(), NULL);
const Wells* wells = wellsManager.c_wells();
}
// There is 3 wells in total in the deck at the 3rd schedule step.
// PROD1 is shut and should therefore not be counted.
// The new well is therefore the secound well.
BOOST_CHECK_EQUAL(2 , wells->number_of_wells);
BOOST_CHECK_EQUAL( wells->name[0] , "INJ1");
BOOST_CHECK_EQUAL( wells->name[1] , "NEW");
BOOST_AUTO_TEST_CASE(WellsEqual) {
const std::string filename = "wells_manager_data.data";
Opm::ParserPtr parser(new Opm::Parser());
Opm::ParseContext parseContext;
Opm::DeckConstPtr deck(parser->parseFile(filename, parseContext));
check_controls_epoch3( wellsManager.c_wells()->ctrls );
}
Opm::EclipseStateConstPtr eclipseState(new Opm::EclipseState(*deck, parseContext));
Opm::GridManager gridManager(eclipseState->getInputGrid());
}
Opm::WellsManager wellsManager0(eclipseState, 0, *gridManager.c_grid(), NULL);
Opm::WellsManager wellsManager1(eclipseState, 1, *gridManager.c_grid(), NULL);
BOOST_CHECK(wells_equal( wellsManager0.c_wells() , wellsManager0.c_wells(),false));
BOOST_CHECK(!wells_equal( wellsManager0.c_wells() , wellsManager1.c_wells(),false));
}
BOOST_AUTO_TEST_CASE(ControlsEqual) {
const std::string filename = "wells_manager_data.data";
Opm::ParseContext parseContext;
Opm::ParserPtr parser(new Opm::Parser());
Opm::DeckConstPtr deck(parser->parseFile(filename, parseContext));
BOOST_AUTO_TEST_CASE(WellsEqual) {
const std::string filename = "wells_manager_data.data";
Opm::ParserPtr parser(new Opm::Parser());
Opm::ParseContext parseContext;
Opm::DeckConstPtr deck(parser->parseFile(filename , parseContext));
Opm::EclipseStateConstPtr eclipseState(new Opm::EclipseState(*deck, parseContext));
Opm::GridManager gridManager(eclipseState->getInputGrid());
Opm::EclipseStateConstPtr eclipseState(new Opm::EclipseState(deck , parseContext));
Opm::GridManager gridManager(deck);
Opm::WellsManager wellsManager0(eclipseState, 0, *gridManager.c_grid(), NULL);
Opm::WellsManager wellsManager1(eclipseState, 1, *gridManager.c_grid(), NULL);
Opm::WellsManager wellsManager0(eclipseState , 0 , *gridManager.c_grid(), NULL);
Opm::WellsManager wellsManager1(eclipseState , 1 , *gridManager.c_grid(), NULL);
BOOST_CHECK(well_controls_equal( wellsManager0.c_wells()->ctrls[0] , wellsManager0.c_wells()->ctrls[0] , false));
BOOST_CHECK(well_controls_equal( wellsManager0.c_wells()->ctrls[1] , wellsManager0.c_wells()->ctrls[1] , false));
BOOST_CHECK(well_controls_equal( wellsManager1.c_wells()->ctrls[0] , wellsManager1.c_wells()->ctrls[0] , false));
BOOST_CHECK(well_controls_equal( wellsManager1.c_wells()->ctrls[1] , wellsManager1.c_wells()->ctrls[1] , false));
BOOST_CHECK( wells_equal( wellsManager0.c_wells() , wellsManager0.c_wells(),false) );
BOOST_CHECK( !wells_equal( wellsManager0.c_wells() , wellsManager1.c_wells(),false) );
}
BOOST_CHECK(!well_controls_equal( wellsManager0.c_wells()->ctrls[0] , wellsManager0.c_wells()->ctrls[1] , false));
BOOST_CHECK(!well_controls_equal( wellsManager0.c_wells()->ctrls[1] , wellsManager0.c_wells()->ctrls[0] , false));
BOOST_CHECK(!well_controls_equal( wellsManager1.c_wells()->ctrls[0] , wellsManager0.c_wells()->ctrls[0] , false));
BOOST_CHECK(!well_controls_equal( wellsManager1.c_wells()->ctrls[1] , wellsManager0.c_wells()->ctrls[1] , false));
}
BOOST_AUTO_TEST_CASE(WellShutOK) {
const std::string filename = "wells_manager_data.data";
Opm::ParserPtr parser(new Opm::Parser());
Opm::ParseContext parseContext;
Opm::DeckConstPtr deck(parser->parseFile(filename, parseContext));
BOOST_AUTO_TEST_CASE(ControlsEqual) {
const std::string filename = "wells_manager_data.data";
Opm::ParseContext parseContext;
Opm::ParserPtr parser(new Opm::Parser());
Opm::DeckConstPtr deck(parser->parseFile(filename , parseContext));
Opm::EclipseStateConstPtr eclipseState(new Opm::EclipseState(*deck, parseContext));
Opm::GridManager gridManager(eclipseState->getInputGrid());
Opm::EclipseStateConstPtr eclipseState(new Opm::EclipseState(deck , parseContext));
Opm::GridManager gridManager(deck);
Opm::WellsManager wellsManager2(eclipseState, 2, *gridManager.c_grid(), NULL);
Opm::WellsManager wellsManager0(eclipseState , 0 , *gridManager.c_grid(), NULL);
Opm::WellsManager wellsManager1(eclipseState , 1 , *gridManager.c_grid(), NULL);
// Shut wells are not added to the deck. i.e number of wells should be 2-1
BOOST_CHECK(wellsManager2.c_wells()->number_of_wells == 1);
BOOST_CHECK( well_controls_equal( wellsManager0.c_wells()->ctrls[0] , wellsManager0.c_wells()->ctrls[0] , false));
BOOST_CHECK( well_controls_equal( wellsManager0.c_wells()->ctrls[1] , wellsManager0.c_wells()->ctrls[1] , false));
BOOST_CHECK( well_controls_equal( wellsManager1.c_wells()->ctrls[0] , wellsManager1.c_wells()->ctrls[0] , false));
BOOST_CHECK( well_controls_equal( wellsManager1.c_wells()->ctrls[1] , wellsManager1.c_wells()->ctrls[1] , false));
//BOOST_CHECK_NO_THROW( Opm::WellsManager wellsManager2(eclipseState , 2 , *gridManager.c_grid(), NULL));
}
BOOST_CHECK( !well_controls_equal( wellsManager0.c_wells()->ctrls[0] , wellsManager0.c_wells()->ctrls[1] , false));
BOOST_CHECK( !well_controls_equal( wellsManager0.c_wells()->ctrls[1] , wellsManager0.c_wells()->ctrls[0] , false));
BOOST_CHECK( !well_controls_equal( wellsManager1.c_wells()->ctrls[0] , wellsManager0.c_wells()->ctrls[0] , false));
BOOST_CHECK( !well_controls_equal( wellsManager1.c_wells()->ctrls[1] , wellsManager0.c_wells()->ctrls[1] , false));
}
BOOST_AUTO_TEST_CASE(WellSTOPOK) {
const std::string filename = "wells_manager_data_wellSTOP.data";
Opm::ParserPtr parser(new Opm::Parser());
Opm::ParseContext parseContext;
Opm::DeckConstPtr deck(parser->parseFile(filename, parseContext));
Opm::EclipseStateConstPtr eclipseState(new Opm::EclipseState(*deck, parseContext));
Opm::GridManager gridManager(eclipseState->getInputGrid());
BOOST_AUTO_TEST_CASE(WellShutOK) {
const std::string filename = "wells_manager_data.data";
Opm::ParserPtr parser(new Opm::Parser());
Opm::ParseContext parseContext;
Opm::DeckConstPtr deck(parser->parseFile(filename , parseContext));
Opm::EclipseStateConstPtr eclipseState(new Opm::EclipseState(deck , parseContext));
Opm::GridManager gridManager(deck);
Opm::WellsManager wellsManager2(eclipseState , 2 , *gridManager.c_grid(), NULL);
// Shut wells are not added to the deck. i.e number of wells should be 2-1
BOOST_CHECK( wellsManager2.c_wells()->number_of_wells == 1);
//BOOST_CHECK_NO_THROW( Opm::WellsManager wellsManager2(eclipseState , 2 , *gridManager.c_grid(), NULL));
}
BOOST_AUTO_TEST_CASE(WellSTOPOK) {
const std::string filename = "wells_manager_data_wellSTOP.data";
Opm::ParserPtr parser(new Opm::Parser());
Opm::ParseContext parseContext;
Opm::DeckConstPtr deck(parser->parseFile(filename , parseContext));
Opm::EclipseStateConstPtr eclipseState(new Opm::EclipseState(deck , parseContext));
Opm::GridManager gridManager(deck);
Opm::WellsManager wellsManager(eclipseState , 0 , *gridManager.c_grid(), NULL);
Opm::WellsManager wellsManager(eclipseState, 0, *gridManager.c_grid(), NULL);
const Wells* wells = wellsManager.c_wells();
const struct WellControls* ctrls0 = wells->ctrls[0];