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/linalg/ParallelIstlInformation.hpp>
28 #include <opm/simulators/wells/RegionAttributeHelpers.hpp>
29 #include <opm/simulators/utils/DeferredLoggingErrorHelpers.hpp>
30 #include <dune/grid/common/gridenums.hh>
35 #include <type_traits>
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();
125 OPM_BEGIN_PARALLEL_TRY_CATCH();
127 const auto& elemEndIt = gridView.template end<0>();
128 for (
auto elemIt = gridView.template begin</*codim=*/0>();
133 const auto& elem = *elemIt;
134 if (elem.partitionType() != Dune::InteriorEntity)
137 elemCtx.updatePrimaryStencil(elem);
138 elemCtx.updatePrimaryIntensiveQuantities(0);
139 const unsigned cellIdx = elemCtx.globalSpaceIndex(0, 0);
140 const auto& intQuants = elemCtx.intensiveQuantities(0, 0);
141 const auto& fs = intQuants.fluidState();
143 const double pv_cell =
144 simulator.model().dofTotalVolume(cellIdx)
145 * intQuants.porosity().value();
148 double hydrocarbon = 1.0;
149 const auto& pu = phaseUsage_;
151 hydrocarbon -= fs.saturation(FluidSystem::waterPhaseIdx).value();
154 const int reg = rmap_.region(cellIdx);
158 const double hydrocarbonPV = pv_cell*hydrocarbon;
159 if (hydrocarbonPV > 0.) {
160 auto& attr = attributes_hpv[reg];
161 attr.pv += hydrocarbonPV;
163 attr.rs += fs.Rs().value() * hydrocarbonPV;
164 attr.rv += fs.Rv().value() * hydrocarbonPV;
167 attr.pressure += fs.pressure(FluidSystem::oilPhaseIdx).value() * hydrocarbonPV;
168 attr.temperature += fs.temperature(FluidSystem::oilPhaseIdx).value() * hydrocarbonPV;
171 attr.pressure += fs.pressure(FluidSystem::gasPhaseIdx).value() * hydrocarbonPV;
172 attr.temperature += fs.temperature(FluidSystem::gasPhaseIdx).value() * hydrocarbonPV;
174 attr.saltConcentration += fs.saltConcentration().value() * hydrocarbonPV;
178 auto& attr = attributes_pv[reg];
181 attr.rs += fs.Rs().value() * pv_cell;
182 attr.rv += fs.Rv().value() * pv_cell;
185 attr.pressure += fs.pressure(FluidSystem::oilPhaseIdx).value() * pv_cell;
186 attr.temperature += fs.temperature(FluidSystem::oilPhaseIdx).value() * pv_cell;
188 attr.pressure += fs.pressure(FluidSystem::gasPhaseIdx).value() * pv_cell;
189 attr.temperature += fs.temperature(FluidSystem::gasPhaseIdx).value() * pv_cell;
192 attr.pressure += fs.pressure(FluidSystem::waterPhaseIdx).value() * pv_cell;
193 attr.temperature += fs.temperature(FluidSystem::waterPhaseIdx).value() * pv_cell;
195 attr.saltConcentration += fs.saltConcentration().value() * pv_cell;
199 OPM_END_PARALLEL_TRY_CATCH(
"SurfaceToReservoirVoidage::defineState() failed: ", simulator.vanguard().grid().comm());
201 for (
const auto& reg : rmap_.activeRegions()) {
202 auto& ra = attr_.attributes(reg);
203 const double hpv_sum = comm.sum(attributes_hpv[reg].pv);
206 const auto& attri_hpv = attributes_hpv[reg];
207 const double p_hpv_sum = comm.sum(attri_hpv.pressure);
208 const double T_hpv_sum = comm.sum(attri_hpv.temperature);
209 const double rs_hpv_sum = comm.sum(attri_hpv.rs);
210 const double rv_hpv_sum = comm.sum(attri_hpv.rv);
211 const double sc_hpv_sum = comm.sum(attri_hpv.saltConcentration);
213 ra.pressure = p_hpv_sum / hpv_sum;
214 ra.temperature = T_hpv_sum / hpv_sum;
215 ra.rs = rs_hpv_sum / hpv_sum;
216 ra.rv = rv_hpv_sum / hpv_sum;
218 ra.saltConcentration = sc_hpv_sum / hpv_sum;
221 const auto& attri_pv = attributes_pv[reg];
222 const double pv_sum = comm.sum(attri_pv.pv);
224 const double p_pv_sum = comm.sum(attri_pv.pressure);
225 const double T_pv_sum = comm.sum(attri_pv.temperature);
226 const double rs_pv_sum = comm.sum(attri_pv.rs);
227 const double rv_pv_sum = comm.sum(attri_pv.rv);
228 const double sc_pv_sum = comm.sum(attri_pv.saltConcentration);
230 ra.pressure = p_pv_sum / pv_sum;
231 ra.temperature = T_pv_sum / pv_sum;
232 ra.rs = rs_pv_sum / pv_sum;
233 ra.rv = rv_pv_sum / pv_sum;
235 ra.saltConcentration = sc_pv_sum / pv_sum;
245 typedef typename RegionMapping<Region>::RegionId
RegionId;
274 template <
class Coeff>
278 const auto& pu = phaseUsage_;
279 const auto& ra = attr_.attributes(r);
280 const double p = ra.pressure;
281 const double T = ra.temperature;
282 const double saltConcentration = ra.saltConcentration;
288 std::fill(& coeff[0], & coeff[0] + phaseUsage_.num_phases, 0.0);
293 const double bw = FluidSystem::waterPvt().inverseFormationVolumeFactor(pvtRegionIdx, T, p, saltConcentration);
295 coeff[iw] = 1.0 / bw;
303 const double detR = 1.0 - (Rs * Rv);
308 const double bo = FluidSystem::oilPvt().inverseFormationVolumeFactor(pvtRegionIdx, T, p, Rs);
309 const double den = bo * detR;
311 coeff[io] += 1.0 / den;
314 coeff[ig] -= ra.rv / den;
321 const double bg = FluidSystem::gasPvt().inverseFormationVolumeFactor(pvtRegionIdx, T, p, Rv);
322 const double den = bg * detR;
324 coeff[ig] += 1.0 / den;
327 coeff[io] -= ra.rs / den;
332 template <
class Coeff>
334 calcInjCoeff(
const RegionId r,
const int pvtRegionIdx, Coeff& coeff)
const
336 const auto& pu = phaseUsage_;
337 const auto& ra = attr_.attributes(r);
338 const double p = ra.pressure;
339 const double T = ra.temperature;
340 const double saltConcentration = ra.saltConcentration;
346 std::fill(& coeff[0], & coeff[0] + phaseUsage_.num_phases, 0.0);
351 const double bw = FluidSystem::waterPvt().inverseFormationVolumeFactor(pvtRegionIdx, T, p, saltConcentration);
353 coeff[iw] = 1.0 / bw;
357 const double bo = FluidSystem::oilPvt().inverseFormationVolumeFactor(pvtRegionIdx, T, p, 0.0);
358 coeff[io] += 1.0 / bo;
362 const double bg = FluidSystem::gasPvt().inverseFormationVolumeFactor(pvtRegionIdx, T, p, 0.0);
363 coeff[ig] += 1.0 / bg;
385 template <
class Rates >
388 Rates& voidage_rates)
const
390 assert(voidage_rates.size() == surface_rates.size());
392 std::fill(voidage_rates.begin(), voidage_rates.end(), 0.0);
394 const auto& pu = phaseUsage_;
395 const auto& ra = attr_.attributes(r);
396 const double p = ra.pressure;
397 const double T = ra.temperature;
398 const double saltConcentration = ra.saltConcentration;
407 const double bw = FluidSystem::waterPvt().inverseFormationVolumeFactor(pvtRegionIdx, T, p, saltConcentration);
409 voidage_rates[iw] = surface_rates[iw] / bw;
415 if (io >= 0 && ig >= 0) {
416 b = surface_rates[ig]/(surface_rates[io]+1.0e-15);
419 double Rs = std::min(a, b);
423 if (io >= 0 && ig >= 0) {
424 b = surface_rates[io]/(surface_rates[ig]+1.0e-15);
427 double Rv = std::min(a, b);
430 const double detR = 1.0 - (Rs * Rv);
435 const double bo = FluidSystem::oilPvt().inverseFormationVolumeFactor(pvtRegionIdx, T, p, Rs);
436 const double den = bo * detR;
438 voidage_rates[io] = surface_rates[io];
441 voidage_rates[io] -= Rv * surface_rates[ig];
444 voidage_rates[io] /= den;
450 const double bg = FluidSystem::gasPvt().inverseFormationVolumeFactor(pvtRegionIdx, T, p, Rv);
451 const double den = bg * detR;
453 voidage_rates[ig] = surface_rates[ig];
456 voidage_rates[ig] -= Rs * surface_rates[io];
459 voidage_rates[ig] /= den;
476 template <
class SolventModule>
480 const auto& ra = attr_.attributes(r);
481 const double p = ra.pressure;
482 const double T = ra.temperature;
483 const auto& solventPvt = SolventModule::solventPvt();
484 const double bs = solventPvt.inverseFormationVolumeFactor(pvtRegionIdx, T, p);
497 const RegionMapping<Region> rmap_;
509 , saltConcentration(0.0)
517 double saltConcentration;
520 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:245
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:478
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:387
void calcCoeff(const RegionId r, const int pvtRegionIdx, Coeff &coeff) const
Compute coefficients for surface-to-reservoir voidage conversion.
Definition: RateConverter.hpp:276
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:395
int oil(const PhaseUsage &pu)
Numerical ID of active oil phase.
Definition: RegionAttributeHelpers.hpp:375
int water(const PhaseUsage &pu)
Numerical ID of active water phase.
Definition: RegionAttributeHelpers.hpp:355
bool water(const PhaseUsage &pu)
Active water predicate.
Definition: RegionAttributeHelpers.hpp:309
bool oil(const PhaseUsage &pu)
Active oil predicate.
Definition: RegionAttributeHelpers.hpp:322
bool gas(const PhaseUsage &pu)
Active gas predicate.
Definition: RegionAttributeHelpers.hpp:335
This file contains a set of helper functions used by VFPProd / VFPInj.
Definition: BlackoilPhases.hpp:26
PhaseUsage phaseUsage(const Phases &phases)
Determine the active phases.
Definition: phaseUsageFromDeck.cpp:33
Definition: BlackoilPhases.hpp:45