diff --git a/dune/istl/paamg/amg.hh b/dune/istl/paamg/amg.hh
index 57b84751829c4a91f1162372388b90979a7182a7..016aceb970b6c22546a52e72dd5a967bf0d37ad5 100644
--- a/dune/istl/paamg/amg.hh
+++ b/dune/istl/paamg/amg.hh
@@ -208,7 +208,8 @@ namespace Dune
* @param pinfo The fine level parallel information.
*/
template<class C>
- void createHierarchies(C& criterion, Operator& matrix,
+ void createHierarchies(C& criterion,
+ const std::shared_ptr<const Operator>& matrixptr,
const PI& pinfo);
/**
* @brief A struct that holds the context of the current level.
@@ -383,12 +384,12 @@ namespace Dune
// TODO: reestablish compile time checks.
//static_assert(static_cast<int>(PI::category)==static_cast<int>(S::category),
// "Matrix and Solver must match in terms of category!");
- createHierarchies(criterion, const_cast<Operator&>(matrix), pinfo);
+ auto matrixptr = stackobject_to_shared_ptr(matrix);
+ createHierarchies(criterion, matrixptr, pinfo);
}
-
template<class M, class X, class S, class PI, class A>
- AMG<M,X,S,PI,A>::AMG(std::shared_ptr<const Operator> matrix,
+ AMG<M,X,S,PI,A>::AMG(std::shared_ptr<const Operator> matrixptr,
const ParameterTree& configuration,
const ParallelInformation& pinfo) :
smoothers_(new Hierarchy<Smoother,A>),
@@ -438,7 +439,7 @@ namespace Dune
verbosity_ = configuration.get("verbosity",2);
criterion.setDebugLevel (verbosity_);
- createHierarchies(criterion, const_cast<Operator&>(*matrix), pinfo);
+ createHierarchies(criterion, matrixptr, pinfo);
}
template<class M, class X, class S, class PI, class A>
@@ -528,11 +529,14 @@ namespace Dune
template<class M, class X, class S, class PI, class A>
template<class C>
- void AMG<M,X,S,PI,A>::createHierarchies(C& criterion, Operator& matrix,
- const PI& pinfo)
+ void AMG<M,X,S,PI,A>::createHierarchies(C& criterion,
+ const std::shared_ptr<const Operator>& matrixptr,
+ const PI& pinfo)
{
Timer watch;
- matrices_.reset(new OperatorHierarchy(matrix, pinfo));
+ matrices_ = std::make_shared<OperatorHierarchy>(
+ std::const_pointer_cast<Operator>(matrixptr),
+ stackobject_to_shared_ptr(const_cast<PI&>(pinfo)));
matrices_->template build<NegateSet<typename PI::OwnerSet> >(criterion);
@@ -648,18 +652,15 @@ namespace Dune
else
// No smoother to make x consistent! Do it by hand
matrices_->parallelInformation().coarsest()->copyOwnerToAll(x,x);
- Range* copy = new Range(b);
if(rhs_)
delete rhs_;
- rhs_ = new Hierarchy<Range,A>(copy);
- Domain* dcopy = new Domain(x);
+ rhs_ = new Hierarchy<Range,A>(std::make_shared<Range>(b));
if(lhs_)
delete lhs_;
- lhs_ = new Hierarchy<Domain,A>(dcopy);
- dcopy = new Domain(x);
+ lhs_ = new Hierarchy<Domain,A>(std::make_shared<Domain>(x));
if(update_)
delete update_;
- update_ = new Hierarchy<Domain,A>(dcopy);
+ update_ = new Hierarchy<Domain,A>(std::make_shared<Domain>(x));
matrices_->coarsenVector(*rhs_);
matrices_->coarsenVector(*lhs_);
matrices_->coarsenVector(*update_);
diff --git a/dune/istl/paamg/fastamg.hh b/dune/istl/paamg/fastamg.hh
index 5b712656d2e6d93fc0cfaa74dbc4b213b2f3d05c..a89b7e6d01203aef0389ad7abf45711076bed335 100644
--- a/dune/istl/paamg/fastamg.hh
+++ b/dune/istl/paamg/fastamg.hh
@@ -167,7 +167,8 @@ namespace Dune
* @param pinfo The fine level parallel information.
*/
template<class C>
- void createHierarchies(C& criterion, Operator& matrix,
+ void createHierarchies(C& criterion,
+ const std::shared_ptr<const Operator>& matrixptr,
const PI& pinfo);
/**
@@ -342,7 +343,8 @@ namespace Dune
// TODO: reestablish compile time checks.
//static_assert(static_cast<int>(PI::category)==static_cast<int>(S::category),
// "Matrix and Solver must match in terms of category!");
- createHierarchies(criterion, const_cast<Operator&>(matrix), pinfo);
+ auto matrixptr = stackobject_to_shared_ptr(matrix);
+ createHierarchies(criterion, matrixptr, pinfo);
}
template<class M, class X, class PI, class A>
@@ -367,11 +369,14 @@ namespace Dune
template<class M, class X, class PI, class A>
template<class C>
- void FastAMG<M,X,PI,A>::createHierarchies(C& criterion, Operator& matrix,
- const PI& pinfo)
+ void FastAMG<M,X,PI,A>::createHierarchies(C& criterion,
+ const std::shared_ptr<const Operator>& matrixptr,
+ const PI& pinfo)
{
Timer watch;
- matrices_.reset(new OperatorHierarchy(matrix, pinfo));
+ matrices_ = std::make_shared<OperatorHierarchy>(
+ std::const_pointer_cast<Operator>(matrixptr),
+ stackobject_to_shared_ptr(const_cast<PI&>(pinfo)));
matrices_->template build<NegateSet<typename PI::OwnerSet> >(criterion);
@@ -483,16 +488,13 @@ namespace Dune
watch1.reset();
// No smoother to make x consistent! Do it by hand
matrices_->parallelInformation().coarsest()->copyOwnerToAll(x,x);
- Range* copy = new Range(b);
if(rhs_)
delete rhs_;
- rhs_ = new Hierarchy<Range,A>(copy);
- Domain* dcopy = new Domain(x);
+ rhs_ = new Hierarchy<Range,A>(std::make_shared<Range>(b));
if(lhs_)
delete lhs_;
- lhs_ = new Hierarchy<Domain,A>(dcopy);
- dcopy = new Domain(x);
- residual_ = new Hierarchy<Domain,A>(dcopy);
+ lhs_ = new Hierarchy<Domain,A>(std::make_shared<Domain>(x));
+ residual_ = new Hierarchy<Domain,A>(std::make_shared<Domain>(x));
matrices_->coarsenVector(*rhs_);
matrices_->coarsenVector(*lhs_);
matrices_->coarsenVector(*residual_);
diff --git a/dune/istl/paamg/hierarchy.hh b/dune/istl/paamg/hierarchy.hh
index 5cff788bca09d0267f488406217eab711deeb60f..e0cedb939bdf2c3a1686f5203a4ed7927d39f746 100644
--- a/dune/istl/paamg/hierarchy.hh
+++ b/dune/istl/paamg/hierarchy.hh
@@ -111,19 +111,27 @@ namespace Dune
typedef typename ConstructionTraits<T>::Arguments Arguments;
+ /**
+ * @brief Construct a new hierarchy.
+ * @param first std::shared_ptr to the first element in the hierarchy.
+ */
+ Hierarchy(const std::shared_ptr<MemberType> & first);
+
/**
* @brief Construct a new hierarchy.
* @param first The first element in the hierarchy.
+ * @todo remove in 2.6
*/
- Hierarchy(MemberType& first);
+ Hierarchy(MemberType& first) DUNE_DEPRECATED_MSG("Use the constructor taking a shared_ptr");
/**
* @brief Construct a new hierarchy.
* @param first Pointer to the first element in the hierarchy.
* @warning Hierarchy will be responsible for the memory
* management of the pointer.
+ * @todo remove in 2.6
*/
- Hierarchy(MemberType* first);
+ Hierarchy(MemberType* first) DUNE_DEPRECATED_MSG("Use the constructor taking a shared_ptr");
/**
* @brief Construct a new empty hierarchy.
@@ -294,12 +302,16 @@ namespace Dune
~Hierarchy();
private:
+ /** @brief fix memory management of the finest element in the hierarchy
+
+ This object is passed in from outside, so that we have to
+ deal with shared ownership.
+ */
+ std::shared_ptr<MemberType> originalFinest_;
/** @brief The finest element in the hierarchy. */
Element* finest_;
/** @brief The coarsest element in the hierarchy. */
Element* coarsest_;
- /** @brief Whether the first element was not allocated by us. */
- Element* nonAllocated_;
/** @brief The allocator for the list elements. */
Allocator allocator_;
/** @brief The number of levels in the hierarchy. */
@@ -357,9 +369,8 @@ namespace Dune
* @param fineMatrix The matrix to coarsen.
* @param pinfo The information about the parallel data decomposition at the first level.
*/
- MatrixHierarchy(const MatrixOperator& fineMatrix,
- const ParallelInformation& pinfo=ParallelInformation());
-
+ MatrixHierarchy(std::shared_ptr<MatrixOperator> fineMatrix,
+ std::shared_ptr<ParallelInformation> pinfo = std::make_shared<ParallelInformation>());
~MatrixHierarchy();
@@ -657,10 +668,10 @@ namespace Dune
}
template<class M, class IS, class A>
- MatrixHierarchy<M,IS,A>::MatrixHierarchy(const MatrixOperator& fineOperator,
- const ParallelInformation& pinfo)
- : matrices_(const_cast<MatrixOperator&>(fineOperator)),
- parallelInformation_(const_cast<ParallelInformation&>(pinfo))
+ MatrixHierarchy<M,IS,A>::MatrixHierarchy(std::shared_ptr<MatrixOperator> fineMatrix,
+ std::shared_ptr<ParallelInformation> pinfo)
+ : matrices_(fineMatrix),
+ parallelInformation_(pinfo)
{
#warning reenable category checks
// static_assert((static_cast<int>(MatrixOperator::category) ==
@@ -1242,17 +1253,28 @@ namespace Dune
template<class T, class A>
Hierarchy<T,A>::Hierarchy()
- : finest_(0), coarsest_(0), nonAllocated_(0), allocator_(), levels_(0)
+ : finest_(0), coarsest_(0), allocator_(), levels_(0)
{}
+ template<class T, class A>
+ Hierarchy<T,A>::Hierarchy(const std::shared_ptr<MemberType> & first)
+ : originalFinest_(first), allocator_()
+ {
+ finest_ = allocator_.allocate(1,0);
+ finest_->element_ = originalFinest_.get();
+ finest_->redistributed_ = nullptr;
+ coarsest_ = finest_;
+ coarsest_->coarser_ = coarsest_->finer_ = nullptr;
+ levels_ = 1;
+ }
+
template<class T, class A>
Hierarchy<T,A>::Hierarchy(MemberType& first)
- : allocator_()
+ : originalFinest_(stackobject_to_shared_ptr(first)), allocator_()
{
finest_ = allocator_.allocate(1,0);
- finest_->element_ = &first;
+ finest_->element_ = originalFinest_.get();
finest_->redistributed_ = nullptr;
- nonAllocated_ = finest_;
coarsest_ = finest_;
coarsest_->coarser_ = coarsest_->finer_ = nullptr;
levels_ = 1;
@@ -1260,12 +1282,11 @@ namespace Dune
template<class T, class A>
Hierarchy<T,A>::Hierarchy(MemberType* first)
- : allocator_()
+ : originalFinest_(first), allocator_()
{
finest_ = allocator_.allocate(1,0);
- finest_->element_ = first;
+ finest_->element_ = originalFinest_.get();
finest_->redistributed_ = nullptr;
- nonAllocated_ = nullptr;
coarsest_ = finest_;
coarsest_->coarser_ = coarsest_->finer_ = nullptr;
levels_ = 1;
@@ -1276,7 +1297,9 @@ namespace Dune
while(coarsest_) {
Element* current = coarsest_;
coarsest_ = coarsest_->finer_;
- if(current != nonAllocated_) {
+ // we changed the internal behaviour
+ // now the finest level is _always_ managedd by a shared_ptr
+ if(current != finest_) {
if(current->redistributed_)
ConstructionTraits<T>::deconstruct(current->redistributed_);
ConstructionTraits<T>::deconstruct(current->element_);
@@ -1289,12 +1312,12 @@ namespace Dune
template<class T, class A>
Hierarchy<T,A>::Hierarchy(const Hierarchy& other)
- : nonAllocated_(), allocator_(other.allocator_),
+ : allocator_(other.allocator_),
levels_(other.levels_)
{
if(!other.finest_)
{
- finest_=coarsest_=nonAllocated_=nullptr;
+ finest_=coarsest_=nullptr;
return;
}
finest_=allocator_.allocate(1,0);
@@ -1339,9 +1362,14 @@ namespace Dune
void Hierarchy<T,A>::addCoarser(Arguments& args)
{
if(!coarsest_) {
+ // we have no levels at all...
assert(!finest_);
+ // allocate into the shared_ptr
+ originalFinest_ =
+ std::shared_ptr<MemberType>(
+ ConstructionTraits<MemberType>::construct(args));
coarsest_ = allocator_.allocate(1,0);
- coarsest_->element_ = ConstructionTraits<MemberType>::construct(args);
+ coarsest_->element_ = originalFinest_.get();
finest_ = coarsest_;
coarsest_->finer_ = nullptr;
}else{
@@ -1360,9 +1388,14 @@ namespace Dune
void Hierarchy<T,A>::addFiner(Arguments& args)
{
if(!finest_) {
+ // we have no levels at all...
assert(!coarsest_);
+ // allocate into the shared_ptr
+ originalFinest_ =
+ std::shared_ptr<MemberType>(
+ ConstructionTraits<MemberType>::construct(args));
finest_ = allocator_.allocate(1,0);
- finest_->element = ConstructionTraits<T>::construct(args);
+ finest_->element = originalFinest_.get();
coarsest_ = finest_;
coarsest_->coarser_ = coarsest_->finer_ = nullptr;
}else{
diff --git a/dune/istl/paamg/test/hierarchytest.cc b/dune/istl/paamg/test/hierarchytest.cc
index 30d6cda2c9d826591b6b407b863c91f7db70546d..3182772df35fb6f4175f01e94a6c5409274a7067 100644
--- a/dune/istl/paamg/test/hierarchytest.cc
+++ b/dune/istl/paamg/test/hierarchytest.cc
@@ -43,8 +43,11 @@ void testHierarchy(int N)
typedef Dune::Amg::Hierarchy<Vector> VHierarchy;
Operator op(mat, pinfo);
- Hierarchy hierarchy(op, pinfo);
- VHierarchy vh(b);
+ Hierarchy hierarchy(
+ Dune::stackobject_to_shared_ptr(op),
+ Dune::stackobject_to_shared_ptr(pinfo));
+ VHierarchy vh(
+ Dune::stackobject_to_shared_ptr(b));
typedef Dune::Amg::CoarsenCriterion<Dune::Amg::SymmetricCriterion<BCRSMat,Dune::Amg::FirstDiagonal> >
Criterion;