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 
36 namespace Opm
37 {
38 
39 class ConvergenceReport;
40 class DeferredLogger;
41 class ParallelWellInfo;
42 class Schedule;
43 class SummaryState;
44 class WellInterfaceGeneric;
45 class WellState;
46 
47 template<class Scalar>
49 {
50 protected:
51  // sparsity pattern for the matrices
52  //[A C^T [x = [ res
53  // B D ] x_well] res_well]
54 
55  // the vector type for the res_well and x_well
56  using VectorBlockWellType = Dune::DynamicVector<Scalar>;
57  using BVectorWell = Dune::BlockVector<VectorBlockWellType>;
58 
59  // the matrix type for the diagonal matrix D
60  using DiagMatrixBlockWellType = Dune::DynamicMatrix<Scalar>;
61  using DiagMatWell = Dune::BCRSMatrix<DiagMatrixBlockWellType>;
62 
63  // the matrix type for the non-diagonal matrix B and C^T
64  using OffDiagMatrixBlockWellType = Dune::DynamicMatrix<Scalar>;
65  using OffDiagMatWell = Dune::BCRSMatrix<OffDiagMatrixBlockWellType>;
66 
67 public:
68 #if HAVE_CUDA || HAVE_OPENCL
70  void getNumBlocks(unsigned int& _nnzs) const;
71 #endif
72 
73 protected:
74  StandardWellGeneric(int Bhp,
75  const WellInterfaceGeneric& baseif);
76 
77  // calculate a relaxation factor to avoid overshoot of total rates
78  static double relaxationFactorRate(const std::vector<double>& primary_variables,
79  const BVectorWell& dwells);
80 
81  // relaxation factor considering only one fraction value
82  static double relaxationFactorFraction(const double old_value,
83  const double dx);
84 
85  double calculateThpFromBhp(const WellState& well_state,
86  const std::vector<double>& rates,
87  const double bhp,
88  DeferredLogger& deferred_logger) const;
89 
90  // checking the convergence of the well control equations
91  void checkConvergenceControlEq(const WellState& well_state,
92  ConvergenceReport& report,
93  DeferredLogger& deferred_logger,
94  const double max_residual_allowed) const;
95 
96  void checkConvergencePolyMW(const std::vector<double>& res,
97  ConvergenceReport& report,
98  const double maxResidualAllowed) const;
99 
100  void computeConnectionPressureDelta();
101 
102  std::optional<double> computeBhpAtThpLimitInj(const std::function<std::vector<double>(const double)>& frates,
103  const SummaryState& summary_state,
104  DeferredLogger& deferred_logger) const;
105  std::optional<double> computeBhpAtThpLimitProdWithAlq(const std::function<std::vector<double>(const double)>& frates,
106  const SummaryState& summary_state,
107  DeferredLogger& deferred_logger,
108  double alq_value) const;
109 
110  void gliftDebug(const std::string &msg,
111  DeferredLogger& deferred_logger) const;
112  bool checkGliftNewtonIterationIdxOk(const int report_step_idx,
113  const int iteration_idx,
114  const Schedule& schedule,
115  DeferredLogger& deferred_logger) const;
116  bool doGasLiftOptimize(const WellState &well_state,
117  const int report_step_idx,
118  const int iteration_idx,
119  const Schedule& schedule,
120  DeferredLogger& deferred_logger) const;
121 
122  // Base interface reference
123  const WellInterfaceGeneric& baseif_;
124 
125  // residuals of the well equations
126  BVectorWell resWell_;
127 
128  // densities of the fluid in each perforation
129  std::vector<double> perf_densities_;
130  // pressure drop between different perforations
131  std::vector<double> perf_pressure_diffs_;
132 
133  // Enable GLIFT debug mode. This will enable output of logging messages.
134  bool glift_debug = false;
135  // Optimize only wells under THP control
136  bool glift_optimize_only_thp_wells = true;
137 
138  // two off-diagonal matrices
139  OffDiagMatWell duneB_;
140  OffDiagMatWell duneC_;
141  // diagonal matrix for the well
142  DiagMatWell invDuneD_;
143 
144  // Wrapper for the parallel application of B for distributed wells
146 
147  // several vector used in the matrix calculation
148  mutable BVectorWell Bx_;
149  mutable BVectorWell invDrw_;
150 
151  double getRho() const { return perf_densities_[0]; }
152 
153 private:
154  int Bhp_; // index of Bhp
155 };
156 
157 }
158 
159 #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:49
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:58
This file contains a set of helper functions used by VFPProd / VFPInj.
Definition: BlackoilPhases.hpp:26