20 #ifndef OPM_POLYMERPROPERTIES_HEADER_INCLUDED
21 #define OPM_POLYMERPROPERTIES_HEADER_INCLUDED
23 #include <opm/parser/eclipse/Deck/Deck.hpp>
24 #include <opm/parser/eclipse/EclipseState/EclipseState.hpp>
25 #include <opm/parser/eclipse/EclipseState/Tables/PlyadsTable.hpp>
26 #include <opm/parser/eclipse/EclipseState/Tables/PlymaxTable.hpp>
27 #include <opm/parser/eclipse/EclipseState/Tables/PlyrockTable.hpp>
28 #include <opm/parser/eclipse/EclipseState/Tables/PlyshlogTable.hpp>
29 #include <opm/parser/eclipse/EclipseState/Tables/PlyviscTable.hpp>
30 #include <opm/parser/eclipse/EclipseState/Tables/TableManager.hpp>
31 #include <opm/parser/eclipse/Units/Dimension.hpp>
32 #include <opm/parser/eclipse/Units/UnitSystem.hpp>
37 #include <opm/common/ErrorMacros.hpp>
50 enum AdsorptionBehaviour { Desorption = 1, NoDesorption = 2 };
71 AdsorptionBehaviour ads_index,
72 const std::vector<double>& c_vals_visc,
73 const std::vector<double>& visc_mult_vals,
74 const std::vector<double>& c_vals_ads,
75 const std::vector<double>& ads_vals,
76 const std::vector<double>& water_vel_vals,
77 const std::vector<double>& shear_vrf_vals
80 mix_param_(mix_param),
81 rock_density_(rock_density),
82 dead_pore_vol_(dead_pore_vol),
83 res_factor_(res_factor),
84 c_max_ads_(c_max_ads),
85 ads_index_(ads_index),
86 c_vals_visc_(c_vals_visc),
87 visc_mult_vals_(visc_mult_vals),
88 c_vals_ads_(c_vals_ads),
90 water_vel_vals_(water_vel_vals),
91 shear_vrf_vals_(shear_vrf_vals)
97 readFromDeck(deck, eclipseState);
100 void set(
double c_max,
103 double dead_pore_vol,
106 AdsorptionBehaviour ads_index,
107 const std::vector<double>& c_vals_visc,
108 const std::vector<double>& visc_mult_vals,
109 const std::vector<double>& c_vals_ads,
110 const std::vector<double>& ads_vals,
111 const std::vector<double>& water_vel_vals,
112 const std::vector<double>& shear_vrf_vals
116 mix_param_ = mix_param;
117 rock_density_ = rock_density;
118 dead_pore_vol_ = dead_pore_vol;
119 res_factor_ = res_factor;
120 c_max_ads_ = c_max_ads;
121 c_vals_visc_ = c_vals_visc;
122 visc_mult_vals_ = visc_mult_vals;
123 c_vals_ads_ = c_vals_ads;
124 ads_vals_ = ads_vals;
125 ads_index_ = ads_index;
126 water_vel_vals_ = water_vel_vals;
127 shear_vrf_vals_ = shear_vrf_vals;
130 void readFromDeck(
const Opm::Deck& deck,
const Opm::EclipseState& eclipseState)
133 const auto& tables = eclipseState.getTableManager();
134 const auto& plymaxTable = tables.getPlymaxTables().getTable<PlymaxTable>(0);
135 const auto& plmixparRecord = deck.getKeyword(
"PLMIXPAR").getRecord(0);
138 assert(plymaxTable.numRows() == 1);
140 c_max_ = plymaxTable.getPolymerConcentrationColumn()[0];
141 mix_param_ = plmixparRecord.getItem(
"TODD_LONGSTAFF").getSIDouble(0);
144 const auto& plyrockTable = tables.getPlyrockTables().getTable<PlyrockTable>(0);
147 assert(plyrockTable.numRows() == 1);
149 dead_pore_vol_ = plyrockTable.getDeadPoreVolumeColumn()[0];
150 res_factor_ = plyrockTable.getResidualResistanceFactorColumn()[0];
151 rock_density_ = plyrockTable.getRockDensityFactorColumn()[0];
152 ads_index_ =
static_cast<AdsorptionBehaviour
>(plyrockTable.getAdsorbtionIndexColumn()[0]);
153 c_max_ads_ = plyrockTable.getMaxAdsorbtionColumn()[0];
156 const auto& plyviscTable = tables.getPlyviscTables().getTable<PlyviscTable>(0);
159 c_vals_visc_ = plyviscTable.getPolymerConcentrationColumn().vectorCopy( );
160 visc_mult_vals_ = plyviscTable.getViscosityMultiplierColumn().vectorCopy( );
163 const auto& plyadsTable = tables.getPlyadsTables().getTable<PlyadsTable>(0);
165 c_vals_ads_ = plyadsTable.getPolymerConcentrationColumn().vectorCopy( );
166 ads_vals_ = plyadsTable.getAdsorbedPolymerColumn().vectorCopy( );
168 has_plyshlog_ = deck.hasKeyword(
"PLYSHLOG");
169 has_shrate_ = deck.hasKeyword(
"SHRATE");
173 const auto& plyshlogTable = tables.getPlyshlogTables().getTable<PlyshlogTable>(0);
175 water_vel_vals_ = plyshlogTable.getWaterVelocityColumn().vectorCopy( );
176 shear_vrf_vals_ = plyshlogTable.getShearMultiplierColumn().vectorCopy( );
179 Opm::UnitSystem unitSystem = deck.getActiveUnitSystem();
182 siFactor = unitSystem.parse(
"1/Time").getSIScaling();
183 const auto& shrateKeyword = deck.getKeyword(
"SHRATE");
184 std::vector<double> shrate_readin = shrateKeyword.getSIDoubleData();
185 if (shrate_readin.size() == 1) {
186 shrate_ = shrate_readin[0];
187 }
else if (shrate_readin.size() == 0) {
190 OPM_THROW(std::logic_error,
"Only NTPVT == 1 is allowed for SHRATE keyword now !\n");
193 siFactor = unitSystem.parse(
"Length/Time").getSIScaling();
196 for (
size_t i = 0; i < water_vel_vals_.size(); ++i) {
197 water_vel_vals_[i] *= siFactor;
201 plyshlog_ref_conc_ = plyshlogTable.getRefPolymerConcentration();
203 if (plyshlogTable.hasRefSalinity()) {
204 has_plyshlog_ref_salinity_ =
true;
205 plyshlog_ref_salinity_ = plyshlogTable.getRefSalinity();
207 has_plyshlog_ref_salinity_ =
false;
210 if (plyshlogTable.hasRefTemperature()) {
211 has_plyshlog_ref_temp_ =
true;
212 plyshlog_ref_temp_ = plyshlogTable.getRefTemperature();
214 has_plyshlog_ref_temp_ =
false;
221 double mixParam()
const;
223 double rockDensity()
const;
225 double deadPoreVol()
const;
227 double resFactor()
const;
229 double cMaxAds()
const;
231 int adsIndex()
const;
263 double shearVrf(
const double velocity)
const;
265 double shearVrfWithDer(
const double velocity,
double& der)
const;
267 double viscMult(
double c)
const;
269 double viscMultWithDer(
double c,
double* der)
const;
271 void simpleAdsorption(
double c,
double& c_ads)
const;
273 void simpleAdsorptionWithDer(
double c,
double& c_ads,
274 double& dc_ads_dc)
const;
276 void adsorption(
double c,
double cmax,
double& c_ads)
const;
278 void adsorptionWithDer(
double c,
double cmax,
279 double& c_ads,
double& dc_ads_dc)
const;
281 void effectiveVisc(
const double c,
const double mu_w,
282 double& mu_w_eff)
const;
284 void effectiveViscWithDer(
const double c,
const double visc
286 ,
double dmu_w_eff_dc)
const;
288 void effectiveInvVisc(
const double c,
const double mu_w,
289 double& inv_mu_w_eff)
const;
291 void effectiveInvViscWithDer(
const double c,
293 double& inv_mu_w_eff,
294 double& dinv_mu_w_eff_dc)
const;
296 void effectiveInvPolyVisc(
const double c,
298 double& inv_mu_p_eff)
const;
300 void effectiveInvPolyViscWithDer(
const double c,
302 double& inv_mu_p_eff,
303 double& d_inv_mu_p_eff_dc)
const;
305 void effectiveRelperm(
const double c,
307 const double* relperm,
308 double& eff_relperm_wat)
const;
310 void effectiveRelpermWithDer (
const double c,
312 const double* relperm,
313 const double* drelperm_ds,
314 double& eff_relperm_wat,
315 double& deff_relperm_wat_ds,
316 double& deff_relperm_wat_dc)
const;
318 void effectiveMobilities(
const double c,
321 const double* relperm,
324 void effectiveMobilitiesWithDer(
const double c,
327 const double* relperm,
328 const double* drelpermds,
331 double& dmobwatdc)
const;
333 void effectiveMobilitiesBoth(
const double c,
336 const double* relperm,
337 const double* drelperm_ds,
341 bool if_with_der)
const;
343 void effectiveTotalMobility(
const double c,
346 const double* relperm,
347 double& totmob)
const;
349 void effectiveTotalMobilityWithDer(
const double c,
352 const double* relperm,
353 const double* drelpermds,
355 double* dtotmob_dsdc)
const;
357 void effectiveTotalMobilityBoth(
const double c,
360 const double* relperm,
361 const double* drelperm_ds,
363 double* dtotmob_dsdc,
364 bool if_with_der)
const;
366 void computeMc(
const double& c,
double& mc)
const;
368 void computeMcWithDer(
const double& c,
double& mc,
369 double& dmc_dc)
const;
371 void computeMcBoth(
const double& c,
double& mc,
372 double& dmc_dc,
bool if_with_der)
const;
375 bool computeShearMultLog(std::vector<double>& water_vel, std::vector<double>& visc_mult, std::vector<double>& shear_mult)
const;
380 double rock_density_;
381 double dead_pore_vol_;
391 AdsorptionBehaviour ads_index_;
392 std::vector<double> c_vals_visc_;
393 std::vector<double> visc_mult_vals_;
394 std::vector<double> c_vals_ads_;
395 std::vector<double> ads_vals_;
396 std::vector<double> water_vel_vals_;
397 std::vector<double> shear_vrf_vals_;
399 double plyshlog_ref_conc_;
400 double plyshlog_ref_salinity_;
401 double plyshlog_ref_temp_;
402 bool has_plyshlog_ref_salinity_;
403 bool has_plyshlog_ref_temp_;
406 void simpleAdsorptionBoth(
double c,
double& c_ads,
407 double& dc_ads_dc,
bool if_with_der)
const;
408 void adsorptionBoth(
double c,
double cmax,
409 double& c_ads,
double& dc_ads_dc,
410 bool if_with_der)
const;
412 void effectiveInvViscBoth(
const double c,
const double mu_w,
413 double& inv_mu_w_eff,
414 double& dinv_mu_w_eff_dc,
bool if_with_der)
const;
416 void effectiveInvPolyViscBoth(
const double c,
418 double& inv_mu_p_eff,
419 double& dinv_mu_p_eff_dc,
420 const bool if_with_der)
const;
422 void effectiveRelpermBoth(
const double c,
424 const double* relperm,
425 const double* drelperm_ds,
426 double& eff_relperm_wat,
427 double& deff_relperm_wat_ds,
428 double& deff_relperm_wat_dc,
429 bool if_with_der)
const;
435 #endif // OPM_POLYMERPROPERTIES_HEADER_INCLUDED
bool computeShearMultLog(std::vector< double > &water_vel, std::vector< double > &visc_mult, std::vector< double > &shear_mult) const
Computing the shear multiplier based on the water velocity/shear rate with PLYSHLOG keyword...
Definition: PolymerProperties.cpp:474
double shrate() const
the value of SHRATE
Definition: PolymerProperties.cpp:116
Definition: PolymerProperties.hpp:43
bool hasShrate() const
indicate whether SHRATE keyword is specified
Definition: PolymerProperties.cpp:110
PolymerProperties(double c_max, double mix_param, double rock_density, double dead_pore_vol, double res_factor, double c_max_ads, AdsorptionBehaviour ads_index, const std::vector< double > &c_vals_visc, const std::vector< double > &visc_mult_vals, const std::vector< double > &c_vals_ads, const std::vector< double > &ads_vals, const std::vector< double > &water_vel_vals, const std::vector< double > &shear_vrf_vals)
Construct from parameters.
Definition: PolymerProperties.hpp:65
const std::vector< double > & shearViscosityReductionFactor() const
the viscosity reduction factor PLYSHLOG table
Definition: PolymerProperties.cpp:75
bool hasPlyshlogRefTemp() const
indicate whether reference temperature is specified in PLYSHLOG
Definition: PolymerProperties.cpp:90
double plyshlogRefTemp() const
the reference temperature in PLYSHLOG
Definition: PolymerProperties.cpp:100
double plyshlogRefSalinity() const
the reference salinity in PLYSHLOG
Definition: PolymerProperties.cpp:95
bool hasPlyshlogRefSalinity() const
indicate wheter reference salinity is specified in PLYSHLOG
Definition: PolymerProperties.cpp:85
double plyshlogRefConc() const
the reference polymer concentration in PLYSHLOG
Definition: PolymerProperties.cpp:80
bool hasPlyshlog() const
indicate whehter PLYSHLOG is specified
Definition: PolymerProperties.cpp:105
const std::vector< double > & shearWaterVelocity() const
the water velocity or water shear rate in PLYSHLOG table
Definition: PolymerProperties.cpp:69