RateConverter.hpp
Go to the documentation of this file.
1/*
2 Copyright 2014, 2015 SINTEF ICT, Applied Mathematics.
3 Copyright 2014, 2015 Statoil ASA.
4 Copyright 2017, IRIS
5
6 This file is part of the Open Porous Media Project (OPM).
7
8 OPM is free software: you can redistribute it and/or modify
9 it under the terms of the GNU General Public License as published by
10 the Free Software Foundation, either version 3 of the License, or
11 (at your option) any later version.
12
13 OPM is distributed in the hope that it will be useful,
14 but WITHOUT ANY WARRANTY; without even the implied warranty of
15 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
16 GNU General Public License for more details.
17
18 You should have received a copy of the GNU General Public License
19 along with OPM. If not, see <http://www.gnu.org/licenses/>.
20*/
21
22#ifndef OPM_RATECONVERTER_HPP_HEADER_INCLUDED
23#define OPM_RATECONVERTER_HPP_HEADER_INCLUDED
24
25#include <opm/grid/utility/RegionMapping.hpp>
26
28
31
32#include <dune/grid/common/gridenums.hh>
33#include <dune/grid/common/rangegenerators.hh>
34
35#include <array>
36#include <cassert>
37#include <unordered_map>
38#include <utility>
39#include <vector>
40
51namespace Opm {
52 namespace RateConverter {
53
69 template <class FluidSystem, class Region>
71 public:
72 using Scalar = typename FluidSystem::Scalar;
79 SurfaceToReservoirVoidage(const Region& region)
80 : rmap_ (region)
81 , attr_ (rmap_, Attributes())
82 {}
83
92 template <typename ElementContext, class Simulator>
93 void defineState(const Simulator& simulator)
94 {
95 // create map from cell to region and set all attributes to
96 // zero
97 for (const auto& reg : rmap_.activeRegions()) {
98 auto& ra = attr_.attributes(reg);
99 ra.pressure = 0.0;
100 ra.temperature = 0.0;
101 ra.rs = 0.0;
102 ra.rv = 0.0;
103 ra.pv = 0.0;
104 ra.saltConcentration = 0.0;
105 ra.rsw = 0.0;
106 ra.rvw = 0.0;
107 }
108
109 // quantities for pore volume average
110 std::unordered_map<RegionId, Attributes> attributes_pv;
111
112 // quantities for hydrocarbon volume average
113 std::unordered_map<RegionId, Attributes> attributes_hpv;
114
115 for (const auto& reg : rmap_.activeRegions()) {
116 attributes_pv.insert({reg, Attributes()});
117 attributes_hpv.insert({reg, Attributes()});
118 }
119
120 ElementContext elemCtx( simulator );
121 const auto& gridView = simulator.gridView();
122 const auto& comm = gridView.comm();
123
125 for (const auto& elem : elements(gridView, Dune::Partitions::interior)) {
126 elemCtx.updatePrimaryStencil(elem);
127 elemCtx.updatePrimaryIntensiveQuantities(/*timeIdx=*/0);
128 const unsigned cellIdx = elemCtx.globalSpaceIndex(/*spaceIdx=*/0, /*timeIdx=*/0);
129 const auto& intQuants = elemCtx.intensiveQuantities(/*spaceIdx=*/0, /*timeIdx=*/0);
130 const auto& fs = intQuants.fluidState();
131 // use pore volume weighted averages.
132 const Scalar pv_cell =
133 simulator.model().dofTotalVolume(simulator.vanguard().gridEquilIdxToGridIdx(cellIdx))
134 * intQuants.porosity().value();
135
136 // only count oil and gas filled parts of the domain
137 Scalar hydrocarbon = 1.0;
138 if (FluidSystem::phaseIsActive(FluidSystem::waterPhaseIdx)) {
139 hydrocarbon -= fs.saturation(FluidSystem::waterPhaseIdx).value();
140 }
141
142 const int reg = rmap_.region(cellIdx);
143 assert(reg >= 0);
144
145 // sum p, rs, rv, and T.
146 const Scalar hydrocarbonPV = pv_cell*hydrocarbon;
147 if (hydrocarbonPV > 0.) {
148 auto& attr = attributes_hpv[reg];
149 attr.pv += hydrocarbonPV;
150 if (FluidSystem::phaseIsActive(FluidSystem::oilPhaseIdx) && FluidSystem::phaseIsActive(FluidSystem::gasPhaseIdx)) {
151 attr.rs += fs.Rs().value() * hydrocarbonPV;
152 attr.rv += fs.Rv().value() * hydrocarbonPV;
153 }
154 if (FluidSystem::phaseIsActive(FluidSystem::oilPhaseIdx)) {
155 attr.pressure += fs.pressure(FluidSystem::oilPhaseIdx).value() * hydrocarbonPV;
156 attr.temperature += fs.temperature(FluidSystem::oilPhaseIdx).value() * hydrocarbonPV;
157 } else {
158 assert(FluidSystem::phaseIsActive(FluidSystem::gasPhaseIdx));
159 attr.pressure += fs.pressure(FluidSystem::gasPhaseIdx).value() * hydrocarbonPV;
160 attr.temperature += fs.temperature(FluidSystem::gasPhaseIdx).value() * hydrocarbonPV;
161 }
162 attr.saltConcentration += fs.saltConcentration().value() * hydrocarbonPV;
163 if (FluidSystem::enableDissolvedGasInWater()) {
164 attr.rsw += fs.Rsw().value() * hydrocarbonPV; // scale with total volume?
165 }
166 if (FluidSystem::enableVaporizedWater()) {
167 attr.rvw += fs.Rvw().value() * hydrocarbonPV; // scale with total volume?
168 }
169 }
170
171 if (pv_cell > 0.) {
172 auto& attr = attributes_pv[reg];
173 attr.pv += pv_cell;
174 if (FluidSystem::phaseIsActive(FluidSystem::oilPhaseIdx) && FluidSystem::phaseIsActive(FluidSystem::gasPhaseIdx)) {
175 attr.rs += fs.Rs().value() * pv_cell;
176 attr.rv += fs.Rv().value() * pv_cell;
177 }
178 if (FluidSystem::phaseIsActive(FluidSystem::oilPhaseIdx)) {
179 attr.pressure += fs.pressure(FluidSystem::oilPhaseIdx).value() * pv_cell;
180 attr.temperature += fs.temperature(FluidSystem::oilPhaseIdx).value() * pv_cell;
181 } else if (FluidSystem::phaseIsActive(FluidSystem::gasPhaseIdx)) {
182 attr.pressure += fs.pressure(FluidSystem::gasPhaseIdx).value() * pv_cell;
183 attr.temperature += fs.temperature(FluidSystem::gasPhaseIdx).value() * pv_cell;
184 } else {
185 assert(FluidSystem::phaseIsActive(FluidSystem::waterPhaseIdx));
186 attr.pressure += fs.pressure(FluidSystem::waterPhaseIdx).value() * pv_cell;
187 attr.temperature += fs.temperature(FluidSystem::waterPhaseIdx).value() * pv_cell;
188 }
189 attr.saltConcentration += fs.saltConcentration().value() * pv_cell;
190 if (FluidSystem::enableDissolvedGasInWater()) {
191 attr.rsw += fs.Rsw().value() * pv_cell;
192 }
193 if (FluidSystem::enableVaporizedWater()) {
194 attr.rvw += fs.Rvw().value() * pv_cell;
195 }
196 }
197 }
198
199 OPM_END_PARALLEL_TRY_CATCH("SurfaceToReservoirVoidage::defineState() failed: ", simulator.vanguard().grid().comm());
200
201 this->sumRates(attributes_hpv,
202 attributes_pv,
203 comm);
204 }
205
212
240 template <class Coeff>
241 void
242 calcCoeff(const RegionId r, const int pvtRegionIdx, Coeff& coeff) const;
243
244 template <class Coeff , class Rates>
245 void
246 calcCoeff(const RegionId r, const int pvtRegionIdx, const Rates& surface_rates, Coeff& coeff) const;
247
248 template <class Coeff>
249 void
250 calcCoeff( const int pvtRegionIdx,
251 const Scalar p,
252 const Scalar rs,
253 const Scalar rv,
254 const Scalar rsw,
255 const Scalar rvw,
256 const Scalar T,
257 const Scalar saltConcentration,
258 Coeff& coeff) const;
259
260 template <class Coeff>
261 void
262 calcInjCoeff(const RegionId r, const int pvtRegionIdx, Coeff& coeff) const;
263
286 template <class Rates>
288 const int pvtRegionIdx,
289 const Rates& surface_rates,
290 Rates& voidage_rates) const;
291
326 template <typename SurfaceRates, typename VoidageRates>
327 void calcReservoirVoidageRates(const int pvtRegionIdx,
328 const Scalar p,
329 const Scalar rs,
330 const Scalar rv,
331 const Scalar rsw,
332 const Scalar rvw,
333 const Scalar T,
334 const Scalar saltConcentration,
335 const SurfaceRates& surface_rates,
336 VoidageRates& voidage_rates) const;
337
338 template <class Rates>
339 std::pair<Scalar, Scalar>
341 const Scalar rvMax,
342 const Rates& surface_rates) const;
343
356 template <class SolventModule>
357 void
358 calcCoeffSolvent(const RegionId r, const int pvtRegionIdx, Scalar& coeff) const
359 {
360 const auto& ra = attr_.attributes(r);
361 const Scalar p = ra.pressure;
362 const Scalar T = ra.temperature;
363 const auto& solventPvt = SolventModule::solventPvt();
364 const Scalar bs = solventPvt.inverseFormationVolumeFactor(pvtRegionIdx, T, p);
365 coeff = 1.0 / bs;
366 }
367
368 private:
372 const RegionMapping<Region> rmap_;
373
377 struct Attributes {
378 Attributes()
379 : data{0.0}
380 , pressure(data[0])
381 , temperature(data[1])
382 , rs(data[2])
383 , rv(data[3])
384 , rsw(data[4])
385 , rvw(data[5])
386 , pv(data[6])
387 , saltConcentration(data[7])
388 {}
389
390 Attributes(const Attributes& rhs)
391 : Attributes()
392 {
393 this->data = rhs.data;
394 }
395
396 Attributes& operator=(const Attributes& rhs)
397 {
398 this->data = rhs.data;
399 return *this;
400 }
401
402 std::array<Scalar,8> data;
403 Scalar& pressure;
404 Scalar& temperature;
405 Scalar& rs;
406 Scalar& rv;
407 Scalar& rsw;
408 Scalar& rvw;
409 Scalar& pv;
410 Scalar& saltConcentration;
411 };
412
413 void sumRates(std::unordered_map<RegionId,Attributes>& attributes_hpv,
414 std::unordered_map<RegionId,Attributes>& attributes_pv,
416
417 RegionAttributeHelpers::RegionAttributes<RegionId, Attributes> attr_;
418 };
419
420 } // namespace RateConverter
421} // namespace Opm
422
423#endif /* OPM_RATECONVERTER_HPP_HEADER_INCLUDED */
#define OPM_END_PARALLEL_TRY_CATCH(prefix, comm)
Catch exception and throw in a parallel try-catch clause.
Definition: DeferredLoggingErrorHelpers.hpp:192
#define OPM_BEGIN_PARALLEL_TRY_CATCH()
Macro to setup the try of a parallel try-catch.
Definition: DeferredLoggingErrorHelpers.hpp:158
Definition: RateConverter.hpp:70
std::pair< Scalar, Scalar > inferDissolvedVaporisedRatio(const Scalar rsMax, const Scalar rvMax, const Rates &surface_rates) const
SurfaceToReservoirVoidage(const Region &region)
Definition: RateConverter.hpp:79
void calcReservoirVoidageRates(const RegionId r, const int pvtRegionIdx, const Rates &surface_rates, Rates &voidage_rates) const
void calcCoeff(const RegionId r, const int pvtRegionIdx, const Rates &surface_rates, Coeff &coeff) const
void calcInjCoeff(const RegionId r, const int pvtRegionIdx, Coeff &coeff) const
void calcCoeff(const RegionId r, const int pvtRegionIdx, Coeff &coeff) const
void calcCoeff(const int pvtRegionIdx, const Scalar p, const Scalar rs, const Scalar rv, const Scalar rsw, const Scalar rvw, const Scalar T, const Scalar saltConcentration, Coeff &coeff) const
void defineState(const Simulator &simulator)
Definition: RateConverter.hpp:93
typename RegionMapping< Region >::RegionId RegionId
Definition: RateConverter.hpp:211
void calcReservoirVoidageRates(const int pvtRegionIdx, const Scalar p, const Scalar rs, const Scalar rv, const Scalar rsw, const Scalar rvw, const Scalar T, const Scalar saltConcentration, const SurfaceRates &surface_rates, VoidageRates &voidage_rates) const
void calcCoeffSolvent(const RegionId r, const int pvtRegionIdx, Scalar &coeff) const
Definition: RateConverter.hpp:358
typename FluidSystem::Scalar Scalar
Definition: RateConverter.hpp:72
Manages the initializing and running of time dependent problems.
Definition: simulator.hh:84
Vanguard & vanguard()
Return a reference to the grid manager of simulation.
Definition: simulator.hh:234
const GridView & gridView() const
Return the grid view for which the simulation is done.
Definition: simulator.hh:246
Model & model()
Return the physical model used in the simulation.
Definition: simulator.hh:252
Dune::Communication< MPIComm > Communication
Definition: ParallelCommunication.hpp:30
Definition: blackoilboundaryratevector.hh:39
int RegionId
Definition: WellConstraints.hpp:37