diff --git a/opm/core/pressure/tpfa/cfs_tpfa_residual.c b/opm/core/pressure/tpfa/cfs_tpfa_residual.c index 001f9a3bf..4c2ca10ed 100644 --- a/opm/core/pressure/tpfa/cfs_tpfa_residual.c +++ b/opm/core/pressure/tpfa/cfs_tpfa_residual.c @@ -5,7 +5,9 @@ #include #include + #include +#include #include #include @@ -806,7 +808,7 @@ welleq_coeff_resv(int np, struct cfs_tpfa_res_data *h, pflux = h->pimpl->flux_work; dpflux_w = pflux + (1 * np); dpflux_c = dpflux_w + (1 * np); - distr = ctrl->distr + (ctrl->current * np); + distr = well_controls_get_current_distr( ctrl ); *res = *w2c = *w2w = 0.0; for (p = 0; p < np; p++) { @@ -833,7 +835,7 @@ welleq_coeff_surfrate(int i, int np, struct cfs_tpfa_res_data *h, pflux = h->pimpl->compflux_p + (i * (1 * np)); dpflux_w = h->pimpl->compflux_deriv_p + (i * (2 * np)); dpflux_c = dpflux_w + (1 * (1 * np)); - distr = ctrl->distr + (ctrl->current * (1 * np)); + distr = well_controls_get_current_distr( ctrl ); *res = *w2c = *w2w = 0.0; for (p = 0; p < np; p++) { @@ -867,14 +869,14 @@ assemble_completion_to_well(int i, int w, int c, int nc, int np, W = wells->W; ctrl = W->ctrls[ w ]; - if (ctrl->current < 0) { + if (well_controls_get_current(ctrl) < 0) { /* Interpreting a negative current control index to mean a shut well */ welleq_coeff_shut(np, h, &res, &w2c, &w2w); } else { - switch (ctrl->type[ ctrl->current ]) { + switch (well_controls_get_current_type(ctrl)) { case BHP : - welleq_coeff_bhp(np, pw - ctrl->target[ ctrl->current ], + welleq_coeff_bhp(np, pw - well_controls_get_current_target( ctrl ), h, &res, &w2c, &w2w); break; @@ -931,7 +933,7 @@ assemble_well_contrib(struct cfs_tpfa_res_wells *wells , for (w = i = 0; w < W->number_of_wells; w++) { pw = wpress[ w ]; - is_open = W->ctrls[w]->current >= 0; + is_open = (well_controls_get_current(W->ctrls[w]) >= 0); for (; i < W->well_connpos[w + 1]; i++, pmobp += np) { @@ -956,10 +958,10 @@ assemble_well_contrib(struct cfs_tpfa_res_wells *wells , } ctrl = W->ctrls[ w ]; - if ((ctrl->current >= 0) && /* OPEN? */ - (ctrl->type[ ctrl->current ] != BHP)) { - - h->F[ nc + w ] -= dt * ctrl->target[ ctrl->current ]; + if ((well_controls_get_current(ctrl) >= 0) && /* OPEN? */ + (well_controls_get_current_type(ctrl) != BHP)) { + + h->F[ nc + w ] -= dt * well_controls_get_current_target(ctrl); } else { is_neumann = 0; diff --git a/opm/core/pressure/tpfa/ifs_tpfa.c b/opm/core/pressure/tpfa/ifs_tpfa.c index e3f4506c2..5130ba0b3 100644 --- a/opm/core/pressure/tpfa/ifs_tpfa.c +++ b/opm/core/pressure/tpfa/ifs_tpfa.c @@ -8,6 +8,7 @@ #include #include +#include #include #include @@ -229,7 +230,7 @@ assemble_bhp_well(int nc, int w, ctrls = W->ctrls[ w ]; wdof = nc + w; - bhp = ctrls->target[ ctrls->current ]; + bhp = well_controls_get_current_target(ctrls); jw = csrmatrix_elm_index(wdof, wdof, h->A); @@ -268,7 +269,7 @@ assemble_rate_well(int nc, int w, ctrls = W->ctrls[ w ]; wdof = nc + w; - resv = ctrls->target[ ctrls->current ]; + resv = well_controls_get_current_target(ctrls); jww = csrmatrix_elm_index(wdof, wdof, h->A); @@ -362,7 +363,7 @@ assemble_well_contrib(int nc , for (w = 0; w < W->number_of_wells; w++) { ctrls = W->ctrls[ w ]; - if (ctrls->current < 0) { + if (well_controls_get_current(ctrls) < 0) { /* Treat this well as a shut well, isolated from the domain. */ @@ -370,9 +371,9 @@ assemble_well_contrib(int nc , } else { - assert (ctrls->current < ctrls->num); + assert (well_controls_get_current(ctrls) < well_controls_get_num(ctrls)); - switch (ctrls->type[ ctrls->current ]) { + switch (well_controls_get_current_type(ctrls)) { case BHP: *all_rate = 0; assemble_bhp_well (nc, w, W, mt, wdp, h); @@ -383,8 +384,9 @@ assemble_well_contrib(int nc , /* Ensure minimal consistency. A PRODUCER should * specify a phase distribution of all ones in the * case of RESV controls. */ + const double * distr = well_controls_get_current_distr( ctrls ); for (p = 0; p < np; p++) { - if (ctrls->distr[np * ctrls->current + p] != 1.0) { + if (distr[p] != 1.0) { *ok = 0; break; } @@ -550,7 +552,7 @@ well_solution(const struct UnstructuredGrid *G , if (soln->well_press != NULL) { /* Extract BHP directly from solution vector for non-shut wells */ for (w = 0; w < F->W->number_of_wells; w++) { - if (F->W->ctrls[w]->current >= 0) { + if (well_controls_get_current(F->W->ctrls[w]) >= 0) { soln->well_press[w] = h->x[G->number_of_cells + w]; } } @@ -563,8 +565,8 @@ well_solution(const struct UnstructuredGrid *G , for (w = i = 0; w < F->W->number_of_wells; w++) { bhp = h->x[G->number_of_cells + w]; - shut = (F->W->ctrls[w]->current < 0); - + shut = (well_controls_get_current(F->W->ctrls[w]) < 0); + for (; i < F->W->well_connpos[ w + 1 ]; i++) { c = F->W->well_cells[ i ]; diff --git a/opm/core/simulator/WellState.hpp b/opm/core/simulator/WellState.hpp index 9849bd680..51d2c3b98 100644 --- a/opm/core/simulator/WellState.hpp +++ b/opm/core/simulator/WellState.hpp @@ -21,6 +21,7 @@ #define OPM_WELLSTATE_HEADER_INCLUDED #include +#include #include namespace Opm @@ -50,21 +51,22 @@ namespace Opm // 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)) { + if ((well_controls_get_current(ctrl) < 0) || // SHUT + (well_controls_get_current_type(ctrl) != BHP)) { 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()[first_cell]; } else { - bhp_[w] = ctrl->target[ctrl->current]; + bhp_[w] = well_controls_get_current_target( ctrl ); } + // 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]; + if ((well_controls_get_current(ctrl) >= 0) && // open well + (well_controls_get_current_type(ctrl) == SURFACE_RATE)) { + const double rate_target = well_controls_get_current_target(ctrl); + const double * distr = well_controls_get_current_distr( ctrl ); 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; + wellrates_[np*w + p] = rate_target * distr[p]; } } } diff --git a/opm/core/utility/miscUtilities.cpp b/opm/core/utility/miscUtilities.cpp index 6ed26f1a8..92dffaef9 100644 --- a/opm/core/utility/miscUtilities.cpp +++ b/opm/core/utility/miscUtilities.cpp @@ -22,6 +22,7 @@ #include #include #include +#include #include #include #include @@ -457,22 +458,25 @@ namespace Opm } src.resize(num_cells); for (int w = 0; w < wells.number_of_wells; ++w) { - const int cur = wells.ctrls[w]->current; - if (wells.ctrls[w]->num != 1) { + const int cur = well_controls_get_current(wells.ctrls[w]); + if (well_controls_get_num(wells.ctrls[w]) != 1) { OPM_MESSAGE("In wellsToSrc(): well has more than one control, all but current control will be ignored."); } - if (wells.ctrls[w]->type[cur] != RESERVOIR_RATE) { + if (well_controls_iget_type(wells.ctrls[w] , cur) != RESERVOIR_RATE) { OPM_THROW(std::runtime_error, "In wellsToSrc(): well is something other than RESERVOIR_RATE."); } if (wells.well_connpos[w+1] - wells.well_connpos[w] != 1) { OPM_THROW(std::runtime_error, "In wellsToSrc(): well has multiple perforations."); } - for (int p = 0; p < np; ++p) { - if (wells.ctrls[w]->distr[np*cur + p] != 1.0) { - OPM_THROW(std::runtime_error, "In wellsToSrc(): well not controlled on total rate."); + { + const double * distr = well_controls_iget_distr( wells.ctrls[w] , cur); + for (int p = 0; p < np; ++p) { + if (distr[p] != 1.0) { + OPM_THROW(std::runtime_error, "In wellsToSrc(): well not controlled on total rate."); + } } } - double flow = wells.ctrls[w]->target[cur]; + double flow = well_controls_iget_target(wells.ctrls[w] , cur); if (wells.type[w] == INJECTOR) { flow *= wells.comp_frac[np*w + 0]; // Obtaining water rate for inflow source. } diff --git a/opm/core/well_controls.h b/opm/core/well_controls.h new file mode 100644 index 000000000..3e0626707 --- /dev/null +++ b/opm/core/well_controls.h @@ -0,0 +1,108 @@ +/* + Copyright 2012 SINTEF ICT, Applied Mathematics. + + This file is part of the Open Porous Media project (OPM). + + OPM is free software: you can redistribute it and/or modify + it under the terms of the GNU General Public License as published by + the Free Software Foundation, either version 3 of the License, or + (at your option) any later version. + + OPM is distributed in the hope that it will be useful, + but WITHOUT ANY WARRANTY; without even the implied warranty of + MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + GNU General Public License for more details. + + You should have received a copy of the GNU General Public License + along with OPM. If not, see . +*/ + +#ifndef OPM_WELL_CONTROLS_H_INCLUDED +#define OPM_WELL_CONTROLS_H_INCLUDED + +#include + + + +#ifdef __cplusplus +extern "C" { +#endif + +enum WellControlType { + BHP, /**< Well constrained by BHP target */ + RESERVOIR_RATE, /**< Well constrained by reservoir volume flow rate */ + SURFACE_RATE /**< Well constrained by surface volume flow rate */ +}; + +struct WellControls; + +bool +well_controls_equal(const struct WellControls *ctrls1, const struct WellControls *ctrls2); + +//int +//well_controls_reserve(int nctrl, int nphases, struct WellControls *ctrl); + +struct WellControls * +well_controls_create(void); + +void +well_controls_destroy(struct WellControls *ctrl); + + +int +well_controls_get_num(const struct WellControls *ctrl); + +int +well_controls_get_cpty(const struct WellControls *ctrl); + +int +well_controls_get_current( const struct WellControls * ctrl); + +void +well_controls_set_current( struct WellControls * ctrl, int current); + +void +well_controls_invert_current( struct WellControls * ctrl ); + +int +well_controls_add_new(enum WellControlType type , double target , const double * distr , struct WellControls * ctrl); + +enum WellControlType +well_controls_iget_type(const struct WellControls * ctrl, int control_index); + +enum WellControlType +well_controls_get_current_type(const struct WellControls * ctrl); + +void +well_controls_iset_type( struct WellControls * ctrls , int control_index , enum WellControlType type); + +void +well_controls_iset_target(struct WellControls * ctrl, int control_index , double target); + +double +well_controls_iget_target(const struct WellControls * ctrl, int control_index); + +double +well_controls_get_current_target(const struct WellControls * ctrl); + +const double * +well_controls_iget_distr(const struct WellControls * ctrl, int control_index); + +void +well_controls_iset_distr(const struct WellControls * ctrl, int control_index, const double * distr); + +const double * +well_controls_get_current_distr(const struct WellControls * ctrl); + +void +well_controls_assert_number_of_phases(struct WellControls * ctrl , int number_of_phases); + +void +well_controls_clear(struct WellControls * ctrl); + + +#ifdef __cplusplus +} +#endif + +#endif /* OPM_WELL_CONTROLS_H_INCLUDED */ diff --git a/opm/core/wells.h b/opm/core/wells.h index c58e06613..dfc9be176 100644 --- a/opm/core/wells.h +++ b/opm/core/wells.h @@ -20,6 +20,8 @@ #ifndef OPM_WELLS_H_INCLUDED #define OPM_WELLS_H_INCLUDED +#include +#include /** * \file @@ -41,75 +43,6 @@ enum WellType { PRODUCER /**< Well is a producer */ }; -/** - * Type of well control equation or inequality constraint. - */ -enum WellControlType { - BHP, /**< Well constrained by BHP target */ - RESERVOIR_RATE, /**< Well constrained by reservoir volume flow rate */ - SURFACE_RATE /**< Well constrained by surface volume flow rate */ -}; - -/** - * Controls for a single well. - * Each control specifies a well rate or bottom-hole pressure. Only - * one control can be active at a time, indicated by current. The - * meaning of each control's target value depends on the control type: - * - * - BHP -> target pressure in Pascal. - * - RESERVOIR_RATE -> target reservoir volume rate in cubic(meter)/second - * - SURFACE_RATE -> target surface volume rate in cubic(meter)/second - * - * The sign convention for RATE targets is as follows: - * - * - (+) Fluid flowing into reservoir, i.e. injecting. - * - (-) Fluid flowing out of reservoir, i.e. producing. - * - * For *_RATE controls, the distribution of phases used for the control - * is also needed. For example, a total rate control should have 1.0 - * for each phase, whereas a control on oil rate should have 1.0 for - * the oil phase and 0.0 for the rest. For BHP controls, this is unused. - * The active control acts as an equality constraint, whereas the - * non-active controls should be interpreted as inequality - * constraints (upper or lower bounds). For instance, a PRODUCER's - * BHP constraint defines a minimum acceptable bottom-hole pressure - * value for the well. - */ -struct WellControls -{ - /** - * Number of controls. - */ - int num; - - /** - * Array of control types. - */ - enum WellControlType *type; - - /** - * Array of control targets. - */ - double *target; - - /** - * Array of rate control distributions, - * Wells::number_of_phases numbers for each control - */ - double *distr; - - /** - * Index of current active control. - */ - int current; - - /** - * Internal management structure. - */ - void *data; -}; - - /** * Data structure aggregating static information about all wells in a scenario. @@ -282,6 +215,8 @@ add_well(enum WellType type , * \param[in,out] W Existing set of well controls. * \return Non-zero (true) if successful and zero (false) otherwise. */ + + int append_well_controls(enum WellControlType type , double target, @@ -327,6 +262,10 @@ destroy_wells(struct Wells *W); struct Wells * clone_wells(const struct Wells *W); +bool +wells_equal(const struct Wells *W1, const struct Wells *W2); + + #ifdef __cplusplus } diff --git a/opm/core/wells/WellsGroup.cpp b/opm/core/wells/WellsGroup.cpp index d65555ed3..6c1059a29 100644 --- a/opm/core/wells/WellsGroup.cpp +++ b/opm/core/wells/WellsGroup.cpp @@ -1,6 +1,5 @@ /* Copyright 2012 SINTEF ICT, Applied Mathematics. - This file is part of the Open Porous Media project (OPM). OPM is free software: you can redistribute it and/or modify @@ -20,6 +19,7 @@ #include "config.h" #include #include +#include #include #include @@ -655,18 +655,19 @@ namespace Opm // Check constraints. bool is_producer = (wells_->type[self_index_] == PRODUCER); - const WellControls& ctrls = *wells_->ctrls[self_index_]; - for (int ctrl_index = 0; ctrl_index < ctrls.num; ++ctrl_index) { - if (ctrl_index == ctrls.current || ctrl_index == group_control_index_) { + const WellControls * ctrls = wells_->ctrls[self_index_]; + for (int ctrl_index = 0; ctrl_index < well_controls_get_num(ctrls); ++ctrl_index) { + if (ctrl_index == well_controls_get_current(ctrls) || ctrl_index == group_control_index_) { // We do not check constraints that either were used // as the active control, or that come from group control. continue; } bool ctrl_violated = false; - switch (ctrls.type[ctrl_index]) { + switch (well_controls_iget_type(ctrls , ctrl_index)) { + case BHP: { const double my_well_bhp = well_bhp[self_index_]; - const double my_target_bhp = ctrls.target[ctrl_index]; + const double my_target_bhp = well_controls_iget_target( ctrls , ctrl_index); ctrl_violated = is_producer ? (my_target_bhp > my_well_bhp) : (my_target_bhp < my_well_bhp); if (ctrl_violated) { @@ -676,12 +677,14 @@ namespace Opm } break; } + case RESERVOIR_RATE: { double my_rate = 0.0; + const double * ctrls_distr = well_controls_iget_distr( ctrls , ctrl_index ); for (int phase = 0; phase < np; ++phase) { - my_rate += ctrls.distr[np*ctrl_index + phase]*well_reservoirrates_phase[np*self_index_ + phase]; + my_rate += ctrls_distr[phase] * well_reservoirrates_phase[np*self_index_ + phase]; } - const double my_rate_target = ctrls.target[ctrl_index]; + const double my_rate_target = well_controls_iget_target(ctrls , ctrl_index); ctrl_violated = std::fabs(my_rate) - std::fabs(my_rate_target)> std::max(std::abs(my_rate), std::abs(my_rate_target))*1e-6; if (ctrl_violated) { std::cout << "RESERVOIR_RATE limit violated for well " << name() << ":\n"; @@ -690,12 +693,14 @@ namespace Opm } break; } + case SURFACE_RATE: { double my_rate = 0.0; + const double * ctrls_distr = well_controls_iget_distr( ctrls , ctrl_index ); for (int phase = 0; phase < np; ++phase) { - my_rate += ctrls.distr[np*ctrl_index + phase]*well_surfacerates_phase[np*self_index_ + phase]; + my_rate += ctrls_distr[phase] * well_surfacerates_phase[np*self_index_ + phase]; } - const double my_rate_target = ctrls.target[ctrl_index]; + const double my_rate_target = well_controls_iget_target(ctrls , ctrl_index); ctrl_violated = std::fabs(my_rate) > std::fabs(my_rate_target); if (ctrl_violated) { std::cout << "SURFACE_RATE limit violated for well " << name() << ":\n"; @@ -742,9 +747,9 @@ namespace Opm void WellNode::shutWell() { if (shut_well_) { - // We set the tilde of the current control + // We set the tilde of the current control // set_current_control(self_index_, -1, wells_); - wells_->ctrls[self_index_]->current = ~ wells_->ctrls[self_index_]->current; + well_controls_invert_current(wells_->ctrls[self_index_]); } else { const double target = 0.0; @@ -753,16 +758,16 @@ namespace Opm if (group_control_index_ < 0) { // The well only had its own controls, no group controls. append_well_controls(SURFACE_RATE, target, distr, self_index_, wells_); - group_control_index_ = wells_->ctrls[self_index_]->num - 1; + group_control_index_ = well_controls_get_num(wells_->ctrls[self_index_]) - 1; } else { // We will now modify the last control, that // "belongs to" the group control. - const int np = wells_->number_of_phases; - wells_->ctrls[self_index_]->type[group_control_index_] = SURFACE_RATE; - wells_->ctrls[self_index_]->target[group_control_index_] = target; - std::copy(distr, distr + np, wells_->ctrls[self_index_]->distr + np * group_control_index_); + + well_controls_iset_type( wells_->ctrls[self_index_] , group_control_index_ , SURFACE_RATE); + well_controls_iset_target( wells_->ctrls[self_index_] , group_control_index_ , target); + well_controls_iset_distr(wells_->ctrls[self_index_] , group_control_index_ , distr); } - wells_->ctrls[self_index_]->current = ~ wells_->ctrls[self_index_]->current; + well_controls_invert_current(wells_->ctrls[self_index_]); } } @@ -808,14 +813,13 @@ namespace Opm if (group_control_index_ < 0) { // The well only had its own controls, no group controls. append_well_controls(wct, target, distr, self_index_, wells_); - group_control_index_ = wells_->ctrls[self_index_]->num - 1; + group_control_index_ = well_controls_get_num(wells_->ctrls[self_index_]) - 1; } else { // We will now modify the last control, that // "belongs to" the group control. - const int np = wells_->number_of_phases; - wells_->ctrls[self_index_]->type[group_control_index_] = wct; - wells_->ctrls[self_index_]->target[group_control_index_] = target; - std::copy(distr, distr + np, wells_->ctrls[self_index_]->distr + np*group_control_index_); + well_controls_iset_type(wells_->ctrls[self_index_] , group_control_index_ , wct); + well_controls_iset_target(wells_->ctrls[self_index_] , group_control_index_ ,target); + well_controls_iset_distr(wells_->ctrls[self_index_] , group_control_index_ , distr); } set_current_control(self_index_, group_control_index_, wells_); } @@ -920,14 +924,13 @@ namespace Opm if (group_control_index_ < 0) { // The well only had its own controls, no group controls. append_well_controls(wct, ntarget, distr, self_index_, wells_); - group_control_index_ = wells_->ctrls[self_index_]->num - 1; + group_control_index_ = well_controls_get_num(wells_->ctrls[self_index_]) - 1; } else { // We will now modify the last control, that // "belongs to" the group control. - const int np = wells_->number_of_phases; - wells_->ctrls[self_index_]->type[group_control_index_] = wct; - wells_->ctrls[self_index_]->target[group_control_index_] = ntarget; - std::copy(distr, distr + np, wells_->ctrls[self_index_]->distr + np*group_control_index_); + well_controls_iset_type(wells_->ctrls[self_index_] , group_control_index_ , wct); + well_controls_iset_target(wells_->ctrls[self_index_] , group_control_index_ , ntarget); + well_controls_iset_distr(wells_->ctrls[self_index_] , group_control_index_ , distr); } set_current_control(self_index_, group_control_index_, wells_); } diff --git a/opm/core/wells/WellsManager.cpp b/opm/core/wells/WellsManager.cpp index 82ac19128..007c8b52c 100644 --- a/opm/core/wells/WellsManager.cpp +++ b/opm/core/wells/WellsManager.cpp @@ -22,6 +22,7 @@ #include #include #include +#include #include #include #include @@ -505,7 +506,7 @@ namespace Opm int ok = 1; int control_pos[5] = { -1, -1, -1, -1, -1 }; if (ok && wci_line.surface_flow_max_rate_ >= 0.0) { - control_pos[InjectionControl::RATE] = w_->ctrls[wix]->num; + control_pos[InjectionControl::RATE] = well_controls_get_num(w_->ctrls[wix]); double distr[3] = { 0.0, 0.0, 0.0 }; if (wci_line.injector_type_ == "WATER") { distr[pu.phase_pos[BlackoilPhases::Aqua]] = 1.0; @@ -521,7 +522,7 @@ namespace Opm distr, wix, w_); } if (ok && wci_line.reservoir_flow_max_rate_ >= 0.0) { - control_pos[InjectionControl::RESV] = w_->ctrls[wix]->num; + control_pos[InjectionControl::RESV] = well_controls_get_num(w_->ctrls[wix]); double distr[3] = { 0.0, 0.0, 0.0 }; if (wci_line.injector_type_ == "WATER") { distr[pu.phase_pos[BlackoilPhases::Aqua]] = 1.0; @@ -537,7 +538,7 @@ namespace Opm distr, wix, w_); } if (ok && wci_line.BHP_limit_ > 0.0) { - control_pos[InjectionControl::BHP] = w_->ctrls[wix]->num; + control_pos[InjectionControl::BHP] = well_controls_get_num(w_->ctrls[wix]); ok = append_well_controls(BHP, wci_line.BHP_limit_, NULL, wix, w_); } @@ -618,7 +619,7 @@ namespace Opm if (!pu.phase_used[BlackoilPhases::Liquid]) { OPM_THROW(std::runtime_error, "Oil phase not active and ORAT control specified."); } - control_pos[ProductionControl::ORAT] = w_->ctrls[wix]->num; + control_pos[ProductionControl::ORAT] = well_controls_get_num(w_->ctrls[wix]); double distr[3] = { 0.0, 0.0, 0.0 }; distr[pu.phase_pos[BlackoilPhases::Liquid]] = 1.0; ok = append_well_controls(SURFACE_RATE, -wcp_line.oil_max_rate_, @@ -628,7 +629,7 @@ namespace Opm if (!pu.phase_used[BlackoilPhases::Aqua]) { OPM_THROW(std::runtime_error, "Water phase not active and WRAT control specified."); } - control_pos[ProductionControl::WRAT] = w_->ctrls[wix]->num; + control_pos[ProductionControl::WRAT] = well_controls_get_num(w_->ctrls[wix]); double distr[3] = { 0.0, 0.0, 0.0 }; distr[pu.phase_pos[BlackoilPhases::Aqua]] = 1.0; ok = append_well_controls(SURFACE_RATE, -wcp_line.water_max_rate_, @@ -638,7 +639,7 @@ namespace Opm if (!pu.phase_used[BlackoilPhases::Vapour]) { OPM_THROW(std::runtime_error, "Gas phase not active and GRAT control specified."); } - control_pos[ProductionControl::GRAT] = w_->ctrls[wix]->num; + control_pos[ProductionControl::GRAT] = well_controls_get_num(w_->ctrls[wix]); double distr[3] = { 0.0, 0.0, 0.0 }; distr[pu.phase_pos[BlackoilPhases::Vapour]] = 1.0; ok = append_well_controls(SURFACE_RATE, -wcp_line.gas_max_rate_, @@ -651,7 +652,7 @@ namespace Opm if (!pu.phase_used[BlackoilPhases::Liquid]) { OPM_THROW(std::runtime_error, "Oil phase not active and LRAT control specified."); } - control_pos[ProductionControl::LRAT] = w_->ctrls[wix]->num; + control_pos[ProductionControl::LRAT] = well_controls_get_num(w_->ctrls[wix]); double distr[3] = { 0.0, 0.0, 0.0 }; distr[pu.phase_pos[BlackoilPhases::Aqua]] = 1.0; distr[pu.phase_pos[BlackoilPhases::Liquid]] = 1.0; @@ -659,13 +660,13 @@ namespace Opm distr, wix, w_); } if (ok && wcp_line.reservoir_flow_max_rate_ >= 0.0) { - control_pos[ProductionControl::RESV] = w_->ctrls[wix]->num; + control_pos[ProductionControl::RESV] = well_controls_get_num(w_->ctrls[wix]); double distr[3] = { 1.0, 1.0, 1.0 }; ok = append_well_controls(RESERVOIR_RATE, -wcp_line.reservoir_flow_max_rate_, distr, wix, w_); } if (ok && wcp_line.BHP_limit_ > 0.0) { - control_pos[ProductionControl::BHP] = w_->ctrls[wix]->num; + control_pos[ProductionControl::BHP] = well_controls_get_num(w_->ctrls[wix]); ok = append_well_controls(BHP, wcp_line.BHP_limit_, NULL, wix, w_); } @@ -754,17 +755,17 @@ namespace Opm } const int index = it->second; if (line.openshutflag_ == "SHUT") { - int& cur_ctrl = w_->ctrls[index]->current; + int cur_ctrl = well_controls_get_current(w_->ctrls[index]); if (cur_ctrl >= 0) { - cur_ctrl = ~cur_ctrl; + well_controls_invert_current(w_->ctrls[index]); } - assert(w_->ctrls[index]->current < 0); + assert(well_controls_get_current(w_->ctrls[index]) < 0); } else if (line.openshutflag_ == "OPEN") { - int& cur_ctrl = w_->ctrls[index]->current; + int cur_ctrl = well_controls_get_current(w_->ctrls[index]); if (cur_ctrl < 0) { - cur_ctrl = ~cur_ctrl; + well_controls_invert_current(w_->ctrls[index]); } - assert(w_->ctrls[index]->current >= 0); + assert(well_controls_get_current(w_->ctrls[index]) >= 0); } else { OPM_THROW(std::runtime_error, "Unknown Open/close keyword: \"" << line.openshutflag_<< "\". Allowed values: OPEN, SHUT."); } diff --git a/opm/core/wells/well_controls.c b/opm/core/wells/well_controls.c new file mode 100644 index 000000000..46eae54da --- /dev/null +++ b/opm/core/wells/well_controls.c @@ -0,0 +1,318 @@ +/*=========================================================================== +// +// File: newwells.c +// +// Created: 2012-02-03 11:28:40+0100 +// +// Authors: Knut-Andreas Lie +// Jostein R. Natvig +// Halvor M. Nilsen +// Atgeirr F. Rasmussen +// Xavier Raynaud +// Bård Skaflestad +// +//==========================================================================*/ + + +/* + Copyright 2012 SINTEF ICT, Applied Mathematics. + Copyright 2012 Statoil ASA. + + This file is part of the Open Porous Media Project (OPM). + + OPM is free software: you can redistribute it and/or modify + it under the terms of the GNU General Public License as published by + the Free Software Foundation, either version 3 of the License, or + (at your option) any later version. + + OPM is distributed in the hope that it will be useful, + but WITHOUT ANY WARRANTY; without even the implied warranty of + MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + GNU General Public License for more details. + + You should have received a copy of the GNU General Public License + along with OPM. If not, see . +*/ + +#include "config.h" + +#include +#include +#include +#include +#include + +/** + * Controls for a single well. + * Each control specifies a well rate or bottom-hole pressure. Only + * one control can be active at a time, indicated by current. The + * meaning of each control's target value depends on the control type: + * + * - BHP -> target pressure in Pascal. + * - RESERVOIR_RATE -> target reservoir volume rate in cubic(meter)/second + * - SURFACE_RATE -> target surface volume rate in cubic(meter)/second + * + * The sign convention for RATE targets is as follows: + * + * - (+) Fluid flowing into reservoir, i.e. injecting. + * - (-) Fluid flowing out of reservoir, i.e. producing. + * + * For *_RATE controls, the distribution of phases used for the control + * is also needed. For example, a total rate control should have 1.0 + * for each phase, whereas a control on oil rate should have 1.0 for + * the oil phase and 0.0 for the rest. For BHP controls, this is unused. + * The active control acts as an equality constraint, whereas the + * non-active controls should be interpreted as inequality + * constraints (upper or lower bounds). For instance, a PRODUCER's + * BHP constraint defines a minimum acceptable bottom-hole pressure + * value for the well. + */ + +struct WellControls +{ + /** + * Number of controls. + */ + int num; + + int number_of_phases; + + /** + * Array of control types. + */ + enum WellControlType *type; + + /** + * Array of control targets. + */ + double *target; + + /** + * Array of rate control distributions, + * number_of_phases numbers for each control + */ + double *distr; + + /** + * Index of current active control. + */ + int current; + + /* + The capacity allocated. + */ + int cpty; +}; + + +/* ---------------------------------------------------------------------- */ +void +well_controls_destroy(struct WellControls *ctrl) +/* ---------------------------------------------------------------------- */ +{ + if (ctrl != NULL) { + free (ctrl->distr); + free (ctrl->target); + free (ctrl->type); + } + + free(ctrl); +} + + +/* ---------------------------------------------------------------------- */ +struct WellControls * +well_controls_create(void) +/* ---------------------------------------------------------------------- */ +{ + struct WellControls *ctrl; + + ctrl = malloc(1 * sizeof *ctrl); + + if (ctrl != NULL) { + /* Initialise empty control set */ + ctrl->num = 0; + ctrl->number_of_phases = 0; + ctrl->type = NULL; + ctrl->target = NULL; + ctrl->distr = NULL; + ctrl->current = -1; + ctrl->cpty = 0; + } + + return ctrl; +} + + +/* ---------------------------------------------------------------------- */ +static int +well_controls_reserve(int nctrl, struct WellControls *ctrl) +/* ---------------------------------------------------------------------- */ +{ + int c, p, ok; + void *type, *target, *distr; + + type = realloc(ctrl->type , nctrl * 1 * sizeof *ctrl->type ); + target = realloc(ctrl->target, nctrl * 1 * sizeof *ctrl->target); + distr = realloc(ctrl->distr , nctrl * ctrl->number_of_phases * sizeof *ctrl->distr ); + + ok = 0; + if (type != NULL) { ctrl->type = type ; ok++; } + if (target != NULL) { ctrl->target = target; ok++; } + if (distr != NULL) { ctrl->distr = distr ; ok++; } + + if (ok == 3) { + for (c = ctrl->cpty; c < nctrl; c++) { + ctrl->type [c] = BHP; + ctrl->target[c] = -1.0; + } + + for (p = ctrl->cpty * ctrl->number_of_phases; p < nctrl * ctrl->number_of_phases; ++p) { + ctrl->distr[ p ] = 0.0; + } + + ctrl->cpty = nctrl; + } + + return ok == 3; +} + + +int well_controls_get_num(const struct WellControls *ctrl) { + return ctrl->num; +} + + +int well_controls_get_cpty(const struct WellControls *ctrl) { + return ctrl->cpty; +} + + +int well_controls_get_current( const struct WellControls * ctrl) { + return ctrl->current; +} + +void +well_controls_set_current( struct WellControls * ctrl, int current) { + ctrl->current = current; +} + +void +well_controls_invert_current( struct WellControls * ctrl ) { + ctrl->current = ~ctrl->current; +} + + +enum WellControlType +well_controls_iget_type(const struct WellControls * ctrl, int control_index) { + return ctrl->type[control_index]; +} + + +enum WellControlType +well_controls_get_current_type(const struct WellControls * ctrl) { + return well_controls_iget_type( ctrl , ctrl->current); +} + + +void +well_controls_iset_type( struct WellControls * ctrls , int control_index , enum WellControlType type) { + ctrls->type[control_index] = type; +} + + +double +well_controls_iget_target(const struct WellControls * ctrl, int control_index) { + return ctrl->target[control_index]; +} + +double +well_controls_get_current_target(const struct WellControls * ctrl) { + return ctrl->target[ctrl->current]; +} + +void +well_controls_iset_target(struct WellControls * ctrl, int control_index , double target) { + ctrl->target[control_index] = target; +} + + +const double * +well_controls_iget_distr(const struct WellControls * ctrl, int control_index) { + int offset = control_index * ctrl->number_of_phases; + return &ctrl->distr[offset]; +} + + +const double * +well_controls_get_current_distr(const struct WellControls * ctrl) { + return well_controls_iget_distr( ctrl , ctrl->current ); +} + + + +void +well_controls_iset_distr(const struct WellControls * ctrl, int control_index, const double * distr) { + int offset = control_index * ctrl->number_of_phases; + for (int p=0; p < ctrl->number_of_phases; p++) + ctrl->distr[offset + p] = distr[p]; +} + + +void +well_controls_assert_number_of_phases(struct WellControls * ctrl , int number_of_phases) { + if (ctrl->num == 0) + ctrl->number_of_phases = number_of_phases; + + assert( ctrl->number_of_phases == number_of_phases ); +} + +void +well_controls_clear(struct WellControls * ctrl) { + ctrl->num = 0; + ctrl->number_of_phases = 0; +} + + + +int +well_controls_add_new(enum WellControlType type , double target , const double * distr , struct WellControls * ctrl) { + if (ctrl->num == ctrl->cpty) { + int new_cpty = 2*ctrl->cpty; + if (new_cpty == ctrl->num) + new_cpty += 1; + + if (!well_controls_reserve( new_cpty , ctrl)) + return 0; + } + + well_controls_iset_type( ctrl , ctrl->num , type); + well_controls_iset_target( ctrl , ctrl->num , target); + + if (distr != NULL) + well_controls_iset_distr( ctrl , ctrl->num , distr); + + ctrl->num += 1; + return 1; +} + + +bool +well_controls_equal(const struct WellControls *ctrls1, const struct WellControls *ctrls2) +/* ---------------------------------------------------------------------- */ +{ + bool are_equal = true; + are_equal = (ctrls1->num == ctrls2->num); + are_equal = are_equal && (ctrls1->number_of_phases == ctrls2->number_of_phases); + if (!are_equal) { + return are_equal; + } + + are_equal = are_equal && (memcmp(ctrls1->type, ctrls2->type, ctrls1->num * sizeof *ctrls1->type ) == 0); + are_equal = are_equal && (memcmp(ctrls1->target, ctrls2->target, ctrls1->num * sizeof *ctrls1->target ) == 0); + are_equal = are_equal && (memcmp(ctrls1->distr, ctrls2->distr, ctrls1->num * ctrls1->number_of_phases * sizeof *ctrls1->distr ) == 0); + are_equal = are_equal && (ctrls1->cpty == ctrls2->cpty); + + return are_equal; +} + diff --git a/opm/core/wells/wells.c b/opm/core/wells/wells.c index c1327f74b..8f2e2b199 100644 --- a/opm/core/wells/wells.c +++ b/opm/core/wells/wells.c @@ -35,15 +35,15 @@ */ #include "config.h" #include +#include #include #include #include #include -struct WellControlMgmt { - int cpty; -}; +#include + struct WellMgmt { int well_cpty; @@ -51,28 +51,6 @@ struct WellMgmt { }; -static void -destroy_ctrl_mgmt(struct WellControlMgmt *m) -{ - free(m); -} - - -static struct WellControlMgmt * -create_ctrl_mgmt(void) -{ - struct WellControlMgmt *m; - - m = malloc(1 * sizeof *m); - - if (m != NULL) { - m->cpty = 0; - } - - return m; -} - - static void destroy_well_mgmt(struct WellMgmt *m) { @@ -96,86 +74,6 @@ create_well_mgmt(void) } -/* ---------------------------------------------------------------------- */ -static void -well_controls_destroy(struct WellControls *ctrl) -/* ---------------------------------------------------------------------- */ -{ - if (ctrl != NULL) { - destroy_ctrl_mgmt(ctrl->data); - free (ctrl->distr); - free (ctrl->target); - free (ctrl->type); - } - - free(ctrl); -} - - -/* ---------------------------------------------------------------------- */ -static struct WellControls * -well_controls_create(void) -/* ---------------------------------------------------------------------- */ -{ - struct WellControls *ctrl; - - ctrl = malloc(1 * sizeof *ctrl); - - if (ctrl != NULL) { - /* Initialise empty control set */ - ctrl->num = 0; - ctrl->type = NULL; - ctrl->target = NULL; - ctrl->distr = NULL; - ctrl->current = -1; - - ctrl->data = create_ctrl_mgmt(); - - if (ctrl->data == NULL) { - well_controls_destroy(ctrl); - ctrl = NULL; - } - } - - return ctrl; -} - - -/* ---------------------------------------------------------------------- */ -static int -well_controls_reserve(int nctrl, int nphases, struct WellControls *ctrl) -/* ---------------------------------------------------------------------- */ -{ - int c, p, ok; - void *type, *target, *distr; - - struct WellControlMgmt *m; - - type = realloc(ctrl->type , nctrl * 1 * sizeof *ctrl->type ); - target = realloc(ctrl->target, nctrl * 1 * sizeof *ctrl->target); - distr = realloc(ctrl->distr , nctrl * nphases * sizeof *ctrl->distr ); - - ok = 0; - if (type != NULL) { ctrl->type = type ; ok++; } - if (target != NULL) { ctrl->target = target; ok++; } - if (distr != NULL) { ctrl->distr = distr ; ok++; } - - if (ok == 3) { - m = ctrl->data; - for (c = m->cpty; c < nctrl; c++) { - ctrl->type [c] = BHP; - ctrl->target[c] = -1.0; - } - - for (p = m->cpty * nphases; p < nctrl * nphases; ++p) { - ctrl->distr[ p ] = 0.0; - } - - m->cpty = nctrl; - } - - return ok == 3; -} /* ---------------------------------------------------------------------- */ @@ -254,7 +152,7 @@ initialise_new_wells(int nwells, struct Wells *W) } for (w = m->well_cpty, ok = 1; ok && (w < nwells); w++) { - W->ctrls[w] = well_controls_create(); + W->ctrls[w] = well_controls_create( ); ok = W->ctrls[w] != NULL; } @@ -529,39 +427,16 @@ append_well_controls(enum WellControlType type, struct Wells *W) /* ---------------------------------------------------------------------- */ { - int ok, alloc, np; struct WellControls *ctrl; - struct WellControlMgmt *m; assert (W != NULL); assert ((0 <= well_index) && (well_index < W->number_of_wells)); ctrl = W->ctrls[well_index]; - np = W->number_of_phases; - assert (ctrl != NULL); - - m = ctrl->data; - ok = ctrl->num < m->cpty; - - if (! ok) { - alloc = alloc_size(ctrl->num, 1, m->cpty); - ok = well_controls_reserve(alloc, np, ctrl); - } - - if (ok) { - ctrl->type [ctrl->num] = type ; - ctrl->target[ctrl->num] = target; - - if (distr != NULL) { - memcpy(ctrl->distr + (ctrl->num * np), distr, - np * sizeof *ctrl->distr); - } - - ctrl->num += 1; - } - - return ok; + + well_controls_assert_number_of_phases( ctrl , W->number_of_phases); + return well_controls_add_new(type , target , distr , ctrl); } @@ -572,12 +447,10 @@ set_current_control(int well_index, int current_control, struct Wells *W) { assert (W != NULL); assert ((0 <= well_index) && (well_index < W->number_of_wells)); - assert (W->ctrls[well_index] != NULL); - - assert (current_control < W->ctrls[well_index]->num); - - W->ctrls[well_index]->current = current_control; + assert (current_control < well_controls_get_num(W->ctrls[well_index])); + + well_controls_set_current(W->ctrls[well_index] , current_control); } @@ -589,12 +462,13 @@ clear_well_controls(int well_index, struct Wells *W) assert (W != NULL); assert ((0 <= well_index) && (well_index < W->number_of_wells)); - if (W->ctrls[well_index] != NULL) { - W->ctrls[well_index]->num = 0; - } + if (W->ctrls[well_index] != NULL) + well_controls_clear( W->ctrls[well_index] ); } + + /* ---------------------------------------------------------------------- */ struct Wells * clone_wells(const struct Wells *W) @@ -607,17 +481,17 @@ clone_wells(const struct Wells *W) enum WellControlType type; const struct WellControls *ctrls; - struct Wells *ret; + struct Wells *newWells; if (W == NULL) { - ret = NULL; + newWells = NULL; } else { np = W->number_of_phases; - ret = create_wells(W->number_of_phases, W->number_of_wells, - W->well_connpos[ W->number_of_wells ]); + newWells = create_wells(W->number_of_phases, W->number_of_wells, + W->well_connpos[ W->number_of_wells ]); - if (ret != NULL) { + if (newWells != NULL) { pos = W->well_connpos[ 0 ]; ok = 1; @@ -629,7 +503,7 @@ clone_wells(const struct Wells *W) comp_frac = W->comp_frac != NULL ? W->comp_frac + w*np : NULL; ok = add_well(W->type[ w ], W->depth_ref[ w ], nperf, - comp_frac, cells, WI, W->name[ w ], ret); + comp_frac, cells, WI, W->name[ w ], newWells); /* Capacity should be sufficient from create_wells() */ assert (ok); @@ -637,34 +511,75 @@ clone_wells(const struct Wells *W) ctrls = W->ctrls[ w ]; if (ok) { - ok = well_controls_reserve(ctrls->num, np, ret->ctrls[ w ]); + for (c = 0; ok && (c < well_controls_get_num(ctrls)); c++) { + type = well_controls_iget_type( ctrls , c ); + target = well_controls_iget_target( ctrls , c ); + distr = well_controls_iget_distr( ctrls , c ); - for (c = 0; ok && (c < ctrls->num); c++) { - type = ctrls->type [ c ]; - target = ctrls->target[ c ]; - distr = & ctrls->distr [ c * np ]; - - ok = append_well_controls(type, target, distr, w, ret); - - /* Capacity *should* be sufficient from - * well_controls_reserve() */ + ok = append_well_controls(type, target, distr, w, newWells); + assert (ok); } } if (ok) { - set_current_control(w, ctrls->current, ret); + set_current_control(w, well_controls_get_current( ctrls) , newWells); } pos = W->well_connpos[w + 1]; } if (! ok) { - destroy_wells(ret); - ret = NULL; + destroy_wells(newWells); + newWells = NULL; } } } - return ret; + return newWells; } + +/* ---------------------------------------------------------------------- */ +bool +wells_equal(const struct Wells *W1, const struct Wells *W2) +/* ---------------------------------------------------------------------- */ +{ + bool are_equal = true; + are_equal = (W1->number_of_wells == W2->number_of_wells); + are_equal = are_equal && (W1->number_of_phases == W2->number_of_phases); + if (!are_equal) { + return are_equal; + } + + for (int i=0; inumber_of_wells; i++) { + are_equal = are_equal && (strcmp(W1->name[i], W2->name[i]) == 0); + are_equal = are_equal && (W1->type[i] == W2->type[i]); + are_equal = are_equal && (W1->depth_ref[i] == W2->depth_ref[i]); + are_equal = are_equal && (well_controls_equal(W1->ctrls[i], W2->ctrls[i])); + } + + + { + struct WellMgmt* mgmt1 = W1->data; + struct WellMgmt* mgmt2 = W2->data; + are_equal = are_equal && (mgmt1->perf_cpty == mgmt2->perf_cpty); + are_equal = are_equal && (mgmt1->well_cpty == mgmt2->well_cpty); + } + + are_equal = are_equal && (memcmp(W1->comp_frac, W2->comp_frac, W1->number_of_wells * W1->number_of_phases * sizeof *W1->comp_frac ) == 0); + are_equal = are_equal && (memcmp(W1->well_connpos, W2->well_connpos, (1 + W1->number_of_wells) * sizeof *W1->well_connpos ) == 0); + if (!are_equal) { + return are_equal; + } + + { + int number_of_perforations = W1->well_connpos[W1->number_of_wells]; + + are_equal = are_equal && (memcmp(W1->well_cells, W2->well_cells, number_of_perforations * sizeof *W1->well_cells ) == 0); + are_equal = are_equal && (memcmp(W1->WI, W2->WI, number_of_perforations * sizeof *W1->WI ) == 0); + } + + return are_equal; +} + +/* ---------------------------------------------------------------------- */ diff --git a/tests/test_wellcontrols.cpp b/tests/test_wellcontrols.cpp new file mode 100644 index 000000000..096b9c6a0 --- /dev/null +++ b/tests/test_wellcontrols.cpp @@ -0,0 +1,105 @@ +/* + Copyright 2014 Statoil. + + This file is part of the Open Porous Media project (OPM). + + OPM is free software: you can redistribute it and/or modify + it under the terms of the GNU General Public License as published by + the Free Software Foundation, either version 3 of the License, or + (at your option) any later version. + + OPM is distributed in the hope that it will be useful, + but WITHOUT ANY WARRANTY; without even the implied warranty of + MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + GNU General Public License for more details. + + You should have received a copy of the GNU General Public License + along with OPM. If not, see . +*/ + +#include + +#if HAVE_DYNAMIC_BOOST_TEST +#define BOOST_TEST_DYN_LINK +#endif + +#define NVERBOSE // Suppress own messages when throw()ing + +#define BOOST_TEST_MODULE WellsModuleTest +#include + +#include +#include + +#include +#include +#include + + +BOOST_AUTO_TEST_CASE(Construction) +{ + struct WellControls * ctrls = well_controls_create(); + + well_controls_set_current( ctrls , 1 ); + BOOST_CHECK_EQUAL( 1 , well_controls_get_current( ctrls )); + well_controls_set_current( ctrls , 2 ); + BOOST_CHECK_EQUAL( 2 , well_controls_get_current( ctrls )); + + well_controls_invert_current( ctrls ); + BOOST_CHECK( well_controls_get_current( ctrls ) < 0 ); + well_controls_invert_current( ctrls ); + BOOST_CHECK_EQUAL( 2 , well_controls_get_current( ctrls )); + + { + enum WellControlType type1 = BHP; + enum WellControlType type2 = SURFACE_RATE; + int num_phases = 3; + double dist1[3] = {0 , 1 , 2}; + double dist2[3] = {10, 11 , 12}; + double target = 77; + + well_controls_assert_number_of_phases( ctrls , num_phases ); + well_controls_add_new( type1 , target , dist1 , ctrls ); + well_controls_add_new( type2 , 2*target , dist2 , ctrls ); + + BOOST_CHECK_EQUAL( target , well_controls_iget_target(ctrls , 0 )); + BOOST_CHECK_EQUAL( type1 , well_controls_iget_type(ctrls , 0 )); + + BOOST_CHECK_EQUAL( 2*target , well_controls_iget_target(ctrls , 1 )); + BOOST_CHECK_EQUAL( type2 , well_controls_iget_type(ctrls , 1 )); + well_controls_set_current( ctrls , 1 ); + BOOST_CHECK_EQUAL( type2 , well_controls_get_current_type( ctrls )); + + BOOST_CHECK_EQUAL( well_controls_iget_target( ctrls , 1 ) , well_controls_get_current_target( ctrls )); + + { + const double * d1 = well_controls_iget_distr( ctrls , 0 ); + const double * d2 = well_controls_iget_distr( ctrls , 1 ); + BOOST_CHECK( memcmp(d1 , dist1 , num_phases * sizeof * d1 ) == 0); + BOOST_CHECK( memcmp(d2 , dist2 , num_phases * sizeof * d2 ) == 0); + } + } + well_controls_iset_target( ctrls , 0 , 123); + BOOST_CHECK_EQUAL( 123 , well_controls_iget_target( ctrls , 0 )); + well_controls_iset_target( ctrls , 1 , 456); + BOOST_CHECK_EQUAL( 456 , well_controls_iget_target( ctrls , 1 )); + + well_controls_iset_type( ctrls , 0 , SURFACE_RATE); + BOOST_CHECK_EQUAL( SURFACE_RATE , well_controls_iget_type( ctrls , 0 )); + well_controls_iset_type( ctrls , 1 , BHP); + BOOST_CHECK_EQUAL( BHP, well_controls_iget_type( ctrls , 1 )); + + + { + double newDist[3] = {77,78,79}; + const double * tmp; + well_controls_iset_distr( ctrls , 0 , newDist ); + tmp = well_controls_iget_distr( ctrls , 0); + BOOST_CHECK( memcmp(tmp , newDist , 3 * sizeof * tmp ) == 0); + } + + + well_controls_destroy( ctrls ); +} + + diff --git a/tests/test_wells.cpp b/tests/test_wells.cpp index b044c79fe..9e76d899d 100644 --- a/tests/test_wells.cpp +++ b/tests/test_wells.cpp @@ -29,6 +29,7 @@ #include #include +#include #include #include @@ -104,20 +105,20 @@ BOOST_AUTO_TEST_CASE(Controls) if (ok1 && ok2) { WellControls* ctrls = W->ctrls[0]; - BOOST_CHECK_EQUAL(ctrls->num , 2); - BOOST_CHECK_EQUAL(ctrls->current, -1); + BOOST_CHECK_EQUAL(well_controls_get_num(ctrls) , 2); + BOOST_CHECK_EQUAL(well_controls_get_current(ctrls), -1); set_current_control(0, 0, W.get()); - BOOST_CHECK_EQUAL(ctrls->current, 0); + BOOST_CHECK_EQUAL(well_controls_get_current(ctrls), 0); set_current_control(0, 1, W.get()); - BOOST_CHECK_EQUAL(ctrls->current, 1); + BOOST_CHECK_EQUAL(well_controls_get_current(ctrls), 1); - BOOST_CHECK_EQUAL(ctrls->type[0], BHP); - BOOST_CHECK_EQUAL(ctrls->type[1], SURFACE_RATE); + BOOST_CHECK_EQUAL(well_controls_iget_type(ctrls , 0) , BHP); + BOOST_CHECK_EQUAL(well_controls_iget_type(ctrls , 1) , SURFACE_RATE); - BOOST_CHECK_EQUAL(ctrls->target[0], 1.0); - BOOST_CHECK_EQUAL(ctrls->target[1], 1.0); + BOOST_CHECK_EQUAL(well_controls_iget_target(ctrls , 0), 1.0); + BOOST_CHECK_EQUAL(well_controls_iget_target(ctrls , 1), 1.0); } } } @@ -188,18 +189,48 @@ BOOST_AUTO_TEST_CASE(Copy) WellControls* c1 = W1->ctrls[w]; WellControls* c2 = W2->ctrls[w]; - BOOST_CHECK_EQUAL(c2->num , c1->num ); - BOOST_CHECK_EQUAL(c2->current, c1->current); + BOOST_CHECK_EQUAL(well_controls_get_num(c2) , well_controls_get_num(c1)); + BOOST_CHECK_EQUAL(well_controls_get_current(c2) , well_controls_get_current(c1)); - for (int c = 0; c < c1->num; ++c) { - BOOST_CHECK_EQUAL(c2->type [c], c1->type [c]); - BOOST_CHECK_EQUAL(c2->target[c], c1->target[c]); + for (int c = 0; c < well_controls_get_num(c1); ++c) { + BOOST_CHECK_EQUAL(well_controls_iget_type(c2, c) , well_controls_iget_type(c1 , c)); + BOOST_CHECK_EQUAL(well_controls_iget_target(c2, c) , well_controls_iget_target(c1 , c)); - for (int p = 0; p < W1->number_of_phases; ++p) { - BOOST_CHECK_EQUAL(c2->distr[c*W1->number_of_phases + p], - c1->distr[c*W1->number_of_phases + p]); + { + const double * dist1 = well_controls_iget_distr(c1 , c ); + const double * dist2 = well_controls_iget_distr(c2 , c ); + + for (int p = 0; p < W1->number_of_phases; ++p) + BOOST_CHECK_EQUAL( dist1[p] , dist2[p]); } } + BOOST_CHECK( well_controls_equal( c1 , c2 )); } } } + +BOOST_AUTO_TEST_CASE(Equals_WellsEqual_ReturnsTrue) { + const int nphases = 2; + const int nwells = 2; + const int nperfs = 2; + + std::shared_ptr W1(create_wells(nphases, nwells, nperfs), + destroy_wells); + std::shared_ptr W2(create_wells(nphases, nwells, nperfs), + destroy_wells); + + BOOST_CHECK(wells_equal(W1.get(), W2.get())); +} + + +BOOST_AUTO_TEST_CASE(Equals_WellsDiffer_ReturnsFalse) { + const int nphases = 2; + const int nperfs = 2; + + std::shared_ptr W1(create_wells(nphases, 2, nperfs), + destroy_wells); + std::shared_ptr W2(create_wells(nphases, 3, nperfs), + destroy_wells); + + BOOST_CHECK(!wells_equal(W1.get(), W2.get())); +} diff --git a/tests/test_wellsmanager.cpp b/tests/test_wellsmanager.cpp new file mode 100644 index 000000000..b3465fac6 --- /dev/null +++ b/tests/test_wellsmanager.cpp @@ -0,0 +1,236 @@ +/* + Copyright 2013 Statoil ASA. + + This file is part of the Open Porous Media project (OPM). + + OPM is free software: you can redistribute it and/or modify + it under the terms of the GNU General Public License as published by + the Free Software Foundation, either version 3 of the License, or + (at your option) any later version. + + OPM is distributed in the hope that it will be useful, + but WITHOUT ANY WARRANTY; without even the implied warranty of + MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + GNU General Public License for more details. + + You should have received a copy of the GNU General Public License + along with OPM. If not, see . +*/ + +#include + +#if HAVE_DYNAMIC_BOOST_TEST +#define BOOST_TEST_DYN_LINK +#endif + +#define NVERBOSE // Suppress own messages when throw()ing + +#define BOOST_TEST_MODULE WellsManagerTests +#include +#include +#include +#include + +#include +#include + + +void wells_static_check(const Wells* wells) { + BOOST_CHECK_EQUAL(2 , wells->number_of_wells); + BOOST_CHECK_EQUAL(3 , wells->number_of_phases); + + BOOST_CHECK_EQUAL("INJ1" , wells->name[0]); + BOOST_CHECK_EQUAL("PROD1" , wells->name[1]); + + /* The mapping from well number into the wells->WI and wells->well_cells arrays. */ + BOOST_CHECK_EQUAL(0 , wells->well_connpos[0]); + BOOST_CHECK_EQUAL(1 , wells->well_connpos[1]); + BOOST_CHECK_EQUAL(2 , wells->well_connpos[2]); + + /* Connection factor */ + BOOST_CHECK_CLOSE(1.2279166666666664e-12 , wells->WI[0] , 0.001); + BOOST_CHECK_CLOSE(1.2279166666666664e-12 , wells->WI[1] , 0.001); + + /* Completed cells */ + BOOST_CHECK_EQUAL(0 , wells->well_cells[0]); + BOOST_CHECK_EQUAL(9 + 2*10 + 2*10*10 , wells->well_cells[1]); +} + + +/* + The number of controls is determined by looking at which elements + have been given explicit - non-default - values in the WCONxxxx + keyword. Is that at all interesting? +*/ + + +void check_controls_epoch0( struct WellControls ** ctrls) { + // The injector + { + const struct WellControls * ctrls0 = ctrls[0]; + BOOST_CHECK_EQUAL( 3 , well_controls_get_num(ctrls0)); // The number of controls for the injector == 3?? + + BOOST_CHECK_EQUAL( SURFACE_RATE , well_controls_iget_type(ctrls0 , 0) ); + BOOST_CHECK_EQUAL( RESERVOIR_RATE , well_controls_iget_type(ctrls0 , 1) ); + BOOST_CHECK_EQUAL( BHP , well_controls_iget_type(ctrls0 , 2) ); + + // The different targets + BOOST_CHECK_EQUAL( 100.0 / 86400 , well_controls_iget_target(ctrls0,0)); + BOOST_CHECK_EQUAL( 200.0 / 86400 , well_controls_iget_target(ctrls0,1)); + BOOST_CHECK_EQUAL( 400 * 100000 , well_controls_iget_target(ctrls0,2)); + + // Which control is active + BOOST_CHECK_EQUAL( 0 , well_controls_get_current(ctrls0) ); + + // The phase distribution in the active target + { + const double * distr = well_controls_iget_distr( ctrls0 , 0 ); + BOOST_CHECK_EQUAL( 0 , distr[0] ); // Water + BOOST_CHECK_EQUAL( 0 , distr[1] ); // Oil + BOOST_CHECK_EQUAL( 1 , distr[2] ); // Gas + } + } + + // The producer + { + const struct WellControls * ctrls1 = ctrls[1]; + BOOST_CHECK_EQUAL( 2 , well_controls_get_num( ctrls1 )); // The number of controls for the producer == 2?? + BOOST_CHECK_EQUAL( SURFACE_RATE , well_controls_iget_type(ctrls1 , 0) ); + BOOST_CHECK_EQUAL( BHP , well_controls_iget_type(ctrls1 , 1) ); + + // The different targets + BOOST_CHECK_EQUAL( -20000.0 / 86400 , well_controls_iget_target(ctrls1,0)); + BOOST_CHECK_EQUAL( 1000 * 100000 , well_controls_iget_target(ctrls1,1)); + + // Which control is active + BOOST_CHECK_EQUAL( 0 , well_controls_get_current(ctrls1)); + + // The phase distribution in the active target + { + const double * distr = well_controls_iget_distr( ctrls1 , 0 ); + BOOST_CHECK_EQUAL( 0 , distr[0] ); // Water + BOOST_CHECK_EQUAL( 1 , distr[1] ); // Oil + BOOST_CHECK_EQUAL( 0 , distr[2] ); // Gas + } + } +} + + + + +void check_controls_epoch1( struct WellControls ** ctrls) { + // The injector + { + const struct WellControls * ctrls0 = ctrls[0]; + BOOST_CHECK_EQUAL( 3 , well_controls_get_num(ctrls0)); // The number of controls for the injector == 3?? + + BOOST_CHECK_EQUAL( SURFACE_RATE , well_controls_iget_type(ctrls0 , 0 )); + BOOST_CHECK_EQUAL( RESERVOIR_RATE , well_controls_iget_type(ctrls0 , 1 )); + BOOST_CHECK_EQUAL( BHP , well_controls_iget_type(ctrls0 , 2 )); + + // The different targets + BOOST_CHECK_CLOSE( 10.0 / 86400 , well_controls_iget_target(ctrls0 , 0) , 0.001); + BOOST_CHECK_CLOSE( 20.0 / 86400 , well_controls_iget_target(ctrls0 , 1) , 0.001); + BOOST_CHECK_CLOSE( 40 * 100000 , well_controls_iget_target(ctrls0 , 2) , 0.001); + + // Which control is active + BOOST_CHECK_EQUAL( 1 , well_controls_get_current(ctrls0)); + + { + const double * distr = well_controls_iget_distr( ctrls0 , 1 ); + BOOST_CHECK_EQUAL( 1 , distr[0] ); // Water + BOOST_CHECK_EQUAL( 0 , distr[1] ); // Oil + BOOST_CHECK_EQUAL( 0 , distr[2] ); // Gas + } + } + + // The producer + { + const struct WellControls * ctrls1 = ctrls[1]; + BOOST_CHECK_EQUAL( 3 , well_controls_get_num(ctrls1)); // The number of controls for the producer - now 3. + BOOST_CHECK_EQUAL( SURFACE_RATE , well_controls_iget_type(ctrls1 , 0) ); + BOOST_CHECK_EQUAL( RESERVOIR_RATE , well_controls_iget_type(ctrls1 , 1) ); + BOOST_CHECK_EQUAL( BHP , well_controls_iget_type(ctrls1 , 2) ); + + // The different targets + BOOST_CHECK_CLOSE( -999.0 / 86400 , well_controls_iget_target(ctrls1 , 0), 0.001); + BOOST_CHECK_CLOSE( -123.0 / 86400 , well_controls_iget_target(ctrls1 , 1), 0.001); + BOOST_CHECK_CLOSE( 100 * 100000 , well_controls_iget_target(ctrls1 , 2), 0.001); + + // Which control is active + BOOST_CHECK_EQUAL( 1 , well_controls_get_current(ctrls1) ); + + { + const double * distr = well_controls_iget_distr( ctrls1 , 1 ); + BOOST_CHECK_EQUAL( 1 , distr[0] ); // Water + BOOST_CHECK_EQUAL( 1 , distr[1] ); // Oil + BOOST_CHECK_EQUAL( 1 , distr[2] ); // Gas + } + } +} + + + + +BOOST_AUTO_TEST_CASE(Constructor_Works) { + Opm::EclipseGridParser Deck("wells_manager_data.data"); + Opm::GridManager gridManager(Deck); + + Deck.setCurrentEpoch(0); + { + Opm::WellsManager wellsManager(Deck, *gridManager.c_grid(), NULL); + const Wells* wells = wellsManager.c_wells(); + wells_static_check( wells ); + check_controls_epoch0( wells->ctrls ); + } + + + Deck.setCurrentEpoch(1); + { + Opm::WellsManager wellsManager(Deck, *gridManager.c_grid(), NULL); + const Wells* wells = wellsManager.c_wells(); + + wells_static_check( wells ); + check_controls_epoch1( wells->ctrls ); + } +} + + + +BOOST_AUTO_TEST_CASE(WellsEqual) { + Opm::EclipseGridParser Deck("wells_manager_data.data"); + Opm::GridManager gridManager(Deck); + + Deck.setCurrentEpoch(0); + Opm::WellsManager wellsManager0(Deck, *gridManager.c_grid(), NULL); + + Deck.setCurrentEpoch(1); + Opm::WellsManager wellsManager1(Deck, *gridManager.c_grid(), NULL); + + BOOST_CHECK( wells_equal( wellsManager0.c_wells() , wellsManager0.c_wells()) ); + BOOST_CHECK( !wells_equal( wellsManager0.c_wells() , wellsManager1.c_wells()) ); +} + + +BOOST_AUTO_TEST_CASE(ControlsEqual) { + Opm::EclipseGridParser Deck("wells_manager_data.data"); + Opm::GridManager gridManager(Deck); + + Deck.setCurrentEpoch(0); + Opm::WellsManager wellsManager0(Deck, *gridManager.c_grid(), NULL); + + Deck.setCurrentEpoch(1); + Opm::WellsManager wellsManager1(Deck, *gridManager.c_grid(), NULL); + + BOOST_CHECK( well_controls_equal( wellsManager0.c_wells()->ctrls[0] , wellsManager0.c_wells()->ctrls[0])); + BOOST_CHECK( well_controls_equal( wellsManager0.c_wells()->ctrls[1] , wellsManager0.c_wells()->ctrls[1])); + BOOST_CHECK( well_controls_equal( wellsManager1.c_wells()->ctrls[0] , wellsManager1.c_wells()->ctrls[0])); + BOOST_CHECK( well_controls_equal( wellsManager1.c_wells()->ctrls[1] , wellsManager1.c_wells()->ctrls[1])); + + BOOST_CHECK( !well_controls_equal( wellsManager0.c_wells()->ctrls[0] , wellsManager0.c_wells()->ctrls[1])); + BOOST_CHECK( !well_controls_equal( wellsManager0.c_wells()->ctrls[1] , wellsManager0.c_wells()->ctrls[0])); + BOOST_CHECK( !well_controls_equal( wellsManager1.c_wells()->ctrls[0] , wellsManager0.c_wells()->ctrls[0])); + BOOST_CHECK( !well_controls_equal( wellsManager1.c_wells()->ctrls[1] , wellsManager0.c_wells()->ctrls[1])); +} + + diff --git a/tests/wells_manager_data.data b/tests/wells_manager_data.data new file mode 100755 index 000000000..4c810afab --- /dev/null +++ b/tests/wells_manager_data.data @@ -0,0 +1,58 @@ +OIL +GAS +WATER + +DIMENS + 10 10 5 / + +GRID + +DXV +10*1000.0 / + +DYV +10*1000.0 / + +DZV +10.0 20.0 30.0 10.0 5.0 / + +SCHEDULE + +WELSPECS + 'INJ1' 'G' 1 1 8335 'GAS' / + 'PROD1' 'G' 10 10 8400 'OIL' / +/ + +COMPDAT + 'INJ1' 1 1 1 1 'OPEN' 1 10.6092 0.5 / + 'PROD1' 10 3 3 3 'OPEN' 0 10.6092 0.5 / +/ + +WCONPROD + 'PROD1' 'OPEN' 'ORAT' 20000 4* 1000 / + / + +WCONINJE + 'INJ1' 'GAS' 'OPEN' 'RATE' 100 200 400 / + / + + +DATES + 1 'FEB' 2000 / +/ + +WCONPROD + 'PROD1' 'OPEN' 'RESV' 999 3* 123 100 / +/ + +WCONINJE + 'INJ1' 'WATER' 'OPEN' 'RESV' 10 20 40 / +/ + + + +TSTEP +14.0 + / + +END