21 #ifndef OPM_PRESSURE_TRANSFER_POLICY_HEADER_INCLUDED
22 #define OPM_PRESSURE_TRANSFER_POLICY_HEADER_INCLUDED
31 template <
class FineOperator,
class CoarseOperator,
class Communication,
bool transpose = false>
36 typedef Communication ParallelInformation;
37 typedef typename FineOperator::domain_type FineVectorType;
41 : communication_(&
const_cast<Communication&
>(comm))
43 , pressure_var_index_(pressure_var_index)
49 using CoarseMatrix =
typename CoarseOperator::matrix_type;
50 const auto& fineLevelMatrix = fineOperator.getmat();
51 coarseLevelMatrix_.reset(
new CoarseMatrix(fineLevelMatrix.N(), fineLevelMatrix.M(), CoarseMatrix::row_wise));
52 auto createIter = coarseLevelMatrix_->createbegin();
54 for (
const auto& row : fineLevelMatrix) {
55 for (
auto col = row.begin(), cend = row.end(); col != cend; ++col) {
56 createIter.insert(col.index());
62 coarseLevelCommunication_.reset(communication_, [](Communication*) {});
64 this->
lhs_.resize(this->coarseLevelMatrix_->M());
65 this->
rhs_.resize(this->coarseLevelMatrix_->N());
66 using OperatorArgs =
typename Dune::Amg::ConstructionTraits<CoarseOperator>::Arguments;
67 #if DUNE_VERSION_NEWER(DUNE_ISTL, 2, 7)
68 OperatorArgs oargs(coarseLevelMatrix_, *coarseLevelCommunication_);
69 this->
operator_ = Dune::Amg::ConstructionTraits<CoarseOperator>::construct(oargs);
71 OperatorArgs oargs(*coarseLevelMatrix_, *coarseLevelCommunication_);
72 this->
operator_.reset(Dune::Amg::ConstructionTraits<CoarseOperator>::construct(oargs));
78 const auto& fineMatrix = fineOperator.getmat();
79 *coarseLevelMatrix_ = 0;
80 auto rowCoarse = coarseLevelMatrix_->begin();
81 for (
auto row = fineMatrix.begin(), rowEnd = fineMatrix.end(); row != rowEnd; ++row, ++rowCoarse) {
82 assert(row.index() == rowCoarse.index());
83 auto entryCoarse = rowCoarse->begin();
84 for (
auto entry = row->begin(), entryEnd = row->end(); entry != entryEnd; ++entry, ++entryCoarse) {
85 assert(entry.index() == entryCoarse.index());
88 const auto& bw = weights_[entry.index()];
89 for (
size_t i = 0; i < bw.size(); ++i) {
90 matrix_el += (*entry)[pressure_var_index_][i] * bw[i];
93 const auto& bw = weights_[row.index()];
94 for (
size_t i = 0; i < bw.size(); ++i) {
95 matrix_el += (*entry)[i][pressure_var_index_] * bw[i];
98 (*entryCoarse) = matrix_el;
101 assert(rowCoarse == coarseLevelMatrix_->end());
109 auto end = fine.end(), begin = fine.begin();
111 for (
auto block = begin; block != end; ++block) {
112 const auto& bw = weights_[block.index()];
115 rhs_el = (*block)[pressure_var_index_];
117 for (
size_t i = 0; i < block->size(); ++i) {
118 rhs_el += (*block)[i] * bw[i];
121 this->
rhs_[block - begin] = rhs_el;
129 auto end = fine.end(), begin = fine.begin();
131 for (
auto block = begin; block != end; ++block) {
133 const auto& bw = weights_[block.index()];
134 for (
size_t i = 0; i < block->size(); ++i) {
135 (*block)[i] = this->
lhs_[block - begin] * bw[i];
138 (*block)[pressure_var_index_] = this->
lhs_[block - begin];
148 const Communication& getCoarseLevelCommunication()
const
150 return *coarseLevelCommunication_;
153 std::size_t getPressureIndex()
const
155 return pressure_var_index_;
158 Communication* communication_;
159 const FineVectorType& weights_;
160 const std::size_t pressure_var_index_;
161 std::shared_ptr<Communication> coarseLevelCommunication_;
162 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: PressureTransferPolicy.hpp:33
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:127
virtual void createCoarseLevelSystem(const FineOperator &fineOperator) override
Algebraically creates the coarse level system.
Definition: PressureTransferPolicy.hpp:47
virtual void calculateCoarseEntries(const FineOperator &fineOperator) override
???.
Definition: PressureTransferPolicy.hpp:76
virtual PressureTransferPolicy * clone() const override
Clone the current object.
Definition: PressureTransferPolicy.hpp:143
This file contains a set of helper functions used by VFPProd / VFPInj.
Definition: BlackoilPhases.hpp:27
Algebraic twolevel methods.