My Project
MultisegmentWellPrimaryVariables.hpp
1/*
2 Copyright 2017 SINTEF Digital, Mathematics and Cybernetics.
3 Copyright 2017 Statoil ASA.
4
5 This file is part of the Open Porous Media project (OPM).
6
7 OPM is free software: you can redistribute it and/or modify
8 it under the terms of the GNU General Public License as published by
9 the Free Software Foundation, either version 3 of the License, or
10 (at your option) any later version.
11
12 OPM is distributed in the hope that it will be useful,
13 but WITHOUT ANY WARRANTY; without even the implied warranty of
14 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
15 GNU General Public License for more details.
16
17 You should have received a copy of the GNU General Public License
18 along with OPM. If not, see <http://www.gnu.org/licenses/>.
19*/
20
21
22#ifndef OPM_MULTISEGMENTWELL_PRIMARY_VARIABLES_HEADER_INCLUDED
23#define OPM_MULTISEGMENTWELL_PRIMARY_VARIABLES_HEADER_INCLUDED
24
25#include <opm/material/densead/Evaluation.hpp>
26
27#include <opm/simulators/wells/MultisegmentWellEquations.hpp>
28
29#include <array>
30#include <vector>
31
32namespace Opm
33{
34
35class DeferredLogger;
36template<class Scalar> class MultisegmentWellGeneric;
37template<class FluidSystem, class Indices, class Scalar> class WellInterfaceIndices;
38class WellState;
39
40template<class FluidSystem, class Indices, class Scalar>
42{
43public:
44 // TODO: for now, not considering the polymer, solvent and so on to simplify the development process.
45
46 // TODO: we need to have order for the primary variables and also the order for the well equations.
47 // sometimes, they are similar, while sometimes, they can have very different forms.
48
49 // Table showing the primary variable indices, depending on what phases are present:
50 //
51 // WOG OG WG WO W/O/G (single phase)
52 // WQTotal 0 0 0 0 0
53 // WFrac 1 -1000 1 1 -1000
54 // GFrac 2 1 -1000 -1000 -1000
55 // Spres 3 2 2 2 1
56
57 static constexpr bool has_water = (Indices::waterSwitchIdx >= 0);
58 static constexpr bool has_gas = (Indices::compositionSwitchIdx >= 0);
59 static constexpr bool has_oil = (Indices::numPhases - has_gas - has_water) > 0;
60
61 // In the implementation, one should use has_wfrac_variable
62 // rather than has_water to check if you should do something
63 // with the variable at the WFrac location, similar for GFrac.
64 static constexpr bool has_wfrac_variable = has_water && Indices::numPhases > 1;
65 static constexpr bool has_gfrac_variable = has_gas && has_oil;
66
67 static constexpr int WQTotal = 0;
68 static constexpr int WFrac = has_wfrac_variable ? 1 : -1000;
69 static constexpr int GFrac = has_gfrac_variable ? has_wfrac_variable + 1 : -1000;
70 static constexpr int SPres = has_wfrac_variable + has_gfrac_variable + 1;
71
72 // the number of well equations TODO: it should have a more general strategy for it
73 static constexpr int numWellEq = Indices::numPhases + 1;
74
75 using EvalWell = DenseAd::Evaluation<double, /*size=*/Indices::numEq + numWellEq>;
76
78 using BVectorWell = typename Equations::BVectorWell;
79
81 : well_(well)
82 {}
83
85 void resize(const int numSegments);
86
88 void init();
89
91 void update(const WellState& well_state, const bool stop_or_zero_rate_target);
92
94 void updateNewton(const BVectorWell& dwells,
95 const double relaxation_factor,
96 const double DFLimit,
97 const bool stop_or_zero_rate_target,
98 const double max_pressure_change);
99
102 const double rho,
103 const bool stop_or_zero_rate_target,
104 WellState& well_state,
105 DeferredLogger& deferred_logger) const;
106
109 EvalWell volumeFractionScaled(const int seg,
110 const int compIdx) const;
111
114 EvalWell surfaceVolumeFraction(const int seg,
115 const int compIdx) const;
116
118 EvalWell getSegmentRateUpwinding(const int seg,
119 const int seg_upwind,
120 const size_t comp_idx) const;
121
123 EvalWell getBhp() const;
124
126 EvalWell getSegmentPressure(const int seg) const;
127
129 EvalWell getSegmentRate(const int seg,
130 const int comp_idx) const;
131
133 EvalWell getQs(const int comp_idx) const;
134
136 EvalWell getWQTotal() const;
137
139 const std::array<EvalWell,numWellEq>& eval(const int idx) const
140 { return evaluation_[idx]; }
141
142private:
144 void processFractions(const int seg);
145
147 EvalWell volumeFraction(const int seg,
148 const unsigned compIdx) const;
149
152 std::vector<std::array<double, numWellEq>> value_;
153
156 std::vector<std::array<EvalWell, numWellEq>> evaluation_;
157
159};
160
161}
162
163#endif // OPM_MULTISEGMENTWELL_PRIMARY_VARIABLES_HEADER_INCLUDED
Definition: DeferredLogger.hpp:57
Definition: MultisegmentWellGeneric.hpp:42
Definition: MultisegmentWellPrimaryVariables.hpp:42
EvalWell getQs(const int comp_idx) const
Returns scaled rate for a component.
Definition: MultisegmentWellPrimaryVariables.cpp:603
void resize(const int numSegments)
Resize values and evaluations.
Definition: MultisegmentWellPrimaryVariables.cpp:46
void copyToWellState(const MultisegmentWellGeneric< Scalar > &mswell, const double rho, const bool stop_or_zero_rate_target, WellState &well_state, DeferredLogger &deferred_logger) const
Copy values to well state.
Definition: MultisegmentWellPrimaryVariables.cpp:208
void updateNewton(const BVectorWell &dwells, const double relaxation_factor, const double DFLimit, const bool stop_or_zero_rate_target, const double max_pressure_change)
Update values from newton update vector.
Definition: MultisegmentWellPrimaryVariables.cpp:155
EvalWell getSegmentPressure(const int seg) const
Get pressure for a segment.
Definition: MultisegmentWellPrimaryVariables.cpp:578
EvalWell getBhp() const
Get bottomhole pressure.
Definition: MultisegmentWellPrimaryVariables.cpp:586
EvalWell getWQTotal() const
Get WQTotal.
Definition: MultisegmentWellPrimaryVariables.cpp:611
EvalWell volumeFractionScaled(const int seg, const int compIdx) const
Returns scaled volume fraction for a component in a segment.
Definition: MultisegmentWellPrimaryVariables.cpp:505
void update(const WellState &well_state, const bool stop_or_zero_rate_target)
Copy values from well state.
Definition: MultisegmentWellPrimaryVariables.cpp:67
EvalWell getSegmentRateUpwinding(const int seg, const int seg_upwind, const size_t comp_idx) const
Returns upwinding rate for a component in a segment.
Definition: MultisegmentWellPrimaryVariables.cpp:538
void init()
Initialize evaluations from values.
Definition: MultisegmentWellPrimaryVariables.cpp:54
EvalWell surfaceVolumeFraction(const int seg, const int compIdx) const
Returns surface volume fraction for a component in a segment.
Definition: MultisegmentWellPrimaryVariables.cpp:522
const std::array< EvalWell, numWellEq > & eval(const int idx) const
Returns a const ref to an evaluation.
Definition: MultisegmentWellPrimaryVariables.hpp:139
EvalWell getSegmentRate(const int seg, const int comp_idx) const
Get rate for a component in a segment.
Definition: MultisegmentWellPrimaryVariables.cpp:594
Definition: WellInterfaceIndices.hpp:33
The state of a set of wells, tailored for use by the fully implicit blackoil simulator.
Definition: WellState.hpp:60
This file contains a set of helper functions used by VFPProd / VFPInj.
Definition: BlackoilPhases.hpp:27