blackoilpolymermodules.hh
Go to the documentation of this file.
1// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
2// vi: set et ts=4 sw=4 sts=4:
3/*
4 This file is part of the Open Porous Media project (OPM).
5
6 OPM is free software: you can redistribute it and/or modify
7 it under the terms of the GNU General Public License as published by
8 the Free Software Foundation, either version 2 of the License, or
9 (at your option) any later version.
10
11 OPM is distributed in the hope that it will be useful,
12 but WITHOUT ANY WARRANTY; without even the implied warranty of
13 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14 GNU General Public License for more details.
15
16 You should have received a copy of the GNU General Public License
17 along with OPM. If not, see <http://www.gnu.org/licenses/>.
18
19 Consult the COPYING file in the top-level source directory of this
20 module for the precise wording of the license and the list of
21 copyright holders.
22*/
28#ifndef EWOMS_BLACK_OIL_POLYMER_MODULE_HH
29#define EWOMS_BLACK_OIL_POLYMER_MODULE_HH
30
31#include "blackoilproperties.hh"
32
35
36#include <opm/common/OpmLog/OpmLog.hpp>
37
38#if HAVE_ECL_INPUT
39#include <opm/input/eclipse/EclipseState/EclipseState.hpp>
40#include <opm/input/eclipse/EclipseState/Tables/PlyadsTable.hpp>
41#include <opm/input/eclipse/EclipseState/Tables/PlymaxTable.hpp>
42#include <opm/input/eclipse/EclipseState/Tables/PlyrockTable.hpp>
43#include <opm/input/eclipse/EclipseState/Tables/PlyshlogTable.hpp>
44#include <opm/input/eclipse/EclipseState/Tables/PlyviscTable.hpp>
45#endif
46
47#include <dune/common/fvector.hh>
48
49#include <cmath>
50#include <stdexcept>
51#include <string>
52
53namespace Opm {
59template <class TypeTag, bool enablePolymerV = getPropValue<TypeTag, Properties::EnablePolymer>()>
61{
74
75 using Toolbox = MathToolbox<Evaluation>;
76
77 using TabulatedFunction = typename BlackOilPolymerParams<Scalar>::TabulatedFunction;
78 using TabulatedTwoDFunction = typename BlackOilPolymerParams<Scalar>::TabulatedTwoDFunction;
79
80 static constexpr unsigned polymerConcentrationIdx = Indices::polymerConcentrationIdx;
81 static constexpr unsigned polymerMoleWeightIdx = Indices::polymerMoleWeightIdx;
82 static constexpr unsigned contiPolymerEqIdx = Indices::contiPolymerEqIdx;
83 static constexpr unsigned contiPolymerMolarWeightEqIdx = Indices::contiPolymerMWEqIdx;
84 static constexpr unsigned waterPhaseIdx = FluidSystem::waterPhaseIdx;
85
86
87 static constexpr unsigned enablePolymer = enablePolymerV;
88 static constexpr bool enablePolymerMolarWeight = getPropValue<TypeTag, Properties::EnablePolymerMW>();
89
90 static constexpr unsigned numEq = getPropValue<TypeTag, Properties::NumEq>();
91 static constexpr unsigned numPhases = FluidSystem::numPhases;
92
93public:
94#if HAVE_ECL_INPUT
98 static void initFromState(const EclipseState& eclState)
99 {
100 // some sanity checks: if polymers are enabled, the POLYMER keyword must be
101 // present, if polymers are disabled the keyword must not be present.
102 if (enablePolymer && !eclState.runspec().phases().active(Phase::POLYMER)) {
103 throw std::runtime_error("Non-trivial polymer treatment requested at compile time, but "
104 "the deck does not contain the POLYMER keyword");
105 }
106 else if (!enablePolymer && eclState.runspec().phases().active(Phase::POLYMER)) {
107 throw std::runtime_error("Polymer treatment disabled at compile time, but the deck "
108 "contains the POLYMER keyword");
109 }
110
111 if (enablePolymerMolarWeight && !eclState.runspec().phases().active(Phase::POLYMW)) {
112 throw std::runtime_error("Polymer molecular weight tracking is enabled at compile time, but "
113 "the deck does not contain the POLYMW keyword");
114 }
115 else if (!enablePolymerMolarWeight && eclState.runspec().phases().active(Phase::POLYMW)) {
116 throw std::runtime_error("Polymer molecular weight tracking is disabled at compile time, but the deck "
117 "contains the POLYMW keyword");
118 }
119
120 if (enablePolymerMolarWeight && !enablePolymer) {
121 throw std::runtime_error("Polymer molecular weight tracking is enabled while polymer treatment "
122 "is disabled at compile time");
123 }
124
125 if (!eclState.runspec().phases().active(Phase::POLYMER))
126 return; // polymer treatment is supposed to be disabled
127
128 const auto& tableManager = eclState.getTableManager();
129
130 unsigned numSatRegions = tableManager.getTabdims().getNumSatTables();
131 params_.setNumSatRegions(numSatRegions);
132
133 // initialize the objects which deal with the PLYROCK keyword
134 const auto& plyrockTables = tableManager.getPlyrockTables();
135 if (!plyrockTables.empty()) {
136 assert(numSatRegions == plyrockTables.size());
137 for (unsigned satRegionIdx = 0; satRegionIdx < numSatRegions; ++ satRegionIdx) {
138 const auto& plyrockTable = plyrockTables.template getTable<PlyrockTable>(satRegionIdx);
139 params_.setPlyrock(satRegionIdx,
140 plyrockTable.getDeadPoreVolumeColumn()[0],
141 plyrockTable.getResidualResistanceFactorColumn()[0],
142 plyrockTable.getRockDensityFactorColumn()[0],
143 static_cast<typename BlackOilPolymerParams<Scalar>::AdsorptionBehaviour>(plyrockTable.getAdsorbtionIndexColumn()[0]),
144 plyrockTable.getMaxAdsorbtionColumn()[0]);
145 }
146 }
147 else {
148 throw std::runtime_error("PLYROCK must be specified in POLYMER runs\n");
149 }
150
151 // initialize the objects which deal with the PLYADS keyword
152 const auto& plyadsTables = tableManager.getPlyadsTables();
153 if (!plyadsTables.empty()) {
154 assert(numSatRegions == plyadsTables.size());
155 for (unsigned satRegionIdx = 0; satRegionIdx < numSatRegions; ++ satRegionIdx) {
156 const auto& plyadsTable = plyadsTables.template getTable<PlyadsTable>(satRegionIdx);
157 // Copy data
158 const auto& c = plyadsTable.getPolymerConcentrationColumn();
159 const auto& ads = plyadsTable.getAdsorbedPolymerColumn();
160 params_.plyadsAdsorbedPolymer_[satRegionIdx].setXYContainers(c, ads);
161 }
162 }
163 else {
164 throw std::runtime_error("PLYADS must be specified in POLYMER runs\n");
165 }
166
167
168 unsigned numPvtRegions = tableManager.getTabdims().getNumPVTTables();
169 params_.plyviscViscosityMultiplierTable_.resize(numPvtRegions);
170
171 // initialize the objects which deal with the PLYVISC keyword
172 const auto& plyviscTables = tableManager.getPlyviscTables();
173 if (!plyviscTables.empty()) {
174 // different viscosity model is used for POLYMW
175 if (enablePolymerMolarWeight) {
176 OpmLog::warning("PLYVISC should not be used in POLYMW runs, "
177 "it will have no effect. A viscosity model based on PLYVMH is used instead.\n");
178 }
179 else {
180 assert(numPvtRegions == plyviscTables.size());
181 for (unsigned pvtRegionIdx = 0; pvtRegionIdx < numPvtRegions; ++ pvtRegionIdx) {
182 const auto& plyadsTable = plyviscTables.template getTable<PlyviscTable>(pvtRegionIdx);
183 // Copy data
184 const auto& c = plyadsTable.getPolymerConcentrationColumn();
185 const auto& visc = plyadsTable.getViscosityMultiplierColumn();
186 params_.plyviscViscosityMultiplierTable_[pvtRegionIdx].setXYContainers(c, visc);
187 }
188 }
189 }
190 else if (!enablePolymerMolarWeight) {
191 throw std::runtime_error("PLYVISC must be specified in POLYMER runs\n");
192 }
193
194 // initialize the objects which deal with the PLYMAX keyword
195 const auto& plymaxTables = tableManager.getPlymaxTables();
196 const unsigned numMixRegions = plymaxTables.size();
197 params_.setNumMixRegions(numMixRegions, enablePolymerMolarWeight);
198 if (!plymaxTables.empty()) {
199 for (unsigned mixRegionIdx = 0; mixRegionIdx < numMixRegions; ++ mixRegionIdx) {
200 const auto& plymaxTable = plymaxTables.template getTable<PlymaxTable>(mixRegionIdx);
201 params_.plymaxMaxConcentration_[mixRegionIdx] = plymaxTable.getPolymerConcentrationColumn()[0];
202 }
203 }
204 else {
205 throw std::runtime_error("PLYMAX must be specified in POLYMER runs\n");
206 }
207
208 if (!eclState.getTableManager().getPlmixparTable().empty()) {
209 if (enablePolymerMolarWeight) {
210 OpmLog::warning("PLMIXPAR should not be used in POLYMW runs, it will have no effect.\n");
211 }
212 else {
213 const auto& plmixparTable = eclState.getTableManager().getPlmixparTable();
214 // initialize the objects which deal with the PLMIXPAR keyword
215 for (unsigned mixRegionIdx = 0; mixRegionIdx < numMixRegions; ++ mixRegionIdx) {
216 params_.plymixparToddLongstaff_[mixRegionIdx] = plmixparTable[mixRegionIdx].todd_langstaff;
217 }
218 }
219 }
220 else if (!enablePolymerMolarWeight) {
221 throw std::runtime_error("PLMIXPAR must be specified in POLYMER runs\n");
222 }
223
224 params_.hasPlyshlog_ = eclState.getTableManager().hasTables("PLYSHLOG");
225 params_.hasShrate_ = eclState.getTableManager().useShrate();
226
227 if ((params_.hasPlyshlog_ || params_.hasShrate_) && enablePolymerMolarWeight) {
228 OpmLog::warning("PLYSHLOG and SHRATE should not be used in POLYMW runs, they will have no effect.\n");
229 }
230
231 if (params_.hasPlyshlog_ && !enablePolymerMolarWeight) {
232 const auto& plyshlogTables = tableManager.getPlyshlogTables();
233 assert(numPvtRegions == plyshlogTables.size());
234 params_.plyshlogShearEffectRefMultiplier_.resize(numPvtRegions);
235 params_.plyshlogShearEffectRefLogVelocity_.resize(numPvtRegions);
236 for (unsigned pvtRegionIdx = 0; pvtRegionIdx < numPvtRegions; ++ pvtRegionIdx) {
237 const auto& plyshlogTable = plyshlogTables.template getTable<PlyshlogTable>(pvtRegionIdx);
238
239 Scalar plyshlogRefPolymerConcentration = plyshlogTable.getRefPolymerConcentration();
240 auto waterVelocity = plyshlogTable.getWaterVelocityColumn().vectorCopy();
241 auto shearMultiplier = plyshlogTable.getShearMultiplierColumn().vectorCopy();
242
243 // do the unit version here for the waterVelocity
244 UnitSystem unitSystem = eclState.getDeckUnitSystem();
245 double siFactor = params_.hasShrate_? unitSystem.parse("1/Time").getSIScaling() : unitSystem.parse("Length/Time").getSIScaling();
246 for (size_t i = 0; i < waterVelocity.size(); ++i) {
247 waterVelocity[i] *= siFactor;
248 // for plyshlog the input must be stored as logarithms
249 // the interpolation is then done the log-space.
250 waterVelocity[i] = std::log(waterVelocity[i]);
251 }
252
253 Scalar refViscMult = params_.plyviscViscosityMultiplierTable_[pvtRegionIdx].eval(plyshlogRefPolymerConcentration, /*extrapolate=*/true);
254 // convert the table using referece conditions
255 for (size_t i = 0; i < waterVelocity.size(); ++i) {
256 shearMultiplier[i] *= refViscMult;
257 shearMultiplier[i] -= 1;
258 shearMultiplier[i] /= (refViscMult - 1);
259 shearMultiplier[i] = shearMultiplier[i];
260 }
261 params_.plyshlogShearEffectRefMultiplier_[pvtRegionIdx].resize(waterVelocity.size());
262 params_.plyshlogShearEffectRefLogVelocity_[pvtRegionIdx].resize(waterVelocity.size());
263
264 for (size_t i = 0; i < waterVelocity.size(); ++i) {
265 params_.plyshlogShearEffectRefMultiplier_[pvtRegionIdx][i] = shearMultiplier[i];
266 params_.plyshlogShearEffectRefLogVelocity_[pvtRegionIdx][i] = waterVelocity[i];
267 }
268 }
269 }
270
271 if (params_.hasShrate_ && !enablePolymerMolarWeight) {
272 if (!params_.hasPlyshlog_) {
273 throw std::runtime_error("PLYSHLOG must be specified if SHRATE is used in POLYMER runs\n");
274 }
275 const auto& shrateTable = eclState.getTableManager().getShrateTable();
276 params_.shrate_.resize(numPvtRegions);
277 for (unsigned pvtRegionIdx = 0; pvtRegionIdx < numPvtRegions; ++ pvtRegionIdx) {
278 if (shrateTable.empty()) {
279 params_.shrate_[pvtRegionIdx] = 4.8; //default;
280 }
281 else if (shrateTable.size() == numPvtRegions) {
282 params_.shrate_[pvtRegionIdx] = shrateTable[pvtRegionIdx].rate;
283 }
284 else {
285 throw std::runtime_error("SHRATE must either have 0 or number of NUMPVT entries\n");
286 }
287 }
288 }
289
290 if constexpr (enablePolymerMolarWeight) {
291 const auto& plyvmhTable = eclState.getTableManager().getPlyvmhTable();
292 if (!plyvmhTable.empty()) {
293 assert(plyvmhTable.size() == numMixRegions);
294 for (size_t regionIdx = 0; regionIdx < numMixRegions; ++regionIdx) {
295 params_.plyvmhCoefficients_[regionIdx].k_mh = plyvmhTable[regionIdx].k_mh;
296 params_.plyvmhCoefficients_[regionIdx].a_mh = plyvmhTable[regionIdx].a_mh;
297 params_.plyvmhCoefficients_[regionIdx].gamma = plyvmhTable[regionIdx].gamma;
298 params_.plyvmhCoefficients_[regionIdx].kappa = plyvmhTable[regionIdx].kappa;
299 }
300 }
301 else {
302 throw std::runtime_error("PLYVMH keyword must be specified in POLYMW rus \n");
303 }
304
305 // handling PLYMWINJ keyword
306 const auto& plymwinjTables = tableManager.getPlymwinjTables();
307 for (const auto& table : plymwinjTables) {
308 const int tableNumber = table.first;
309 const auto& plymwinjtable = table.second;
310 const std::vector<double>& throughput = plymwinjtable.getThroughputs();
311 const std::vector<double>& watervelocity = plymwinjtable.getVelocities();
312 const std::vector<std::vector<double>>& molecularweight = plymwinjtable.getMoleWeights();
313 TabulatedTwoDFunction tablefunc(throughput, watervelocity, molecularweight, true, false);
314 params_.plymwinjTables_[tableNumber] = std::move(tablefunc);
315 }
316
317 // handling SKPRWAT keyword
318 const auto& skprwatTables = tableManager.getSkprwatTables();
319 for (const auto& table : skprwatTables) {
320 const int tableNumber = table.first;
321 const auto& skprwattable = table.second;
322 const std::vector<double>& throughput = skprwattable.getThroughputs();
323 const std::vector<double>& watervelocity = skprwattable.getVelocities();
324 const std::vector<std::vector<double>>& skinpressure = skprwattable.getSkinPressures();
325 TabulatedTwoDFunction tablefunc(throughput, watervelocity, skinpressure, true, false);
326 params_.skprwatTables_[tableNumber] = std::move(tablefunc);
327 }
328
329 // handling SKPRPOLY keyword
330 const auto& skprpolyTables = tableManager.getSkprpolyTables();
331 for (const auto& table : skprpolyTables) {
332 const int tableNumber = table.first;
333 const auto& skprpolytable = table.second;
334 const std::vector<double>& throughput = skprpolytable.getThroughputs();
335 const std::vector<double>& watervelocity = skprpolytable.getVelocities();
336 const std::vector<std::vector<double>>& skinpressure = skprpolytable.getSkinPressures();
337 const double refPolymerConcentration = skprpolytable.referenceConcentration();
339 {refPolymerConcentration,
340 TabulatedTwoDFunction(throughput, watervelocity, skinpressure, true, false)};
341 params_.skprpolyTables_[tableNumber] = std::move(tablefunc);
342 }
343 }
344 }
345#endif
346
350 static TabulatedTwoDFunction& getPlymwinjTable(const int tableNumber)
351 {
352 const auto iterTable = params_.plymwinjTables_.find(tableNumber);
353 if (iterTable != params_.plymwinjTables_.end()) {
354 return iterTable->second;
355 }
356 else {
357 throw std::runtime_error(" the PLYMWINJ table " + std::to_string(tableNumber) + " does not exist\n");
358 }
359 }
360
364 static TabulatedTwoDFunction& getSkprwatTable(const int tableNumber)
365 {
366 const auto iterTable = params_.skprwatTables_.find(tableNumber);
367 if (iterTable != params_.skprwatTables_.end()) {
368 return iterTable->second;
369 }
370 else {
371 throw std::runtime_error(" the SKPRWAT table " + std::to_string(tableNumber) + " does not exist\n");
372 }
373 }
374
379 getSkprpolyTable(const int tableNumber)
380 {
381 const auto iterTable = params_.skprpolyTables_.find(tableNumber);
382 if (iterTable != params_.skprpolyTables_.end()) {
383 return iterTable->second;
384 }
385 else {
386 throw std::runtime_error(" the SKPRPOLY table " + std::to_string(tableNumber) + " does not exist\n");
387 }
388 }
389
393 static void registerParameters()
394 {
395 if constexpr (enablePolymer)
397 }
398
402 static void registerOutputModules(Model& model,
403 Simulator& simulator)
404 {
405 if constexpr (enablePolymer)
406 model.addOutputModule(new VtkBlackOilPolymerModule<TypeTag>(simulator));
407 }
408
409 static bool primaryVarApplies(unsigned pvIdx)
410 {
411 if constexpr (enablePolymer) {
412 if constexpr (enablePolymerMolarWeight)
413 return pvIdx == polymerConcentrationIdx || pvIdx == polymerMoleWeightIdx;
414 else
415 return pvIdx == polymerConcentrationIdx;
416 }
417 else
418 return false;
419 }
420
421 static std::string primaryVarName(unsigned pvIdx)
422 {
423 assert(primaryVarApplies(pvIdx));
424
425 if (pvIdx == polymerConcentrationIdx) {
426 return "polymer_waterconcentration";
427 }
428 else {
429 return "polymer_molecularweight";
430 }
431 }
432
433 static Scalar primaryVarWeight([[maybe_unused]] unsigned pvIdx)
434 {
435 assert(primaryVarApplies(pvIdx));
436
437 // TODO: it may be beneficial to chose this differently.
438 return static_cast<Scalar>(1.0);
439 }
440
441 static bool eqApplies(unsigned eqIdx)
442 {
443 if constexpr (enablePolymer) {
444 if constexpr (enablePolymerMolarWeight)
445 return eqIdx == contiPolymerEqIdx || eqIdx == contiPolymerMolarWeightEqIdx;
446 else
447 return eqIdx == contiPolymerEqIdx;
448 }
449 else
450 return false;
451 }
452
453 static std::string eqName(unsigned eqIdx)
454 {
455 assert(eqApplies(eqIdx));
456
457 if (eqIdx == contiPolymerEqIdx)
458 return "conti^polymer";
459 else
460 return "conti^polymer_molecularweight";
461 }
462
463 static Scalar eqWeight([[maybe_unused]] unsigned eqIdx)
464 {
465 assert(eqApplies(eqIdx));
466
467 // TODO: it may be beneficial to chose this differently.
468 return static_cast<Scalar>(1.0);
469 }
470
471 // must be called after water storage is computed
472 template <class LhsEval>
473 static void addStorage(Dune::FieldVector<LhsEval, numEq>& storage,
474 const IntensiveQuantities& intQuants)
475 {
476 if constexpr (enablePolymer) {
477 const auto& fs = intQuants.fluidState();
478
479 LhsEval surfaceVolumeWater =
480 Toolbox::template decay<LhsEval>(fs.saturation(waterPhaseIdx))
481 * Toolbox::template decay<LhsEval>(fs.invB(waterPhaseIdx))
482 * Toolbox::template decay<LhsEval>(intQuants.porosity());
483
484 // avoid singular matrix if no water is present.
485 surfaceVolumeWater = max(surfaceVolumeWater, 1e-10);
486
487 // polymer in water phase
488 const LhsEval massPolymer = surfaceVolumeWater
489 * Toolbox::template decay<LhsEval>(intQuants.polymerConcentration())
490 * (1.0 - Toolbox::template decay<LhsEval>(intQuants.polymerDeadPoreVolume()));
491
492 // polymer in solid phase
493 const LhsEval adsorptionPolymer =
494 Toolbox::template decay<LhsEval>(1.0 - intQuants.porosity())
495 * Toolbox::template decay<LhsEval>(intQuants.polymerRockDensity())
496 * Toolbox::template decay<LhsEval>(intQuants.polymerAdsorption());
497
498 LhsEval accumulationPolymer = massPolymer + adsorptionPolymer;
499
500 storage[contiPolymerEqIdx] += accumulationPolymer;
501
502 // tracking the polymer molecular weight
503 if constexpr (enablePolymerMolarWeight) {
504 accumulationPolymer = max(accumulationPolymer, 1e-10);
505
506 storage[contiPolymerMolarWeightEqIdx] += accumulationPolymer
507 * Toolbox::template decay<LhsEval> (intQuants.polymerMoleWeight());
508 }
509 }
510 }
511
512 static void computeFlux([[maybe_unused]] RateVector& flux,
513 [[maybe_unused]] const ElementContext& elemCtx,
514 [[maybe_unused]] unsigned scvfIdx,
515 [[maybe_unused]] unsigned timeIdx)
516
517 {
518 if constexpr (enablePolymer) {
519 const auto& extQuants = elemCtx.extensiveQuantities(scvfIdx, timeIdx);
520
521 const unsigned upIdx = extQuants.upstreamIndex(FluidSystem::waterPhaseIdx);
522 const unsigned inIdx = extQuants.interiorIndex();
523 const auto& up = elemCtx.intensiveQuantities(upIdx, timeIdx);
524 const unsigned contiWaterEqIdx = Indices::conti0EqIdx + Indices::canonicalToActiveComponentIndex(FluidSystem::waterCompIdx);
525
526 if (upIdx == inIdx) {
527 flux[contiPolymerEqIdx] =
528 extQuants.volumeFlux(waterPhaseIdx)
529 *up.fluidState().invB(waterPhaseIdx)
530 *up.polymerViscosityCorrection()
531 /extQuants.polymerShearFactor()
532 *up.polymerConcentration();
533
534 // modify water
535 flux[contiWaterEqIdx] /=
536 extQuants.waterShearFactor();
537 }
538 else {
539 flux[contiPolymerEqIdx] =
540 extQuants.volumeFlux(waterPhaseIdx)
541 *decay<Scalar>(up.fluidState().invB(waterPhaseIdx))
542 *decay<Scalar>(up.polymerViscosityCorrection())
543 /decay<Scalar>(extQuants.polymerShearFactor())
544 *decay<Scalar>(up.polymerConcentration());
545
546 // modify water
547 flux[contiWaterEqIdx] /=
548 decay<Scalar>(extQuants.waterShearFactor());
549 }
550
551 // flux related to transport of polymer molecular weight
552 if constexpr (enablePolymerMolarWeight) {
553 if (upIdx == inIdx)
554 flux[contiPolymerMolarWeightEqIdx] =
555 flux[contiPolymerEqIdx]*up.polymerMoleWeight();
556 else
557 flux[contiPolymerMolarWeightEqIdx] =
558 flux[contiPolymerEqIdx]*decay<Scalar>(up.polymerMoleWeight());
559 }
560 }
561 }
562
566 static Scalar computeUpdateError(const PrimaryVariables&,
567 const EqVector&)
568 {
569 // do not consider consider the cange of polymer primary variables for
570 // convergence
571 // TODO: maybe this should be changed
572 return static_cast<Scalar>(0.0);
573 }
574
575 template <class DofEntity>
576 static void serializeEntity(const Model& model, std::ostream& outstream, const DofEntity& dof)
577 {
578 if constexpr (enablePolymer) {
579 unsigned dofIdx = model.dofMapper().index(dof);
580 const PrimaryVariables& priVars = model.solution(/*timeIdx=*/0)[dofIdx];
581 outstream << priVars[polymerConcentrationIdx];
582 outstream << priVars[polymerMoleWeightIdx];
583 }
584 }
585
586 template <class DofEntity>
587 static void deserializeEntity(Model& model, std::istream& instream, const DofEntity& dof)
588 {
589 if constexpr (enablePolymer) {
590 unsigned dofIdx = model.dofMapper().index(dof);
591 PrimaryVariables& priVars0 = model.solution(/*timeIdx=*/0)[dofIdx];
592 PrimaryVariables& priVars1 = model.solution(/*timeIdx=*/1)[dofIdx];
593
594 instream >> priVars0[polymerConcentrationIdx];
595 instream >> priVars0[polymerMoleWeightIdx];
596
597 // set the primary variables for the beginning of the current time step.
598 priVars1[polymerConcentrationIdx] = priVars0[polymerConcentrationIdx];
599 priVars1[polymerMoleWeightIdx] = priVars0[polymerMoleWeightIdx];
600 }
601 }
602
603 static const Scalar plyrockDeadPoreVolume(const ElementContext& elemCtx,
604 unsigned scvIdx,
605 unsigned timeIdx)
606 {
607 unsigned satnumRegionIdx = elemCtx.problem().satnumRegionIndex(elemCtx, scvIdx, timeIdx);
608 return params_.plyrockDeadPoreVolume_[satnumRegionIdx];
609 }
610
611 static const Scalar plyrockResidualResistanceFactor(const ElementContext& elemCtx,
612 unsigned scvIdx,
613 unsigned timeIdx)
614 {
615 unsigned satnumRegionIdx = elemCtx.problem().satnumRegionIndex(elemCtx, scvIdx, timeIdx);
616 return params_.plyrockResidualResistanceFactor_[satnumRegionIdx];
617 }
618
619 static const Scalar plyrockRockDensityFactor(const ElementContext& elemCtx,
620 unsigned scvIdx,
621 unsigned timeIdx)
622 {
623 unsigned satnumRegionIdx = elemCtx.problem().satnumRegionIndex(elemCtx, scvIdx, timeIdx);
624 return params_.plyrockRockDensityFactor_[satnumRegionIdx];
625 }
626
627 static const Scalar plyrockAdsorbtionIndex(const ElementContext& elemCtx,
628 unsigned scvIdx,
629 unsigned timeIdx)
630 {
631 unsigned satnumRegionIdx = elemCtx.problem().satnumRegionIndex(elemCtx, scvIdx, timeIdx);
632 return params_.plyrockAdsorbtionIndex_[satnumRegionIdx];
633 }
634
635 static const Scalar plyrockMaxAdsorbtion(const ElementContext& elemCtx,
636 unsigned scvIdx,
637 unsigned timeIdx)
638 {
639 unsigned satnumRegionIdx = elemCtx.problem().satnumRegionIndex(elemCtx, scvIdx, timeIdx);
640 return params_.plyrockMaxAdsorbtion_[satnumRegionIdx];
641 }
642
643 static const TabulatedFunction& plyadsAdsorbedPolymer(const ElementContext& elemCtx,
644 unsigned scvIdx,
645 unsigned timeIdx)
646 {
647 unsigned satnumRegionIdx = elemCtx.problem().satnumRegionIndex(elemCtx, scvIdx, timeIdx);
648 return params_.plyadsAdsorbedPolymer_[satnumRegionIdx];
649 }
650
651 static const TabulatedFunction& plyviscViscosityMultiplierTable(const ElementContext& elemCtx,
652 unsigned scvIdx,
653 unsigned timeIdx)
654 {
655 unsigned pvtnumRegionIdx = elemCtx.problem().pvtRegionIndex(elemCtx, scvIdx, timeIdx);
656 return params_.plyviscViscosityMultiplierTable_[pvtnumRegionIdx];
657 }
658
659 static const TabulatedFunction& plyviscViscosityMultiplierTable(unsigned pvtnumRegionIdx)
660 {
661 return params_.plyviscViscosityMultiplierTable_[pvtnumRegionIdx];
662 }
663
664 static const Scalar plymaxMaxConcentration(const ElementContext& elemCtx,
665 unsigned scvIdx,
666 unsigned timeIdx)
667 {
668 unsigned polymerMixRegionIdx = elemCtx.problem().plmixnumRegionIndex(elemCtx, scvIdx, timeIdx);
669 return params_.plymaxMaxConcentration_[polymerMixRegionIdx];
670 }
671
672 static const Scalar plymixparToddLongstaff(const ElementContext& elemCtx,
673 unsigned scvIdx,
674 unsigned timeIdx)
675 {
676 unsigned polymerMixRegionIdx = elemCtx.problem().plmixnumRegionIndex(elemCtx, scvIdx, timeIdx);
677 return params_.plymixparToddLongstaff_[polymerMixRegionIdx];
678 }
679
681 plyvmhCoefficients(const ElementContext& elemCtx,
682 const unsigned scvIdx,
683 const unsigned timeIdx)
684 {
685 const unsigned polymerMixRegionIdx = elemCtx.problem().plmixnumRegionIndex(elemCtx, scvIdx, timeIdx);
686 return params_.plyvmhCoefficients_[polymerMixRegionIdx];
687 }
688
689 static bool hasPlyshlog()
690 {
691 return params_.hasPlyshlog_;
692 }
693
694 static bool hasShrate()
695 {
696 return params_.hasShrate_;
697 }
698
699 static const Scalar shrate(unsigned pvtnumRegionIdx)
700 {
701 return params_.shrate_[pvtnumRegionIdx];
702 }
703
710 template <class Evaluation>
711 static Evaluation computeShearFactor(const Evaluation& polymerConcentration,
712 unsigned pvtnumRegionIdx,
713 const Evaluation& v0)
714 {
715 using ToolboxLocal = MathToolbox<Evaluation>;
716
717 const auto& viscosityMultiplierTable = params_.plyviscViscosityMultiplierTable_[pvtnumRegionIdx];
718 Scalar viscosityMultiplier = viscosityMultiplierTable.eval(scalarValue(polymerConcentration), /*extrapolate=*/true);
719
720 const Scalar eps = 1e-14;
721 // return 1.0 if the polymer has no effect on the water.
722 if (std::abs((viscosityMultiplier - 1.0)) < eps)
723 return ToolboxLocal::createConstant(v0, 1.0);
724
725 const std::vector<Scalar>& shearEffectRefLogVelocity = params_.plyshlogShearEffectRefLogVelocity_[pvtnumRegionIdx];
726 auto v0AbsLog = log(abs(v0));
727 // return 1.0 if the velocity /sharte is smaller than the first velocity entry.
728 if (v0AbsLog < shearEffectRefLogVelocity[0])
729 return ToolboxLocal::createConstant(v0, 1.0);
730
731 // compute shear factor from input
732 // Z = (1 + (P - 1) * M(v)) / P
733 // where M(v) is computed from user input
734 // and P = viscosityMultiplier
735 const std::vector<Scalar>& shearEffectRefMultiplier = params_.plyshlogShearEffectRefMultiplier_[pvtnumRegionIdx];
736 size_t numTableEntries = shearEffectRefLogVelocity.size();
737 assert(shearEffectRefMultiplier.size() == numTableEntries);
738
739 std::vector<Scalar> shearEffectMultiplier(numTableEntries, 1.0);
740 for (size_t i = 0; i < numTableEntries; ++i) {
741 shearEffectMultiplier[i] = (1.0 + (viscosityMultiplier - 1.0)*shearEffectRefMultiplier[i]) / viscosityMultiplier;
742 shearEffectMultiplier[i] = log(shearEffectMultiplier[i]);
743 }
744 // store the logarithmic velocity and logarithmic multipliers in a table for easy look up and
745 // linear interpolation in the logarithmic space.
746 TabulatedFunction logShearEffectMultiplier = TabulatedFunction(numTableEntries, shearEffectRefLogVelocity, shearEffectMultiplier, /*bool sortInputs =*/ false);
747
748 // Find sheared velocity (v) that satisfies
749 // F = log(v) + log (Z) - log(v0) = 0;
750
751 // Set up the function
752 // u = log(v)
753 auto F = [&logShearEffectMultiplier, &v0AbsLog](const Evaluation& u) {
754 return u + logShearEffectMultiplier.eval(u, true) - v0AbsLog;
755 };
756 // and its derivative
757 auto dF = [&logShearEffectMultiplier](const Evaluation& u) {
758 return 1 + logShearEffectMultiplier.evalDerivative(u, true);
759 };
760
761 // Solve F = 0 using Newton
762 // Use log(v0) as initial value for u
763 auto u = v0AbsLog;
764 bool converged = false;
765 // TODO make this into parameters
766 for (int i = 0; i < 20; ++i) {
767 auto f = F(u);
768 auto df = dF(u);
769 u -= f/df;
770 if (std::abs(scalarValue(f)) < 1e-12) {
771 converged = true;
772 break;
773 }
774 }
775 if (!converged) {
776 throw std::runtime_error("Not able to compute shear velocity. \n");
777 }
778
779 // return the shear factor
780 return exp(logShearEffectMultiplier.eval(u, /*extrapolate=*/true));
781
782 }
783
784 const Scalar molarMass() const
785 {
786 return 0.25; // kg/mol
787 }
788
789private:
790 static BlackOilPolymerParams<Scalar> params_;
791};
792
793
794template <class TypeTag, bool enablePolymerV>
795BlackOilPolymerParams<typename BlackOilPolymerModule<TypeTag, enablePolymerV>::Scalar>
796BlackOilPolymerModule<TypeTag, enablePolymerV>::params_;
797
805template <class TypeTag, bool enablePolymerV = getPropValue<TypeTag, Properties::EnablePolymer>()>
807{
809
817
819
820 enum { numPhases = getPropValue<TypeTag, Properties::NumPhases>() };
821 static constexpr int polymerConcentrationIdx = Indices::polymerConcentrationIdx;
822 static constexpr int waterPhaseIdx = FluidSystem::waterPhaseIdx;
823 static constexpr bool enablePolymerMolarWeight = getPropValue<TypeTag, Properties::EnablePolymerMW>();
824 static constexpr int polymerMoleWeightIdx = Indices::polymerMoleWeightIdx;
825
826
827public:
828
834 void polymerPropertiesUpdate_(const ElementContext& elemCtx,
835 unsigned dofIdx,
836 unsigned timeIdx)
837 {
838 const auto linearizationType = elemCtx.linearizationType();
839 const PrimaryVariables& priVars = elemCtx.primaryVars(dofIdx, timeIdx);
840 polymerConcentration_ = priVars.makeEvaluation(polymerConcentrationIdx, timeIdx, linearizationType);
841 if constexpr (enablePolymerMolarWeight) {
842 polymerMoleWeight_ = priVars.makeEvaluation(polymerMoleWeightIdx, timeIdx, linearizationType);
843 }
844
845 // permeability reduction due to polymer
846 const Scalar& maxAdsorbtion = PolymerModule::plyrockMaxAdsorbtion(elemCtx, dofIdx, timeIdx);
847 const auto& plyadsAdsorbedPolymer = PolymerModule::plyadsAdsorbedPolymer(elemCtx, dofIdx, timeIdx);
848 polymerAdsorption_ = plyadsAdsorbedPolymer.eval(polymerConcentration_, /*extrapolate=*/true);
850 const Scalar& maxPolymerAdsorption = elemCtx.problem().maxPolymerAdsorption(elemCtx, dofIdx, timeIdx);
851 polymerAdsorption_ = std::max(Evaluation(maxPolymerAdsorption) , polymerAdsorption_);
852 }
853
854 // compute resitanceFactor
855 const Scalar& residualResistanceFactor = PolymerModule::plyrockResidualResistanceFactor(elemCtx, dofIdx, timeIdx);
856 const Evaluation resistanceFactor = 1.0 + (residualResistanceFactor - 1.0) * polymerAdsorption_ / maxAdsorbtion;
857
858 // compute effective viscosities
859 if constexpr (!enablePolymerMolarWeight) {
860 const Scalar cmax = PolymerModule::plymaxMaxConcentration(elemCtx, dofIdx, timeIdx);
861 const auto& fs = asImp_().fluidState_;
862 const Evaluation& muWater = fs.viscosity(waterPhaseIdx);
863 const auto& viscosityMultiplier = PolymerModule::plyviscViscosityMultiplierTable(elemCtx, dofIdx, timeIdx);
864 const Evaluation viscosityMixture = viscosityMultiplier.eval(polymerConcentration_, /*extrapolate=*/true) * muWater;
865
866 // Do the Todd-Longstaff mixing
867 const Scalar plymixparToddLongstaff = PolymerModule::plymixparToddLongstaff(elemCtx, dofIdx, timeIdx);
868 const Evaluation viscosityPolymer = viscosityMultiplier.eval(cmax, /*extrapolate=*/true) * muWater;
869 const Evaluation viscosityPolymerEffective = pow(viscosityMixture, plymixparToddLongstaff) * pow(viscosityPolymer, 1.0 - plymixparToddLongstaff);
870 const Evaluation viscosityWaterEffective = pow(viscosityMixture, plymixparToddLongstaff) * pow(muWater, 1.0 - plymixparToddLongstaff);
871
872 const Evaluation cbar = polymerConcentration_ / cmax;
873 // waterViscosity / effectiveWaterViscosity
874 waterViscosityCorrection_ = muWater * ((1.0 - cbar) / viscosityWaterEffective + cbar / viscosityPolymerEffective);
875 // effectiveWaterViscosity / effectivePolymerViscosity
876 polymerViscosityCorrection_ = (muWater / waterViscosityCorrection_) / viscosityPolymerEffective;
877 }
878 else { // based on PLYVMH
879 const auto& plyvmhCoefficients = PolymerModule::plyvmhCoefficients(elemCtx, dofIdx, timeIdx);
880 const Scalar k_mh = plyvmhCoefficients.k_mh;
881 const Scalar a_mh = plyvmhCoefficients.a_mh;
882 const Scalar gamma = plyvmhCoefficients.gamma;
883 const Scalar kappa = plyvmhCoefficients.kappa;
884
885 // viscosity model based on Mark-Houwink equation and Huggins equation
886 // 1000 is a emperical constant, most likely related to unit conversion
887 const Evaluation intrinsicViscosity = k_mh * pow(polymerMoleWeight_ * 1000., a_mh);
888 const Evaluation x = polymerConcentration_ * intrinsicViscosity;
889 waterViscosityCorrection_ = 1.0 / (1.0 + gamma * (x + kappa * x * x));
891 }
892
893 // adjust water mobility
894 asImp_().mobility_[waterPhaseIdx] *= waterViscosityCorrection_ / resistanceFactor;
895
896 // update rock properties
899 }
900
901 const Evaluation& polymerConcentration() const
902 { return polymerConcentration_; }
903
904 const Evaluation& polymerMoleWeight() const
905 {
906 if constexpr (enablePolymerMolarWeight)
907 return polymerMoleWeight_;
908 else
909 throw std::logic_error("polymerMoleWeight() is called but polymer milecular weight is disabled");
910 }
911
912 const Scalar& polymerDeadPoreVolume() const
913 { return polymerDeadPoreVolume_; }
914
915 const Evaluation& polymerAdsorption() const
916 { return polymerAdsorption_; }
917
918 const Scalar& polymerRockDensity() const
919 { return polymerRockDensity_; }
920
921 // effectiveWaterViscosity / effectivePolymerViscosity
922 const Evaluation& polymerViscosityCorrection() const
924
925 // waterViscosity / effectiveWaterViscosity
926 const Evaluation& waterViscosityCorrection() const
927 { return waterViscosityCorrection_; }
928
929
930protected:
931 Implementation& asImp_()
932 { return *static_cast<Implementation*>(this); }
933
935 // polymer molecular weight
942
943
944};
945
946template <class TypeTag>
948{
952
953public:
954 void polymerPropertiesUpdate_(const ElementContext&,
955 unsigned,
956 unsigned)
957 { }
958
959 const Evaluation& polymerMoleWeight() const
960 { throw std::logic_error("polymerMoleWeight() called but polymer molecular weight is disabled"); }
961
962 const Evaluation& polymerConcentration() const
963 { throw std::runtime_error("polymerConcentration() called but polymers are disabled"); }
964
965 const Evaluation& polymerDeadPoreVolume() const
966 { throw std::runtime_error("polymerDeadPoreVolume() called but polymers are disabled"); }
967
968 const Evaluation& polymerAdsorption() const
969 { throw std::runtime_error("polymerAdsorption() called but polymers are disabled"); }
970
971 const Evaluation& polymerRockDensity() const
972 { throw std::runtime_error("polymerRockDensity() called but polymers are disabled"); }
973
974 const Evaluation& polymerViscosityCorrection() const
975 { throw std::runtime_error("polymerViscosityCorrection() called but polymers are disabled"); }
976
977 const Evaluation& waterViscosityCorrection() const
978 { throw std::runtime_error("waterViscosityCorrection() called but polymers are disabled"); }
979};
980
981
989template <class TypeTag, bool enablePolymerV = getPropValue<TypeTag, Properties::EnablePolymer>()>
991{
993
1001
1002 static constexpr unsigned gasPhaseIdx = FluidSystem::gasPhaseIdx;
1003 static constexpr int dimWorld = GridView::dimensionworld;
1004 static constexpr unsigned waterPhaseIdx = FluidSystem::waterPhaseIdx;
1005
1006 using Toolbox = MathToolbox<Evaluation>;
1008 using DimVector = Dune::FieldVector<Scalar, dimWorld>;
1009 using DimEvalVector = Dune::FieldVector<Evaluation, dimWorld>;
1010
1011public:
1018 template <class Dummy = bool> // we need to make this method a template to avoid
1019 // compiler errors if it is not instantiated!
1020 void updateShearMultipliersPerm(const ElementContext&,
1021 unsigned,
1022 unsigned)
1023 {
1024 throw std::runtime_error("The extension of the blackoil model for polymers is not yet "
1025 "implemented for problems specified using permeabilities.");
1026 }
1027
1034 template <class Dummy = bool> // we need to make this method a template to avoid
1035 // compiler errors if it is not instantiated!
1036 void updateShearMultipliers(const ElementContext& elemCtx,
1037 unsigned scvfIdx,
1038 unsigned timeIdx)
1039 {
1040
1041 waterShearFactor_ = 1.0;
1042 polymerShearFactor_ = 1.0;
1043
1045 return;
1046
1047 const ExtensiveQuantities& extQuants = asImp_();
1048 unsigned upIdx = extQuants.upstreamIndex(waterPhaseIdx);
1049 unsigned interiorDofIdx = extQuants.interiorIndex();
1050 unsigned exteriorDofIdx = extQuants.exteriorIndex();
1051 const auto& up = elemCtx.intensiveQuantities(upIdx, timeIdx);
1052 const auto& intQuantsIn = elemCtx.intensiveQuantities(interiorDofIdx, timeIdx);
1053 const auto& intQuantsEx = elemCtx.intensiveQuantities(exteriorDofIdx, timeIdx);
1054
1055 // compute water velocity from flux
1056 Evaluation poroAvg = intQuantsIn.porosity()*0.5 + intQuantsEx.porosity()*0.5;
1057 unsigned pvtnumRegionIdx = elemCtx.problem().pvtRegionIndex(elemCtx, scvfIdx, timeIdx);
1058 const Evaluation& Sw = up.fluidState().saturation(waterPhaseIdx);
1059 unsigned cellIdx = elemCtx.globalSpaceIndex(scvfIdx, timeIdx);
1060 const auto& materialLawManager = elemCtx.problem().materialLawManager();
1061 const auto& scaledDrainageInfo =
1062 materialLawManager->oilWaterScaledEpsInfoDrainage(cellIdx);
1063 const Scalar& Swcr = scaledDrainageInfo.Swcr;
1064
1065 // guard against zero porosity and no mobile water
1066 Evaluation denom = max(poroAvg * (Sw - Swcr), 1e-12);
1067 Evaluation waterVolumeVelocity = extQuants.volumeFlux(waterPhaseIdx) / denom;
1068
1069 // if shrate is specified. Compute shrate based on the water velocity
1071 const Evaluation& relWater = up.relativePermeability(waterPhaseIdx);
1072 Scalar trans = elemCtx.problem().transmissibility(elemCtx, interiorDofIdx, exteriorDofIdx);
1073 if (trans > 0.0) {
1074 Scalar faceArea = elemCtx.stencil(timeIdx).interiorFace(scvfIdx).area();
1075 auto dist = elemCtx.pos(interiorDofIdx, timeIdx) - elemCtx.pos(exteriorDofIdx, timeIdx);
1076 // compute permeability from transmissibility.
1077 Scalar absPerm = trans / faceArea * dist.two_norm();
1078 waterVolumeVelocity *=
1079 PolymerModule::shrate(pvtnumRegionIdx)*sqrt(poroAvg*Sw / (relWater*absPerm));
1080 assert(isfinite(waterVolumeVelocity));
1081 }
1082 }
1083
1084 // compute share factors for water and polymer
1085 waterShearFactor_ =
1086 PolymerModule::computeShearFactor(up.polymerConcentration(),
1087 pvtnumRegionIdx,
1088 waterVolumeVelocity);
1089 polymerShearFactor_ =
1090 PolymerModule::computeShearFactor(up.polymerConcentration(),
1091 pvtnumRegionIdx,
1092 waterVolumeVelocity*up.polymerViscosityCorrection());
1093
1094 }
1095
1096 const Evaluation& polymerShearFactor() const
1097 { return polymerShearFactor_; }
1098
1099 const Evaluation& waterShearFactor() const
1100 { return waterShearFactor_; }
1101
1102
1103private:
1104 Implementation& asImp_()
1105 { return *static_cast<Implementation*>(this); }
1106
1107 Evaluation polymerShearFactor_;
1108 Evaluation waterShearFactor_;
1109
1110};
1111
1112template <class TypeTag>
1114{
1117
1118public:
1119 void updateShearMultipliers(const ElementContext&,
1120 unsigned,
1121 unsigned)
1122 { }
1123
1124 void updateShearMultipliersPerm(const ElementContext&,
1125 unsigned,
1126 unsigned)
1127 { }
1128
1129 const Evaluation& polymerShearFactor() const
1130 { throw std::runtime_error("polymerShearFactor() called but polymers are disabled"); }
1131
1132 const Evaluation& waterShearFactor() const
1133 { throw std::runtime_error("waterShearFactor() called but polymers are disabled"); }
1134};
1135
1136
1137} // namespace Opm
1138
1139#endif
Contains the parameters required to extend the black-oil model by polymer.
Declares the properties required by the black oil model.
const Evaluation & polymerShearFactor() const
Definition: blackoilpolymermodules.hh:1129
void updateShearMultipliersPerm(const ElementContext &, unsigned, unsigned)
Definition: blackoilpolymermodules.hh:1124
void updateShearMultipliers(const ElementContext &, unsigned, unsigned)
Definition: blackoilpolymermodules.hh:1119
const Evaluation & waterShearFactor() const
Definition: blackoilpolymermodules.hh:1132
Provides the polymer specific extensive quantities to the generic black-oil module's extensive quanti...
Definition: blackoilpolymermodules.hh:991
void updateShearMultipliers(const ElementContext &elemCtx, unsigned scvfIdx, unsigned timeIdx)
Method which calculates the shear factor based on flow velocity.
Definition: blackoilpolymermodules.hh:1036
const Evaluation & waterShearFactor() const
Definition: blackoilpolymermodules.hh:1099
void updateShearMultipliersPerm(const ElementContext &, unsigned, unsigned)
Method which calculates the shear factor based on flow velocity.
Definition: blackoilpolymermodules.hh:1020
const Evaluation & polymerShearFactor() const
Definition: blackoilpolymermodules.hh:1096
const Evaluation & polymerConcentration() const
Definition: blackoilpolymermodules.hh:962
const Evaluation & waterViscosityCorrection() const
Definition: blackoilpolymermodules.hh:977
const Evaluation & polymerRockDensity() const
Definition: blackoilpolymermodules.hh:971
const Evaluation & polymerDeadPoreVolume() const
Definition: blackoilpolymermodules.hh:965
const Evaluation & polymerMoleWeight() const
Definition: blackoilpolymermodules.hh:959
void polymerPropertiesUpdate_(const ElementContext &, unsigned, unsigned)
Definition: blackoilpolymermodules.hh:954
const Evaluation & polymerAdsorption() const
Definition: blackoilpolymermodules.hh:968
const Evaluation & polymerViscosityCorrection() const
Definition: blackoilpolymermodules.hh:974
Provides the volumetric quantities required for the equations needed by the polymers extension of the...
Definition: blackoilpolymermodules.hh:807
const Evaluation & polymerAdsorption() const
Definition: blackoilpolymermodules.hh:915
Scalar polymerRockDensity_
Definition: blackoilpolymermodules.hh:938
const Evaluation & polymerMoleWeight() const
Definition: blackoilpolymermodules.hh:904
const Evaluation & polymerViscosityCorrection() const
Definition: blackoilpolymermodules.hh:922
const Scalar & polymerRockDensity() const
Definition: blackoilpolymermodules.hh:918
void polymerPropertiesUpdate_(const ElementContext &elemCtx, unsigned dofIdx, unsigned timeIdx)
Update the intensive properties needed to handle polymers from the primary variables.
Definition: blackoilpolymermodules.hh:834
Evaluation waterViscosityCorrection_
Definition: blackoilpolymermodules.hh:941
Scalar polymerDeadPoreVolume_
Definition: blackoilpolymermodules.hh:937
Evaluation polymerAdsorption_
Definition: blackoilpolymermodules.hh:939
Evaluation polymerMoleWeight_
Definition: blackoilpolymermodules.hh:936
const Evaluation & waterViscosityCorrection() const
Definition: blackoilpolymermodules.hh:926
Evaluation polymerConcentration_
Definition: blackoilpolymermodules.hh:934
const Evaluation & polymerConcentration() const
Definition: blackoilpolymermodules.hh:901
Implementation & asImp_()
Definition: blackoilpolymermodules.hh:931
const Scalar & polymerDeadPoreVolume() const
Definition: blackoilpolymermodules.hh:912
Evaluation polymerViscosityCorrection_
Definition: blackoilpolymermodules.hh:940
Contains the high level supplements required to extend the black oil model by polymer.
Definition: blackoilpolymermodules.hh:61
static const Scalar plyrockAdsorbtionIndex(const ElementContext &elemCtx, unsigned scvIdx, unsigned timeIdx)
Definition: blackoilpolymermodules.hh:627
static const TabulatedFunction & plyviscViscosityMultiplierTable(const ElementContext &elemCtx, unsigned scvIdx, unsigned timeIdx)
Definition: blackoilpolymermodules.hh:651
static std::string primaryVarName(unsigned pvIdx)
Definition: blackoilpolymermodules.hh:421
static const TabulatedFunction & plyadsAdsorbedPolymer(const ElementContext &elemCtx, unsigned scvIdx, unsigned timeIdx)
Definition: blackoilpolymermodules.hh:643
const Scalar molarMass() const
Definition: blackoilpolymermodules.hh:784
static TabulatedTwoDFunction & getPlymwinjTable(const int tableNumber)
get the PLYMWINJ table
Definition: blackoilpolymermodules.hh:350
static void serializeEntity(const Model &model, std::ostream &outstream, const DofEntity &dof)
Definition: blackoilpolymermodules.hh:576
static Scalar eqWeight(unsigned eqIdx)
Definition: blackoilpolymermodules.hh:463
static void registerParameters()
Register all run-time parameters for the black-oil polymer module.
Definition: blackoilpolymermodules.hh:393
static bool eqApplies(unsigned eqIdx)
Definition: blackoilpolymermodules.hh:441
static Scalar primaryVarWeight(unsigned pvIdx)
Definition: blackoilpolymermodules.hh:433
static BlackOilPolymerParams< Scalar >::SkprpolyTable & getSkprpolyTable(const int tableNumber)
get the SKPRPOLY table
Definition: blackoilpolymermodules.hh:379
static bool hasShrate()
Definition: blackoilpolymermodules.hh:694
static const TabulatedFunction & plyviscViscosityMultiplierTable(unsigned pvtnumRegionIdx)
Definition: blackoilpolymermodules.hh:659
static const Scalar plyrockResidualResistanceFactor(const ElementContext &elemCtx, unsigned scvIdx, unsigned timeIdx)
Definition: blackoilpolymermodules.hh:611
static void registerOutputModules(Model &model, Simulator &simulator)
Register all polymer specific VTK and ECL output modules.
Definition: blackoilpolymermodules.hh:402
static Scalar computeUpdateError(const PrimaryVariables &, const EqVector &)
Return how much a Newton-Raphson update is considered an error.
Definition: blackoilpolymermodules.hh:566
static const BlackOilPolymerParams< Scalar >::PlyvmhCoefficients & plyvmhCoefficients(const ElementContext &elemCtx, const unsigned scvIdx, const unsigned timeIdx)
Definition: blackoilpolymermodules.hh:681
static const Scalar shrate(unsigned pvtnumRegionIdx)
Definition: blackoilpolymermodules.hh:699
static bool primaryVarApplies(unsigned pvIdx)
Definition: blackoilpolymermodules.hh:409
static void addStorage(Dune::FieldVector< LhsEval, numEq > &storage, const IntensiveQuantities &intQuants)
Definition: blackoilpolymermodules.hh:473
static const Scalar plyrockRockDensityFactor(const ElementContext &elemCtx, unsigned scvIdx, unsigned timeIdx)
Definition: blackoilpolymermodules.hh:619
static const Scalar plymaxMaxConcentration(const ElementContext &elemCtx, unsigned scvIdx, unsigned timeIdx)
Definition: blackoilpolymermodules.hh:664
static bool hasPlyshlog()
Definition: blackoilpolymermodules.hh:689
static std::string eqName(unsigned eqIdx)
Definition: blackoilpolymermodules.hh:453
static const Scalar plyrockDeadPoreVolume(const ElementContext &elemCtx, unsigned scvIdx, unsigned timeIdx)
Definition: blackoilpolymermodules.hh:603
static const Scalar plymixparToddLongstaff(const ElementContext &elemCtx, unsigned scvIdx, unsigned timeIdx)
Definition: blackoilpolymermodules.hh:672
static TabulatedTwoDFunction & getSkprwatTable(const int tableNumber)
get the SKPRWAT table
Definition: blackoilpolymermodules.hh:364
static void deserializeEntity(Model &model, std::istream &instream, const DofEntity &dof)
Definition: blackoilpolymermodules.hh:587
static const Scalar plyrockMaxAdsorbtion(const ElementContext &elemCtx, unsigned scvIdx, unsigned timeIdx)
Definition: blackoilpolymermodules.hh:635
static void computeFlux(RateVector &flux, const ElementContext &elemCtx, unsigned scvfIdx, unsigned timeIdx)
Definition: blackoilpolymermodules.hh:512
static Evaluation computeShearFactor(const Evaluation &polymerConcentration, unsigned pvtnumRegionIdx, const Evaluation &v0)
Computes the shear factor.
Definition: blackoilpolymermodules.hh:711
VTK output module for the black oil model's polymer related quantities.
Definition: vtkblackoilpolymermodule.hh:71
Definition: blackoilboundaryratevector.hh:37
typename Properties::Detail::GetPropImpl< TypeTag, Property >::type::type GetPropType
get the type alias defined in the property (equivalent to old macro GET_PROP_TYPE(....
Definition: propertysystem.hh:235
Definition: blackoilpolymerparams.hh:99
Definition: blackoilpolymerparams.hh:106
Struct holding the parameters for the BlackOilPolymerModule class.
Definition: blackoilpolymerparams.hh:41
AdsorptionBehaviour
Definition: blackoilpolymerparams.hh:45
IntervalTabulated2DFunction< Scalar > TabulatedTwoDFunction
Definition: blackoilpolymerparams.hh:43
Tabulated1DFunction< Scalar > TabulatedFunction
Definition: blackoilpolymerparams.hh:42