22 #ifndef OPM_WELLHELPERS_HEADER_INCLUDED
23 #define OPM_WELLHELPERS_HEADER_INCLUDED
25 #include <opm/core/wells.h>
40 double rateToCompare(
const std::vector<double>& well_phase_flow_rate,
46 for (
int phase = 0; phase < num_phases; ++phase) {
49 rate += well_phase_flow_rate[well*num_phases + phase] * distr[phase];
55 bool constraintBroken(
const std::vector<double>& bhp,
56 const std::vector<double>& thp,
57 const std::vector<double>& well_phase_flow_rate,
60 const WellType& well_type,
61 const WellControls* wc,
64 const WellControlType ctrl_type = well_controls_iget_type(wc, ctrl_index);
65 const double target = well_controls_iget_target(wc, ctrl_index);
66 const double* distr = well_controls_iget_distr(wc, ctrl_index);
75 broken = bhp[well] > target;
79 broken = thp[well] > target;
84 broken = rateToCompare(well_phase_flow_rate,
85 well, num_phases, distr) > target;
95 broken = bhp[well] < target;
99 broken = thp[well] < target;
107 broken = rateToCompare(well_phase_flow_rate,
108 well, num_phases, distr) < target;
115 OPM_THROW(std::logic_error,
"Can only handle INJECTOR and PRODUCER wells.");
133 double computeHydrostaticCorrection(
const Wells& wells,
const int w,
double vfp_ref_depth,
134 const double& rho,
const double gravity) {
135 if ( wells.well_connpos[w] == wells.well_connpos[w+1] )
143 const double well_ref_depth = wells.depth_ref[w];
144 const double dh = vfp_ref_depth - well_ref_depth;
145 const double dp = rho*gravity*dh;
151 double computeHydrostaticCorrection(
const double well_ref_depth,
const double vfp_ref_depth,
152 const double rho,
const double gravity) {
153 const double dh = vfp_ref_depth - well_ref_depth;
154 const double dp = rho * gravity * dh;
159 template <
class Vector>
161 Vector computeHydrostaticCorrection(
const Wells& wells,
const Vector vfp_ref_depth,
162 const Vector& well_perforation_densities,
const double gravity) {
163 const int nw = wells.number_of_wells;
164 Vector retval = Vector::Zero(nw);
166 #pragma omp parallel for schedule(static)
167 for (
int i=0; i<nw; ++i) {
168 const int perf = wells.well_connpos[i];
169 retval[i] = computeHydrostaticCorrection(wells, i, vfp_ref_depth[i], well_perforation_densities[perf], gravity);
176 std::vector<double> computeHydrostaticCorrection(
const Wells& wells,
const std::vector<double>& vfp_ref_depth,
177 const std::vector<double>& well_perforation_densities,
const double gravity) {
178 const int nw = wells.number_of_wells;
179 std::vector<double> retval(nw,0.0);
181 #pragma omp parallel for schedule(static)
182 for (
int i=0; i<nw; ++i) {
183 const int perf = wells.well_connpos[i];
184 retval[i] = computeHydrostaticCorrection(wells, i, vfp_ref_depth[i], well_perforation_densities[perf], gravity);