24 #ifndef OPM_BLACKOILWELLMODEL_HEADER_INCLUDED 25 #define OPM_BLACKOILWELLMODEL_HEADER_INCLUDED 27 #include <opm/common/OpmLog/OpmLog.hpp> 29 #include <opm/common/utility/platform_dependent/disable_warnings.h> 30 #include <opm/common/utility/platform_dependent/reenable_warnings.h> 35 #include <opm/parser/eclipse/EclipseState/Schedule/Schedule.hpp> 37 #include <opm/core/wells.h> 38 #include <opm/core/wells/DynamicListEconLimited.hpp> 39 #include <opm/core/wells/WellCollection.hpp> 40 #include <opm/core/simulator/SimulatorReport.hpp> 41 #include <opm/autodiff/VFPProperties.hpp> 42 #include <opm/autodiff/WellHelpers.hpp> 43 #include <opm/autodiff/BlackoilModelEnums.hpp> 44 #include <opm/autodiff/WellDensitySegmented.hpp> 45 #include <opm/autodiff/BlackoilPropsAdFromDeck.hpp> 46 #include <opm/autodiff/BlackoilDetails.hpp> 47 #include <opm/autodiff/BlackoilModelParameters.hpp> 48 #include <opm/autodiff/WellStateFullyImplicitBlackoil.hpp> 50 #include <opm/autodiff/WellInterface.hpp> 51 #include <opm/autodiff/StandardWell.hpp> 52 #include <opm/autodiff/MultisegmentWell.hpp> 53 #include<dune/common/fmatrix.hh> 54 #include<dune/istl/bcrsmatrix.hh> 55 #include<dune/istl/matrixmatrix.hh> 57 #include <opm/material/densead/Math.hpp> 59 #include <opm/simulators/WellSwitchingLogger.hpp> 65 template<
typename TypeTag>
72 typedef typename GET_PROP_TYPE(TypeTag, Grid) Grid;
73 typedef typename GET_PROP_TYPE(TypeTag, FluidSystem) FluidSystem;
74 typedef typename GET_PROP_TYPE(TypeTag, ElementContext) ElementContext;
75 typedef typename GET_PROP_TYPE(TypeTag, Indices) BlackoilIndices;
76 typedef typename GET_PROP_TYPE(TypeTag, Simulator) Simulator;
77 typedef typename GET_PROP_TYPE(TypeTag, Scalar) Scalar;
79 static const int numEq = BlackoilIndices::numEq;
80 static const int solventSaturationIdx = BlackoilIndices::solventSaturationIdx;
84 typedef Dune::FieldVector<Scalar, numEq > VectorBlockType;
85 typedef Dune::BlockVector<VectorBlockType> BVector;
87 typedef Ewoms::BlackOilPolymerModule<TypeTag> PolymerModule;
91 SurfaceToReservoirVoidage<BlackoilPropsAdFromDeck::FluidSystem, std::vector<int> >;
95 WellCollection* well_collection,
96 const std::vector< const Well* >& wells_ecl,
99 const bool terminal_output,
100 const int current_index,
101 const std::vector<int>& pvt_region_idx);
103 void init(
const PhaseUsage phase_usage_arg,
104 const std::vector<bool>& active_arg,
105 const double gravity_arg,
106 const std::vector<double>& depth_arg,
110 void setVFPProperties(
const VFPProperties* vfp_properties_arg);
113 SimulatorReport assemble(Simulator& ebosSimulator,
114 const int iterationIdx,
119 void apply( BVector& r)
const;
122 void apply(
const BVector& x, BVector& Ax)
const;
125 void applyScaleAdd(
const Scalar alpha,
const BVector& x, BVector& Ax)
const;
129 void recoverWellSolutionAndUpdateWellState(
const BVector& x,
WellState& well_state)
const;
131 int numWells()
const;
136 void setWellsActive(
const bool wells_active);
141 bool getWellConvergence(
const Simulator& ebosSimulator,
142 const std::vector<Scalar>& B_avg)
const;
146 const int current_step,
147 const Wells* wells_struct,
149 DynamicListEconLimited& list_econ_limited)
const;
151 WellCollection* wellCollection()
const;
157 const std::vector< const Well* > wells_ecl_;
162 const int number_of_wells_;
164 const int number_of_phases_;
168 using WellInterfacePtr = std::unique_ptr<WellInterface<TypeTag> >;
173 std::vector<WellInterfacePtr > well_container_;
178 static std::vector<WellInterfacePtr > createWellContainer(
const Wells* wells,
179 const std::vector<const Well*>& wells_ecl,
180 const bool use_multisegment_well,
185 WellCollection* well_collection_;
187 bool terminal_output_;
190 int current_timeIdx_;
192 PhaseUsage phase_usage_;
193 std::vector<bool> active_;
195 const std::vector<int>& pvt_region_idx_;
198 int number_of_cells_;
203 mutable BVector scaleAddRes_;
205 void updateWellControls(
WellState& xw)
const;
207 void updateGroupControls(
WellState& well_state)
const;
210 void updatePrimaryVariables(
const WellState& well_state)
const;
212 void setupCompressedToCartesian(
const int* global_cell,
int number_of_cells, std::map<int,int>& cartesian_to_compressed )
const;
214 void computeRepRadiusPerfLength(
const Grid& grid);
217 void computeAverageFormationFactor(
const Simulator& ebosSimulator,
218 std::vector<double>& B_avg)
const;
220 void applyVREPGroupControl(
WellState& well_state)
const;
222 void computeWellVoidageRates(
const WellState& well_state,
223 std::vector<double>& well_voidage_rates,
224 std::vector<double>& voidage_conversion_coeffs)
const;
228 void computeWellPotentials(
const Simulator& ebosSimulator,
230 std::vector<double>& well_potentials)
const;
232 const std::vector<double>& wellPerfEfficiencyFactors()
const;
234 void calculateEfficiencyFactors();
244 SimulatorReport solveWellEq(Simulator& ebosSimulator,
248 void initPrimaryVariablesEvaluation()
const;
251 int numComponents()
const 253 if (numPhases() == 2) {
256 int numComp = FluidSystem::numComponents;
264 int numPhases()
const;
266 int flowPhaseToEbosPhaseIdx(
const int phaseIdx )
const;
268 void resetWellControlFromState(
const WellState& xw)
const;
270 void assembleWellEq(Simulator& ebosSimulator,
273 bool only_wells)
const;
277 void prepareTimeStep(
const Simulator& ebos_simulator,
280 void prepareGroupControl(
const Simulator& ebos_simulator,
288 #include "BlackoilWellModel_impl.hpp" A thin wrapper class that holds one VFPProdProperties and one VFPInjProperties object.
Definition: VFPProperties.hpp:37
bool localWellsActive() const
return true if wells are available on this process
Definition: BlackoilWellModel_impl.hpp:378
Class for handling the blackoil well model.
Definition: BlackoilWellModel.hpp:66
bool wellsActive() const
return true if wells are available in the reservoir
Definition: BlackoilWellModel_impl.hpp:354
a struct to collect information about the convergence checking
Definition: WellInterface.hpp:117
This file contains a set of helper functions used by VFPProd / VFPInj.
Definition: AdditionalObjectDeleter.hpp:22
Facility for converting component rates at surface conditions to phase (voidage) rates at reservoir c...
Solver parameters for the BlackoilModel.
Definition: BlackoilModelParameters.hpp:31
The state of a set of wells, tailored for use by the fully implicit blackoil simulator.
Definition: WellStateFullyImplicitBlackoil.hpp:44
void calculateExplicitQuantities(const Simulator &ebosSimulator, const WellState &xw) const
Calculating the explict quantities used in the well calculation.
Definition: BlackoilWellModel_impl.hpp:540
void updateListEconLimited(const Schedule &schedule, const int current_step, const Wells *wells_struct, const WellState &well_state, DynamicListEconLimited &list_econ_limited) const
upate the dynamic lists related to economic limits
Definition: BlackoilWellModel_impl.hpp:579