22#ifndef OPM_MULTISEGMENTWELL_PRIMARY_VARIABLES_HEADER_INCLUDED
23#define OPM_MULTISEGMENTWELL_PRIMARY_VARIABLES_HEADER_INCLUDED
25#include <opm/material/densead/Evaluation.hpp>
27#include <opm/simulators/wells/MultisegmentWellEquations.hpp>
28#include <opm/input/eclipse/Schedule/SummaryState.hpp>
29#include <opm/input/eclipse/Units/Units.hpp>
43template<
class Flu
idSystem,
class Indices>
60 static constexpr bool has_water = (Indices::waterSwitchIdx >= 0);
61 static constexpr bool has_gas = (Indices::compositionSwitchIdx >= 0);
62 static constexpr bool has_oil = (Indices::numPhases - has_gas - has_water) > 0;
67 static constexpr bool has_wfrac_variable = has_water && Indices::numPhases > 1;
68 static constexpr bool has_gfrac_variable = has_gas && has_oil;
70 static constexpr int WQTotal = 0;
71 static constexpr int WFrac = has_wfrac_variable ? 1 : -1000;
72 static constexpr int GFrac = has_gfrac_variable ? has_wfrac_variable + 1 : -1000;
73 static constexpr int SPres = has_wfrac_variable + has_gfrac_variable + 1;
76 static constexpr int numWellEq = Indices::numPhases + 1;
78 using Scalar =
typename FluidSystem::Scalar;
79 using EvalWell = DenseAd::Evaluation<Scalar, Indices::numEq + numWellEq>;
82 using BVectorWell =
typename Equations::BVectorWell;
89 void resize(
const int numSegments);
96 const bool stop_or_zero_rate_target);
100 const Scalar relaxation_factor,
101 const Scalar DFLimit,
102 const bool stop_or_zero_rate_target,
103 const Scalar max_pressure_change);
109 const SummaryState& summary_state,
115 const int compIdx)
const;
120 const int compIdx)
const;
124 const int seg_upwind,
125 const std::size_t comp_idx)
const;
135 const int comp_idx)
const;
138 EvalWell
getQs(
const int comp_idx)
const;
144 const std::array<EvalWell, numWellEq>&
eval(
const int idx)
const
145 {
return evaluation_[idx]; }
148 const std::array<Scalar, numWellEq>&
value(
const int idx)
const
149 {
return value_[idx]; }
152 void setValue(
const int idx,
const std::array<Scalar, numWellEq>& val)
153 { value_[idx] = val; }
160 void processFractions(
const int seg);
163 EvalWell volumeFraction(
const int seg,
164 const unsigned compIdx)
const;
168 std::vector<std::array<Scalar, numWellEq>> value_;
172 std::vector<std::array<EvalWell, numWellEq>> evaluation_;
176 static constexpr double bhp_lower_limit = 1. * unit::barsa - 1. * unit::Pascal;
177 static constexpr double seg_pres_lower_limit = 0.;
Definition DeferredLogger.hpp:57
Definition MultisegmentWellPrimaryVariables.hpp:39
Definition MultisegmentWellPrimaryVariables.hpp:45
void resize(const int numSegments)
Resize values and evaluations.
Definition MultisegmentWellPrimaryVariables.cpp:49
EvalWell volumeFractionScaled(const int seg, const int compIdx) const
Returns scaled volume fraction for a component in a segment.
Definition MultisegmentWellPrimaryVariables.cpp:522
void outputLowLimitPressureSegments(DeferredLogger &deferred_logger) const
output the segments with pressure close to lower pressure limit for debugging purpose
Definition MultisegmentWellPrimaryVariables.cpp:636
EvalWell getSegmentRate(const int seg, const int comp_idx) const
Get rate for a component in a segment.
Definition MultisegmentWellPrimaryVariables.cpp:611
void init()
Initialize evaluations from values.
Definition MultisegmentWellPrimaryVariables.cpp:57
EvalWell getWQTotal() const
Get WQTotal.
Definition MultisegmentWellPrimaryVariables.cpp:628
void copyToWellState(const MultisegmentWellGeneric< Scalar > &mswell, const Scalar rho, WellState< Scalar > &well_state, const SummaryState &summary_state, DeferredLogger &deferred_logger) const
Copy values to well state.
Definition MultisegmentWellPrimaryVariables.cpp:215
const std::array< EvalWell, numWellEq > & eval(const int idx) const
Returns a const ref to an array of evaluations.
Definition MultisegmentWellPrimaryVariables.hpp:144
EvalWell getBhp() const
Get bottomhole pressure.
Definition MultisegmentWellPrimaryVariables.cpp:603
EvalWell surfaceVolumeFraction(const int seg, const int compIdx) const
Returns surface volume fraction for a component in a segment.
Definition MultisegmentWellPrimaryVariables.cpp:539
void setValue(const int idx, const std::array< Scalar, numWellEq > &val)
Set a value array. Note that this does not also set the corresponding evaluation.
Definition MultisegmentWellPrimaryVariables.hpp:152
EvalWell getQs(const int comp_idx) const
Returns scaled rate for a component.
Definition MultisegmentWellPrimaryVariables.cpp:620
void update(const WellState< Scalar > &well_state, const bool stop_or_zero_rate_target)
Copy values from well state.
Definition MultisegmentWellPrimaryVariables.cpp:70
EvalWell getSegmentPressure(const int seg) const
Get pressure for a segment.
Definition MultisegmentWellPrimaryVariables.cpp:595
void updateNewton(const BVectorWell &dwells, const Scalar relaxation_factor, const Scalar DFLimit, const bool stop_or_zero_rate_target, const Scalar max_pressure_change)
Update values from newton update vector.
Definition MultisegmentWellPrimaryVariables.cpp:159
EvalWell getSegmentRateUpwinding(const int seg, const int seg_upwind, const std::size_t comp_idx) const
Returns upwinding rate for a component in a segment.
Definition MultisegmentWellPrimaryVariables.cpp:555
const std::array< Scalar, numWellEq > & value(const int idx) const
Returns a value array.
Definition MultisegmentWellPrimaryVariables.hpp:148
Definition WellInterfaceIndices.hpp:34
The state of a set of wells, tailored for use by the fully implicit blackoil simulator.
Definition WellState.hpp:62
This file contains a set of helper functions used by VFPProd / VFPInj.
Definition blackoilboundaryratevector.hh:37