00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022 #ifndef OPM_MULTISEGMENTWELLS_HEADER_INCLUDED
00023 #define OPM_MULTISEGMENTWELLS_HEADER_INCLUDED
00024
00025 #include <dune/common/parallel/mpihelper.hh>
00026
00027 #include <opm/common/utility/platform_dependent/disable_warnings.h>
00028 #include <Eigen/Eigen>
00029 #include <Eigen/Sparse>
00030 #include <opm/common/utility/platform_dependent/reenable_warnings.h>
00031
00032 #include <cassert>
00033
00034 #include <opm/core/props/BlackoilPhases.hpp>
00035 #include <opm/core/wells/WellCollection.hpp>
00036
00037 #include <opm/autodiff/AutoDiffBlock.hpp>
00038 #include <opm/autodiff/AutoDiffHelpers.hpp>
00039 #include <opm/autodiff/BlackoilModelEnums.hpp>
00040 #include <opm/autodiff/BlackoilPropsAdFromDeck.hpp>
00041 #include <opm/autodiff/LinearisedBlackoilResidual.hpp>
00042 #include <opm/autodiff/WellHelpers.hpp>
00043 #include <opm/autodiff/VFPProperties.hpp>
00044
00045 #include <opm/autodiff/WellMultiSegment.hpp>
00046 #include <opm/autodiff/WellDensitySegmented.hpp>
00047 #include <opm/simulators/WellSwitchingLogger.hpp>
00048
00049
00050
00051 namespace Opm {
00052
00053
00055 class MultisegmentWells {
00056 public:
00057
00058 using ADB = AutoDiffBlock<double>;
00059 using Vector = ADB::V;
00060
00061
00062 struct MultisegmentWellOps {
00063 explicit MultisegmentWellOps(const std::vector<WellMultiSegmentConstPtr>& wells_ms);
00064 Eigen::SparseMatrix<double> w2p;
00065 Eigen::SparseMatrix<double> p2w;
00066 Eigen::SparseMatrix<double> w2s;
00067 Eigen::SparseMatrix<double> s2w;
00068 Eigen::SparseMatrix<double> s2p;
00069 Eigen::SparseMatrix<double> p2s;
00070 Eigen::SparseMatrix<double> s2s_inlets;
00071 Eigen::SparseMatrix<double> s2s_outlet;
00072 Eigen::SparseMatrix<double> topseg2w;
00073 AutoDiffMatrix eliminate_topseg;
00074 std::vector<int> well_cells;
00075 Vector conn_trans_factors;
00076 bool has_multisegment_wells;
00077 };
00078
00079
00080
00081 using DataBlock = Eigen::Array<double,
00082 Eigen::Dynamic,
00083 Eigen::Dynamic,
00084 Eigen::RowMajor>;
00085 using Communication =
00086 Dune::CollectiveCommunication<typename Dune::MPIHelper
00087 ::MPICommunicator>;
00088
00089
00090
00091
00092 MultisegmentWells(const Wells* wells_arg,
00093 WellCollection* well_collection,
00094 const std::vector< const Well* >& wells_ecl,
00095 const int time_step);
00096
00097 std::vector<WellMultiSegmentConstPtr> createMSWellVector(const Wells* wells_arg,
00098 const std::vector< const Well* >& wells_ecl,
00099 const int time_step);
00100
00101 void init(const BlackoilPropsAdFromDeck* fluid_arg,
00102 const std::vector<bool>* active_arg,
00103 const std::vector<PhasePresence>* pc_arg,
00104 const VFPProperties* vfp_properties_arg,
00105 const double gravity_arg,
00106 const Vector& depth_arg);
00107
00108 const std::vector<WellMultiSegmentConstPtr>& msWells() const;
00109 const MultisegmentWellOps& wellOps() const;
00110
00111 const Wells& wells() const;
00112
00113 const Wells* wellsPointer() const;
00114
00115 int numPhases() const { return num_phases_; };
00116
00117 int numWells() const { return wells_multisegment_.size(); }
00118
00119 int numSegment() const { return nseg_total_; };
00120
00121 int numPerf() const { return nperf_total_; };
00122
00123 bool wellsActive() const { return wells_active_; };
00124
00125 void setWellsActive(const bool wells_active) { wells_active_ = wells_active; };
00126
00127 bool localWellsActive() const { return ! wells_multisegment_.empty(); };
00128
00129 int numWellVars() const { return (num_phases_ + 1) * nseg_total_; };
00130
00131 template <class ReservoirResidualQuant, class SolutionState>
00132 void extractWellPerfProperties(const SolutionState& state,
00133 const std::vector<ReservoirResidualQuant>& rq,
00134 std::vector<ADB>& mob_perfcells,
00135 std::vector<ADB>& b_perfcells) const;
00136
00137 Vector& wellPerforationCellPressureDiffs() { return well_perforation_cell_pressure_diffs_; };
00138
00139 Vector& wellSegmentPerforationDepthDiffs() { return well_segment_perforation_depth_diffs_; };
00140
00141 const Vector& wellPerforationCellDensities() const { return well_perforation_cell_densities_; };
00142 Vector& wellPerforationCellDensities() { return well_perforation_cell_densities_; };
00143
00144 const std::vector<Vector>& segmentCompSurfVolumeInitial() const { return segment_comp_surf_volume_initial_; };
00145 std::vector<Vector>& segmentCompSurfVolumeInitial() { return segment_comp_surf_volume_initial_; };
00146
00147 const std::vector<ADB>& segmentCompSurfVolumeCurrent() const { return segment_comp_surf_volume_current_; };
00148
00149 const std::vector<int>& topWellSegments() const { return top_well_segments_; };
00150 std::vector<int>& topWellSegments() { return top_well_segments_; };
00151
00152 Vector& segVDt() { return segvdt_; };
00153
00154 const Vector& wellPerforationDensities() const { return well_perforation_densities_; };
00155 Vector& wellPerforationDensities() { return well_perforation_densities_; };
00156
00157 const Vector& wellPerforationPressureDiffs() const { return well_perforation_pressure_diffs_; };
00158 Vector& wellPerforationPressureDiffs() { return well_perforation_pressure_diffs_; };
00159
00160 template <class WellState>
00161 void
00162 updateWellState(const Vector& dwells,
00163 const double dpmaxrel,
00164 WellState& well_state) const;
00165
00166
00167
00168 template <class SolutionState>
00169 void
00170 computeWellFlux(const SolutionState& state,
00171 const std::vector<ADB>& mob_perfcells,
00172 const std::vector<ADB>& b_perfcells,
00173 Vector& aliveWells,
00174 std::vector<ADB>& cq_s) const;
00175
00176 template <class SolutionState, class WellState>
00177 void updatePerfPhaseRatesAndPressures(const std::vector<ADB>& cq_s,
00178 const SolutionState& state,
00179 WellState& xw) const;
00180
00181
00182
00183
00184 template <class SolutionState>
00185 void
00186 computeSegmentFluidProperties(const SolutionState& state);
00187
00188 void
00189 computeSegmentPressuresDelta(const double grav);
00190
00191 template <class SolutionState>
00192 void
00193 addWellFluxEq(const std::vector<ADB>& cq_s,
00194 const SolutionState& state,
00195 LinearisedBlackoilResidual& residual);
00196
00197 template <class SolutionState, class WellState>
00198 void
00199 addWellControlEq(const SolutionState& state,
00200 const WellState& xw,
00201 const Vector& aliveWells,
00202 LinearisedBlackoilResidual& residual);
00203
00204 template <class WellState>
00205 void
00206 updateWellControls(WellState& xw) const;
00207
00208
00209
00210 void
00211 variableStateWellIndices(std::vector<int>& indices,
00212 int& next) const;
00213
00214 template <class SolutionState>
00215 void
00216 variableStateExtractWellsVars(const std::vector<int>& indices,
00217 std::vector<ADB>& vars,
00218 SolutionState& state) const;
00219
00220 std::vector<int>
00221 variableWellStateIndices() const;
00222
00223 template <class WellState>
00224 void
00225 variableWellStateInitials(const WellState& xw,
00226 std::vector<Vector>& vars0) const;
00227
00228 template <class SolutionState, class WellState>
00229 void computeWellConnectionPressures(const SolutionState& state,
00230 const WellState& xw,
00231 const std::vector<ADB>& kr_adb,
00232 const std::vector<ADB>& fluid_density);
00233
00234 WellCollection* wellCollection() const;
00235
00236 void calculateEfficiencyFactors();
00237
00238 const Vector& wellPerfEfficiencyFactors() const;
00239
00240
00241 protected:
00242
00243 bool wells_active_;
00244 std::vector<WellMultiSegmentConstPtr> wells_multisegment_;
00245 MultisegmentWellOps wops_ms_;
00246
00247 WellCollection* well_collection_;
00248
00249
00250
00251
00252 Vector well_perforation_efficiency_factors_;
00253
00254 const int num_phases_;
00255 int nseg_total_;
00256 int nperf_total_;
00257
00258
00259
00260
00261
00262 const Wells* wells_;
00263
00264 const BlackoilPropsAdFromDeck* fluid_;
00265 const std::vector<bool>* active_;
00266 const std::vector<PhasePresence>* phase_condition_;
00267 const VFPProperties* vfp_properties_;
00268 double gravity_;
00269
00270
00271 Vector perf_cell_depth_;
00272
00273
00274
00275
00276
00277
00278 Vector well_perforation_cell_pressure_diffs_;
00279
00280
00281 ADB well_segment_perforation_pressure_diffs_;
00282
00283
00284 Vector well_segment_perforation_depth_diffs_;
00285
00286
00287 Vector perf_cell_depth_diffs_;
00288
00289
00290
00291
00292 Vector well_perforation_cell_densities_;
00293
00294
00295
00296 ADB well_segment_densities_;
00297
00298
00299
00300
00301 ADB well_segment_pressures_delta_;
00302
00303
00304
00305 std::vector<Vector> segment_comp_surf_volume_initial_;
00306
00307
00308 std::vector<ADB> segment_comp_surf_volume_current_;
00309
00310
00311 ADB segment_mass_flow_rates_;
00312
00313
00314
00315
00316 ADB segment_viscosities_;
00317
00318 std::vector<int> top_well_segments_;
00319
00320
00321
00322 Vector segvdt_;
00323
00324
00325
00326
00327
00328 Vector well_perforation_densities_;
00329 Vector well_perforation_pressure_diffs_;
00330
00331 };
00332
00333 }
00334
00335 #include "MultisegmentWells_impl.hpp"
00336
00337 #endif // OPM_MULTISEGMENTWELLS_HEADER_INCLUDED