diff --git a/opm/core/simulator/WellState.hpp b/opm/core/simulator/WellState.hpp index 3915ed0e0..b8b7a31d0 100644 --- a/opm/core/simulator/WellState.hpp +++ b/opm/core/simulator/WellState.hpp @@ -31,31 +31,48 @@ namespace Opm { public: /// Allocate and initialize if wells is non-null. + /// Also tries to give useful initial values to the bhp() and + /// wellRates() fields, depending on controls. The + /// perfRates() field is filled with zero, and perfPress() + /// with -1e100. template void init(const Wells* wells, const State& state) { if (wells) { const int nw = wells->number_of_wells; + const int np = wells->number_of_phases; bhp_.resize(nw); - // Initialize bhp to be target pressure - // if bhp-controlled well, otherwise set - // to pressure in first perforation cell. + wellrates_.resize(nw * np, 0.0); for (int w = 0; w < nw; ++w) { const WellControls* ctrl = wells->ctrls[w]; - + // Initialize bhp to be target pressure if + // bhp-controlled well, otherwise set to a little + // above or below (depending on if the well is an + // injector or producer) pressure in first perforation + // cell. if ((ctrl->current < 0) || // SHUT (ctrl->type[ctrl->current] != BHP)) { - const int cell = wells->well_cells[wells->well_connpos[w]]; + const int first_cell = wells->well_cells[wells->well_connpos[w]]; const double safety_factor = (wells->type[w] == INJECTOR) ? 1.01 : 0.99; - bhp_[w] = safety_factor*state.pressure()[cell]; - } - else { + bhp_[w] = safety_factor*state.pressure()[first_cell]; + } else { bhp_[w] = ctrl->target[ctrl->current]; } + // Initialize well rates to match controls if type is SURFACE_RATE + if ((ctrl->current >= 0) || // open well + (ctrl->type[ctrl->current] != SURFACE_RATE)) { + const double rate_target = ctrl->target[ctrl->current]; + for (int p = 0; p < np; ++p) { + const double phase_distr = ctrl->distr[np * ctrl->current + p]; + wellrates_[np*w + p] = rate_target * phase_distr; + } + } } + // The perforation rates and perforation pressures are + // not expected to be consistent with bhp_ and wellrates_ + // after init(). perfrates_.resize(wells->well_connpos[nw], 0.0); perfpress_.resize(wells->well_connpos[nw], -1e100); - wellrates_.resize(wells->well_connpos[nw] * wells->number_of_phases, 0.0); } }