20 #ifndef OPM_BLACKOIMULTISEGMENTLMODEL_IMPL_HEADER_INCLUDED
21 #define OPM_BLACKOIMULTISEGMENTLMODEL_IMPL_HEADER_INCLUDED
23 #include <opm/autodiff/BlackoilMultiSegmentModel.hpp>
25 #include <opm/autodiff/AutoDiffBlock.hpp>
26 #include <opm/autodiff/AutoDiffHelpers.hpp>
27 #include <opm/autodiff/GridHelpers.hpp>
28 #include <opm/autodiff/BlackoilPropsAdFromDeck.hpp>
29 #include <opm/autodiff/GeoProps.hpp>
30 #include <opm/autodiff/WellDensitySegmented.hpp>
31 #include <opm/autodiff/VFPProperties.hpp>
32 #include <opm/autodiff/VFPProdProperties.hpp>
33 #include <opm/autodiff/VFPInjProperties.hpp>
35 #include <opm/core/grid.h>
36 #include <opm/core/linalg/LinearSolverInterface.hpp>
37 #include <opm/core/linalg/ParallelIstlInformation.hpp>
38 #include <opm/core/props/rock/RockCompressibility.hpp>
39 #include <opm/common/ErrorMacros.hpp>
40 #include <opm/common/Exceptions.hpp>
41 #include <opm/parser/eclipse/Units/Units.hpp>
42 #include <opm/core/well_controls.h>
43 #include <opm/core/utility/parameters/ParameterGroup.hpp>
44 #include <opm/parser/eclipse/EclipseState/EclipseState.hpp>
64 const RockCompressibility* rock_comp_props,
67 std::shared_ptr< const EclipseState > eclState,
68 const bool has_disgas,
69 const bool has_vapoil,
70 const bool terminal_output)
71 :
Base(param, grid, fluid, geo, rock_comp_props, well_model, linsolver,
72 eclState, has_disgas, has_vapoil, terminal_output)
84 const ReservoirState& reservoir_state,
85 const WellState& well_state)
88 pvdt_ = geo_.poreVolume() / dt;
90 updatePrimalVariableFromState(reservoir_state);
93 const int nw = wellsMultiSegment().size();
95 if ( !msWellOps().has_multisegment_wells ) {
96 wellModel().segVDt() = V::Zero(nw);
100 const int nseg_total = well_state.numSegments();
101 std::vector<double> segment_volume;
102 segment_volume.reserve(nseg_total);
103 for (
int w = 0; w < nw; ++w) {
104 WellMultiSegmentConstPtr well = wellsMultiSegment()[w];
105 const std::vector<double>& segment_volume_well = well->segmentVolume();
106 segment_volume.insert(segment_volume.end(), segment_volume_well.begin(), segment_volume_well.end());
108 assert(
int(segment_volume.size()) == nseg_total);
109 wellModel().segVDt() = Eigen::Map<V>(segment_volume.data(), nseg_total) / dt;
116 template <
class Gr
id>
120 Base::makeConstantState(state);
129 template <
class Gr
id>
133 WellState& well_state,
134 const bool initial_assembly)
136 using namespace Opm::AutoDiffGrid;
151 wellModel().updateWellControls(well_state);
154 if (asImpl().wellModel().wellCollection()->groupControlActive()) {
156 Base::applyVREPGroupControl(reservoir_state, well_state);
158 asImpl().wellModel().wellCollection()->updateWellTargets(well_state.wellRates());
162 SolutionState state = asImpl().variableState(reservoir_state, well_state);
164 if (initial_assembly) {
167 asImpl().makeConstantState(state0);
170 asImpl().computeAccum(state0, 0);
171 wellModel().computeSegmentFluidProperties(state0);
172 const int np = numPhases();
173 assert(np ==
int(wellModel().segmentCompSurfVolumeInitial().size()));
174 for (
int phase = 0; phase < np; ++phase) {
175 wellModel().segmentCompSurfVolumeInitial()[phase] = wellModel().segmentCompSurfVolumeCurrent()[phase].value();
178 asImpl().computeWellConnectionPressures(state0, well_state);
191 asImpl().assembleMassBalanceEq(state);
194 if ( ! wellsActive() ) {
195 SimulatorReport report;
199 wellModel().computeSegmentFluidProperties(state);
201 const double gravity = detail::getGravity(geo_.gravity(), UgGridHelpers::dimensions(grid_));
202 wellModel().computeSegmentPressuresDelta(gravity);
204 std::vector<ADB> mob_perfcells;
205 std::vector<ADB> b_perfcells;
206 SimulatorReport report;
207 wellModel().extractWellPerfProperties(state, sd_.rq, mob_perfcells, b_perfcells);
208 if (param_.solve_welleq_initially_ && initial_assembly) {
210 report = asImpl().solveWellEq(mob_perfcells, b_perfcells, reservoir_state, state, well_state);
216 std::vector<ADB> cq_s;
217 wellModel().computeWellFlux(state, mob_perfcells, b_perfcells, aliveWells, cq_s);
218 wellModel().updatePerfPhaseRatesAndPressures(cq_s, state, well_state);
219 wellModel().addWellFluxEq(cq_s, state, residual_);
220 asImpl().addWellContributionToMassBalanceEq(cq_s, state, well_state);
221 wellModel().addWellControlEq(state, well_state, aliveWells, residual_);
229 template <
class Gr
id>
232 const std::vector<ADB>& b_perfcells,
233 const ReservoirState& reservoir_state,
234 SolutionState& state,
235 WellState& well_state)
237 SimulatorReport report = Base::solveWellEq(mob_perfcells, b_perfcells, reservoir_state, state, well_state);
239 if (report.converged) {
242 const int np = numPhases();
243 const int nseg_total =well_state.numSegments();
247 ADB::V new_segp = Eigen::Map<ADB::V>(well_state.segPress().data(), nseg_total);
250 std::vector<ADB::M> old_segp_derivs = state.segp.derivative();
251 state.segp =
ADB::function(std::move(new_segp), std::move(old_segp_derivs));
257 const DataBlock segrates = Eigen::Map<const DataBlock>(well_state.segPhaseRates().data(), nseg_total, np).transpose();
258 ADB::V new_segqs = Eigen::Map<const V>(segrates.data(), nseg_total * np);
259 std::vector<ADB::M> old_segqs_derivs = state.segqs.derivative();
260 state.segqs =
ADB::function(std::move(new_segqs), std::move(old_segqs_derivs));
265 asImpl().computeWellConnectionPressures(state, well_state);
274 template <
class Gr
id>
276 BlackoilMultiSegmentModel<Grid>::
277 computeWellConnectionPressures(
const SolutionState& state,
278 const WellState& well_state)
280 const int np = numPhases();
281 const std::vector<ADB> kr_adb = Base::computeRelPerm(state);
282 std::vector<ADB> fluid_density(np,
ADB::null());
284 for (
int phaseIdx = 0; phaseIdx < np; ++phaseIdx) {
285 const int canonicalPhaseIdx = canph_[phaseIdx];
286 fluid_density[phaseIdx] = fluidDensity(canonicalPhaseIdx, sd_.rq[phaseIdx].b, state.rs, state.rv);
288 wellModel().computeWellConnectionPressures(state, well_state, kr_adb, fluid_density);
293 #endif // OPM_BLACKOILMODELBASE_IMPL_HEADER_INCLUDED
BlackoilMultiSegmentModel(const typename Base::ModelParameters ¶m, const Grid &grid, const BlackoilPropsAdFromDeck &fluid, const DerivedGeology &geo, const RockCompressibility *rock_comp_props, const MultisegmentWells &well_model, const NewtonIterationBlackoilInterface &linsolver, std::shared_ptr< const EclipseState > eclState, const bool has_disgas, const bool has_vapoil, const bool terminal_output)
Construct the model.
Definition: BlackoilMultiSegmentModel_impl.hpp:60
static AutoDiffBlock null()
Construct an empty AutoDiffBlock.
Definition: AutoDiffBlock.hpp:104
Definition: BlackoilMultiSegmentModel.hpp:35
void prepareStep(const SimulatorTimerInterface &timer, const ReservoirState &reservoir_state, const WellState &well_state)
Called once before each time step.
Definition: BlackoilMultiSegmentModel_impl.hpp:83
virtual double currentStepLength() const =0
Current step length.
SimulatorReport assemble(const ReservoirState &reservoir_state, WellState &well_state, const bool initial_assembly)
Assemble the residual and Jacobian of the nonlinear system.
Definition: BlackoilMultiSegmentModel_impl.hpp:132
A model implementation for three-phase black oil with support for multi-segment wells.
Definition: BlackoilMultiSegmentModel.hpp:55
static AutoDiffBlock constant(V &&val)
Create an AutoDiffBlock representing a constant.
Definition: AutoDiffBlock.hpp:111
Interface class for SimulatorTimer objects, to be improved.
Definition: SimulatorTimerInterface.hpp:35
Class containing static geological properties that are derived from grid and petrophysical properties...
Definition: GeoProps.hpp:59
static AutoDiffBlock function(V &&val, std::vector< M > &&jac)
Create an AutoDiffBlock by directly specifying values and jacobians.
Definition: AutoDiffBlock.hpp:184
Interface class for (linear) solvers for the fully implicit black-oil system.
Definition: NewtonIterationBlackoilInterface.hpp:31
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