My Project
StandardWellGeneric.hpp
1/*
2 Copyright 2017 SINTEF Digital, Mathematics and Cybernetics.
3 Copyright 2017 Statoil ASA.
4 Copyright 2016 - 2017 IRIS AS.
5
6 This file is part of the Open Porous Media project (OPM).
7
8 OPM is free software: you can redistribute it and/or modify
9 it under the terms of the GNU General Public License as published by
10 the Free Software Foundation, either version 3 of the License, or
11 (at your option) any later version.
12
13 OPM is distributed in the hope that it will be useful,
14 but WITHOUT ANY WARRANTY; without even the implied warranty of
15 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
16 GNU General Public License for more details.
17
18 You should have received a copy of the GNU General Public License
19 along with OPM. If not, see <http://www.gnu.org/licenses/>.
20*/
21
22
23#ifndef OPM_STANDARDWELL_GENERIC_HEADER_INCLUDED
24#define OPM_STANDARDWELL_GENERIC_HEADER_INCLUDED
25
26#include <dune/common/dynmatrix.hh>
27#include <dune/common/dynvector.hh>
28#include <dune/istl/bcrsmatrix.hh>
29#include <dune/istl/bvector.hh>
30
31#include <opm/simulators/wells/WellHelpers.hpp>
32
33#include <optional>
34#include <vector>
35#include <functional>
36
37namespace Opm
38{
39
40class ConvergenceReport;
41class DeferredLogger;
42class ParallelWellInfo;
43class Schedule;
44class SummaryState;
45class WellInterfaceGeneric;
46class WellState;
47
48template<class Scalar>
50{
51protected:
52 // sparsity pattern for the matrices
53 //[A C^T [x = [ res
54 // B D ] x_well] res_well]
55
56 // the vector type for the res_well and x_well
57 using VectorBlockWellType = Dune::DynamicVector<Scalar>;
58 using BVectorWell = Dune::BlockVector<VectorBlockWellType>;
59
60 // the matrix type for the diagonal matrix D
61 using DiagMatrixBlockWellType = Dune::DynamicMatrix<Scalar>;
62 using DiagMatWell = Dune::BCRSMatrix<DiagMatrixBlockWellType>;
63
64 // the matrix type for the non-diagonal matrix B and C^T
65 using OffDiagMatrixBlockWellType = Dune::DynamicMatrix<Scalar>;
66 using OffDiagMatWell = Dune::BCRSMatrix<OffDiagMatrixBlockWellType>;
67
68public:
70 void getNumBlocks(unsigned int& _nnzs) const;
71
72protected:
73 StandardWellGeneric(int Bhp,
74 const WellInterfaceGeneric& baseif);
75
76 // calculate a relaxation factor to avoid overshoot of total rates
77 static double relaxationFactorRate(const std::vector<double>& primary_variables,
78 const BVectorWell& dwells);
79
80 // relaxation factor considering only one fraction value
81 static double relaxationFactorFraction(const double old_value,
82 const double dx);
83
84 double calculateThpFromBhp(const WellState& well_state,
85 const std::vector<double>& rates,
86 const double bhp,
87 DeferredLogger& deferred_logger) const;
88
89 // checking the convergence of the well control equations
90 void checkConvergenceControlEq(const WellState& well_state,
91 ConvergenceReport& report,
92 DeferredLogger& deferred_logger,
93 const double max_residual_allowed) const;
94
95 void checkConvergencePolyMW(const std::vector<double>& res,
96 ConvergenceReport& report,
97 const double maxResidualAllowed) const;
98
99 void computeConnectionPressureDelta();
100
101 std::optional<double> computeBhpAtThpLimitInj(const std::function<std::vector<double>(const double)>& frates,
102 const SummaryState& summary_state,
103 DeferredLogger& deferred_logger) const;
104 std::optional<double> computeBhpAtThpLimitProdWithAlq(const std::function<std::vector<double>(const double)>& frates,
105 const SummaryState& summary_state,
106 DeferredLogger& deferred_logger,
107 double maxPerfPress,
108 double alq_value) const;
109
110 // Base interface reference
111 const WellInterfaceGeneric& baseif_;
112
113 // residuals of the well equations
114 BVectorWell resWell_;
115
116 // densities of the fluid in each perforation
117 std::vector<double> perf_densities_;
118 // pressure drop between different perforations
119 std::vector<double> perf_pressure_diffs_;
120
121 // two off-diagonal matrices
122 OffDiagMatWell duneB_;
123 OffDiagMatWell duneC_;
124 // diagonal matrix for the well
125 DiagMatWell invDuneD_;
126 DiagMatWell duneD_;
127
128 // Wrapper for the parallel application of B for distributed wells
130
131 // several vector used in the matrix calculation
132 mutable BVectorWell Bx_;
133 mutable BVectorWell invDrw_;
134
135 double getRho() const { return perf_densities_[0]; }
136
137private:
138 int Bhp_; // index of Bhp
139};
140
141}
142
143#endif // OPM_STANDARDWELL_GENERIC_HEADER_INCLUDED
Represents the convergence status of the whole simulator, to make it possible to query and store the ...
Definition: ConvergenceReport.hpp:36
Definition: DeferredLogger.hpp:57
Definition: StandardWellGeneric.hpp:50
void getNumBlocks(unsigned int &_nnzs) const
get the number of blocks of the C and B matrices, used to allocate memory in a WellContributions obje...
Definition: StandardWellGeneric.cpp:524
Definition: WellInterfaceGeneric.hpp:51
The state of a set of wells, tailored for use by the fully implicit blackoil simulator.
Definition: WellState.hpp:56
A wrapper around the B matrix for distributed wells.
Definition: WellHelpers.hpp:49
This file contains a set of helper functions used by VFPProd / VFPInj.
Definition: BlackoilPhases.hpp:27