22#ifndef OPM_WELLOPERATORS_HEADER_INCLUDED
23#define OPM_WELLOPERATORS_HEADER_INCLUDED
25#include <dune/istl/operators.hh>
26#include <opm/simulators/linalg/matrixblock.hh>
48template <
class X,
class Y>
52 using PressureMatrix = Dune::BCRSMatrix<Opm::MatrixBlock<double, 1, 1>>;
53 virtual void addWellPressureEquations(PressureMatrix& jacobian,
const X& weights,
const bool use_well_weights)
const = 0;
54 virtual void addWellPressureEquationsStruct(PressureMatrix& jacobian)
const = 0;
55 virtual int getNumberOfExtraEquations()
const = 0;
58template <
class WellModel,
class X,
class Y>
63 using field_type =
typename Base::field_type;
64 using PressureMatrix =
typename Base::PressureMatrix;
73 void apply(
const X& x, Y& y)
const override
79 virtual void applyscaleadd(field_type alpha,
const X& x, Y& y)
const override
81 wellMod_.applyScaleAdd(alpha, x, y);
89 Dune::SolverCategory::Category
category()
const override
91 return Dune::SolverCategory::sequential;
93 void addWellPressureEquations(PressureMatrix& jacobian,
const X& weights,
const bool use_well_weights)
const override
95 wellMod_.addWellPressureEquations(jacobian, weights, use_well_weights);
97 void addWellPressureEquationsStruct(PressureMatrix& jacobian)
const override
99 wellMod_.addWellPressureEquationsStruct(jacobian);
101 int getNumberOfExtraEquations()
const override
103 return wellMod_.numLocalWellsEnd();
107 const WellModel& wellMod_;
121template<
class M,
class X,
class Y,
bool overlapping >
125 typedef M matrix_type;
126 typedef X domain_type;
127 typedef Y range_type;
128 typedef typename X::field_type field_type;
129 using PressureMatrix = Dune::BCRSMatrix<Opm::MatrixBlock<double, 1, 1>>;
131 typedef Dune::OwnerOverlapCopyCommunication<int,int> communication_type;
133 typedef Dune::CollectiveCommunication< int > communication_type;
136 Dune::SolverCategory::Category category()
const override
139 Dune::SolverCategory::overlapping : Dune::SolverCategory::sequential;
145 const std::shared_ptr< communication_type >& comm = std::shared_ptr< communication_type >())
146 : A_( A ), wellOper_( wellOper ), comm_(comm)
150 virtual void apply(
const X& x, Y& y )
const override
155 wellOper_.apply(x, y );
164 virtual void applyscaleadd (field_type alpha,
const X& x, Y& y)
const override
169 wellOper_.applyscaleadd( alpha, x, y );
177 virtual const matrix_type& getmat()
const override {
return A_; }
179 void addWellPressureEquations(PressureMatrix& jacobian,
const X& weights,
const bool use_well_weights)
const
181 wellOper_.addWellPressureEquations(jacobian, weights, use_well_weights);
183 void addWellPressureEquationsStruct(PressureMatrix& jacobian)
const
185 wellOper_.addWellPressureEquationsStruct(jacobian);
187 int getNumberOfExtraEquations()
const
189 return wellOper_.getNumberOfExtraEquations();
193 const matrix_type& A_ ;
195 std::shared_ptr< communication_type > comm_;
207template<
class M,
class X,
class Y,
bool overlapping >
211 typedef M matrix_type;
212 typedef X domain_type;
213 typedef Y range_type;
214 typedef typename X::field_type field_type;
215 using PressureMatrix = Dune::BCRSMatrix<Opm::MatrixBlock<double, 1, 1>>;
217 typedef Dune::OwnerOverlapCopyCommunication<int,int> communication_type;
219 typedef Dune::CollectiveCommunication< int > communication_type;
223 Dune::SolverCategory::Category category()
const override
226 Dune::SolverCategory::overlapping : Dune::SolverCategory::sequential;
232 const size_t interiorSize )
233 : A_( A ), wellOper_( wellOper ), interiorSize_(interiorSize)
236 virtual void apply(
const X& x, Y& y )
const override
238 for (
auto row = A_.begin(); row.index() < interiorSize_; ++row)
241 auto endc = (*row).end();
242 for (
auto col = (*row).begin(); col != endc; ++col)
243 (*col).umv(x[col.index()], y[row.index()]);
247 wellOper_.apply(x, y );
249 ghostLastProject( y );
253 virtual void applyscaleadd (field_type alpha,
const X& x, Y& y)
const override
255 for (
auto row = A_.begin(); row.index() < interiorSize_; ++row)
257 auto endc = (*row).end();
258 for (
auto col = (*row).begin(); col != endc; ++col)
259 (*col).usmv(alpha, x[col.index()], y[row.index()]);
262 wellOper_.applyscaleadd( alpha, x, y );
264 ghostLastProject( y );
267 virtual const matrix_type& getmat()
const override {
return A_; }
269 void addWellPressureEquations(PressureMatrix& jacobian,
const X& weights,
const bool use_well_weights)
const
271 wellOper_.addWellPressureEquations(jacobian, weights, use_well_weights);
273 void addWellPressureEquationsStruct(PressureMatrix& jacobian)
const
275 wellOper_.addWellPressureEquationsStruct(jacobian);
277 int getNumberOfExtraEquations()
const
279 return wellOper_.getNumberOfExtraEquations();
283 void ghostLastProject(Y& y)
const
285 size_t end = y.size();
286 for (
size_t i = interiorSize_; i < end; ++i)
290 const matrix_type& A_ ;
292 size_t interiorSize_;
Definition: WellOperators.hpp:60
void apply(const X &x, Y &y) const override
apply operator to x: The input vector is consistent and the output must also be consistent on the in...
Definition: WellOperators.hpp:73
virtual void applyscaleadd(field_type alpha, const X &x, Y &y) const override
apply operator to x, scale and add:
Definition: WellOperators.hpp:79
Dune::SolverCategory::Category category() const override
Category for operator.
Definition: WellOperators.hpp:89
Adapter to combine a matrix and another linear operator into a combined linear operator.
Definition: WellOperators.hpp:209
WellModelGhostLastMatrixAdapter(const M &A, const Opm::LinearOperatorExtra< X, Y > &wellOper, const size_t interiorSize)
constructor: just store a reference to a matrix
Definition: WellOperators.hpp:230
Adapter to combine a matrix and another linear operator into a combined linear operator.
Definition: WellOperators.hpp:123
WellModelMatrixAdapter(const M &A, const Opm::LinearOperatorExtra< X, Y > &wellOper, const std::shared_ptr< communication_type > &comm=std::shared_ptr< communication_type >())
constructor: just store a reference to a matrix
Definition: WellOperators.hpp:143
This file contains a set of helper functions used by VFPProd / VFPInj.
Definition: BlackoilPhases.hpp:27