OilfieldToolsNet/opm-simulators
4
0
Fork 0
mirror of https://github.com/OPM/opm-simulators.git synced 2025-02-25 18:55:30 -06:00
Code Issues Packages Projects Releases Wiki Activity

fix errors for debugging.

Browse Source
This commit is contained in:
Liu Ming
2013-12-12 21:58:25 +08:00
parent 25160d019a
commit 2eaf24decf
5 changed files with 22 additions and 25 deletions
Download Patch File Download Diff File Expand all files Collapse all files

BIN
opm/polymer/fullyimplicit/.AutoDiffBlock.hpp.swp
View File

Binary file not shown.

24
opm/polymer/fullyimplicit/FullyImplicitTwophasePolymerSolver.cpp
View File

@@ -228,7 +228,7 @@ typedef Eigen::Array<double,
ADB ADB
FullyImplicitTwophasePolymerSolver:: FullyImplicitTwophasePolymerSolver::
computeCmax(const ADB& c) computeCmax(const ADB& c) const
{ {
const int nc = c.value().size(); const int nc = c.value().size();
V cmax(nc); V cmax(nc);
@@ -258,8 +258,8 @@ typedef Eigen::Array<double,
const std::vector<ADB> kr = computeRelPerm(state); const std::vector<ADB> kr = computeRelPerm(state);
const ADB cmax = computeCmax(state.concentration); const ADB cmax = computeCmax(state.concentration);
const ADB ads = adsorption(state.concentration, cmax); const ADB ads = polymer_props_ad_.adsorption(state.concentration, cmax);
const ADB krw_eff = polymer_props_ad_.effectiveRelPerm(c, cmax, kr[0], state.saturation[0]); const ADB krw_eff = polymer_props_ad_.effectiveRelPerm(state.concentration, cmax, kr[0], state.saturation[0]);
const ADB mc = computeMc(state); const ADB mc = computeMc(state);
const std::vector<ADB> mflux = computeMassFlux(trans, mc, kr[0], krw_eff, state); const std::vector<ADB> mflux = computeMassFlux(trans, mc, kr[0], krw_eff, state);
@@ -274,7 +274,7 @@ typedef Eigen::Array<double,
residual_[2] = pvdt * (state.saturation[0] * state.concentration residual_[2] = pvdt * (state.saturation[0] * state.concentration
- old_state.saturation[0] * old_state.concentration) - old_state.saturation[0] * old_state.concentration)
+ pvdt * rho_r * (1. - phi) / phi * ads + pvdt * rho_r * (1. - phi) / phi * ads
+ ops_.div * mflux[3] - srouce[3]; + ops_.div * mflux[3] - source[3];
} }
@@ -289,7 +289,7 @@ typedef Eigen::Array<double,
std::vector<ADB> mflux; std::vector<ADB> mflux;
ADB inv_wat_eff_vis = polymer_props_ad_.effectiveInvWaterVisc(state.concentration, mus); ADB inv_wat_eff_vis = polymer_props_ad_.effectiveInvWaterVisc(state.concentration, mus);
ADB wat_mob = krw_eff * inv_wat_eff_vis; ADB wat_mob = krw_eff * inv_wat_eff_vis;
ADB oil_mob = kr[1] / V::Constant(kr[1].size(), 1, mus[1]); ADB oil_mob = kro / V::Constant(kro.size(), 1, mus[1]);
ADB poly_mob = mc * krw_eff * inv_wat_eff_vis; ADB poly_mob = mc * krw_eff * inv_wat_eff_vis;
@@ -325,11 +325,10 @@ typedef Eigen::Array<double,
insrc.push_back(*it); insrc.push_back(*it);
outsrc.push_back(0.0); outsrc.push_back(0.0);
} else { } else {
outsrc.emplace_back(0); outsrc.push_back(0);
insrc.emplace_back(0); insrc.push_back(0);
} }
} }
const V source = Eigen::Map<const V>(& src[0], grid_.number_of_cells);
const V outSrc = Eigen::Map<const V>(& outsrc[0], grid_.number_of_cells); const V outSrc = Eigen::Map<const V>(& outsrc[0], grid_.number_of_cells);
const V inSrc = Eigen::Map<const V>(& insrc[0], grid_.number_of_cells); const V inSrc = Eigen::Map<const V>(& insrc[0], grid_.number_of_cells);
const V polyin = Eigen::Map<const V>(& polymer_inflow_c[0], grid_.number_of_cells); const V polyin = Eigen::Map<const V>(& polymer_inflow_c[0], grid_.number_of_cells);
@@ -338,7 +337,6 @@ typedef Eigen::Array<double,
// compute the in-fracflow. // compute the in-fracflow.
V zero = V::Zero(grid_.number_of_cells); V zero = V::Zero(grid_.number_of_cells);
V one = V::Ones(grid_.number_of_cells); V one = V::Ones(grid_.number_of_cells);
return f_out * outSrc + f_in * inSrc ;
std::vector<ADB> source; std::vector<ADB> source;
//water source //water source
@@ -346,7 +344,9 @@ typedef Eigen::Array<double,
//oil source //oil source
source.push_back(f[1] * outSrc + zero * inSrc); source.push_back(f[1] * outSrc + zero * inSrc);
//polymer source //polymer source
source.push_back(f[0] * outSrc * c + one * inSrc * polyin) source.push_back(f[0] * outSrc * c + one * inSrc * polyin);
return source;
} }
@@ -359,8 +359,8 @@ typedef Eigen::Array<double,
{ {
const double* mus = fluid_.viscosity(); const double* mus = fluid_.viscosity();
ADB inv_wat_eff_vis = polymer_props_ad_.effectiveInvWaterVisc(c, mus); ADB inv_wat_eff_vis = polymer_props_ad_.effectiveInvWaterVisc(c, mus);
ADB wat_mob = kr[0] * inv_wat_eff_vis; ADB wat_mob = krw_eff * inv_wat_eff_vis;
ADB oil_mob = kr[1] / V::Constant(kr[1].size(), 1, mus[1]); ADB oil_mob = kro / V::Constant(kro.size(), 1, mus[1]);
ADB total_mob = wat_mob + oil_mob; ADB total_mob = wat_mob + oil_mob;
std::vector<ADB> fracflow; std::vector<ADB> fracflow;

7
opm/polymer/fullyimplicit/FullyImplicitTwophasePolymerSolver.hpp
View File

@@ -90,11 +90,6 @@ namespace Opm {
computeFracFlow(const ADB& kro, computeFracFlow(const ADB& kro,
const ADB& krw_eff, const ADB& krw_eff,
const ADB& c) const; const ADB& c) const;
ADB
computePolymerMassFlux(const V& trans,
const ADB& mc,
const std::vector<ADB>& kr,
const SolutionState& state) const;
double double
residualNorm() const; residualNorm() const;
ADB ADB
@@ -103,6 +98,8 @@ namespace Opm {
const std::vector<double>& polymer_inflow_c, const std::vector<double>& polymer_inflow_c,
const SolutionState& state) const; const SolutionState& state) const;
ADB
computeCmax(const ADB& c) const;
ADB ADB
computeMc(const SolutionState& state) const; computeMc(const SolutionState& state) const;
ADB ADB

12
opm/polymer/fullyimplicit/PolymerPropsAd.cpp
View File

@@ -326,14 +326,14 @@ namespace Opm {
const int nc = c.size(); const int nc = c.size();
V one = V::Ones(nc); V one = V::Ones(nc);
V ads = adsorption(c); V ads = adsorption(c, cmax_cells);
double max_ads = polymer_props_.cMaxAds(); double max_ads = polymer_props_.cMaxAds();
double res_factor = polymer_props_.resFactor(); double res_factor = polymer_props_.resFactor();
double factor = (res_factor -1.) / max_ads; double factor = (res_factor -1.) / max_ads;
V rk = one + factor * ads; V rk = one + factor * ads;
V krw_eff = krw / rk; V krw_eff = krw / rk;
return eff_relperm; return krw_eff;
} }
@@ -347,20 +347,20 @@ namespace Opm {
V one = V::Ones(nc); V one = V::Ones(nc);
ADB ads = adsorption(c); ADB ads = adsorption(c, cmax_cells);
V krw_eff = effectiveRelPerm(c.value(), cmax_cells.value(), krw.value()); V krw_eff = effectiveRelPerm(c.value(), cmax_cells.value(), krw.value());
double max_ads = polymer_props_.cMaxAds(); double max_ads = polymer_props_.cMaxAds();
double res_factor = polymer_props_.resFactor(); double res_factor = polymer_props_.resFactor();
double factor = (res_factor - 1.) / max_ads; double factor = (res_factor - 1.) / max_ads;
ADB rk = one + ads * factor; ADB rk = one + ads * factor;
ADB::M dkrw_ds = krw.derivative() / rk.derivative(); ADB dkrw_ds = krw / rk.value();
ADB::M dkrw_dc = -krw.value() / (rk.value * rk.value()) * ads.derivative() * factor; ADB dkrw_dc = -krw.value() / (rk.value() * rk.value()) * ads * factor;
const int num_blocks = c.numBlocks(); const int num_blocks = c.numBlocks();
std::vector<ADB::M> jacs(num_blocks); std::vector<ADB::M> jacs(num_blocks);
for (int block = 0; block < num_blocks; ++block) { for (int block = 0; block < num_blocks; ++block) {
jac[block] = dkrw_ds * sw.derivative()[block] + dkrw_dc * c.derivative()[block]; jacs[block] = dkrw_ds.derivative()[block] * sw.derivative()[block] + dkrw_dc.derivative()[block] * c.derivative()[block];
} }
return ADB::function(krw_eff, jacs); return ADB::function(krw_eff, jacs);

4
opm/polymer/fullyimplicit/PolymerPropsAd.hpp
View File

@@ -83,10 +83,10 @@ namespace Opm {
polymerWaterVelocityRatio(const ADB& c) const; polymerWaterVelocityRatio(const ADB& c) const;
V V
adsorption(const V& c) const; adsorption(const V& c, const V& cmax_cells) const;
ADB ADB
adsorption(const ADB& c) const; adsorption(const ADB& c, const ADB& cmax_cells) const;
V V
effectiveRelPerm(const V& c, const V& cmax_cells, const V& relperm) const; effectiveRelPerm(const V& c, const V& cmax_cells, const V& relperm) const;