00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022 #ifndef OPM_WELLINTERFACE_HEADER_INCLUDED
00023 #define OPM_WELLINTERFACE_HEADER_INCLUDED
00024
00025 #include <opm/common/OpmLog/OpmLog.hpp>
00026
00027
00028 #include <opm/parser/eclipse/EclipseState/Schedule/Well.hpp>
00029 #include <opm/core/wells.h>
00030 #include <opm/core/well_controls.h>
00031 #include <opm/core/props/BlackoilPhases.hpp>
00032 #include <opm/core/wells/WellsManager.hpp>
00033
00034 #include <opm/autodiff/VFPProperties.hpp>
00035 #include <opm/autodiff/VFPInjProperties.hpp>
00036 #include <opm/autodiff/VFPProdProperties.hpp>
00037 #include <opm/autodiff/WellHelpers.hpp>
00038 #include <opm/autodiff/WellStateFullyImplicitBlackoil.hpp>
00039 #include <opm/autodiff/BlackoilModelParameters.hpp>
00040
00041 #include <opm/simulators/WellSwitchingLogger.hpp>
00042
00043 #include<dune/common/fmatrix.hh>
00044 #include<dune/istl/bcrsmatrix.hh>
00045 #include<dune/istl/matrixmatrix.hh>
00046
00047 #include <opm/material/densead/Math.hpp>
00048 #include <opm/material/densead/Evaluation.hpp>
00049
00050 #include <string>
00051 #include <memory>
00052 #include <vector>
00053 #include <cassert>
00054
00055 namespace Opm
00056 {
00057
00058
00059 template<typename TypeTag>
00060 class WellInterface
00061 {
00062 public:
00063
00064 using WellState = WellStateFullyImplicitBlackoil;
00065
00066 typedef BlackoilModelParameters ModelParameters;
00067 typedef typename GET_PROP_TYPE(TypeTag, Grid) Grid;
00068 typedef typename GET_PROP_TYPE(TypeTag, Simulator) Simulator;
00069 typedef typename GET_PROP_TYPE(TypeTag, FluidSystem) FluidSystem;
00070 typedef typename GET_PROP_TYPE(TypeTag, Indices) BlackoilIndices;
00071 typedef typename GET_PROP_TYPE(TypeTag, IntensiveQuantities) IntensiveQuantities;
00072 typedef typename GET_PROP_TYPE(TypeTag, MaterialLaw) MaterialLaw;
00073
00074 static const int numEq = BlackoilIndices::numEq;
00075 typedef double Scalar;
00076
00077 typedef Dune::FieldVector<Scalar, numEq > VectorBlockType;
00078 typedef Dune::FieldMatrix<Scalar, numEq, numEq > MatrixBlockType;
00079 typedef Dune::BCRSMatrix <MatrixBlockType> Mat;
00080 typedef Dune::BlockVector<VectorBlockType> BVector;
00081 typedef DenseAd::Evaluation<double, numEq> Eval;
00082
00083 typedef Ewoms::BlackOilPolymerModule<TypeTag> PolymerModule;
00084
00085 static const bool has_solvent = GET_PROP_VALUE(TypeTag, EnableSolvent);
00086 static const bool has_polymer = GET_PROP_VALUE(TypeTag, EnablePolymer);
00087
00089 WellInterface(const Well* well, const int time_step, const Wells* wells, const ModelParameters& param);
00090
00092 virtual ~WellInterface() {}
00093
00095 const std::string& name() const;
00096
00098 const std::vector<int>& cells() {return well_cells_; }
00099
00101 WellType wellType() const;
00102
00104 WellControls* wellControls() const;
00105
00106 void setVFPProperties(const VFPProperties* vfp_properties_arg);
00107
00108 virtual void init(const PhaseUsage* phase_usage_arg,
00109 const std::vector<bool>* active_arg,
00110 const std::vector<double>& depth_arg,
00111 const double gravity_arg,
00112 const int num_cells);
00113
00114 virtual void initPrimaryVariablesEvaluation() const = 0;
00115
00117 struct ConvergenceReport {
00118 struct ProblemWell {
00119 std::string well_name;
00120 std::string phase_name;
00121 };
00122 bool converged = true;
00123 bool nan_residual_found = false;
00124 std::vector<ProblemWell> nan_residual_wells;
00125
00126 bool too_large_residual_found = false;
00127 std::vector<ProblemWell> too_large_residual_wells;
00128
00129 ConvergenceReport& operator+=(const ConvergenceReport& rhs) {
00130 converged = converged && rhs.converged;
00131 nan_residual_found = nan_residual_found || rhs.nan_residual_found;
00132 if (rhs.nan_residual_found) {
00133 for (const ProblemWell& well : rhs.nan_residual_wells) {
00134 nan_residual_wells.push_back(well);
00135 }
00136 }
00137 too_large_residual_found = too_large_residual_found || rhs.too_large_residual_found;
00138 if (rhs.too_large_residual_found) {
00139 for (const ProblemWell& well : rhs.too_large_residual_wells) {
00140 too_large_residual_wells.push_back(well);
00141 }
00142 }
00143 return *this;
00144 }
00145 };
00146
00147 virtual ConvergenceReport getWellConvergence(const std::vector<double>& B_avg) const = 0;
00148
00149 virtual void solveEqAndUpdateWellState(WellState& well_state) = 0;
00150
00151 virtual void assembleWellEq(Simulator& ebosSimulator,
00152 const double dt,
00153 WellState& well_state,
00154 bool only_wells) = 0;
00155
00156 void updateListEconLimited(const WellState& well_state,
00157 DynamicListEconLimited& list_econ_limited) const;
00158
00159 void setWellEfficiencyFactor(const double efficiency_factor);
00160
00161 void computeRepRadiusPerfLength(const Grid& grid, const std::map<int, int>& cartesian_to_compressed);
00162
00165 virtual void recoverWellSolutionAndUpdateWellState(const BVector& x,
00166 WellState& well_state) const = 0;
00167
00169 virtual void apply(const BVector& x, BVector& Ax) const = 0;
00170
00172 virtual void apply(BVector& r) const = 0;
00173
00174
00175 virtual void computeWellPotentials(const Simulator& ebosSimulator,
00176 const WellState& well_state,
00177 std::vector<double>& well_potentials) = 0;
00178
00179 virtual void updateWellStateWithTarget(const int current,
00180 WellState& xw) const = 0;
00181
00182 void updateWellControl(WellState& xw,
00183 wellhelpers::WellSwitchingLogger& logger) const;
00184
00185 virtual void updatePrimaryVariables(const WellState& well_state) const = 0;
00186
00187 virtual void calculateExplicitQuantities(const Simulator& ebosSimulator,
00188 const WellState& xw) = 0;
00189
00190 protected:
00191
00192
00193 static const int INVALIDCONNECTION = -100000;
00194
00195 const Well* well_ecl_;
00196
00197 const int current_step_;
00198
00199
00200 int index_of_well_;
00201
00202
00203 const ModelParameters& param_;
00204
00205
00206
00207 enum WellType well_type_;
00208
00209
00210 int number_of_phases_;
00211
00212
00213
00214 std::vector<double> comp_frac_;
00215
00216
00217 struct WellControls* well_controls_;
00218
00219
00220 int number_of_perforations_;
00221
00222
00223
00224 int first_perf_;
00225
00226
00227 std::vector<double> well_index_;
00228
00229
00230 std::vector<double> perf_depth_;
00231
00232
00233 double ref_depth_;
00234
00235 double well_efficiency_factor_;
00236
00237
00238 std::vector<int> well_cells_;
00239
00240
00241 std::vector<int> saturation_table_number_;
00242
00243
00244 std::vector<double> perf_rep_radius_;
00245
00246
00247 std::vector<double> perf_length_;
00248
00249
00250 std::vector<double> bore_diameters_;
00251
00252 const PhaseUsage* phase_usage_;
00253
00254 bool getAllowCrossFlow() const;
00255
00256 const std::vector<bool>* active_;
00257
00258 const VFPProperties* vfp_properties_;
00259
00260 double gravity_;
00261
00262 const std::vector<bool>& active() const;
00263
00264 const PhaseUsage& phaseUsage() const;
00265
00266 int flowPhaseToEbosCompIdx( const int phaseIdx ) const;
00267
00268 int flowPhaseToEbosPhaseIdx( const int phaseIdx ) const;
00269
00270
00271 int numComponents() const;
00272
00273 double wsolvent() const;
00274
00275 double wpolymer() const;
00276
00277 bool checkRateEconLimits(const WellEconProductionLimits& econ_production_limits,
00278 const WellState& well_state) const;
00279
00280 bool wellHasTHPConstraints() const;
00281
00282
00283 const std::vector<double>& compFrac() const;
00284
00285 double mostStrictBhpFromBhpLimits() const;
00286
00287
00288
00289
00290
00291
00292
00293
00294 using RatioCheckTuple = std::tuple<bool, bool, int, double>;
00295
00296 RatioCheckTuple checkMaxWaterCutLimit(const WellEconProductionLimits& econ_production_limits,
00297 const WellState& well_state) const;
00298
00299 RatioCheckTuple checkRatioEconLimits(const WellEconProductionLimits& econ_production_limits,
00300 const WellState& well_state) const;
00301
00302 };
00303
00304 }
00305
00306 #include "WellInterface_impl.hpp"
00307
00308 #endif // OPM_WELLINTERFACE_HEADER_INCLUDED