22 #ifndef OPM_MULTISEGMENTWELLS_HEADER_INCLUDED
23 #define OPM_MULTISEGMENTWELLS_HEADER_INCLUDED
25 #include <dune/common/parallel/mpihelper.hh>
27 #include <opm/common/utility/platform_dependent/disable_warnings.h>
28 #include <Eigen/Eigen>
29 #include <Eigen/Sparse>
30 #include <opm/common/utility/platform_dependent/reenable_warnings.h>
34 #include <opm/core/props/BlackoilPhases.hpp>
35 #include <opm/core/wells/WellCollection.hpp>
37 #include <opm/autodiff/AutoDiffBlock.hpp>
38 #include <opm/autodiff/AutoDiffHelpers.hpp>
39 #include <opm/autodiff/BlackoilModelEnums.hpp>
40 #include <opm/autodiff/BlackoilPropsAdFromDeck.hpp>
41 #include <opm/autodiff/LinearisedBlackoilResidual.hpp>
42 #include <opm/autodiff/WellHelpers.hpp>
43 #include <opm/autodiff/VFPProperties.hpp>
45 #include <opm/autodiff/WellMultiSegment.hpp>
46 #include <opm/autodiff/WellDensitySegmented.hpp>
47 #include <opm/simulators/WellSwitchingLogger.hpp>
64 Eigen::SparseMatrix<double> w2p;
65 Eigen::SparseMatrix<double> p2w;
66 Eigen::SparseMatrix<double> w2s;
67 Eigen::SparseMatrix<double> s2w;
68 Eigen::SparseMatrix<double> s2p;
69 Eigen::SparseMatrix<double> p2s;
70 Eigen::SparseMatrix<double> s2s_inlets;
71 Eigen::SparseMatrix<double> s2s_outlet;
72 Eigen::SparseMatrix<double> topseg2w;
74 std::vector<int> well_cells;
75 Vector conn_trans_factors;
76 bool has_multisegment_wells;
81 using DataBlock = Eigen::Array<double,
86 Dune::CollectiveCommunication<
typename Dune::MPIHelper
93 WellCollection* well_collection,
94 const std::vector< const Well* >& wells_ecl,
97 std::vector<WellMultiSegmentConstPtr> createMSWellVector(
const Wells* wells_arg,
98 const std::vector< const Well* >& wells_ecl,
102 const std::vector<bool>* active_arg,
103 const std::vector<PhasePresence>* pc_arg,
105 const double gravity_arg,
106 const Vector& depth_arg);
108 const std::vector<WellMultiSegmentConstPtr>& msWells()
const;
111 const Wells& wells()
const;
113 const Wells* wellsPointer()
const;
115 int numPhases()
const {
return num_phases_; };
117 int numWells()
const {
return wells_multisegment_.size(); }
119 int numSegment()
const {
return nseg_total_; };
121 int numPerf()
const {
return nperf_total_; };
123 bool wellsActive()
const {
return wells_active_; };
125 void setWellsActive(
const bool wells_active) { wells_active_ = wells_active; };
127 bool localWellsActive()
const {
return ! wells_multisegment_.empty(); };
129 int numWellVars()
const {
return (num_phases_ + 1) * nseg_total_; };
131 template <
class ReservoirRes
idualQuant,
class SolutionState>
132 void extractWellPerfProperties(
const SolutionState& state,
133 const std::vector<ReservoirResidualQuant>& rq,
134 std::vector<ADB>& mob_perfcells,
135 std::vector<ADB>& b_perfcells)
const;
137 Vector& wellPerforationCellPressureDiffs() {
return well_perforation_cell_pressure_diffs_; };
139 Vector& wellSegmentPerforationDepthDiffs() {
return well_segment_perforation_depth_diffs_; };
141 const Vector& wellPerforationCellDensities()
const {
return well_perforation_cell_densities_; };
142 Vector& wellPerforationCellDensities() {
return well_perforation_cell_densities_; };
144 const std::vector<Vector>& segmentCompSurfVolumeInitial()
const {
return segment_comp_surf_volume_initial_; };
145 std::vector<Vector>& segmentCompSurfVolumeInitial() {
return segment_comp_surf_volume_initial_; };
147 const std::vector<ADB>& segmentCompSurfVolumeCurrent()
const {
return segment_comp_surf_volume_current_; };
149 const std::vector<int>& topWellSegments()
const {
return top_well_segments_; };
150 std::vector<int>& topWellSegments() {
return top_well_segments_; };
152 Vector& segVDt() {
return segvdt_; };
154 const Vector& wellPerforationDensities()
const {
return well_perforation_densities_; };
155 Vector& wellPerforationDensities() {
return well_perforation_densities_; };
157 const Vector& wellPerforationPressureDiffs()
const {
return well_perforation_pressure_diffs_; };
158 Vector& wellPerforationPressureDiffs() {
return well_perforation_pressure_diffs_; };
160 template <
class WellState>
162 updateWellState(
const Vector& dwells,
163 const double dpmaxrel,
164 WellState& well_state)
const;
168 template <
class SolutionState>
170 computeWellFlux(
const SolutionState& state,
171 const std::vector<ADB>& mob_perfcells,
172 const std::vector<ADB>& b_perfcells,
174 std::vector<ADB>& cq_s)
const;
176 template <
class SolutionState,
class WellState>
177 void updatePerfPhaseRatesAndPressures(
const std::vector<ADB>& cq_s,
178 const SolutionState& state,
179 WellState& xw)
const;
184 template <
class SolutionState>
186 computeSegmentFluidProperties(
const SolutionState& state);
189 computeSegmentPressuresDelta(
const double grav);
191 template <
class SolutionState>
193 addWellFluxEq(
const std::vector<ADB>& cq_s,
194 const SolutionState& state,
195 LinearisedBlackoilResidual& residual);
197 template <
class SolutionState,
class WellState>
199 addWellControlEq(
const SolutionState& state,
201 const Vector& aliveWells,
202 LinearisedBlackoilResidual& residual);
204 template <
class WellState>
206 updateWellControls(WellState& xw)
const;
211 variableStateWellIndices(std::vector<int>& indices,
214 template <
class SolutionState>
216 variableStateExtractWellsVars(
const std::vector<int>& indices,
217 std::vector<ADB>& vars,
218 SolutionState& state)
const;
221 variableWellStateIndices()
const;
223 template <
class WellState>
225 variableWellStateInitials(
const WellState& xw,
226 std::vector<Vector>& vars0)
const;
228 template <
class SolutionState,
class WellState>
229 void computeWellConnectionPressures(
const SolutionState& state,
231 const std::vector<ADB>& kr_adb,
232 const std::vector<ADB>& fluid_density);
234 WellCollection* wellCollection()
const;
236 void calculateEfficiencyFactors();
238 const Vector& wellPerfEfficiencyFactors()
const;
244 std::vector<WellMultiSegmentConstPtr> wells_multisegment_;
245 MultisegmentWellOps wops_ms_;
247 WellCollection* well_collection_;
252 Vector well_perforation_efficiency_factors_;
254 const int num_phases_;
264 const BlackoilPropsAdFromDeck* fluid_;
265 const std::vector<bool>* active_;
266 const std::vector<PhasePresence>* phase_condition_;
267 const VFPProperties* vfp_properties_;
271 Vector perf_cell_depth_;
278 Vector well_perforation_cell_pressure_diffs_;
281 ADB well_segment_perforation_pressure_diffs_;
284 Vector well_segment_perforation_depth_diffs_;
287 Vector perf_cell_depth_diffs_;
292 Vector well_perforation_cell_densities_;
296 ADB well_segment_densities_;
301 ADB well_segment_pressures_delta_;
305 std::vector<Vector> segment_comp_surf_volume_initial_;
308 std::vector<ADB> segment_comp_surf_volume_current_;
311 ADB segment_mass_flow_rates_;
316 ADB segment_viscosities_;
318 std::vector<int> top_well_segments_;
328 Vector well_perforation_densities_;
329 Vector well_perforation_pressure_diffs_;
335 #include "MultisegmentWells_impl.hpp"
337 #endif // OPM_MULTISEGMENTWELLS_HEADER_INCLUDED
A thin wrapper class that holds one VFPProdProperties and one VFPInjProperties object.
Definition: VFPProperties.hpp:37
AutoDiffMatrix is a wrapper class that optimizes matrix operations.
Definition: AutoDiffMatrix.hpp:43
Definition: MultisegmentWells.hpp:62
A class for forward-mode automatic differentiation with vector values and sparse jacobian matrices...
Definition: AutoDiffBlock.hpp:95
This class implements the AD-adapted fluid interface for three-phase black-oil.
Definition: BlackoilPropsAdFromDeck.hpp:61
Eigen::Array< double, Eigen::Dynamic, 1 > V
Underlying type for values.
Definition: AutoDiffBlock.hpp:99
Class for handling the multi-segment well model.
Definition: MultisegmentWells.hpp:55