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

Merge pull request #4164 from bska/state-indep-resv-calc

Browse Source 
Make Reservoir Voidage Rate Calculator Independent
This commit is contained in:
Bård Skaflestad 2022-11-10 13:06:18 +01:00 committed by GitHub
commit 622321a8e6
No known key found for this signature in database
GPG Key ID: 4AEE18F83AFDEB23
1 changed files with 155 additions and 77 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

232
opm/simulators/wells/RateConverter.hpp
View File

@ -24,10 +24,12 @@
#include <opm/core/props/BlackoilPhases.hpp>
#include <opm/grid/utility/RegionMapping.hpp>
#include <opm/simulators/wells/RegionAttributeHelpers.hpp>
#include <opm/simulators/utils/DeferredLoggingErrorHelpers.hpp>
#include <opm/simulators/wells/RegionAttributeHelpers.hpp>
#include <dune/grid/common/gridenums.hh>
#include <dune/grid/common/rangegenerators.hh>
#include <algorithm>
#include <cmath>
#include <memory>
@ -36,6 +38,7 @@
#include <unordered_map>
#include <utility>
#include <vector>
/**
* \file
* Facility for converting component rates at surface conditions to
@ -56,13 +59,13 @@ namespace Opm {
* The conversion uses fluid properties evaluated at average
* hydrocarbon pressure in regions or field.
*
* \tparam FluidSystem Fluid system class. Expected to be a BlackOilFluidSystem
* \tparam FluidSystem Fluid system class. Expected to be a
* BlackOilFluidSystem
*
* \tparam Region Type of a forward region mapping. Expected
* to provide indexed access through \code operator[]()
* \endcode as well as inner types \c value_type, \c
* size_type, and \c const_iterator. Typically \code
* std::vector<int> \endcode.
* \tparam Region Type of a forward region mapping. Expected to
* provide indexed access through \code operator[]() \endcode as
* well as inner types \c value_type, \c size_type, and \c
* const_iterator. Typically \code std::vector<int> \endcode.
*/
template <class FluidSystem, class Region>
class SurfaceToReservoirVoidage {
@ -70,18 +73,15 @@ namespace Opm {
/**
* Constructor.
*
* \param[in] region Forward region mapping. Often
* corresponds to the "FIPNUM" mapping of an ECLIPSE input
* deck.
* \param[in] region Forward region mapping. Often corresponds
* to the "FIPNUM" mapping of an ECLIPSE input deck.
*/
SurfaceToReservoirVoidage(const PhaseUsage& phaseUsage,
const Region& region)
const Region& region)
: phaseUsage_(phaseUsage)
, rmap_ (region)
, attr_ (rmap_, Attributes())
{
}
, rmap_ (region)
, attr_ (rmap_, Attributes())
{}
/**
* Compute pore volume averaged hydrocarbon state pressure, rs and rv.
@ -94,9 +94,8 @@ namespace Opm {
template <typename ElementContext, class EbosSimulator>
void defineState(const EbosSimulator& simulator)
{
// create map from cell to region
// and set all attributes to zero
// create map from cell to region and set all attributes to
// zero
for (const auto& reg : rmap_.activeRegions()) {
auto& ra = attr_.attributes(reg);
ra.pressure = 0.0;
@ -355,103 +354,134 @@ namespace Opm {
}
}
/**
* Convert surface volume flow rates to reservoir voidage flow
* rates.
*
* State dependent version. Client must call \code
* defineState() \endcode prior to invoking this member
* function.
*
* \tparam Rates Type representing contiguous collection of
* surface flow rates. Must support direct indexing through
* \code operator[]() \endcode.
*
* \param[in] r Zero based fluid-in-place region index.
*
* \param[in] pvtRegionIdx Zero based PVT region index.
*
* \param[in] surface_rates surface volume flow rates for all
* active phases.
*
* \param[out] voidage_rates reservoir volume flow rates for all
* active phases.
*/
template <class Rates>
void calcReservoirVoidageRates(const RegionId r,
const int pvtRegionIdx,
const Rates& surface_rates,
Rates& voidage_rates) const
{
const auto& ra = this->attr_.attributes(r);
this->calcReservoirVoidageRates(pvtRegionIdx,
ra.pressure, ra.rs, ra.rv,
ra.temperature,
ra.saltConcentration,
surface_rates,
voidage_rates);
}
/**
* Converting surface volume rates to reservoir voidage rates
* Convert surface volume flow rates to reservoir voidage flow
* rates.
*
* \tparam Rates Type representing contiguous collection
* of surface-to-reservoir conversion coefficients. Must
* support direct indexing through \code operator[]()
* \endcode.
* State independent version.
*
* \tparam Rates Type representing contiguous collection of
* surface flow rates. Must support direct indexing through
* \code operator[]() \endcode.
*
* \param[in] r Fluid-in-place region of the well
* \param[in] pvtRegionIdx PVT region of the well
* \param[in] surface_rates surface voluem rates for
* all active phases
* \param[in] pvtRegionIdx PVT region.
*
* \param[out] voidage_rates reservoir volume rates for
* all active phases
* \param[in] p Fluid pressure.
*
* \param[in] rs Dissolved gas/oil ratio.
*
* \param[in] rv Vaporised oil/gas ratio.
*
* \param[in] T Temperature. Unused in non-thermal simulation
* runs.
*
* \param[in] saltConcentration Salt concentration. Unused in
* simulation runs without salt precipitation.
*
* \param[in] surface_rates Surface volume flow rates for all
* active phases.
*
* \param[out] voidage_rates Reservoir volume flow rates for all
* active phases.
*/
template <class Rates >
void
calcReservoirVoidageRates(const RegionId r, const int pvtRegionIdx, const Rates& surface_rates,
Rates& voidage_rates) const
template <class Rates>
void calcReservoirVoidageRates(const int pvtRegionIdx,
const double p,
const double rs,
const double rv,
const double T,
const double saltConcentration,
const Rates& surface_rates,
Rates& voidage_rates) const
{
assert(voidage_rates.size() == surface_rates.size());
std::fill(voidage_rates.begin(), voidage_rates.end(), 0.0);
const auto& pu = phaseUsage_;
const auto& ra = attr_.attributes(r);
const double p = ra.pressure;
const double T = ra.temperature;
const double saltConcentration = ra.saltConcentration;
const auto& pu = this->phaseUsage_;
const auto iw = RegionAttributeHelpers::PhasePos::water(pu);
const auto io = RegionAttributeHelpers::PhasePos::oil (pu);
const auto ig = RegionAttributeHelpers::PhasePos::gas (pu);
const int iw = RegionAttributeHelpers::PhasePos::water(pu);
const int io = RegionAttributeHelpers::PhasePos::oil (pu);
const int ig = RegionAttributeHelpers::PhasePos::gas (pu);
const auto [Rs, Rv] = this->
dissolvedVaporisedRatio(io, ig, rs, rv, surface_rates);
if (RegionAttributeHelpers::PhaseUsed::water(pu)) {
// q[w]_r = q[w]_s / bw
const double bw = FluidSystem::waterPvt().inverseFormationVolumeFactor(pvtRegionIdx, T, p, saltConcentration);
const auto bw = FluidSystem::waterPvt()
.inverseFormationVolumeFactor(pvtRegionIdx, T, p,
saltConcentration);
voidage_rates[iw] = surface_rates[iw] / bw;
}
// Use average Rs and Rv:
auto a = ra.rs;
auto b = a;
if (io >= 0 && ig >= 0) {
b = surface_rates[ig]/(surface_rates[io]+1.0e-15);
}
double Rs = std::min(a, b);
a = ra.rv;
b = a;
if (io >= 0 && ig >= 0) {
b = surface_rates[io]/(surface_rates[ig]+1.0e-15);
}
double Rv = std::min(a, b);
// Determinant of 'R' matrix
const double detR = 1.0 - (Rs * Rv);
const auto detR = 1.0 - (Rs * Rv);
if (RegionAttributeHelpers::PhaseUsed::oil(pu)) {
// q[o]_r = 1/(bo * (1 - rs*rv)) * (q[o]_s - rv*q[g]_s)
const double bo = FluidSystem::oilPvt().inverseFormationVolumeFactor(pvtRegionIdx, T, p, Rs);
const double den = bo * detR;
voidage_rates[io] = surface_rates[io];
if (RegionAttributeHelpers::PhaseUsed::gas(pu)) {
voidage_rates[io] -= Rv * surface_rates[ig];
}
voidage_rates[io] /= den;
const auto bo = FluidSystem::oilPvt()
.inverseFormationVolumeFactor(pvtRegionIdx, T, p, Rs);
voidage_rates[io] /= bo * detR;
}
if (RegionAttributeHelpers::PhaseUsed::gas(pu)) {
// q[g]_r = 1/(bg * (1 - rs*rv)) * (q[g]_s - rs*q[o]_s)
const double bg = FluidSystem::gasPvt().inverseFormationVolumeFactor(pvtRegionIdx, T, p, Rv, 0.0 /*=Rvw*/);
const double den = bg * detR;
voidage_rates[ig] = surface_rates[ig];
if (RegionAttributeHelpers::PhaseUsed::oil(pu)) {
voidage_rates[ig] -= Rs * surface_rates[io];
}
voidage_rates[ig] /= den;
const auto bg = FluidSystem::gasPvt()
.inverseFormationVolumeFactor(pvtRegionIdx, T, p,
Rv, 0.0 /*=Rvw*/);
voidage_rates[ig] /= bg * detR;
}
}
/**
* Compute coefficients for surface-to-reservoir voidage
* conversion for solvent.
@ -510,7 +540,55 @@ namespace Opm {
RegionAttributeHelpers::RegionAttributes<RegionId, Attributes> attr_;
template <typename Rates>
std::pair<double, double>
dissolvedVaporisedRatio(const int io,
const int ig,
const double rs,
const double rv,
const Rates& surface_rates) const
{
if ((io < 0) || (ig < 0)) {
return { rs, rv };
}
#define BURN_RESV_BRIDGES 0
#if !BURN_RESV_BRIDGES
auto b = rs;
if (io >= 0 && ig >= 0) {
b = surface_rates[ig] / (surface_rates[io] + 1.0e-15);
}
const double Rs = std::min(rs, b);
b = rv;
if (io >= 0 && ig >= 0) {
b = surface_rates[io] / (surface_rates[ig] + 1.0e-15);
}
const double Rv = std::min(rv, b);
return { Rs, Rv };
#else // BURN_RESV_BRIDGES
auto eps = std::copysign(1.0e-15, surface_rates[io]);
const auto Rs = surface_rates[ig] / (surface_rates[io] + eps);
eps = std::copysign(1.0e-15, surface_rates[ig]);
const auto Rv = surface_rates[io] / (surface_rates[ig] + eps);
return {
std::clamp(static_cast<double>(Rs), 0.0, rs),
std::clamp(static_cast<double>(Rv), 0.0, rv)
};
#endif // BURN_RESV_BRIDGES
#undef BURN_RESV_BRIDGES
}
};
} // namespace RateConverter
} // namespace Opm