From 139bfb7e625ae2604a7546cbcdedb88d11bcaffa Mon Sep 17 00:00:00 2001 From: Liu Ming Date: Fri, 10 Jan 2014 14:15:24 +0800 Subject: [PATCH] refactor well controls use new well controls gourp. --- ...FullyImplicitCompressiblePolymerSolver.cpp | 21 +++++++++++-------- .../FullyImplicitTwophasePolymerSolver.cpp | 19 ++++++++++------- opm/polymer/fullyimplicit/PolymerPropsAd.cpp | 9 +++++--- opm/polymer/fullyimplicit/PolymerPropsAd.hpp | 6 +++++- 4 files changed, 34 insertions(+), 21 deletions(-) diff --git a/opm/polymer/fullyimplicit/FullyImplicitCompressiblePolymerSolver.cpp b/opm/polymer/fullyimplicit/FullyImplicitCompressiblePolymerSolver.cpp index 672a84077..96ce6bcdc 100644 --- a/opm/polymer/fullyimplicit/FullyImplicitCompressiblePolymerSolver.cpp +++ b/opm/polymer/fullyimplicit/FullyImplicitCompressiblePolymerSolver.cpp @@ -31,7 +31,7 @@ #include #include #include - +#include #include #include #include @@ -590,7 +590,7 @@ namespace { const ADB perf_mob = producer.select(subset(cell_mob, well_cells), perf_mob_injector); const ADB perf_flux = perf_mob * (nkgradp_well); // No gravity term for perforations. - well_perf_rates[phase] = (perf_flux*perf_b); + well_perf_rates[phase] = (perf_flux * perf_b); const ADB well_rates = wops_.p2w * well_perf_rates[phase]; well_rates_all += superset(well_rates, Span(nw, 1, phase*nw), nw*np); @@ -617,15 +617,18 @@ namespace { M rate_distr(nw, np*nw); for (int w = 0; w < nw; ++w) { const WellControls* wc = wells_.ctrls[w]; - if (wc->type[wc->current] == BHP) { - bhp_targets[w] = wc->target[wc->current]; + if (well_controls_get_current_type(wc) == BHP) { + bhp_targets[w] = well_controls_get_current_target(wc); rate_targets[w] = -1e100; - } else if (wc->type[wc->current] == SURFACE_RATE) { + } else if (well_controls_get_current_type(wc) == SURFACE_RATE) { bhp_targets[w] = -1e100; - rate_targets[w] = wc->target[wc->current]; - for (int phase = 0; phase < np; ++phase) { - rate_distr.insert(w, phase*nw + w) = wc->distr[phase]; - } + rate_targets[w] = well_controls_get_current_target(wc); + { + const double* distr = well_controls_get_current_distr(wc); + for (int phase = 0; phase < np; ++phase) { + rate_distr.insert(w, phase*nw + w) = distr[phase]; + } + } } else { OPM_THROW(std::runtime_error, "Can only handle BHP and SURFACE_RATE type controls."); } diff --git a/opm/polymer/fullyimplicit/FullyImplicitTwophasePolymerSolver.cpp b/opm/polymer/fullyimplicit/FullyImplicitTwophasePolymerSolver.cpp index d869b3133..98a78deee 100644 --- a/opm/polymer/fullyimplicit/FullyImplicitTwophasePolymerSolver.cpp +++ b/opm/polymer/fullyimplicit/FullyImplicitTwophasePolymerSolver.cpp @@ -15,7 +15,7 @@ #include #include #include - +#include #include #include #include @@ -499,15 +499,18 @@ namespace { M rate_distr(nw, np*nw); for (int w = 0; w < nw; ++w) { const WellControls* wc = wells_.ctrls[w]; - if (wc->type[wc->current] == BHP) { - bhp_targets[w] = wc->target[wc->current]; + if (well_controls_get_current_type(wc) == BHP) { + bhp_targets[w] = well_controls_get_current_target(wc); rate_targets[w] = -1e100; - } else if (wc->type[wc->current] == SURFACE_RATE) { + } else if (well_controls_get_current_type(wc) == SURFACE_RATE) { bhp_targets[w] = -1e100; - rate_targets[w] = wc->target[wc->current]; - for (int phase = 0; phase < np; ++phase) { - rate_distr.insert(w, phase*nw + w) = wc->distr[phase]; - } + rate_targets[w] = well_controls_get_current_target(wc); + { + const double* distr = well_controls_get_current_distr(wc); + for (int phase = 0; phase < np; ++phase) { + rate_distr.insert(w, phase*nw + w) = distr[phase]; + } + } } else { OPM_THROW(std::runtime_error, "Can only handle BHP type controls."); } diff --git a/opm/polymer/fullyimplicit/PolymerPropsAd.cpp b/opm/polymer/fullyimplicit/PolymerPropsAd.cpp index 37ce0cfb0..89b657f17 100644 --- a/opm/polymer/fullyimplicit/PolymerPropsAd.cpp +++ b/opm/polymer/fullyimplicit/PolymerPropsAd.cpp @@ -1,4 +1,3 @@ - #include #include #include @@ -170,6 +169,8 @@ namespace Opm { } + + double PolymerPropsAd::deadPoreVol() const { @@ -360,8 +361,10 @@ namespace Opm { double res_factor = polymer_props_.resFactor(); double factor = (res_factor - 1.) / max_ads; ADB rk = one + ads * factor; - ADB dkrw_ds = krw / rk.value(); - ADB dkrw_dc = -krw.value() / (rk.value() * rk.value()) * ads * factor; +// ADB dkrw_ds = krw / rk.value(); + ADB dkrw_ds = krw / rk; +// ADB dkrw_dc = -krw.value() / (rk.value() * rk.value()) * ads * factor; + ADB dkrw_dc = -factor * krw / (rk * rk) * ads ; const int num_blocks = c.numBlocks(); std::vector jacs(num_blocks); diff --git a/opm/polymer/fullyimplicit/PolymerPropsAd.hpp b/opm/polymer/fullyimplicit/PolymerPropsAd.hpp index d0d72caaa..f3396bbc5 100644 --- a/opm/polymer/fullyimplicit/PolymerPropsAd.hpp +++ b/opm/polymer/fullyimplicit/PolymerPropsAd.hpp @@ -65,8 +65,10 @@ namespace Opm { */ double rockDensity() const; double deadPoreVol() const; - typedef AutoDiffBlock ADB; + + typedef AutoDiffBlock ADB; typedef ADB::V V; + PolymerPropsAd(const PolymerProperties& polymer_props); ~PolymerPropsAd(); @@ -94,6 +96,8 @@ namespace Opm { ADB effectiveRelPerm(const ADB& c, const ADB& cmax_cells, const ADB& krw, const ADB& sw) const; + + private: const PolymerProperties& polymer_props_; };