22#ifndef OPM_WELLOPERATORS_HEADER_INCLUDED
23#define OPM_WELLOPERATORS_HEADER_INCLUDED
25#include <dune/common/parallel/communication.hh>
26#include <dune/istl/operators.hh>
27#include <dune/istl/bcrsmatrix.hh>
29#include <opm/common/TimingMacros.hpp>
31#include <opm/simulators/linalg/matrixblock.hh>
32#include <dune/common/shared_ptr.hh>
33#include <dune/istl/paamg/smoother.hh>
55template <
class X,
class Y>
59 using field_type =
typename X::field_type;
60 using PressureMatrix = Dune::BCRSMatrix<MatrixBlock<field_type, 1, 1>>;
61 virtual void addWellPressureEquations(PressureMatrix& jacobian,
64 virtual void addWellPressureEquationsStruct(PressureMatrix& jacobian)
const = 0;
65 virtual int getNumberOfExtraEquations()
const = 0;
68template <
class WellModel,
class X,
class Y>
73 using field_type =
typename Base::field_type;
74 using PressureMatrix =
typename Base::PressureMatrix;
84 void apply(
const X& x, Y&
y)
const override
94 wellMod_.applyScaleAdd(alpha, x,
y);
102 Dune::SolverCategory::Category
category()
const override
104 return Dune::SolverCategory::sequential;
107 void addWellPressureEquations(PressureMatrix& jacobian,
115 void addWellPressureEquationsStruct(PressureMatrix& jacobian)
const override
118 wellMod_.addWellPressureEquationsStruct(jacobian);
121 int getNumberOfExtraEquations()
const override
123 return wellMod_.numLocalWellsEnd();
127 const WellModel& wellMod_;
131template <
class WellModel,
class X,
class Y>
137 using field_type =
typename WBase::field_type;
138 using PressureMatrix =
typename WBase::PressureMatrix;
140 void setDomainIndex(
int index) { domainIndex_ = index; }
142 void apply(
const X& x, Y&
y)
const override
145 this->wellMod_.applyDomain(x,
y, domainIndex_);
148 void applyscaleadd(field_type alpha,
const X& x, Y&
y)
const override
151 this->wellMod_.applyScaleAddDomain(alpha, x,
y, domainIndex_);
154 void addWellPressureEquations(PressureMatrix& jacobian,
159 this->wellMod_.addWellPressureEquationsDomain(jacobian, weights,
use_well_weights, domainIndex_);
163 int domainIndex_ = -1;
177template<
class M,
class X,
class Y,
bool overlapping >
181 using matrix_type = M;
182 using domain_type = X;
183 using range_type = Y;
184 using field_type =
typename X::field_type;
185 using PressureMatrix = Dune::BCRSMatrix<MatrixBlock<field_type, 1, 1>>;
187 using communication_type = Dune::OwnerOverlapCopyCommunication<int,int>;
189 using communication_type = Dune::Communication<int>;
192 Dune::SolverCategory::Category category()
const override
195 Dune::SolverCategory::overlapping : Dune::SolverCategory::sequential;
201 const std::shared_ptr<communication_type>& comm = {})
202 : A_(
A ), wellOper_(
wellOper ), comm_(comm)
205 void apply(
const X& x, Y&
y )
const override
211 wellOper_.apply(x,
y);
221 void applyscaleadd (field_type alpha,
const X& x, Y&
y)
const override
224 A_.usmv(alpha, x,
y);
227 wellOper_.applyscaleadd(alpha, x,
y);
236 const matrix_type& getmat()
const override {
return A_; }
238 void addWellPressureEquations(PressureMatrix& jacobian,
246 void addWellPressureEquationsStruct(PressureMatrix& jacobian)
const
249 wellOper_.addWellPressureEquationsStruct(jacobian);
252 int getNumberOfExtraEquations()
const
254 return wellOper_.getNumberOfExtraEquations();
258 const matrix_type& A_ ;
260 std::shared_ptr<communication_type> comm_;
271template<
class M,
class X,
class Y,
bool overlapping >
275 using matrix_type = M;
276 using domain_type = X;
277 using range_type = Y;
278 using field_type =
typename X::field_type;
279 using PressureMatrix = Dune::BCRSMatrix<MatrixBlock<field_type, 1, 1>>;
281 using communication_type = Dune::OwnerOverlapCopyCommunication<int,int>;
283 using communication_type = Dune::Communication<int>;
286 Dune::SolverCategory::Category category()
const override
289 Dune::SolverCategory::overlapping : Dune::SolverCategory::sequential;
299 void apply(
const X& x, Y&
y)
const override
302 for (
auto row = A_.begin();
row.index() < interiorSize_; ++
row)
305 auto endc = (*row).end();
307 (*col).umv(x[
col.index()],
y[
row.index()]);
311 wellOper_.apply(x,
y);
317 void applyscaleadd (field_type alpha,
const X& x, Y&
y)
const override
320 for (
auto row = A_.begin();
row.index() < interiorSize_; ++
row)
322 auto endc = (*row).end();
324 (*col).usmv(alpha, x[
col.index()],
y[
row.index()]);
327 wellOper_.applyscaleadd(alpha, x,
y);
332 const matrix_type& getmat()
const override {
return A_; }
334 void addWellPressureEquations(PressureMatrix& jacobian,
342 void addWellPressureEquationsStruct(PressureMatrix& jacobian)
const
345 wellOper_.addWellPressureEquationsStruct(jacobian);
348 int getNumberOfExtraEquations()
const
350 return wellOper_.getNumberOfExtraEquations();
354 void ghostLastProject(Y&
y)
const
356 std::size_t end =
y.size();
357 for (std::size_t i = interiorSize_; i < end; ++i)
361 const matrix_type& A_ ;
363 std::size_t interiorSize_;
374template<
class M,
class X,
class Y,
class C>
378 typedef M matrix_type;
379 typedef X domain_type;
380 typedef Y range_type;
381 typedef typename X::field_type field_type;
384 typedef C communication_type;
386 Dune::SolverCategory::Category category()
const override
388 return Dune::SolverCategory::overlapping;
393 const communication_type& comm)
396 interiorSize_ = setInteriorSize(comm_);
400 const communication_type& comm)
401 : A_(
A ), comm_(comm)
403 interiorSize_ = setInteriorSize(comm_);
406 virtual void apply(
const X& x, Y&
y )
const override
408 for (
auto row = A_->begin();
row.index() < interiorSize_; ++
row)
411 auto endc = (*row).end();
413 (*col).umv(x[
col.index()],
y[
row.index()]);
416 ghostLastProject(
y );
420 virtual void applyscaleadd (field_type alpha,
const X& x, Y&
y)
const override
422 for (
auto row = A_->begin();
row.index() < interiorSize_; ++
row)
424 auto endc = (*row).end();
426 (*col).usmv(alpha, x[
col.index()],
y[
row.index()]);
429 ghostLastProject(
y );
432 virtual const matrix_type& getmat()
const override {
return *A_; }
434 size_t getInteriorSize()
const {
return interiorSize_;}
437 void ghostLastProject(Y&
y)
const
439 size_t end =
y.size();
440 for (
size_t i = interiorSize_; i < end; ++i)
444 size_t setInteriorSize(
const communication_type& comm)
const
452 if (
idx->local().attribute()==1) {
454 auto loc =
idx->local().local();
463 const std::shared_ptr<const matrix_type> A_ ;
464 const communication_type& comm_;
465 size_t interiorSize_;
473template<
class M,
class X,
class Y,
class C>
474class ConstructionTraits<
Opm::GhostLastMatrixAdapter<M,X,Y,C> >
477 typedef ParallelOperatorArgs<M,C> Arguments;
479 static inline std::shared_ptr<Opm::GhostLastMatrixAdapter<M,X,Y,C>> construct(
const Arguments& args)
481 return std::make_shared<Opm::GhostLastMatrixAdapter<M,X,Y,C>>
482 (args.matrix_, args.comm_);
Definition WellOperators.hpp:133
Dune linear operator that assumes ghost rows are ordered after interior rows.
Definition WellOperators.hpp:376
GhostLastMatrixAdapter(const M &A, const communication_type &comm)
constructor: just store a reference to a matrix
Definition WellOperators.hpp:392
Definition WellOperators.hpp:70
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:84
void applyscaleadd(field_type alpha, const X &x, Y &y) const override
apply operator to x, scale and add:
Definition WellOperators.hpp:91
Dune::SolverCategory::Category category() const override
Category for operator.
Definition WellOperators.hpp:102
Adapter to combine a matrix and another linear operator into a combined linear operator.
Definition WellOperators.hpp:273
WellModelGhostLastMatrixAdapter(const M &A, const LinearOperatorExtra< X, Y > &wellOper, const std::size_t interiorSize)
constructor: just store a reference to a matrix
Definition WellOperators.hpp:293
Adapter to combine a matrix and another linear operator into a combined linear operator.
Definition WellOperators.hpp:179
WellModelMatrixAdapter(const M &A, const LinearOperatorExtra< X, Y > &wellOper, const std::shared_ptr< communication_type > &comm={})
constructor: just store a reference to a matrix
Definition WellOperators.hpp:199
This file contains a set of helper functions used by VFPProd / VFPInj.
Definition blackoilboundaryratevector.hh:37
constexpr auto getPropValue()
get the value data member of a property
Definition propertysystem.hh:242