23 #ifndef OPM_STANDARDWELL_GENERIC_HEADER_INCLUDED
24 #define OPM_STANDARDWELL_GENERIC_HEADER_INCLUDED
26 #include <dune/common/dynmatrix.hh>
27 #include <dune/common/dynvector.hh>
28 #include <dune/istl/bcrsmatrix.hh>
29 #include <dune/istl/bvector.hh>
31 #include <opm/simulators/wells/WellHelpers.hpp>
39 class ConvergenceReport;
41 class ParallelWellInfo;
44 class WellInterfaceGeneric;
47 template<
class Scalar>
56 using VectorBlockWellType = Dune::DynamicVector<Scalar>;
57 using BVectorWell = Dune::BlockVector<VectorBlockWellType>;
60 using DiagMatrixBlockWellType = Dune::DynamicMatrix<Scalar>;
61 using DiagMatWell = Dune::BCRSMatrix<DiagMatrixBlockWellType>;
64 using OffDiagMatrixBlockWellType = Dune::DynamicMatrix<Scalar>;
65 using OffDiagMatWell = Dune::BCRSMatrix<OffDiagMatrixBlockWellType>;
68 #if HAVE_CUDA || HAVE_OPENCL
70 void getNumBlocks(
unsigned int& _nnzs)
const;
78 static double relaxationFactorRate(
const std::vector<double>& primary_variables,
79 const BVectorWell& dwells);
82 static double relaxationFactorFraction(
const double old_value,
85 double calculateThpFromBhp(
const WellState& well_state,
86 const std::vector<double>& rates,
91 void checkConvergenceControlEq(
const WellState& well_state,
94 const double max_residual_allowed)
const;
96 void checkConvergencePolyMW(
const std::vector<double>& res,
98 const double maxResidualAllowed)
const;
100 void computeConnectionPressureDelta();
102 std::optional<double> computeBhpAtThpLimitInj(
const std::function<std::vector<double>(
const double)>& frates,
103 const SummaryState& summary_state,
105 std::optional<double> computeBhpAtThpLimitProdWithAlq(
const std::function<std::vector<double>(
const double)>& frates,
106 const SummaryState& summary_state,
109 double alq_value)
const;
115 BVectorWell resWell_;
118 std::vector<double> perf_densities_;
120 std::vector<double> perf_pressure_diffs_;
123 OffDiagMatWell duneB_;
124 OffDiagMatWell duneC_;
126 DiagMatWell invDuneD_;
132 mutable BVectorWell Bx_;
133 mutable BVectorWell invDrw_;
135 double getRho()
const {
return perf_densities_[0]; }
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:27