22#ifndef OPM_RATECONVERTER_HPP_HEADER_INCLUDED
23#define OPM_RATECONVERTER_HPP_HEADER_INCLUDED
25#include <opm/core/props/BlackoilPhases.hpp>
26#include <opm/grid/utility/RegionMapping.hpp>
27#include <opm/simulators/wells/RegionAttributeHelpers.hpp>
28#include <opm/simulators/utils/DeferredLoggingErrorHelpers.hpp>
29#include <dune/grid/common/gridenums.hh>
30#include <dune/grid/common/rangegenerators.hh>
36#include <unordered_map>
50 namespace RateConverter {
67 template <
class Flu
idSystem,
class Region>
81 , attr_ (rmap_, Attributes())
94 template <
typename ElementContext,
class EbosSimulator>
100 for (
const auto& reg : rmap_.activeRegions()) {
101 auto& ra = attr_.attributes(reg);
103 ra.temperature = 0.0;
107 ra.saltConcentration = 0.0;
112 std::unordered_map<RegionId, Attributes> attributes_pv;
115 std::unordered_map<RegionId, Attributes> attributes_hpv;
117 for (
const auto& reg : rmap_.activeRegions()) {
118 attributes_pv.insert({reg, Attributes()});
119 attributes_hpv.insert({reg, Attributes()});
122 ElementContext elemCtx( simulator );
123 const auto& gridView = simulator.gridView();
124 const auto& comm = gridView.comm();
126 OPM_BEGIN_PARALLEL_TRY_CATCH();
127 for (
const auto& elem : elements(gridView, Dune::Partitions::interior)) {
128 elemCtx.updatePrimaryStencil(elem);
129 elemCtx.updatePrimaryIntensiveQuantities(0);
130 const unsigned cellIdx = elemCtx.globalSpaceIndex(0, 0);
131 const auto& intQuants = elemCtx.intensiveQuantities(0, 0);
132 const auto& fs = intQuants.fluidState();
134 const double pv_cell =
135 simulator.model().dofTotalVolume(cellIdx)
136 * intQuants.porosity().value();
139 double hydrocarbon = 1.0;
140 const auto& pu = phaseUsage_;
142 hydrocarbon -= fs.saturation(FluidSystem::waterPhaseIdx).value();
145 const int reg = rmap_.region(cellIdx);
149 const double hydrocarbonPV = pv_cell*hydrocarbon;
150 if (hydrocarbonPV > 0.) {
151 auto& attr = attributes_hpv[reg];
152 attr.pv += hydrocarbonPV;
154 attr.rs += fs.Rs().value() * hydrocarbonPV;
155 attr.rv += fs.Rv().value() * hydrocarbonPV;
158 attr.pressure += fs.pressure(FluidSystem::oilPhaseIdx).value() * hydrocarbonPV;
159 attr.temperature += fs.temperature(FluidSystem::oilPhaseIdx).value() * hydrocarbonPV;
162 attr.pressure += fs.pressure(FluidSystem::gasPhaseIdx).value() * hydrocarbonPV;
163 attr.temperature += fs.temperature(FluidSystem::gasPhaseIdx).value() * hydrocarbonPV;
165 attr.saltConcentration += fs.saltConcentration().value() * hydrocarbonPV;
169 auto& attr = attributes_pv[reg];
172 attr.rs += fs.Rs().value() * pv_cell;
173 attr.rv += fs.Rv().value() * pv_cell;
176 attr.pressure += fs.pressure(FluidSystem::oilPhaseIdx).value() * pv_cell;
177 attr.temperature += fs.temperature(FluidSystem::oilPhaseIdx).value() * pv_cell;
179 attr.pressure += fs.pressure(FluidSystem::gasPhaseIdx).value() * pv_cell;
180 attr.temperature += fs.temperature(FluidSystem::gasPhaseIdx).value() * pv_cell;
183 attr.pressure += fs.pressure(FluidSystem::waterPhaseIdx).value() * pv_cell;
184 attr.temperature += fs.temperature(FluidSystem::waterPhaseIdx).value() * pv_cell;
186 attr.saltConcentration += fs.saltConcentration().value() * pv_cell;
190 OPM_END_PARALLEL_TRY_CATCH(
"SurfaceToReservoirVoidage::defineState() failed: ", simulator.vanguard().grid().comm());
192 for (
const auto& reg : rmap_.activeRegions()) {
193 auto& ra = attr_.attributes(reg);
194 const double hpv_sum = comm.sum(attributes_hpv[reg].pv);
197 const auto& attri_hpv = attributes_hpv[reg];
198 const double p_hpv_sum = comm.sum(attri_hpv.pressure);
199 const double T_hpv_sum = comm.sum(attri_hpv.temperature);
200 const double rs_hpv_sum = comm.sum(attri_hpv.rs);
201 const double rv_hpv_sum = comm.sum(attri_hpv.rv);
202 const double sc_hpv_sum = comm.sum(attri_hpv.saltConcentration);
204 ra.pressure = p_hpv_sum / hpv_sum;
205 ra.temperature = T_hpv_sum / hpv_sum;
206 ra.rs = rs_hpv_sum / hpv_sum;
207 ra.rv = rv_hpv_sum / hpv_sum;
209 ra.saltConcentration = sc_hpv_sum / hpv_sum;
212 const auto& attri_pv = attributes_pv[reg];
213 const double pv_sum = comm.sum(attri_pv.pv);
215 const double p_pv_sum = comm.sum(attri_pv.pressure);
216 const double T_pv_sum = comm.sum(attri_pv.temperature);
217 const double rs_pv_sum = comm.sum(attri_pv.rs);
218 const double rv_pv_sum = comm.sum(attri_pv.rv);
219 const double sc_pv_sum = comm.sum(attri_pv.saltConcentration);
221 ra.pressure = p_pv_sum / pv_sum;
222 ra.temperature = T_pv_sum / pv_sum;
223 ra.rs = rs_pv_sum / pv_sum;
224 ra.rv = rv_pv_sum / pv_sum;
226 ra.saltConcentration = sc_pv_sum / pv_sum;
236 typedef typename RegionMapping<Region>::RegionId
RegionId;
265 template <
class Coeff>
269 const auto& pu = phaseUsage_;
270 const auto& ra = attr_.attributes(r);
271 const double p = ra.pressure;
272 const double T = ra.temperature;
273 const double saltConcentration = ra.saltConcentration;
279 std::fill(& coeff[0], & coeff[0] + phaseUsage_.num_phases, 0.0);
284 const double bw = FluidSystem::waterPvt().inverseFormationVolumeFactor(pvtRegionIdx, T, p, saltConcentration);
286 coeff[iw] = 1.0 / bw;
294 const double detR = 1.0 - (Rs * Rv);
299 const double bo = FluidSystem::oilPvt().inverseFormationVolumeFactor(pvtRegionIdx, T, p, Rs);
300 const double den = bo * detR;
302 coeff[io] += 1.0 / den;
305 coeff[ig] -= ra.rv / den;
312 const double bg = FluidSystem::gasPvt().inverseFormationVolumeFactor(pvtRegionIdx, T, p, Rv, 0.0 );
313 const double den = bg * detR;
315 coeff[ig] += 1.0 / den;
318 coeff[io] -= ra.rs / den;
323 template <
class Coeff>
325 calcInjCoeff(
const RegionId r,
const int pvtRegionIdx, Coeff& coeff)
const
327 const auto& pu = phaseUsage_;
328 const auto& ra = attr_.attributes(r);
329 const double p = ra.pressure;
330 const double T = ra.temperature;
331 const double saltConcentration = ra.saltConcentration;
337 std::fill(& coeff[0], & coeff[0] + phaseUsage_.num_phases, 0.0);
342 const double bw = FluidSystem::waterPvt().inverseFormationVolumeFactor(pvtRegionIdx, T, p, saltConcentration);
344 coeff[iw] = 1.0 / bw;
348 const double bo = FluidSystem::oilPvt().inverseFormationVolumeFactor(pvtRegionIdx, T, p, 0.0);
349 coeff[io] += 1.0 / bo;
353 const double bg = FluidSystem::gasPvt().inverseFormationVolumeFactor(pvtRegionIdx, T, p, 0.0, 0.0);
354 coeff[ig] += 1.0 / bg;
376 template <
class Rates >
379 Rates& voidage_rates)
const
381 assert(voidage_rates.size() == surface_rates.size());
383 std::fill(voidage_rates.begin(), voidage_rates.end(), 0.0);
385 const auto& pu = phaseUsage_;
386 const auto& ra = attr_.attributes(r);
387 const double p = ra.pressure;
388 const double T = ra.temperature;
389 const double saltConcentration = ra.saltConcentration;
398 const double bw = FluidSystem::waterPvt().inverseFormationVolumeFactor(pvtRegionIdx, T, p, saltConcentration);
400 voidage_rates[iw] = surface_rates[iw] / bw;
406 if (io >= 0 && ig >= 0) {
407 b = surface_rates[ig]/(surface_rates[io]+1.0e-15);
410 double Rs = std::min(a, b);
414 if (io >= 0 && ig >= 0) {
415 b = surface_rates[io]/(surface_rates[ig]+1.0e-15);
418 double Rv = std::min(a, b);
421 const double detR = 1.0 - (Rs * Rv);
426 const double bo = FluidSystem::oilPvt().inverseFormationVolumeFactor(pvtRegionIdx, T, p, Rs);
427 const double den = bo * detR;
429 voidage_rates[io] = surface_rates[io];
432 voidage_rates[io] -= Rv * surface_rates[ig];
435 voidage_rates[io] /= den;
441 const double bg = FluidSystem::gasPvt().inverseFormationVolumeFactor(pvtRegionIdx, T, p, Rv, 0.0 );
442 const double den = bg * detR;
444 voidage_rates[ig] = surface_rates[ig];
447 voidage_rates[ig] -= Rs * surface_rates[io];
450 voidage_rates[ig] /= den;
467 template <
class SolventModule>
471 const auto& ra = attr_.attributes(r);
472 const double p = ra.pressure;
473 const double T = ra.temperature;
474 const auto& solventPvt = SolventModule::solventPvt();
475 const double bs = solventPvt.inverseFormationVolumeFactor(pvtRegionIdx, T, p);
488 const RegionMapping<Region> rmap_;
500 , saltConcentration(0.0)
508 double saltConcentration;
511 RegionAttributeHelpers::RegionAttributes<RegionId, Attributes> attr_;
Convert component rates at surface conditions to phase (voidage) rates at reservoir conditions.
Definition: RateConverter.hpp:68
RegionMapping< Region >::RegionId RegionId
Region identifier.
Definition: RateConverter.hpp:236
SurfaceToReservoirVoidage(const PhaseUsage &phaseUsage, const Region ®ion)
Constructor.
Definition: RateConverter.hpp:77
void calcCoeffSolvent(const RegionId r, const int pvtRegionIdx, double &coeff) const
Compute coefficients for surface-to-reservoir voidage conversion for solvent.
Definition: RateConverter.hpp:469
void calcReservoirVoidageRates(const RegionId r, const int pvtRegionIdx, const Rates &surface_rates, Rates &voidage_rates) const
Converting surface volume rates to reservoir voidage rates.
Definition: RateConverter.hpp:378
void calcCoeff(const RegionId r, const int pvtRegionIdx, Coeff &coeff) const
Compute coefficients for surface-to-reservoir voidage conversion.
Definition: RateConverter.hpp:267
void defineState(const EbosSimulator &simulator)
Compute pore volume averaged hydrocarbon state pressure, rs and rv.
Definition: RateConverter.hpp:95
int gas(const PhaseUsage &pu)
Numerical ID of active gas phase.
Definition: RegionAttributeHelpers.hpp:394
int oil(const PhaseUsage &pu)
Numerical ID of active oil phase.
Definition: RegionAttributeHelpers.hpp:374
int water(const PhaseUsage &pu)
Numerical ID of active water phase.
Definition: RegionAttributeHelpers.hpp:354
bool water(const PhaseUsage &pu)
Active water predicate.
Definition: RegionAttributeHelpers.hpp:308
bool oil(const PhaseUsage &pu)
Active oil predicate.
Definition: RegionAttributeHelpers.hpp:321
bool gas(const PhaseUsage &pu)
Active gas predicate.
Definition: RegionAttributeHelpers.hpp:334
This file contains a set of helper functions used by VFPProd / VFPInj.
Definition: BlackoilPhases.hpp:27
PhaseUsage phaseUsage(const Phases &phases)
Determine the active phases.
Definition: phaseUsageFromDeck.cpp:37
Definition: BlackoilPhases.hpp:46