59 using PressureMatrix = Dune::BCRSMatrix<Opm::MatrixBlock<double, 1, 1>>;
60 virtual void addWellPressureEquations(PressureMatrix& jacobian,
const X& weights,
const bool use_well_weights)
const = 0;
61 virtual void addWellPressureEquationsStruct(PressureMatrix& jacobian)
const = 0;
62 virtual int getNumberOfExtraEquations()
const = 0;
70 using field_type =
typename Base::field_type;
71 using PressureMatrix =
typename Base::PressureMatrix;
80 void apply(
const X& x, Y& y)
const override
87 virtual void applyscaleadd(field_type alpha,
const X& x, Y& y)
const override
90 wellMod_.applyScaleAdd(alpha, x, y);
98 Dune::SolverCategory::Category
category()
const override
100 return Dune::SolverCategory::sequential;
102 void addWellPressureEquations(PressureMatrix& jacobian,
const X& weights,
const bool use_well_weights)
const override
104 OPM_TIMEBLOCK(addWellPressureEquations);
105 wellMod_.addWellPressureEquations(jacobian, weights, use_well_weights);
107 void addWellPressureEquationsStruct(PressureMatrix& jacobian)
const override
109 OPM_TIMEBLOCK(addWellPressureEquationsStruct);
110 wellMod_.addWellPressureEquationsStruct(jacobian);
112 int getNumberOfExtraEquations()
const override
114 return wellMod_.numLocalWellsEnd();
118 const WellModel& wellMod_;
136 typedef M matrix_type;
137 typedef X domain_type;
138 typedef Y range_type;
139 typedef typename X::field_type field_type;
140 using PressureMatrix = Dune::BCRSMatrix<Opm::MatrixBlock<double, 1, 1>>;
142 typedef Dune::OwnerOverlapCopyCommunication<int,int> communication_type;
144 typedef Dune::CollectiveCommunication< int > communication_type;
147 Dune::SolverCategory::Category category()
const override
150 Dune::SolverCategory::overlapping : Dune::SolverCategory::sequential;
156 const std::shared_ptr< communication_type >& comm = std::shared_ptr< communication_type >())
157 : A_( A ), wellOper_( wellOper ), comm_(comm)
161 virtual void apply(
const X& x, Y& y )
const override
163 OPM_TIMEBLOCK(apply);
167 wellOper_.apply(x, y );
176 virtual void applyscaleadd (field_type alpha,
const X& x, Y& y)
const override
178 OPM_TIMEBLOCK(applyscaleadd);
182 wellOper_.applyscaleadd( alpha, x, y );
190 virtual const matrix_type& getmat()
const override {
return A_; }
192 void addWellPressureEquations(PressureMatrix& jacobian,
const X& weights,
const bool use_well_weights)
const
194 OPM_TIMEBLOCK(addWellPressureEquations);
195 wellOper_.addWellPressureEquations(jacobian, weights, use_well_weights);
197 void addWellPressureEquationsStruct(PressureMatrix& jacobian)
const
199 OPM_TIMEBLOCK(addWellPressureEquations);
200 wellOper_.addWellPressureEquationsStruct(jacobian);
202 int getNumberOfExtraEquations()
const
204 return wellOper_.getNumberOfExtraEquations();
208 const matrix_type& A_ ;
210 std::shared_ptr< communication_type > comm_;
226 typedef M matrix_type;
227 typedef X domain_type;
228 typedef Y range_type;
229 typedef typename X::field_type field_type;
230 using PressureMatrix = Dune::BCRSMatrix<Opm::MatrixBlock<double, 1, 1>>;
232 typedef Dune::OwnerOverlapCopyCommunication<int,int> communication_type;
234 typedef Dune::CollectiveCommunication< int > communication_type;
238 Dune::SolverCategory::Category category()
const override
241 Dune::SolverCategory::overlapping : Dune::SolverCategory::sequential;
247 const std::size_t interiorSize )
248 : A_( A ), wellOper_( wellOper ), interiorSize_(interiorSize)
251 virtual void apply(
const X& x, Y& y )
const override
253 OPM_TIMEBLOCK(apply);
254 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).umv(x[col.index()], y[row.index()]);
263 wellOper_.apply(x, y );
265 ghostLastProject( y );
269 virtual void applyscaleadd (field_type alpha,
const X& x, Y& y)
const override
271 OPM_TIMEBLOCK(applyscaleadd);
272 for (
auto row = A_.begin(); row.index() < interiorSize_; ++row)
274 auto endc = (*row).end();
275 for (
auto col = (*row).begin(); col != endc; ++col)
276 (*col).usmv(alpha, x[col.index()], y[row.index()]);
279 wellOper_.applyscaleadd( alpha, x, y );
281 ghostLastProject( y );
284 virtual const matrix_type& getmat()
const override {
return A_; }
286 void addWellPressureEquations(PressureMatrix& jacobian,
const X& weights,
const bool use_well_weights)
const
288 OPM_TIMEBLOCK(addWellPressureEquations);
289 wellOper_.addWellPressureEquations(jacobian, weights, use_well_weights);
291 void addWellPressureEquationsStruct(PressureMatrix& jacobian)
const
293 OPM_TIMEBLOCK(addWellPressureEquationsStruct);
294 wellOper_.addWellPressureEquationsStruct(jacobian);
296 int getNumberOfExtraEquations()
const
298 return wellOper_.getNumberOfExtraEquations();
302 void ghostLastProject(Y& y)
const
304 std::size_t end = y.size();
305 for (std::size_t i = interiorSize_; i < end; ++i)
309 const matrix_type& A_ ;
311 std::size_t interiorSize_;
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:80
virtual void applyscaleadd(field_type alpha, const X &x, Y &y) const override
apply operator to x, scale and add:
Definition WellOperators.hpp:87
WellModelGhostLastMatrixAdapter(const M &A, const Opm::LinearOperatorExtra< X, Y > &wellOper, const std::size_t interiorSize)
constructor: just store a reference to a matrix
Definition WellOperators.hpp:245
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:154