22 #ifndef OPM_WELLOPERATORS_HEADER_INCLUDED
23 #define OPM_WELLOPERATORS_HEADER_INCLUDED
25 #include <dune/istl/operators.hh>
48 template <
class WellModel,
class X,
class Y>
52 using Base = Dune::LinearOperator<X, Y>;
53 using field_type =
typename Base::field_type;
63 void apply (
const X& x, Y& y)
const override
65 wellMod_.apply(x, y );
69 virtual void applyscaleadd (field_type alpha,
const X& x, Y& y)
const override
71 wellMod_.applyScaleAdd( alpha, x, y );
79 Dune::SolverCategory::Category
category()
const override
81 return Dune::SolverCategory::sequential;
84 const WellModel& wellMod_;
98 template<
class M,
class X,
class Y,
bool overlapping >
102 typedef M matrix_type;
103 typedef X domain_type;
104 typedef Y range_type;
105 typedef typename X::field_type field_type;
108 typedef Dune::OwnerOverlapCopyCommunication<int,int> communication_type;
110 typedef Dune::CollectiveCommunication< int > communication_type;
113 Dune::SolverCategory::Category category()
const override
116 Dune::SolverCategory::overlapping : Dune::SolverCategory::sequential;
121 const Dune::LinearOperator<X, Y>& wellOper,
122 const std::shared_ptr< communication_type >& comm = std::shared_ptr< communication_type >())
123 : A_( A ), wellOper_( wellOper ), comm_(comm)
127 virtual void apply(
const X& x, Y& y )
const override
132 wellOper_.apply(x, y );
141 virtual void applyscaleadd (field_type alpha,
const X& x, Y& y)
const override
146 wellOper_.applyscaleadd( alpha, x, y );
154 virtual const matrix_type& getmat()
const override {
return A_; }
157 const matrix_type& A_ ;
158 const Dune::LinearOperator<X, Y>& wellOper_;
159 std::shared_ptr< communication_type > comm_;
171 template<
class M,
class X,
class Y,
bool overlapping >
175 typedef M matrix_type;
176 typedef X domain_type;
177 typedef Y range_type;
178 typedef typename X::field_type field_type;
181 typedef Dune::OwnerOverlapCopyCommunication<int,int> communication_type;
183 typedef Dune::CollectiveCommunication< int > communication_type;
187 Dune::SolverCategory::Category category()
const override
190 Dune::SolverCategory::overlapping : Dune::SolverCategory::sequential;
195 const Dune::LinearOperator<X, Y>& wellOper,
196 const size_t interiorSize )
197 : A_( A ), wellOper_( wellOper ), interiorSize_(interiorSize)
200 virtual void apply(
const X& x, Y& y )
const override
202 for (
auto row = A_.begin(); row.index() < interiorSize_; ++row)
205 auto endc = (*row).end();
206 for (
auto col = (*row).begin(); col != endc; ++col)
207 (*col).umv(x[col.index()], y[row.index()]);
211 wellOper_.apply(x, y );
213 ghostLastProject( y );
217 virtual void applyscaleadd (field_type alpha,
const X& x, Y& y)
const override
219 for (
auto row = A_.begin(); row.index() < interiorSize_; ++row)
221 auto endc = (*row).end();
222 for (
auto col = (*row).begin(); col != endc; ++col)
223 (*col).usmv(alpha, x[col.index()], y[row.index()]);
226 wellOper_.applyscaleadd( alpha, x, y );
228 ghostLastProject( y );
231 virtual const matrix_type& getmat()
const override {
return A_; }
234 void ghostLastProject(Y& y)
const
236 size_t end = y.size();
237 for (
size_t i = interiorSize_; i < end; ++i)
241 const matrix_type& A_ ;
242 const Dune::LinearOperator<X, Y>& wellOper_;
243 size_t interiorSize_;
Linear operator wrapper for well model.
Definition: WellOperators.hpp:50
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:63
virtual void applyscaleadd(field_type alpha, const X &x, Y &y) const override
apply operator to x, scale and add:
Definition: WellOperators.hpp:69
Dune::SolverCategory::Category category() const override
Category for operator.
Definition: WellOperators.hpp:79
Adapter to combine a matrix and another linear operator into a combined linear operator.
Definition: WellOperators.hpp:173
WellModelGhostLastMatrixAdapter(const M &A, const Dune::LinearOperator< X, Y > &wellOper, const size_t interiorSize)
constructor: just store a reference to a matrix
Definition: WellOperators.hpp:194
Adapter to combine a matrix and another linear operator into a combined linear operator.
Definition: WellOperators.hpp:100
WellModelMatrixAdapter(const M &A, const Dune::LinearOperator< 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:120
This file contains a set of helper functions used by VFPProd / VFPInj.
Definition: BlackoilPhases.hpp:27