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>
28#include <opm/simulators/linalg/WellOperators.hpp>
32namespace Opm {
namespace Details {
33 template<
class Scalar>
using PressureMatrixType = Dune::BCRSMatrix<MatrixBlock<Scalar, 1, 1>>;
34 template<
class Scalar>
using PressureVectorType = Dune::BlockVector<Dune::FieldVector<Scalar, 1>>;
35 template<
class Scalar>
using SeqCoarseOperatorType = Dune::MatrixAdapter<PressureMatrixType<Scalar>,
36 PressureVectorType<Scalar>,
37 PressureVectorType<Scalar>>;
38 template<
class Scalar,
class Comm>
39 using ParCoarseOperatorType
41 PressureVectorType<Scalar>,
42 PressureVectorType<Scalar>,
44 template<
class Scalar,
class Comm>
45 using CoarseOperatorType = std::conditional_t<std::is_same<Comm, Dune::Amg::SequentialInformation>::value,
46 SeqCoarseOperatorType<Scalar>,
47 ParCoarseOperatorType<Scalar,Comm>>;
52template <
class FineOperator,
class Communication,
class Scalar,
bool transpose = false>
57 using CoarseOperator =
typename Details::CoarseOperatorType<Scalar,Communication>;
59 using ParallelInformation = Communication;
60 using FineVectorType =
typename FineOperator::domain_type;
64 const FineVectorType& weights,
66 int pressure_var_index)
67 : communication_(&
const_cast<Communication&
>(comm))
69 , pressure_var_index_(pressure_var_index)
75 using CoarseMatrix =
typename CoarseOperator::matrix_type;
76 const auto& fineLevelMatrix = fineOperator.getmat();
77 coarseLevelMatrix_.reset(
new CoarseMatrix(fineLevelMatrix.N(), fineLevelMatrix.M(), CoarseMatrix::row_wise));
78 auto createIter = coarseLevelMatrix_->createbegin();
80 for (
const auto& row : fineLevelMatrix) {
81 for (
auto col = row.begin(), cend = row.end(); col != cend; ++col) {
82 createIter.insert(col.index());
88 coarseLevelCommunication_.reset(communication_, [](Communication*) {});
90 this->
lhs_.resize(this->coarseLevelMatrix_->M());
91 this->
rhs_.resize(this->coarseLevelMatrix_->N());
92 using OperatorArgs =
typename Dune::Amg::ConstructionTraits<CoarseOperator>::Arguments;
93 OperatorArgs oargs(coarseLevelMatrix_, *coarseLevelCommunication_);
94 this->
operator_ = Dune::Amg::ConstructionTraits<CoarseOperator>::construct(oargs);
99 const auto& fineMatrix = fineOperator.getmat();
100 *coarseLevelMatrix_ = 0;
101 auto rowCoarse = coarseLevelMatrix_->begin();
102 for (
auto row = fineMatrix.begin(), rowEnd = fineMatrix.end(); row != rowEnd; ++row, ++rowCoarse) {
103 assert(row.index() == rowCoarse.index());
104 auto entryCoarse = rowCoarse->begin();
105 for (
auto entry = row->begin(), entryEnd = row->end(); entry != entryEnd; ++entry, ++entryCoarse) {
106 assert(entry.index() == entryCoarse.index());
107 Scalar matrix_el = 0;
109 const auto& bw = weights_[entry.index()];
110 for (std::size_t i = 0; i < bw.size(); ++i) {
111 matrix_el += (*entry)[pressure_var_index_][i] * bw[i];
114 const auto& bw = weights_[row.index()];
115 for (std::size_t i = 0; i < bw.size(); ++i) {
116 matrix_el += (*entry)[i][pressure_var_index_] * bw[i];
119 (*entryCoarse) = matrix_el;
122 assert(rowCoarse == coarseLevelMatrix_->end());
130 auto end = fine.end(), begin = fine.begin();
132 for (
auto block = begin; block != end; ++block) {
133 const auto& bw = weights_[block.index()];
136 rhs_el = (*block)[pressure_var_index_];
138 for (std::size_t i = 0; i < block->size(); ++i) {
139 rhs_el += (*block)[i] * bw[i];
142 this->
rhs_[block - begin] = rhs_el;
150 auto end = fine.end(), begin = fine.begin();
152 for (
auto block = begin; block != end; ++block) {
154 const auto& bw = weights_[block.index()];
155 for (std::size_t i = 0; i < block->size(); ++i) {
156 (*block)[i] = this->
lhs_[block - begin] * bw[i];
159 (*block)[pressure_var_index_] = this->
lhs_[block - begin];
169 const Communication& getCoarseLevelCommunication()
const
171 return *coarseLevelCommunication_;
174 std::size_t getPressureIndex()
const
176 return pressure_var_index_;
180 Communication* communication_;
181 const FineVectorType& weights_;
182 const std::size_t pressure_var_index_;
183 std::shared_ptr<Communication> coarseLevelCommunication_;
184 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::domain_type FineDomainType
The type of the domain of the fine level operator.
Definition twolevelmethodcpr.hh:58
CoarseDomainType lhs_
Definition twolevelmethodcpr.hh:141
CoarseRangeType rhs_
Definition twolevelmethodcpr.hh:139
FineOperatorType::range_type FineRangeType
The type of the range of the fine level operator.
Definition twolevelmethodcpr.hh:54
std::shared_ptr< CoarseOperatorType > operator_
Definition twolevelmethodcpr.hh:143
Dune linear operator that assumes ghost rows are ordered after interior rows.
Definition WellOperators.hpp:376
Definition PressureTransferPolicy.hpp:55
void moveToFineLevel(typename ParentType::FineDomainType &fine) override
Updates the fine level linear system after the correction of the coarse levels system.
Definition PressureTransferPolicy.hpp:148
PressureTransferPolicy * clone() const override
Clone the current object.
Definition PressureTransferPolicy.hpp:164
void calculateCoarseEntries(const FineOperator &fineOperator) override
???.
Definition PressureTransferPolicy.hpp:97
void createCoarseLevelSystem(const FineOperator &fineOperator) override
Algebraically creates the coarse level system.
Definition PressureTransferPolicy.hpp:73
Definition PropertyTree.hpp:37
This file contains a set of helper functions used by VFPProd / VFPInj.
Definition blackoilboundaryratevector.hh:37
Algebraic twolevel methods.