polymerUtilities uses now PolymerProperties to compute effective mobilities.

This commit is contained in:
Xavier Raynaud
2012-04-23 11:49:05 +02:00
parent 4794348740
commit e867c53722
5 changed files with 109 additions and 46 deletions

View File

@@ -593,10 +593,10 @@ main(int argc, char** argv)
// Solve pressure.
if (use_gravity) {
computeTotalMobilityOmega(*props, polyprop, allcells, state.saturation(), state.concentration(),
computeTotalMobilityOmega(*props, polyprop, allcells, state.saturation(), state.concentration(), state.maxconcentration(),
totmob, omega);
} else {
computeTotalMobility(*props, polyprop, allcells, state.saturation(), state.concentration(),
computeTotalMobility(*props, polyprop, allcells, state.saturation(), state.concentration(), state.maxconcentration(),
totmob);
}
std::vector<double> empty_vector_for_wells;
@@ -681,7 +681,7 @@ main(int argc, char** argv)
Opm::computeSaturatedVol(porevol, state.saturation(), satvol);
polymass = Opm::computePolymerMass(porevol, state.saturation(), state.concentration(), polyprop.deadPoreVol());
polymass_adsorbed = Opm::computePolymerAdsorbed(*props, polyprop, porevol, state.maxconcentration());
Opm::computeInjectedProduced(*props, polyprop, state.saturation(), state.concentration(),
Opm::computeInjectedProduced(*props, polyprop, state.saturation(), state.concentration(), state.maxconcentration(),
src, simtimer.currentStepLength(), inflow_c,
injected, produced, polyinj, polyprod);
tot_injected[0] += injected[0];

View File

@@ -279,6 +279,53 @@ namespace Opm
}
}
void PolymerProperties::effectiveTotalMobility(const double c,
const double cmax,
const double* visc,
const double* relperm,
double& totmob) const
{
std::vector<double> dummy(2);
double dummy2[4];
effectiveTotalMobilityBoth(c, cmax, visc, relperm, dummy2,
totmob, dummy, false);
}
void PolymerProperties::effectiveTotalMobilityWithDer(const double c,
const double cmax,
const double* visc,
const double* relperm,
const double* drelperm_ds,
double& totmob,
std::vector<double>& dtotmob_dsdc) const
{
effectiveTotalMobilityBoth(c, cmax, visc, relperm, drelperm_ds,
totmob, dtotmob_dsdc, true);
}
void PolymerProperties::effectiveTotalMobilityBoth(const double c,
const double cmax,
const double* visc,
const double* relperm,
const double* drelperm_ds,
double& totmob,
std::vector<double>& dtotmob_dsdc,
bool if_with_der) const
{
std::vector<double> mob(2);
std::vector<double> dmob_ds(2);
double dmobwat_dc;
effectiveMobilitiesBoth(c, cmax, visc, relperm, drelperm_ds,
mob, dmob_ds, dmobwat_dc, if_with_der);
totmob = mob[0] + mob[1];
if (if_with_der) {
dtotmob_dsdc[0] = dmob_ds[0] + dmob_ds[2]; //derivative with respect to s
dtotmob_dsdc[1] = dmobwat_dc; //derivative with respect to c
} else {
dtotmob_dsdc.clear();
}
}
void PolymerProperties::computeMc(const double& c, double& mc) const
{
double dummy;

View File

@@ -145,20 +145,20 @@ namespace Opm
void simpleAdsorption(double c, double& c_ads) const;
void simpleAdsorptionWithDer(double c, double& c_ads,
void simpleAdsorptionWithDer(double c, double& c_ads,
double& dc_ads_dc) const;
void adsorption(double c, double cmax, double& c_ads) const;
void adsorptionWithDer(double c, double cmax,
void adsorptionWithDer(double c, double cmax,
double& c_ads, double& dc_ads_dc) const;
void effectiveVisc(const double c, const double* visc, double* visc_eff) const;
void effectiveInvVisc(const double c, const double* visc,
void effectiveInvVisc(const double c, const double* visc,
double* inv_visc_eff) const;
void effectiveViscWithDer(const double c, const double* visc,
void effectiveViscWithDer(const double c, const double* visc,
double* visc_eff, double* dvisc_eff_dc) const;
void effectiveRelperm(const double c,
@@ -185,26 +185,49 @@ namespace Opm
const double* visc,
const double* relperm,
const double* drelpermds,
std::vector<double>& mob,
std::vector<double>& dmobds,
double& dmobwatdc) const;
std::vector<double>& mob,
std::vector<double>& dmob_ds,
double& dmobwatdc) const;
void effectiveMobilitiesBoth(const double c,
const double cmax,
const double* visc,
const double* relperm,
const double* drelperm_ds,
std::vector<double>& mob,
std::vector<double>& mob,
std::vector<double>& dmob_ds,
double& dmobwat_dc,
bool if_with_der) const;
void effectiveTotalMobility(const double c,
const double cmax,
const double* visc,
const double* relperm,
double& totmob) const;
void effectiveTotalMobilityWithDer(const double c,
const double cmax,
const double* visc,
const double* relperm,
const double* drelpermds,
double& totmob,
std::vector<double>& dtotmob_dsdc) const;
void effectiveTotalMobilityBoth(const double c,
const double cmax,
const double* visc,
const double* relperm,
const double* drelperm_ds,
double& totmob,
std::vector<double>& dtotmob_dsdc,
bool if_with_der) const;
void computeMc(const double& c, double& mc) const;
void computeMcWithDer(const double& c, double& mc,
double& dmc_dc) const;
void computeMcBoth(const double& c, double& mc,
void computeMcBoth(const double& c, double& mc,
double& dmc_dc, bool if_with_der) const;
private:
@@ -219,10 +242,10 @@ namespace Opm
std::vector<double> visc_mult_vals_;
std::vector<double> c_vals_ads_;
std::vector<double> ads_vals_;
void simpleAdsorptionBoth(double c, double& c_ads,
void simpleAdsorptionBoth(double c, double& c_ads,
double& dc_ads_dc, bool if_with_der) const;
void adsorptionBoth(double c, double cmax,
double& c_ads, double& dc_ads_dc,
void adsorptionBoth(double c, double cmax,
double& c_ads, double& dc_ads_dc,
bool if_with_der) const;
void effectiveViscBoth(const double c, const double* visc, double* visc_eff,
double* dvisc_eff_dc, bool if_with_der) const;

View File

@@ -34,22 +34,18 @@ namespace Opm
const std::vector<int>& cells,
const std::vector<double>& s,
const std::vector<double>& c,
const std::vector<double>& cmax,
std::vector<double>& totmob)
{
int num_cells = cells.size();
int num_phases = props.numPhases();
totmob.resize(num_cells);
ASSERT(int(s.size()) == num_cells*num_phases);
std::vector<double> kr(num_cells*num_phases);
std::vector<double> kr(2*num_cells);
props.relperm(num_cells, &s[0], &cells[0], &kr[0], 0);
const double* mu = props.viscosity();
double inv_mu_eff[2] = { 0.0 };
const double* visc = props.viscosity();
for (int cell = 0; cell < num_cells; ++cell) {
totmob[cell] = 0;
polyprops.effectiveInvVisc(c[cell], mu, inv_mu_eff);
for (int phase = 0; phase < num_phases; ++phase) {
totmob[cell] += kr[num_phases*cell + phase]*inv_mu_eff[phase];
}
double* kr_cell = &kr[2*cell];
polyprops.effectiveTotalMobility(c[cell], cmax[cell], visc, kr_cell,
totmob[cell]);
}
}
@@ -70,6 +66,7 @@ namespace Opm
const std::vector<int>& cells,
const std::vector<double>& s,
const std::vector<double>& c,
const std::vector<double>& cmax,
std::vector<double>& totmob,
std::vector<double>& omega)
{
@@ -80,21 +77,16 @@ namespace Opm
ASSERT(int(s.size()) == num_cells*num_phases);
std::vector<double> kr(num_cells*num_phases);
props.relperm(num_cells, &s[0], &cells[0], &kr[0], 0);
const double* mu = props.viscosity();
double inv_mu_eff[2] = { 0.0 };
const double* visc = props.viscosity();
const double* rho = props.density();
std::vector<double> mob(num_phases); // here we assume num_phases=2
for (int cell = 0; cell < num_cells; ++cell) {
totmob[cell] = 0.0;
omega[cell] = 0.0;
polyprops.effectiveInvVisc(c[cell], mu, inv_mu_eff);
for (int phase = 0; phase < num_phases; ++phase) {
totmob[cell] += kr[num_phases*cell + phase]*inv_mu_eff[phase];
}
// Must finish computing totmob before we can use it.
for (int phase = 0; phase < num_phases; ++phase) {
omega[cell] += rho[phase]*(kr[num_phases*cell + phase]*inv_mu_eff[phase])/totmob[cell];
}
}
double* kr_cell = &kr[2*cell];
polyprops.effectiveMobilities(c[cell], cmax[cell], visc, kr_cell,
mob);
totmob[cell] = mob[0] + mob[1];
omega[cell] = rho[0]*mob[0]/totmob[cell] + rho[1]*mob[1]/totmob[cell];
}
}
@@ -120,6 +112,7 @@ namespace Opm
const Opm::PolymerProperties& polyprops,
const std::vector<double>& s,
const std::vector<double>& c,
const std::vector<double>& cmax,
const std::vector<double>& src,
const double dt,
const double inj_c,
@@ -137,8 +130,8 @@ namespace Opm
std::fill(produced, produced + np, 0.0);
polyinj = 0.0;
polyprod = 0.0;
std::vector<double> inv_eff_visc(np);
const double* visc = props.viscosity();
std::vector<double> kr_cell(np);
std::vector<double> mob(np);
for (int cell = 0; cell < num_cells; ++cell) {
if (src[cell] > 0.0) {
@@ -147,13 +140,10 @@ namespace Opm
} else if (src[cell] < 0.0) {
const double flux = -src[cell]*dt;
const double* sat = &s[np*cell];
props.relperm(1, sat, &cell, &mob[0], 0);
polyprops.effectiveInvVisc(c[cell], visc, &inv_eff_visc[0]);
double totmob = 0.0;
for (int p = 0; p < np; ++p) {
mob[p] *= inv_eff_visc[p];
totmob += mob[p];
}
props.relperm(1, sat, &cell, &kr_cell[0], 0);
polyprops.effectiveMobilities(c[cell], cmax[cell], visc,
&kr_cell[0], mob);
double totmob = mob[0] + mob[1];
for (int p = 0; p < np; ++p) {
produced[p] += (mob[p]/totmob)*flux;
}

View File

@@ -42,6 +42,7 @@ namespace Opm
const std::vector<int>& cells,
const std::vector<double>& s,
const std::vector<double>& c,
const std::vector<double>& cmax,
std::vector<double>& totmob);
/// @brief Computes total mobility and omega for a set of s/c values.
@@ -58,6 +59,7 @@ namespace Opm
const std::vector<int>& cells,
const std::vector<double>& s,
const std::vector<double>& c,
const std::vector<double>& cmax,
std::vector<double>& totmob,
std::vector<double>& omega);
@@ -83,6 +85,7 @@ namespace Opm
const Opm::PolymerProperties& polyprops,
const std::vector<double>& s,
const std::vector<double>& c,
const std::vector<double>& cmax,
const std::vector<double>& src,
const double dt,
const double inj_c,