[cleanup][istl] Adjust to shared_ptr change in dune-istl >= 2.7.

This commit is contained in:
Robert Kloefkorn 2019-05-07 16:59:00 +02:00
parent c4bcd9e467
commit 6c77fae891
6 changed files with 129 additions and 68 deletions

View File

@ -127,7 +127,7 @@ scaleMatrixDRS(const Operator& op, std::size_t pressureEqnIndex, const Vector& w
for (auto i = matrix->begin(); i != endi; ++i) { for (auto i = matrix->begin(); i != endi; ++i) {
const BlockVector& bw = weights[i.index()]; const BlockVector& bw = weights[i.index()];
const auto endj = (*i).end(); const auto endj = (*i).end();
for (auto j = (*i).begin(); j != endj; ++j) { for (auto j = (*i).begin(); j != endj; ++j) {
Block& block = *j; Block& block = *j;
BlockVector& bvec = block[pressureEqnIndex]; BlockVector& bvec = block[pressureEqnIndex];
// should introduce limits which also change the weights // should introduce limits which also change the weights
@ -387,7 +387,11 @@ private:
cargs.setMatrix(op.getmat()); cargs.setMatrix(op.getmat());
cargs.setComm(comm); cargs.setComm(comm);
cargs.setArgs(args); cargs.setArgs(args);
#if DUNE_VERSION_NEWER(DUNE_ISTL, 2, 7)
smoother_ = Dune::Amg::ConstructionTraits<Smoother>::construct(cargs);
#else
smoother_.reset(Dune::Amg::ConstructionTraits<Smoother>::construct(cargs)); smoother_.reset(Dune::Amg::ConstructionTraits<Smoother>::construct(cargs));
#endif
} }
} }
@ -484,7 +488,7 @@ private:
Dune::LoopSolver<X> solver(const_cast<typename AMGType::Operator&>(op_), *sp, *prec, Dune::LoopSolver<X> solver(const_cast<typename AMGType::Operator&>(op_), *sp, *prec,
tolerance, maxit, verbosity); tolerance, maxit, verbosity);
solver.apply(x,b,res); solver.apply(x,b,res);
#else #else
if ( !amg_ ) if ( !amg_ )
{ {
Dune::LoopSolver<X> solver(const_cast<typename AMGType::Operator&>(op_), *sp, Dune::LoopSolver<X> solver(const_cast<typename AMGType::Operator&>(op_), *sp,
@ -529,7 +533,7 @@ private:
const CPRParameter* param_; const CPRParameter* param_;
X x_; X x_;
std::unique_ptr<AMGType> amg_; std::unique_ptr<AMGType> amg_;
std::unique_ptr<Smoother> smoother_; std::shared_ptr<Smoother> smoother_;
const typename AMGType::Operator& op_; const typename AMGType::Operator& op_;
Criterion crit_; Criterion crit_;
const Communication& comm_; const Communication& comm_;
@ -578,15 +582,18 @@ private:
}; };
template<class Smoother, class Operator, class Communication> template<class Smoother, class Operator, class Communication>
Smoother* constructSmoother(const Operator& op, std::shared_ptr< Smoother >
const typename Dune::Amg::SmootherTraits<Smoother>::Arguments& smargs, constructSmoother(const Operator& op,
const Communication& comm) const typename Dune::Amg::SmootherTraits<Smoother>::Arguments& smargs,
const Communication& comm)
{ {
typename Dune::Amg::ConstructionTraits<Smoother>::Arguments args; typename Dune::Amg::ConstructionTraits<Smoother>::Arguments args;
args.setMatrix(op.getmat()); args.setMatrix(op.getmat());
args.setComm(comm); args.setComm(comm);
args.setArgs(smargs); args.setArgs(smargs);
return Dune::Amg::ConstructionTraits<Smoother>::construct(args); // for DUNE < 2.7 ConstructionTraits<Smoother>::construct returns a raw
// pointer, therefore std::make_shared cannot be used here
return std::shared_ptr< Smoother > (Dune::Amg::ConstructionTraits<Smoother>::construct(args));
} }
template<class G, class C, class S> template<class G, class C, class S>
@ -755,7 +762,11 @@ public:
using CommunicationArgs = typename Dune::Amg::ConstructionTraits<Communication>::Arguments; using CommunicationArgs = typename Dune::Amg::ConstructionTraits<Communication>::Arguments;
CommunicationArgs commArgs(communication_->communicator(), communication_->getSolverCategory()); CommunicationArgs commArgs(communication_->communicator(), communication_->getSolverCategory());
#if DUNE_VERSION_NEWER(DUNE_ISTL, 2, 7)
coarseLevelCommunication_ = Dune::Amg::ConstructionTraits<Communication>::construct(commArgs);
#else
coarseLevelCommunication_.reset(Dune::Amg::ConstructionTraits<Communication>::construct(commArgs)); coarseLevelCommunication_.reset(Dune::Amg::ConstructionTraits<Communication>::construct(commArgs));
#endif
using Iterator = typename std::vector<bool>::iterator; using Iterator = typename std::vector<bool>::iterator;
using std::get; using std::get;
auto visitedMap = get(Dune::Amg::VertexVisitedTag(), *(get<1>(graphs))); auto visitedMap = get(Dune::Amg::VertexVisitedTag(), *(get<1>(graphs)));
@ -828,7 +839,11 @@ public:
this->rhs_.resize(this->coarseLevelMatrix_->N()); this->rhs_.resize(this->coarseLevelMatrix_->N());
using OperatorArgs = typename Dune::Amg::ConstructionTraits<CoarseOperator>::Arguments; using OperatorArgs = typename Dune::Amg::ConstructionTraits<CoarseOperator>::Arguments;
OperatorArgs oargs(*coarseLevelMatrix_, *coarseLevelCommunication_); OperatorArgs oargs(*coarseLevelMatrix_, *coarseLevelCommunication_);
#if DUNE_VERSION_NEWER(DUNE_ISTL, 2, 7)
this->operator_ = Dune::Amg::ConstructionTraits<CoarseOperator>::construct(oargs);
#else
this->operator_.reset(Dune::Amg::ConstructionTraits<CoarseOperator>::construct(oargs)); this->operator_.reset(Dune::Amg::ConstructionTraits<CoarseOperator>::construct(oargs));
#endif
} }
void calculateCoarseEntriesWithAggregatesMap(const Operator& fineOperator) void calculateCoarseEntriesWithAggregatesMap(const Operator& fineOperator)
@ -1032,13 +1047,13 @@ public:
weights_(weights), weights_(weights),
scaledMatrix_(Detail::scaleMatrixDRS(fineOperator, COMPONENT_INDEX, weights, param)), scaledMatrix_(Detail::scaleMatrixDRS(fineOperator, COMPONENT_INDEX, weights, param)),
scaledMatrixOperator_(Detail::createOperator(fineOperator, *scaledMatrix_, comm)), scaledMatrixOperator_(Detail::createOperator(fineOperator, *scaledMatrix_, comm)),
smoother_(Detail::constructSmoother<Smoother>(*scaledMatrixOperator_, smoother_( Detail::constructSmoother<Smoother>(*scaledMatrixOperator_, smargs, comm)),
smargs, comm)),
levelTransferPolicy_(criterion, comm, param.cpr_pressure_aggregation_), levelTransferPolicy_(criterion, comm, param.cpr_pressure_aggregation_),
coarseSolverPolicy_(&param, smargs, criterion), coarseSolverPolicy_(&param, smargs, criterion),
twoLevelMethod_(*scaledMatrixOperator_, smoother_, twoLevelMethod_(*scaledMatrixOperator_, smoother_,
levelTransferPolicy_, coarseSolverPolicy_, 0, 1) levelTransferPolicy_, coarseSolverPolicy_, 0, 1)
{} {
}
void pre(typename TwoLevelMethod::FineDomainType& x, void pre(typename TwoLevelMethod::FineDomainType& x,
typename TwoLevelMethod::FineRangeType& b) override typename TwoLevelMethod::FineRangeType& b) override

View File

@ -141,12 +141,12 @@ namespace Opm
const Communication& comm) const Communication& comm)
{ {
*scaledMatrix_ = *Detail::scaleMatrixDRS(fineOperator, COMPONENT_INDEX, weights_, param_); *scaledMatrix_ = *Detail::scaleMatrixDRS(fineOperator, COMPONENT_INDEX, weights_, param_);
smoother_.reset(Detail::constructSmoother<Smoother>(*scaledMatrixOperator_, smargs, comm)); smoother_ = Detail::constructSmoother<Smoother>(*scaledMatrixOperator_, smargs, comm);
twoLevelMethod_.updatePreconditioner(*scaledMatrixOperator_, twoLevelMethod_.updatePreconditioner(*scaledMatrixOperator_,
smoother_, smoother_,
coarseSolverPolicy_); coarseSolverPolicy_);
} }
void pre(typename TwoLevelMethod::FineDomainType& x, void pre(typename TwoLevelMethod::FineDomainType& x,
typename TwoLevelMethod::FineRangeType& b) override typename TwoLevelMethod::FineRangeType& b) override
{ {

View File

@ -176,7 +176,7 @@ public:
protected: protected:
const matrix_type& A_ ; const matrix_type& A_ ;
const matrix_type& A_for_precond_ ; const matrix_type& A_for_precond_ ;
const WellModel& wellMod_; const WellModel& wellMod_;
std::shared_ptr< communication_type > comm_; std::shared_ptr< communication_type > comm_;
}; };
@ -246,7 +246,7 @@ protected:
void scaleSystem() void scaleSystem()
{ {
const bool matrix_cont_added = EWOMS_GET_PARAM(TypeTag, bool, MatrixAddWellContributions); const bool matrix_cont_added = EWOMS_GET_PARAM(TypeTag, bool, MatrixAddWellContributions);
if (matrix_cont_added) { if (matrix_cont_added) {
bool form_cpr = true; bool form_cpr = true;
if (parameters_.system_strategy_ == "quasiimpes") { if (parameters_.system_strategy_ == "quasiimpes") {
@ -663,7 +663,7 @@ protected:
// conservation equations, ignoring all other terms. // conservation equations, ignoring all other terms.
Vector getStorageWeights() const Vector getStorageWeights() const
{ {
Vector weights(rhs_->size()); Vector weights(rhs_->size());
BlockVector rhs(0.0); BlockVector rhs(0.0);
rhs[pressureVarIndex] = 1.0; rhs[pressureVarIndex] = 1.0;
int index = 0; int index = 0;
@ -830,7 +830,7 @@ protected:
} }
} }
} }
static void multBlocksVector(Vector& ebosResid_cp, const MatrixBlockType& leftTrans) static void multBlocksVector(Vector& ebosResid_cp, const MatrixBlockType& leftTrans)
{ {
for (auto& bvec : ebosResid_cp) { for (auto& bvec : ebosResid_cp) {

View File

@ -127,19 +127,30 @@ struct ConstructionTraits<Opm::ParallelOverlappingILU0<Matrix,Domain,Range,Paral
{ {
typedef Opm::ParallelOverlappingILU0<Matrix,Domain,Range,ParallelInfo> T; typedef Opm::ParallelOverlappingILU0<Matrix,Domain,Range,ParallelInfo> T;
typedef DefaultParallelConstructionArgs<T,ParallelInfo> Arguments; typedef DefaultParallelConstructionArgs<T,ParallelInfo> Arguments;
static inline Opm::ParallelOverlappingILU0<Matrix,Domain,Range,ParallelInfo>* construct(Arguments& args)
#if DUNE_VERSION_NEWER(DUNE_ISTL, 2, 7)
typedef std::shared_ptr< T > ParallelOverlappingILU0Pointer;
#else
typedef T* ParallelOverlappingILU0Pointer;
#endif
static inline ParallelOverlappingILU0Pointer construct(Arguments& args)
{ {
return new T(args.getMatrix(), return ParallelOverlappingILU0Pointer(
args.getComm(), new T(args.getMatrix(),
args.getArgs().getN(), args.getComm(),
args.getArgs().relaxationFactor, args.getArgs().getN(),
args.getArgs().getMilu()); args.getArgs().relaxationFactor,
args.getArgs().getMilu()) );
} }
#if ! DUNE_VERSION_NEWER(DUNE_ISTL, 2, 7)
// this method is not needed anymore in 2.7 since std::shared_ptr is used
static inline void deconstruct(T* bp) static inline void deconstruct(T* bp)
{ {
delete bp; delete bp;
} }
#endif
}; };
@ -400,7 +411,7 @@ namespace Opm
newRow[ordering[col.index()]] = *col; newRow[ordering[col.index()]] = *col;
} }
} }
// call decomposition on pattern // call decomposition on pattern
switch ( milu ) switch ( milu )
{ {
case MILU_VARIANT::MILU_1: case MILU_VARIANT::MILU_1:

View File

@ -55,16 +55,33 @@ struct ConstructionTraits<Opm::ParallelRestrictedOverlappingSchwarz<Range,
typedef ConstructionTraits<SeqPreconditioner> SeqConstructionTraits; typedef ConstructionTraits<SeqPreconditioner> SeqConstructionTraits;
/// \brief Construct a parallel restricted overlapping schwarz preconditioner. /// \brief Construct a parallel restricted overlapping schwarz preconditioner.
static inline Opm::ParallelRestrictedOverlappingSchwarz<Range, #if DUNE_VERSION_NEWER(DUNE_ISTL, 2, 7)
Domain, typedef std::shared_ptr< Opm::ParallelRestrictedOverlappingSchwarz<Range,
ParallelInfo, Domain,
SeqPreconditioner>* ParallelInfo,
SeqPreconditioner> > ParallelRestrictedOverlappingSchwarzPointer;
#else
typedef Opm::ParallelRestrictedOverlappingSchwarz<Range,
Domain,
ParallelInfo,
SeqPreconditioner>* ParallelRestrictedOverlappingSchwarzPointer;
#endif
static inline ParallelRestrictedOverlappingSchwarzPointer
construct(Arguments& args) construct(Arguments& args)
{ {
return new Opm::ParallelRestrictedOverlappingSchwarz return
<Range,Domain,ParallelInfo,SeqPreconditioner>(*SeqConstructionTraits #if DUNE_VERSION_NEWER(DUNE_ISTL, 2, 7)
::construct(args), std::make_shared(
args.getComm()); #endif
new Opm::ParallelRestrictedOverlappingSchwarz
<Range,Domain,ParallelInfo,SeqPreconditioner>(*SeqConstructionTraits ::construct(args),
args.getComm())
#if DUNE_VERSION_NEWER(DUNE_ISTL, 2, 7)
);
#else
;
#endif
} }
/// \brief Deconstruct and free a parallel restricted overlapping schwarz preconditioner. /// \brief Deconstruct and free a parallel restricted overlapping schwarz preconditioner.

View File

@ -290,7 +290,7 @@ namespace Dune
/** /**
* @brief Update the coarse solver and the hierarchies. * @brief Update the coarse solver and the hierarchies.
*/ */
template<class C> template<class C>
void updateSolver(C& criterion, Operator& /* matrix */, const PI& pinfo); void updateSolver(C& criterion, Operator& /* matrix */, const PI& pinfo);
/** /**
@ -306,9 +306,26 @@ namespace Dune
* @param matrix The fine level matrix operator. * @param matrix The fine level matrix operator.
* @param pinfo The fine level parallel information. * @param pinfo The fine level parallel information.
*/ */
#if DUNE_VERSION_NEWER( DUNE_ISTL, 2, 7 )
template<class C>
void createHierarchies(C& criterion, Operator& matrix,
const PI& pinfo)
{
// create shared_ptr with empty deleter
std::shared_ptr< Operator > op( &matrix, []( Operator* ) {});
std::shared_ptr< PI > pifo( const_cast< PI* > (&pinfo), []( PI * ) {});
createHierarchies( criterion, op, pifo);
}
template<class C>
void createHierarchies(C& criterion, std::shared_ptr< Operator > matrix,
std::shared_ptr< PI > pinfo );
#else
template<class C> template<class C>
void createHierarchies(C& criterion, Operator& matrix, void createHierarchies(C& criterion, Operator& matrix,
const PI& pinfo); const PI& pinfo);
#endif
void setupCoarseSolver(); void setupCoarseSolver();
@ -397,11 +414,11 @@ namespace Dune
/** @brief The solver of the coarsest level. */ /** @brief The solver of the coarsest level. */
std::shared_ptr<CoarseSolver> solver_; std::shared_ptr<CoarseSolver> solver_;
/** @brief The right hand side of our problem. */ /** @brief The right hand side of our problem. */
Hierarchy<Range,A>* rhs_; std::shared_ptr< Hierarchy<Range,A> > rhs_;
/** @brief The left approximate solution of our problem. */ /** @brief The left approximate solution of our problem. */
Hierarchy<Domain,A>* lhs_; std::shared_ptr< Hierarchy<Domain,A> > lhs_;
/** @brief The total update for the outer solver. */ /** @brief The total update for the outer solver. */
Hierarchy<Domain,A>* update_; std::shared_ptr< Hierarchy<Domain,A> > update_;
/** @brief The type of the scalar product for the coarse solver. */ /** @brief The type of the scalar product for the coarse solver. */
using ScalarProduct = Dune::ScalarProduct<X>; using ScalarProduct = Dune::ScalarProduct<X>;
/** @brief Scalar product on the coarse level. */ /** @brief Scalar product on the coarse level. */
@ -440,11 +457,11 @@ namespace Dune
verbosity_(amg.verbosity_) verbosity_(amg.verbosity_)
{ {
if(amg.rhs_) if(amg.rhs_)
rhs_=new Hierarchy<Range,A>(*amg.rhs_); rhs_.reset( new Hierarchy<Range,A>(*amg.rhs_) );
if(amg.lhs_) if(amg.lhs_)
lhs_=new Hierarchy<Domain,A>(*amg.lhs_); lhs_.reset( new Hierarchy<Domain,A>(*amg.lhs_) );
if(amg.update_) if(amg.update_)
update_=new Hierarchy<Domain,A>(*amg.update_); update_.reset( new Hierarchy<Domain,A>(*amg.update_) );
} }
template<class M, class X, class S, class PI, class A> template<class M, class X, class S, class PI, class A>
@ -495,7 +512,6 @@ namespace Dune
createHierarchies(criterion, const_cast<Operator&>(matrix), pinfo); createHierarchies(criterion, const_cast<Operator&>(matrix), pinfo);
} }
template<class M, class X, class S, class PI, class A> template<class M, class X, class S, class PI, class A>
AMGCPR<M,X,S,PI,A>::~AMGCPR() AMGCPR<M,X,S,PI,A>::~AMGCPR()
{ {
@ -505,15 +521,6 @@ namespace Dune
if(coarseSmoother_) if(coarseSmoother_)
coarseSmoother_.reset(); coarseSmoother_.reset();
} }
if(lhs_)
delete lhs_;
lhs_=nullptr;
if(update_)
delete update_;
update_=nullptr;
if(rhs_)
delete rhs_;
rhs_=nullptr;
} }
template<class M, class X, class S, class PI, class A> template<class M, class X, class S, class PI, class A>
@ -539,8 +546,12 @@ namespace Dune
template<class M, class X, class S, class PI, class A> template<class M, class X, class S, class PI, class A>
template<class C> template<class C>
void AMGCPR<M,X,S,PI,A>::createHierarchies(C& criterion, Operator& matrix, #if DUNE_VERSION_NEWER( DUNE_ISTL, 2, 7)
const PI& pinfo) void AMGCPR<M,X,S,PI,A>::createHierarchies(C& criterion, std::shared_ptr< Operator > matrix,
std::shared_ptr< PI > pinfo )
#else
void AMGCPR<M,X,S,PI,A>::createHierarchies(C& criterion, Operator& matrix, const PI& pinfo )
#endif
{ {
Timer watch; Timer watch;
matrices_.reset(new OperatorHierarchy(matrix, pinfo)); matrices_.reset(new OperatorHierarchy(matrix, pinfo));
@ -581,7 +592,11 @@ namespace Dune
cargs.setComm(*matrices_->parallelInformation().coarsest()); cargs.setComm(*matrices_->parallelInformation().coarsest());
} }
#if DUNE_VERSION_NEWER(DUNE_ISTL, 2, 7)
coarseSmoother_ = ConstructionTraits<Smoother>::construct(cargs);
#else
coarseSmoother_.reset(ConstructionTraits<Smoother>::construct(cargs)); coarseSmoother_.reset(ConstructionTraits<Smoother>::construct(cargs));
#endif
#if DUNE_VERSION_NEWER(DUNE_ISTL, 2, 6) #if DUNE_VERSION_NEWER(DUNE_ISTL, 2, 6)
scalarProduct_ = createScalarProduct<X>(cargs.getComm(),category()); scalarProduct_ = createScalarProduct<X>(cargs.getComm(),category());
@ -704,18 +719,24 @@ namespace Dune
else else
// No smoother to make x consistent! Do it by hand // No smoother to make x consistent! Do it by hand
matrices_->parallelInformation().coarsest()->copyOwnerToAll(x,x); matrices_->parallelInformation().coarsest()->copyOwnerToAll(x,x);
Range* copy = new Range(b);
if(rhs_)
delete rhs_; #if DUNE_VERSION_NEWER( DUNE_ISTL, 2, 7)
rhs_ = new Hierarchy<Range,A>(copy); typedef std::shared_ptr< Range > RangePtr ;
Domain* dcopy = new Domain(x); typedef std::shared_ptr< Domain > DomainPtr;
if(lhs_) #else
delete lhs_; typedef Range* RangePtr;
lhs_ = new Hierarchy<Domain,A>(dcopy); typedef Domain* DomainPtr;
dcopy = new Domain(x); #endif
if(update_)
delete update_; // Hierarchy takes ownership of pointers
update_ = new Hierarchy<Domain,A>(dcopy); RangePtr copy( new Range(b) );
rhs_.reset( new Hierarchy<Range,A>(copy) );
DomainPtr dcopy( new Domain(x) );
lhs_.reset( new Hierarchy<Domain,A>(dcopy) );
DomainPtr dcopy2( new Domain(x) );
update_.reset( new Hierarchy<Domain,A>(dcopy2) );
matrices_->coarsenVector(*rhs_); matrices_->coarsenVector(*rhs_);
matrices_->coarsenVector(*lhs_); matrices_->coarsenVector(*lhs_);
matrices_->coarsenVector(*update_); matrices_->coarsenVector(*update_);
@ -1024,14 +1045,11 @@ namespace Dune
smoother->post(*lhs); smoother->post(*lhs);
} }
//delete &(*lhs_->finest()); //delete &(*lhs_->finest());
delete lhs_; lhs_.reset();
lhs_=nullptr;
//delete &(*update_->finest()); //delete &(*update_->finest());
delete update_; update_.reset();
update_=nullptr;
//delete &(*rhs_->finest()); //delete &(*rhs_->finest());
delete rhs_; rhs_.reset();
rhs_=nullptr;
} }
template<class M, class X, class S, class PI, class A> template<class M, class X, class S, class PI, class A>