21#ifndef OPM_PRESSURE_TRANSFER_POLICY_HEADER_INCLUDED
22#define OPM_PRESSURE_TRANSFER_POLICY_HEADER_INCLUDED
26#include <opm/simulators/linalg/PropertyTree.hpp>
27#include <opm/simulators/linalg/matrixblock.hh>
36 using PressureMatrixType = Dune::BCRSMatrix<Opm::MatrixBlock<double, 1, 1>>;
37 using PressureVectorType = Dune::BlockVector<Dune::FieldVector<double, 1>>;
38 using SeqCoarseOperatorType = Dune::MatrixAdapter<PressureMatrixType, PressureVectorType, PressureVectorType>;
40 using ParCoarseOperatorType
41 = Dune::OverlappingSchwarzOperator<PressureMatrixType, PressureVectorType, PressureVectorType, Comm>;
43 using CoarseOperatorType = std::conditional_t<std::is_same<Comm, Dune::Amg::SequentialInformation>::value,
44 SeqCoarseOperatorType,
45 ParCoarseOperatorType<Comm>>;
50template <
class FineOperator,
class Communication,
bool transpose = false>
55 typedef typename Details::CoarseOperatorType<Communication> CoarseOperator;
57 typedef Communication ParallelInformation;
58 typedef typename FineOperator::domain_type FineVectorType;
62 const FineVectorType& weights,
65 : communication_(&
const_cast<Communication&
>(comm))
73 using CoarseMatrix =
typename CoarseOperator::matrix_type;
76 auto createIter = coarseLevelMatrix_->createbegin();
86 coarseLevelCommunication_.reset(communication_, [](Communication*) {});
88 this->
lhs_.resize(this->coarseLevelMatrix_->M());
89 this->
rhs_.resize(this->coarseLevelMatrix_->N());
90 using OperatorArgs =
typename Dune::Amg::ConstructionTraits<CoarseOperator>::Arguments;
92 this->
operator_ = Dune::Amg::ConstructionTraits<CoarseOperator>::construct(
oargs);
98 *coarseLevelMatrix_ = 0;
99 auto rowCoarse = coarseLevelMatrix_->begin();
107 const auto&
bw = weights_[
entry.index()];
108 for (std::size_t
i = 0;
i <
bw.size(); ++
i) {
112 const auto&
bw = weights_[
row.index()];
113 for (std::size_t
i = 0;
i <
bw.size(); ++
i) {
128 auto end =
fine.end(), begin =
fine.begin();
131 const auto&
bw = weights_[
block.index()];
134 rhs_el = (*block)[pressure_var_index_];
136 for (std::size_t i = 0; i < block->size(); ++i) {
137 rhs_el += (*block)[i] * bw[i];
140 this->
rhs_[block - begin] = rhs_el;
148 auto end =
fine.end(), begin =
fine.begin();
152 const auto&
bw = weights_[
block.index()];
153 for (std::size_t
i = 0;
i <
block->size(); ++
i) {
157 (*block)[pressure_var_index_] = this->
lhs_[
block - begin];
167 const Communication& getCoarseLevelCommunication()
const
169 return *coarseLevelCommunication_;
172 std::size_t getPressureIndex()
const
174 return pressure_var_index_;
178 Communication* communication_;
179 const FineVectorType& weights_;
180 const std::size_t pressure_var_index_;
181 std::shared_ptr<Communication> coarseLevelCommunication_;
182 std::shared_ptr<typename CoarseOperator::matrix_type> coarseLevelMatrix_;
Abstract base class for transfer between levels and creation of the coarse level system.
Definition twolevelmethodcpr.hh:44
FineOperatorType::range_type FineRangeType
The type of the range of the fine level operator.
Definition twolevelmethodcpr.hh:54
CoarseDomainType lhs_
The coarse level lhs.
Definition twolevelmethodcpr.hh:141
CoarseRangeType rhs_
The coarse level rhs.
Definition twolevelmethodcpr.hh:139
FineOperatorType::domain_type FineDomainType
The type of the domain of the fine level operator.
Definition twolevelmethodcpr.hh:58
std::shared_ptr< CoarseOperatorType > operator_
the coarse level linear operator.
Definition twolevelmethodcpr.hh:143
Definition AquiferInterface.hpp:35
Definition PressureTransferPolicy.hpp:53
virtual PressureTransferPolicy * clone() const override
Clone the current object.
Definition PressureTransferPolicy.hpp:162
virtual void createCoarseLevelSystem(const FineOperator &fineOperator) override
Algebraically creates the coarse level system.
Definition PressureTransferPolicy.hpp:71
virtual void calculateCoarseEntries(const FineOperator &fineOperator) override
???.
Definition PressureTransferPolicy.hpp:95
virtual void moveToFineLevel(typename ParentType::FineDomainType &fine) override
Updates the fine level linear system after the correction of the coarse levels system.
Definition PressureTransferPolicy.hpp:146
Definition PropertyTree.hpp:37
This file contains a set of helper functions used by VFPProd / VFPInj.
Definition BlackoilPhases.hpp:27
Algebraic twolevel methods.