Commit 084ba1ca authored by Müller, Felix's avatar Müller, Felix Committed by Praetorius, Simon
Browse files

Add a Manager class for GridTransfer with automatic registration of DOFVectors

parent 5abe55d5
...@@ -31,6 +31,7 @@ install(FILES ...@@ -31,6 +31,7 @@ install(FILES
GridFunctionOperator.hpp GridFunctionOperator.hpp
GridFunctions.hpp GridFunctions.hpp
GridTransfer.hpp GridTransfer.hpp
GridTransferManager.hpp
Initfile.hpp Initfile.hpp
InitfileParser.hpp InitfileParser.hpp
LinearAlgebra.hpp LinearAlgebra.hpp
......
#pragma once #pragma once
#include <vector> #include <list>
#include <amdis/linear_algebra/DOFVectorInterface.hpp> #include <amdis/linear_algebra/DOFVectorInterface.hpp>
#include <amdis/Output.hpp>
namespace AMDiS namespace AMDiS
{ {
class GridTransferInterface
{
public:
virtual ~GridTransferInterface() = default;
virtual void attach(DOFVectorInterface*) = 0;
virtual void detach(DOFVectorInterface*) = 0;
virtual bool preAdapt() = 0;
virtual bool adapt() = 0;
virtual void postAdapt() = 0;
};
template <class Grid> template <class Grid>
class GridTransfer class GridTransfer
: public GridTransferInterface
{ {
using Self = GridTransfer; using Self = GridTransfer;
public: public:
GridTransfer(Grid& grid) void bind(Grid& grid)
: grid_(&grid) {
{} grid_ = &grid;
}
void unbind()
{
grid_ = nullptr;
}
/// Attach a data container to the grid transfer, that gets interpolated during grid change /// Attach a data container to the grid transfer, that gets interpolated during grid change
void attach(DOFVectorInterface* vec) virtual void attach(DOFVectorInterface* vec) override
{ {
data_.push_back(vec); data_.push_back(vec);
} }
virtual void detach(DOFVectorInterface* vec) override
{
auto it = std::find(data_.begin(), data_.end(), vec);
if (it != data_.end())
data_.erase(it);
else
warning("DOFVector to detach not found");
}
/// Prepare the grid and the data for the adaption /// Prepare the grid and the data for the adaption
bool preAdapt() virtual bool preAdapt() override
{ {
assert(grid_ != nullptr);
mightCoarsen_ = grid_->preAdapt(); // any element might be coarsened in adapt() mightCoarsen_ = grid_->preAdapt(); // any element might be coarsened in adapt()
for (auto* vec : data_) for (auto* vec : this->data_)
vec->preAdapt(mightCoarsen_); vec->preAdapt(mightCoarsen_);
return mightCoarsen_; return mightCoarsen_;
} }
/// do the grid adaption /// do the grid adaption
bool adapt() virtual bool adapt() override
{ {
assert(grid_ != nullptr);
refined_ = grid_->adapt(); // returns true if a least one entity was refined refined_ = grid_->adapt(); // returns true if a least one entity was refined
return refined_; return refined_;
} }
// Perform data adaption to the new grid // Perform data adaption to the new grid
void postAdapt() virtual void postAdapt() override
{ {
assert(grid_ != nullptr);
if (mightCoarsen_ || refined_) { if (mightCoarsen_ || refined_) {
for (auto* vec : data_) for (auto* vec : this->data_)
vec->postAdapt(refined_); vec->postAdapt(refined_);
} }
grid_->postAdapt(); grid_->postAdapt();
} }
protected: private:
Grid* grid_; Grid* grid_ = nullptr;
std::vector<DOFVectorInterface*> data_; std::list<DOFVectorInterface*> data_;
bool mightCoarsen_ = false; bool mightCoarsen_ = false;
bool refined_ = false; bool refined_ = false;
}; };
......
#pragma once
#include <cstdint>
#include <map>
#include <memory>
#include <amdis/GridTransfer.hpp>
#include <amdis/utility/ConcurrentCache.hpp>
namespace AMDiS
{
/// Static administration class for automatic handling of DOFVectors during grid adaption
class GridTransferManager
{
template <class T1, class T2>
friend class DOFVectorBase;
using Self = GridTransferManager;
using Key = std::uintptr_t;
using Data = std::unique_ptr<GridTransferInterface>;
private:
// Constructors. Since this is a static class you cannot call any constructor.
GridTransferManager() = delete;
public:
/**
* \brief Adapts the grid according to previously set marks and performs data transfer on
* all DOFVectors associated to the given grid.
*
* This function retrieves the grid's stored GridTransfer and calls its adaptation methods.
* During this process the grid will change according to previously set grid element markings.
* All DOFVectors associated to the given grid will be interpolated to the new grid according to
* their DataTransfer member object.
*
* Takes a grid as argument. Marking of the grid needs to be performed prior to calling this.
* Returns true if the grid changed during adaptation.
**/
template <class Grid>
static bool adapt(Grid& grid)
{
auto gridTransfer = GridTransferManager::gridTransfer(grid);
bool adapted = gridTransfer.preAdapt();
adapted |= gridTransfer.adapt();
gridTransfer.postAdapt();
return adapted;
}
/// Returns the GridTransfer corresponding to a Grid, to be used during the adaptation cycle
template <class Grid>
static GridTransfer<Grid>& gridTransfer(Grid& grid)
{
GridTransfer<Grid>* gridTransfer
= dynamic_cast<GridTransfer<Grid>*>(Self::gridTransferInterface(grid));
gridTransfer->bind(grid);
return *gridTransfer;
}
private:
// DOFVector registration methods. They are called automaticly by the DOFVectors.
template <class Vec>
static void attach(Vec& vec)
{
auto const& grid = vec.basis().gridView().grid();
Self::gridTransferInterface(grid)->attach(&vec);
}
template <class Vec>
static void detach(Vec& vec)
{
auto const& grid = vec.basis().gridView().grid();
Self::gridTransferInterface(grid)->detach(&vec);
}
// Returns the GridTransferInterface class,
// used for attaching and detaching DOFVectors to a GridTransfer
template <class Grid>
static GridTransferInterface* gridTransferInterface(Grid const &grid)
{
Key key = Key(&grid);
auto& gridTransferInterfaceUniquePtr = GridTransferCache::get(key, [&](Key const&)
{
return std::make_unique<GridTransfer<Grid>>();
});
return gridTransferInterfaceUniquePtr.get();
}
private:
// Associate to each Grid (represented as address) a GridTransfer
using GridTransferCache = ConcurrentCache< Key, Data, StaticLockedPolicy, std::map<Key, Data> >;
};
} // end namespace AMDiS
...@@ -18,7 +18,6 @@ ...@@ -18,7 +18,6 @@
#include <amdis/DirichletBC.hpp> #include <amdis/DirichletBC.hpp>
//#include <amdis/Estimator.hpp> //#include <amdis/Estimator.hpp>
#include <amdis/Flag.hpp> #include <amdis/Flag.hpp>
#include <amdis/GridTransfer.hpp>
#include <amdis/Initfile.hpp> #include <amdis/Initfile.hpp>
#include <amdis/LinearAlgebra.hpp> #include <amdis/LinearAlgebra.hpp>
#include <amdis/LinearSolvers.hpp> #include <amdis/LinearSolvers.hpp>
...@@ -347,12 +346,6 @@ namespace AMDiS ...@@ -347,12 +346,6 @@ namespace AMDiS
setGrid(Dune::stackobject_to_shared_ptr(grid)); setGrid(Dune::stackobject_to_shared_ptr(grid));
} }
void addAdaptionData(DOFVectorInterface* vec)
{
assert(bool(gridTransfer_));
gridTransfer_->attach(vec);
}
void addMarker(std::shared_ptr<Marker<Grid>> const& marker) void addMarker(std::shared_ptr<Marker<Grid>> const& marker)
{ {
marker_.push_back(marker); marker_.push_back(marker);
...@@ -384,7 +377,6 @@ namespace AMDiS ...@@ -384,7 +377,6 @@ namespace AMDiS
void adoptGrid(std::shared_ptr<Grid> const& grid) void adoptGrid(std::shared_ptr<Grid> const& grid)
{ {
grid_ = grid; grid_ = grid;
gridTransfer_ = std::make_shared<GridTransfer<Grid>>(*grid_);
Parameters::get(name_ + "->mesh", gridName_); Parameters::get(name_ + "->mesh", gridName_);
} }
...@@ -421,9 +413,6 @@ namespace AMDiS ...@@ -421,9 +413,6 @@ namespace AMDiS
/// Grid of this problem. /// Grid of this problem.
std::shared_ptr<Grid> grid_; std::shared_ptr<Grid> grid_;
/// Handling the adaption of the grid
std::shared_ptr<GridTransfer<Grid>> gridTransfer_;
/// Name of the grid /// Name of the grid
std::string gridName_ = "mesh"; std::string gridName_ = "mesh";
......
...@@ -10,6 +10,7 @@ ...@@ -10,6 +10,7 @@
#include <amdis/AdaptInfo.hpp> #include <amdis/AdaptInfo.hpp>
#include <amdis/FileWriter.hpp> #include <amdis/FileWriter.hpp>
#include <amdis/GridFunctionOperator.hpp> #include <amdis/GridFunctionOperator.hpp>
#include <amdis/GridTransferManager.hpp>
#include <amdis/LocalAssembler.hpp> #include <amdis/LocalAssembler.hpp>
#include <amdis/common/Loops.hpp> #include <amdis/common/Loops.hpp>
...@@ -37,7 +38,6 @@ void ProblemStat<Traits>::initialize( ...@@ -37,7 +38,6 @@ void ProblemStat<Traits>::initialize(
adoptFlag.isSet(INIT_SYSTEM) || adoptFlag.isSet(INIT_SYSTEM) ||
adoptFlag.isSet(INIT_FE_SPACE))) { adoptFlag.isSet(INIT_FE_SPACE))) {
grid_ = adoptProblem->grid_; grid_ = adoptProblem->grid_;
gridTransfer_ = adoptProblem->gridTransfer_;
} }
} }
...@@ -123,7 +123,6 @@ void ProblemStat<Traits>::createGrid() ...@@ -123,7 +123,6 @@ void ProblemStat<Traits>::createGrid()
{ {
Parameters::get(name_ + "->mesh", gridName_); Parameters::get(name_ + "->mesh", gridName_);
grid_ = MeshCreator<Grid>::create(gridName_); grid_ = MeshCreator<Grid>::create(gridName_);
gridTransfer_ = std::make_shared<GridTransfer<Grid>>(*grid_);
msg("Create grid:"); msg("Create grid:");
msg("#elements = {}" , grid_->size(0)); msg("#elements = {}" , grid_->size(0));
...@@ -175,10 +174,6 @@ void ProblemStat<Traits>::createMatricesAndVectors() ...@@ -175,10 +174,6 @@ void ProblemStat<Traits>::createMatricesAndVectors()
solution_ = std::make_shared<SystemVector>(*globalBasis_, INTERPOLATE); solution_ = std::make_shared<SystemVector>(*globalBasis_, INTERPOLATE);
rhs_ = std::make_shared<SystemVector>(*globalBasis_, NO_OPERATION); rhs_ = std::make_shared<SystemVector>(*globalBasis_, NO_OPERATION);
assert(bool(gridTransfer_));
gridTransfer_->attach(solution_.get());
gridTransfer_->attach(rhs_.get());
auto localView = globalBasis_->localView(); auto localView = globalBasis_->localView();
AMDiS::forEachNode_(localView.tree(), [&,this](auto const& node, auto treePath) AMDiS::forEachNode_(localView.tree(), [&,this](auto const& node, auto treePath)
{ {
...@@ -333,11 +328,8 @@ adaptGrid(AdaptInfo& adaptInfo) ...@@ -333,11 +328,8 @@ adaptGrid(AdaptInfo& adaptInfo)
{ {
Dune::Timer t; Dune::Timer t;
bool adapted = gridTransfer_->preAdapt(); bool adapted = GridTransferManager::adapt(*grid_);
adapted |= gridTransfer_->adapt(); // NOTE: |= does not short-circuit and works with bools as ||= globalBasis_->update(gridView());
if (adapted)
globalBasis_->update(gridView());
gridTransfer_->postAdapt();
msg("adaptGrid needed {} seconds", t.elapsed()); msg("adaptGrid needed {} seconds", t.elapsed());
return adapted ? MESH_ADAPTED : Flag(0); return adapted ? MESH_ADAPTED : Flag(0);
......
#pragma once #pragma once
#include <cmath> #include <cmath>
#include <utility>
#include <dune/functions/functionspacebases/sizeinfo.hh> #include <dune/functions/functionspacebases/sizeinfo.hh>
#include <amdis/DataTransfer.hpp> #include <amdis/DataTransfer.hpp>
#include <amdis/GridTransferManager.hpp>
#include <amdis/LocalAssemblerList.hpp> #include <amdis/LocalAssemblerList.hpp>
#include <amdis/common/Math.hpp> #include <amdis/common/Math.hpp>
#include <amdis/common/ScalarTypes.hpp> #include <amdis/common/ScalarTypes.hpp>
...@@ -52,18 +54,47 @@ namespace AMDiS ...@@ -52,18 +54,47 @@ namespace AMDiS
, dataTransfer_(DataTransferFactory::create(op, basis)) , dataTransfer_(DataTransferFactory::create(op, basis))
{ {
compress(); compress();
GridTransferManager::attach(*this);
operators_.init(basis); operators_.init(basis);
} }
DOFVectorBase(Self const&) = default; /// Copy constructor
DOFVectorBase(Self&&) = default; DOFVectorBase(Self const& that)
: basis_(that.basis_)
, backend_(that.backend_)
, elementVector_(that.elementVector_)
, operators_(that.operators_)
, dataTransfer_(that.dataTransfer_)
{
GridTransferManager::attach(*this);
}
/// Move constructor
DOFVectorBase(Self&& that)
: basis_(std::move(that.basis_))
, backend_(std::move(that.backend_))
, elementVector_(std::move(that.elementVector_))
, operators_(std::move(that.operators_))
, dataTransfer_(std::move(that.dataTransfer_))
{
GridTransferManager::attach(*this);
}
/// Destructor
virtual ~DOFVectorBase() override
{
GridTransferManager::detach(*this);
}
/// Copy assignment operator /// Copy assignment operator
Self& operator=(Self const& that) Self& operator=(Self const& that)
{ {
GridTransferManager::detach(*this);
basis_ = that.basis_; basis_ = that.basis_;
backend_.resize(that.size()); backend_.resize(that.size());
backend_ = that.backend_; backend_ = that.backend_;
dataTransfer_ = that.dataTransfer_;
GridTransferManager::attach(*this);
return *this; return *this;
} }
......
...@@ -7,6 +7,7 @@ ...@@ -7,6 +7,7 @@
#include <dune/functions/functionspacebases/powerbasis.hh> #include <dune/functions/functionspacebases/powerbasis.hh>
#include <dune/functions/functionspacebases/lagrangebasis.hh> #include <dune/functions/functionspacebases/lagrangebasis.hh>
#include <amdis/GridTransferManager.hpp>
#include <amdis/LinearAlgebra.hpp> #include <amdis/LinearAlgebra.hpp>
#include "Tests.hpp" #include "Tests.hpp"
...@@ -41,26 +42,41 @@ int main(int argc, char** argv) ...@@ -41,26 +42,41 @@ int main(int argc, char** argv)
Dune::FieldVector<double, 2> L; L = 1.0; Dune::FieldVector<double, 2> L; L = 1.0;
auto s = Dune::filledArray<2>(1); auto s = Dune::filledArray<2>(1);
Dune::YaspGrid<2> grid(L, s); Dune::YaspGrid<2> grid(L, s);
auto const& gridView = grid.leafGridView();
// create basis // create basis
auto basis = makeBasis(grid.leafGridView(), auto basis = makeBasis(gridView,
composite(power<2>(lagrange<2>(), flatInterleaved()), lagrange<1>(), flatLexicographic())); composite(power<2>(lagrange<2>(), flatInterleaved()), lagrange<1>(), flatLexicographic()));
using Basis = decltype(basis); using Basis = decltype(basis);
DOFVector<Basis> vec1(basis);
test_dofvector(basis, vec1);
DOFVector<Basis, float> vec2(basis); {
test_dofvector(basis, vec2); DOFVector<Basis> vec1(basis);
test_dofvector(basis, vec1);
DOFVector<Basis, float> vec2(basis);
test_dofvector(basis, vec2);
auto vec3 = makeDOFVector(basis); auto vec3 = makeDOFVector(basis);
test_dofvector(basis, vec3); test_dofvector(basis, vec3);
auto vec4 = makeDOFVector<float>(basis); auto vec4 = makeDOFVector<float>(basis);
test_dofvector(basis, vec4); test_dofvector(basis, vec4);
#if DUNE_HAVE_CXX_CLASS_TEMPLATE_ARGUMENT_DEDUCTION #if DUNE_HAVE_CXX_CLASS_TEMPLATE_ARGUMENT_DEDUCTION
DOFVector vec5(basis); DOFVector vec5(basis);
test_dofvector(basis, vec5); test_dofvector(basis, vec5);
#endif #endif
} }
\ No newline at end of file
// test GridTransferManager registration
{
DOFVector<Basis> vec1(basis);
test_dofvector(basis, vec1);
for (auto const& e : elements(gridView))
grid.mark(1, e);
GridTransferManager::adapt(grid);
AMDIS_TEST_EQ(vec1.size(), basis.dimension());
}
report_errors();
}
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment