Changed raw pointer to shared_ptr in DOFVector and DOFMatrix

examples/CMakeLists.txt

@@ -29,6 +29,6 @@ add_dependencies(examples
   cahn_hilliard.2d)
 
 if (dune-foamgrid_FOUND)
-  add_amdis_executable(NAME surface.2d SOURCES surface.cc DIM 2 DOW 3)
+  add_amdis_executable(NAME surface.2d SOURCES surface.cc DIM 2 DOW 3 ALBERTA_GRID)
   add_dependencies(examples surface.2d)
 endif ()
\ No newline at end of file

examples/surface.cc

+#include "config.h"
 #include <amdis/AMDiS.hpp>
 #include <amdis/Integrate.hpp>
 #include <amdis/LocalOperators.hpp>
@@ -14,7 +15,7 @@ using namespace AMDiS;
 using namespace Dune::Indices;
 
 using Grid = Dune::FoamGrid<2,3>;
-using Basis = LagrangeBasis<typename Grid::LeafGridView, 1>;
+using Basis = LagrangeBasis<Grid, 1>;
 
 struct UnitRadius
 {
@@ -22,6 +23,9 @@ struct UnitRadius
   static double eval(T&&) { return 1.0; }
 };
 
+#if !HAVE_ALBERTA
+#error "Need Alberta!!"
+#endif
 
 // solve the equation -laplace(u) - k^2 u = f on the sphere
 int main(int argc, char** argv)

src/amdis/DirichletBC.inc.hpp

@@ -35,7 +35,7 @@ fillBoundaryCondition(Mat& matrix, Sol& solution, Rhs& rhs, RN const& rowNode, R
 
   Dune::Hybrid::ifElse(std::is_same<RangeType_t<RN>, R>{},
     [&](auto id) {
-      auto subBasis = Dune::Functions::subspaceBasis(rhs.basis(), rowTreePath);
+      auto subBasis = Dune::Functions::subspaceBasis(*rhs.basis(), rowTreePath);
       Dune::Functions::interpolate(subBasis, rhs, values_, dirichletNodes_);
     });
 
@@ -45,7 +45,7 @@ fillBoundaryCondition(Mat& matrix, Sol& solution, Rhs& rhs, RN const& rowNode, R
 
   Dune::Hybrid::ifElse(std::is_same<RangeType_t<CN>, R>{},
     [&](auto id) {
-      auto subBasis = Dune::Functions::subspaceBasis(solution.basis(), colTreePath);
+      auto subBasis = Dune::Functions::subspaceBasis(*solution.basis(), colTreePath);
       Dune::Functions::interpolate(subBasis, solution, values_, dirichletNodes_);
     });
 }

src/amdis/FileWriter.hpp

@@ -107,7 +107,7 @@ namespace AMDiS
   protected:
     GridView const& gridView() const
     {
-      return discreteFct_.basis().gridView();
+      return discreteFct_.basis()->gridView();
     }
 
   private:

src/amdis/PeriodicBC.inc.hpp

@@ -241,7 +241,7 @@ fillBoundaryCondition(Mat& matrix, Sol& solution, Rhs& rhs, RN const& rowNode, R
     rhs[a.second] += rhs[a.first];
     solution[a.second] = solution[a.first];
   }
-  Dune::Functions::interpolate(rhs.basis(), rhs, [](auto const&) { return 0.0; }, periodicNodes_);
+  Dune::Functions::interpolate(*rhs.basis(), rhs, [](auto const&) { return 0.0; }, periodicNodes_);
 }
 
 } // end namespace AMDiS

src/amdis/ProblemInstat.inc.hpp

@@ -45,7 +45,7 @@ void ProblemInstat<Traits>::createUhOld()
   if (oldSolution_)
     warning("oldSolution already created\n");
   else // create oldSolution
-    oldSolution_.reset(new SystemVector(*problemStat_->globalBasis(), DataTransferOperation::INTERPOLATE));
+    oldSolution_.reset(new SystemVector(problemStat_->globalBasis(), DataTransferOperation::INTERPOLATE));
 }
 
 

src/amdis/ProblemStat.inc.hpp

@@ -174,9 +174,9 @@ void ProblemStat<Traits>::initGlobalBasis()
 template <class Traits>
 void ProblemStat<Traits>::createMatricesAndVectors()
 {
-  systemMatrix_ = std::make_shared<SystemMatrix>(*globalBasis_, *globalBasis_);
-  solution_ = std::make_shared<SystemVector>(*globalBasis_, DataTransferOperation::INTERPOLATE);
-  rhs_ = std::make_shared<SystemVector>(*globalBasis_, DataTransferOperation::NO_OPERATION);
+  systemMatrix_ = std::make_shared<SystemMatrix>(globalBasis_, globalBasis_);
+  solution_ = std::make_shared<SystemVector>(globalBasis_, DataTransferOperation::INTERPOLATE);
+  rhs_ = std::make_shared<SystemVector>(globalBasis_, DataTransferOperation::NO_OPERATION);
 
   auto localView = globalBasis_->localView();
   for_each_node(localView.tree(), [&,this](auto const& node, auto treePath)

src/amdis/common/TypeTraits.hpp

@@ -31,11 +31,29 @@ namespace AMDiS
       using type = remove_cvref_t<T>;
     };
 
+    template <class T>
+    struct UnderlyingType<T*>
+    {
+      using type = remove_cvref_t<T>;
+    };
+
     template <class T>
     struct UnderlyingType<std::reference_wrapper<T>>
     {
       using type = remove_cvref_t<T>;
     };
+
+    template <class T>
+    struct UnderlyingType<std::shared_ptr<T>>
+    {
+      using type = remove_cvref_t<T>;
+    };
+
+    template <class T>
+    struct UnderlyingType<std::unique_ptr<T>>
+    {
+      using type = remove_cvref_t<T>;
+    };
   }
 
   template <class T>

src/amdis/gridfunctions/DOFVectorView.hpp

@@ -25,11 +25,6 @@ namespace AMDiS
     using TreePath = TP;
 
   public:
-    /// Constructor forwards to the treePath constructor, with empty TreePath
-    DOFVectorView(DOFVector<GB,VT>& dofVector)
-      : DOFVectorView{dofVector, Dune::TypeTree::hybridTreePath()}
-    {}
-
     /// Constructor. Stores a pointer to the mutable `dofvector`.
     template <class PreTreePath>
     DOFVectorView(DOFVector<GB,VT>& dofVector, PreTreePath const& preTreePath)
@@ -37,6 +32,11 @@ namespace AMDiS
       , mutableDofVector_(&dofVector)
     {}
 
+    /// Constructor forwards to the treePath constructor, with empty TreePath
+    DOFVectorView(DOFVector<GB,VT>& dofVector)
+      : DOFVectorView(dofVector, Dune::TypeTree::hybridTreePath())
+    {}
+
   public:
     /// \brief Interpolation of GridFunction to DOFVector, assuming that there is no
     /// reference to this DOFVector in the expression.
@@ -50,7 +50,7 @@ namespace AMDiS
     template <class Expr, class Tag = tag::average>
     void interpolate_noalias(Expr&& expr, Tag strategy = {})
     {
-      auto const& basis = this->basis();
+      auto const& basis = *this->basis();
       auto const& treePath = this->treePath();
 
       auto&& gridFct = makeGridFunction(FWD(expr), basis.gridView());

src/amdis/gridfunctions/DiscreteFunction.hpp

@@ -66,17 +66,17 @@ namespace AMDiS
     class LocalFunction;
 
   public:
-    /// Constructor forwards to the treePath constructor, with empty TreePath
-    DiscreteFunction(DOFVector<GB,VT> const& dofVector)
-      : DiscreteFunction{dofVector, Dune::TypeTree::hybridTreePath()}
-    {}
-
     /// Constructor. Stores a pointer to the dofVector and a copy of the treePath.
     DiscreteFunction(DOFVector<GB,VT> const& dofVector, TP const& treePath)
       : dofVector_(&dofVector)
       , treePath_(treePath)
-      , entitySet_(dofVector.basis().gridView())
-      , nodeToRangeEntry_(Dune::Functions::makeDefaultNodeToRangeMap(dofVector.basis(), treePath))
+      , entitySet_(dofVector.basis()->gridView())
+      , nodeToRangeEntry_(Dune::Functions::makeDefaultNodeToRangeMap(*dofVector.basis(), treePath))
+    {}
+
+    /// Constructor forwards to the treePath constructor, with empty TreePath
+    DiscreteFunction(DOFVector<GB,VT> const& dofVector)
+      : DiscreteFunction(dofVector, Dune::TypeTree::hybridTreePath())
     {}
 
     /// \brief Evaluate DiscreteFunction in global coordinates. NOTE: expensive
@@ -96,7 +96,7 @@ namespace AMDiS
 
   public:
     /// \brief Return global basis bound to the DOFVector
-    GlobalBasis const& basis() const
+    std::shared_ptr<GlobalBasis const> basis() const
     {
       return dofVector_->basis();
     }

src/amdis/gridfunctions/DiscreteFunction.inc.hpp

@@ -23,15 +23,17 @@ private:
   using Geometry = typename Element::Geometry;
 
 public:
+  /// Constructor. Stores a copy of the DiscreteFunction.
   LocalFunction(DiscreteFunction const& globalFunction)
     : globalFunction_(globalFunction)
-    , localView_(globalFunction_.basis().localView())
+    , localView_(globalFunction_.basis()->localView())
     , subTree_(&child(localView_.tree(), globalFunction_.treePath()))
   {}
 
+  /// Copy constructor.
   LocalFunction(LocalFunction const& other)
     : globalFunction_(other.globalFunction_)
-    , localView_(globalFunction_.basis().localView())
+    , localView_(globalFunction_.basis()->localView())
     , subTree_(&child(localView_.tree(), globalFunction_.treePath()))
   {}
 
@@ -100,15 +102,17 @@ private:
   using Geometry = typename Element::Geometry;
 
 public:
+  /// Constructor. Stores a copy of the DiscreteFunction.
   GradientLocalFunction(DiscreteFunction const& globalFunction)
     : globalFunction_(globalFunction)
-    , localView_(globalFunction_.basis().localView())
+    , localView_(globalFunction_.basis()->localView())
     , subTree_(&child(localView_.tree(), globalFunction_.treePath()))
   {}
 
+  /// Copy constructor.
   GradientLocalFunction(GradientLocalFunction const& other)
     : globalFunction_(other.globalFunction_)
-    , localView_(globalFunction_.basis().localView())
+    , localView_(globalFunction_.basis()->localView())
     , subTree_(&child(localView_.tree(), globalFunction_.treePath()))
   {}
 
@@ -260,7 +264,7 @@ operator()(Domain const& x) const
   using Grid = typename GlobalBasis::GridView::Grid;
   using IS = typename GlobalBasis::GridView::IndexSet;
 
-  auto const& gv = this->basis().gridView();
+  auto const& gv = this->basis()->gridView();
   Dune::HierarchicSearch<Grid,IS> hsearch{gv.grid(), gv.indexSet()};
 
   auto element = hsearch.findEntity(x);

src/amdis/linearalgebra/DOFMatrixBase.hpp

@@ -42,27 +42,33 @@ namespace AMDiS
     using ElementMatrix = Dune::DynamicMatrix<double>;
 
   public:
-    /// Constructor. Constructs new BaseVector.
-    DOFMatrixBase(RowBasis const& rowBasis, ColBasis const& colBasis)
-      : rowBasis_(&rowBasis)
-      , colBasis_(&colBasis)
+    /// Constructor. Stores the shared_ptr to the bases.
+    DOFMatrixBase(std::shared_ptr<RowBasis> rowBasis, std::shared_ptr<ColBasis> colBasis)
+      : rowBasis_(rowBasis)
+      , colBasis_(colBasis)
     {
-      operators_.init(rowBasis, colBasis);
+      operators_.init(*rowBasis_, *colBasis_);
     }
 
+    /// Constructor. Wraps the reference into a non-destroying shared_ptr or moves the basis into a new shared_ptr.
+    template <class RB_, class CB_>
+    DOFMatrixBase(RB_&& rowBasis, CB_&& colBasis)
+      : DOFMatrixBase(Dune::wrap_or_move(FWD(rowBasis)), Dune::wrap_or_move(FWD(colBasis)))
+    {}
+
     DOFMatrixBase(DOFMatrixBase const&) = delete;
     DOFMatrixBase(DOFMatrixBase&&) = delete;
 
     /// Return the row-basis \ref rowBasis of the matrix
-    RowBasis const& rowBasis() const
+    std::shared_ptr<RowBasis const> rowBasis() const
     {
-      return *rowBasis_;
+      return rowBasis_;
     }
 
     /// Return the col-basis \ref colBasis of the matrix
-    ColBasis const& colBasis() const
+    std::shared_ptr<ColBasis const> colBasis() const
     {
-      return *colBasis_;
+      return colBasis_;
     }
 
     /// Return the data-matrix
@@ -152,8 +158,8 @@ namespace AMDiS
 
   protected:
     /// The finite element space / basis associated with the rows/columns
-    RowBasis const* rowBasis_;
-    ColBasis const* colBasis_;
+    std::shared_ptr<RowBasis> rowBasis_;
+    std::shared_ptr<ColBasis> colBasis_;
 
     /// Data backend
     Backend backend_;

src/amdis/linearalgebra/DOFVectorBase.hpp

@@ -10,12 +10,35 @@
 #include <amdis/GridTransferManager.hpp>
 #include <amdis/OperatorList.hpp>
 #include <amdis/common/Math.hpp>
+#include <amdis/common/TypeTraits.hpp>
 #include <amdis/linearalgebra/DOFVectorInterface.hpp>
 #include <amdis/typetree/MultiIndex.hpp>
 #include <amdis/typetree/TreePath.hpp>
 
 namespace AMDiS
 {
+  /// The container that stores a data-vector and a corresponding basis
+  template <class GlobalBasis, class ValueType = double>
+  class DOFVector;
+
+  /// \brief Create a DOFVector from a basis.
+  /**
+   * This generator function accepts the basis as reference, temporary, or
+   * shared_ptr. Internally the reference is wrapped into a non-destroying
+   * shared_ptr and the temporary is moved into a new shared_ptr.
+   *
+   * The DataTransferOperation controls what is done during grid changes with the
+   * DOFVector. The default is interpolation of the data to the new grid. See
+   * \ref DataTransferOperation for more options.
+   **/
+  template <class ValueType = double, class GlobalBasis>
+  DOFVector<Underlying_t<GlobalBasis>, ValueType>
+  makeDOFVector(GlobalBasis&& basis, DataTransferOperation op = DataTransferOperation::INTERPOLATE)
+  {
+    return {FWD(basis), op};
+  }
+
+
   /// \brief The basic container that stores a base vector and a corresponding basis
   /**
    * Basis implementation of DOFVector, i.e. a vector storing all the
@@ -56,16 +79,22 @@ namespace AMDiS
     using DataTransferFactory = AMDiS::DataTransferFactory<Self>;
 
   public:
-    /// Constructor. Constructs new BaseVector.
-    DOFVectorBase(GlobalBasis& basis, DataTransferOperation op = DataTransferOperation::INTERPOLATE)
-      : basis_(&basis)
-      , dataTransfer_(DataTransferFactory::create(op, basis))
+    /// Constructor. Stores the shared_ptr of the basis and creates a new DataTransfer.
+    DOFVectorBase(std::shared_ptr<GlobalBasis> basis, DataTransferOperation op)
+      : basis_(basis)
+      , dataTransfer_(DataTransferFactory::create(op, *basis_))
     {
       compress();
       attachToGridTransfer();
-      operators_.init(basis);
+      operators_.init(*basis_);
     }
 
+    /// Constructor. Wraps the reference into a non-destroying shared_ptr or moves the basis into a new shared_ptr.
+    template <class GB_>
+    DOFVectorBase(GB_&& basis, DataTransferOperation op)
+      : DOFVectorBase(Dune::wrap_or_move(FWD(basis)), op)
+    {}
+
     /// Copy constructor
     DOFVectorBase(Self const& that)
       : basis_(that.basis_)
@@ -119,11 +148,17 @@ namespace AMDiS
     }
 
     /// Return the basis \ref basis_ associated to the vector
-    GlobalBasis const& basis() const
+    std::shared_ptr<GlobalBasis const> basis() const
     {
-      return *basis_;
+      return basis_;
     }
 
+    /// Return the (mutable) basis \ref basis_ associated to the vector
+    // std::shared_ptr<GlobalBasis> basis()
+    // {
+    //   return basis_;
+    // }
+
     /// Return the data-vector
     BaseVector const& vector() const
     {
@@ -240,7 +275,7 @@ namespace AMDiS
 
   private:
     /// The finite element space / basis associated with the data vector
-    GlobalBasis* basis_;
+    std::shared_ptr<GlobalBasis> basis_;
 
     /// Data backend
     Backend backend_;

src/amdis/linearalgebra/eigen/DOFVector.hpp

@@ -85,26 +85,25 @@ namespace AMDiS
   };
 
 
-  template <class BasisType, class ValueType = double>
-  class DOFVector : public DOFVectorBase<BasisType, EigenVector<ValueType>>
+  template <class GlobalBasis, class ValueType>
+  class DOFVector : public DOFVectorBase<GlobalBasis, EigenVector<ValueType>>
   {
-    using Super = DOFVectorBase<BasisType, EigenVector<ValueType>>;
+    using Super = DOFVectorBase<GlobalBasis, EigenVector<ValueType>>;
 
   public:
-    DOFVector(BasisType& basis, DataTransferOperation op = DataTransferOperation::INTERPOLATE)
+    DOFVector(std::shared_ptr<GlobalBasis> basis, DataTransferOperation op = DataTransferOperation::INTERPOLATE)
+      : Super(std::move(basis), op)
+    {}
+
+    DOFVector(GlobalBasis& basis, DataTransferOperation op = DataTransferOperation::INTERPOLATE)
       : Super(basis, op)
     {}
 
+    DOFVector(GlobalBasis&& basis, DataTransferOperation op = DataTransferOperation::INTERPOLATE)
+      : Super(std::move(basis), op)
+    {}
+
     using Super::operator=;
   };
 
-
-  /// Constructor a dofvector from given basis and name
-  template <class ValueType = double, class Basis>
-  DOFVector<Basis, ValueType>
-  makeDOFVector(Basis& basis, DataTransferOperation op = DataTransferOperation::INTERPOLATE)
-  {
-    return {basis, op};
-  }
-
 } // end namespace AMDiS

src/amdis/linearalgebra/istl/DOFVector.hpp

@@ -94,26 +94,25 @@ namespace AMDiS
 	};
 
 
-  template <class BasisType, class ValueType = double>
-  class DOFVector : public DOFVectorBase<BasisType, IstlVector<ValueType>>
+  template <class GlobalBasis, class ValueType>
+  class DOFVector : public DOFVectorBase<GlobalBasis, IstlVector<ValueType>>
   {
-    using Super = DOFVectorBase<BasisType, IstlVector<ValueType>>;
+    using Super = DOFVectorBase<GlobalBasis, IstlVector<ValueType>>;
 
   public:
-    DOFVector(BasisType& basis, DataTransferOperation op = DataTransferOperation::INTERPOLATE)
+    DOFVector(std::shared_ptr<GlobalBasis> basis, DataTransferOperation op = DataTransferOperation::INTERPOLATE)
+      : Super(std::move(basis), op)
+    {}
+
+    DOFVector(GlobalBasis& basis, DataTransferOperation op = DataTransferOperation::INTERPOLATE)
       : Super(basis, op)
     {}
 
+    DOFVector(GlobalBasis&& basis, DataTransferOperation op = DataTransferOperation::INTERPOLATE)
+      : Super(std::move(basis), op)
+    {}
+
     using Super::operator=;
   };
 
-
-  /// Constructor a dofvector from given basis and name
-  template <class ValueType = double, class Basis>
-  DOFVector<Basis, ValueType>
-  makeDOFVector(Basis& basis, DataTransferOperation op = DataTransferOperation::INTERPOLATE)
-  {
-    return {basis, op};
-  }
-
 } // end namespace AMDiS

src/amdis/linearalgebra/mtl/DOFVector.hpp

@@ -85,26 +85,25 @@ namespace AMDiS
   };
 
 
-  template <class BasisType, class ValueType = double>
-  class DOFVector : public DOFVectorBase<BasisType, MtlVector<ValueType>>
+  template <class GlobalBasis, class ValueType>
+  class DOFVector : public DOFVectorBase<GlobalBasis, MtlVector<ValueType>>
   {
-    using Super = DOFVectorBase<BasisType, MtlVector<ValueType>>;
+    using Super = DOFVectorBase<GlobalBasis, MtlVector<ValueType>>;
 
   public:
-    DOFVector(BasisType& basis, DataTransferOperation op = DataTransferOperation::INTERPOLATE)
+    DOFVector(std::shared_ptr<GlobalBasis> basis, DataTransferOperation op = DataTransferOperation::INTERPOLATE)
+      : Super(std::move(basis), op)
+    {}
+
+    DOFVector(GlobalBasis& basis, DataTransferOperation op = DataTransferOperation::INTERPOLATE)
       : Super(basis, op)
     {}
 
+    DOFVector(GlobalBasis&& basis, DataTransferOperation op = DataTransferOperation::INTERPOLATE)
+      : Super(std::move(basis), op)
+    {}
+
     using Super::operator=;
   };
 
-
-  /// Constructor a dofvector from given basis and name
-  template <class ValueType = double, class Basis>
-  DOFVector<Basis, ValueType>
-  makeDOFVector(Basis& basis, DataTransferOperation op = DataTransferOperation::INTERPOLATE)
-  {
-    return {basis, op};
-  }
-
 } // end namespace AMDiS

test/ExpressionsTest.cpp

@@ -94,7 +94,7 @@ int main(int argc, char** argv)
 
   // integration of expressions
 
-  auto gv = u.basis().gridView();
+  auto gv = u.basis()->gridView();
 
   DUNE_UNUSED auto int1 = integrate(op1, gv, 5);
   DUNE_UNUSED auto int2 = integrate(op2, gv, 5);

test/IntegrateTest.cpp

@@ -24,7 +24,7 @@ int main(int argc, char** argv)
   prob.initialize(INIT_ALL);
 
   auto u = prob.solution(0);
-  auto gv = u.basis().gridView();
+  auto gv = u.basis()->gridView();
 
   u << 1.0;
   auto i1 = integrate(1.0, gv);