Commit 9379904a authored by Praetorius, Simon's avatar Praetorius, Simon
Browse files

rebased on develop and cleaned up

parent ea53bf87
#pragma once
// TODO: Cleanup of copied comments
#include <string>
#include <dune/grid/common/grid.hh>
......@@ -10,7 +10,7 @@
namespace AMDiS {
/**
* \ingroup Adaption
*
......@@ -20,60 +20,56 @@ namespace AMDiS {
template <class Traits>
class Marker
{
protected:
using GlobalBasis = typename Traits::GlobalBasis;
using GridView = typename GlobalBasis::GridView;
using Grid = typename GridView::Grid;
using IndexSet = typename GridView::IndexSet;
using Element = typename GridView::template Codim<0>::Entity;
using EstType = std::vector<double>;
using Estimates = std::vector<double>;
public:
/// Constructor.
Marker() {}
/// Constructor.
Marker(std::string name_, std::string component_, const EstType& est_, std::shared_ptr<Grid> grid_)
: name(name_),
component(component_),
grid(grid_),
est(est_),
maximumMarking(false),
p(2),
info(10),
maxRefineLevel(-1),
minRefineLevel(-1),
refineAllowed(true),
coarsenAllowed(false)
Marker(std::string name, std::string component, Estimates const& est,
std::shared_ptr<Grid> const& grid)
: name_(name)
, component_(component)
, grid_(grid)
, est_(est)
, maximumMarking_(false)
{
Parameters::get(name + "->p", p);
Parameters::get(name + "->info", info);
Parameters::get(name + "->max refinement level", maxRefineLevel);
Parameters::get(name + "->min refinement level", minRefineLevel);
Parameters::get(name_ + "->p", p_);
Parameters::get(name_ + "->info", info_);
Parameters::get(name_ + "->max refinement level", maxRefineLevel_);
Parameters::get(name_ + "->min refinement level", minRefineLevel_);
}
/// destructor
virtual ~Marker() {}
/// Marks element with newMark. If \ref maximumMarking is set, the element
/// Marks element with newMark. If \ref maximumMarking_ is set, the element
/// is marked only if the new mark is bigger than the old one. The return
/// value specifies whether the element has been marked, or not.
inline void mark(const Element& elem, int newMark)
void mark(Element const& elem, int newMark)
{
int oldMark = grid->getMark(elem);
int oldMark = grid_->getMark(elem);
if (!maximumMarking || (newMark > oldMark)) {
bool (marked) = grid->mark(newMark, elem);
if (!maximumMarking_ || (newMark > oldMark)) {
bool marked = grid_->mark(newMark, elem);
if (marked) {
if (oldMark > 0) {
elMarkRefine--;
elMarkRefine_--;
} else if (oldMark < 0) {
elMarkCoarsen--;
elMarkCoarsen_--;
}
if (newMark > 0) {
elMarkRefine++;
elMarkRefine_++;
} else if (newMark < 0) {
elMarkCoarsen++;
elMarkCoarsen_++;
}
} else {
msg("Marking failed");
......@@ -88,88 +84,88 @@ namespace AMDiS {
virtual void finishMarking(AdaptInfo& adaptInfo);
/// Marks one element.
virtual void markElement(AdaptInfo& adaptInfo, const Element& elem);
virtual void markElement(AdaptInfo& adaptInfo, Element const& elem);
/// Marking of the mesh.
virtual Flag markGrid(AdaptInfo& adaptInfo);
/// Sets \ref maximumMarking.
inline void setMaximumMarking(bool maxMark)
void setMaximumMarking(bool maxMark)
{
maximumMarking = maxMark;
maximumMarking_ = maxMark;
}
inline int getElMarkRefine()
int getElMarkRefine()
{
return elMarkRefine;
return elMarkRefine_;
}
inline int getElMarkCoarsen()
int getElMarkCoarsen()
{
return elMarkCoarsen;
return elMarkCoarsen_;
}
/// Returns \ref name of the Marker
inline std::string getName() const
std::string getName() const
{
return name;
return name_;
}
/// Creates a scalar marker depending on the strategy set in parameters.
static std::shared_ptr<Marker<Traits> > createMarker(std::string name, std::string component_, const EstType& est_, std::shared_ptr<Grid> grid_);
static std::shared_ptr<Marker<Traits> > createMarker(std::string name, std::string component, Estimates const& est, std::shared_ptr<Grid> const& grid);
protected:
/// Name of the scalar marker.
std::string name;
std::string name_;
/// Problem component to work on
std::string component;
std::string component_;
/// Pointer to the grid
std::shared_ptr<Grid> grid;
std::shared_ptr<Grid> grid_;
/// Pointer to the local error estimates
EstType est;
Estimates est_;
/// estimation sum
double estSum;
double estSum_ = 0.0;
/// estmation maximum
double estMax;
double estMax_ = 0.0;
/// If true, elements are marked only if the new value is bigger than
/// the current marking.
bool maximumMarking;
bool maximumMarking_;
/// Lower limit for error estimation, from which an element is marked for
/// refinement
double markRLimit;
double markRLimit_;
/// Upper limit for error estimation, from which an element is marked for
/// coarsening
double markCLimit;
double markCLimit_;
/// power in estimator norm
double p;
/// power in estimator norm TODO: is this info necessary in marker?
double p_ = 2.0;
/// Info level.
int info;
int info_ = 10;
/// Counter for elements marked for refinement
int elMarkRefine;
int elMarkRefine_ = 0;
/// Counter for elements marked for coarsening
int elMarkCoarsen;
int elMarkCoarsen_ = 0;
/// Maximal level of all elements.
int maxRefineLevel;
int maxRefineLevel_ = -1;
/// Minimal level of all elements.
int minRefineLevel;
int minRefineLevel_ = -1;
bool refineAllowed;
bool refineAllowed_ = true;
bool coarsenAllowed;
bool coarsenAllowed_ = false;
};
......@@ -180,25 +176,25 @@ namespace AMDiS {
* Global refinement.
*/
template <class Traits>
class GRMarker : public Marker<Traits>
class GRMarker
: public Marker<Traits>
{
using GlobalBasis = typename Traits::GlobalBasis;
using GridView = typename GlobalBasis::GridView;
using Grid = typename GridView::Grid;
using IndexSet = typename GridView::IndexSet;
using Element = typename GridView::template Codim<0>::Entity;
using EstType = std::vector<double>;
using Super = Marker<Traits>;
using Grid = typename Super::Grid;
using Element = typename Super::Element;
using Estimates = typename Super::Estimates;
public:
/// Constructor.
GRMarker(std::string name_, std::string component_, const EstType& est_, std::shared_ptr<Grid> grid_)
: Marker<Traits>(name_, component_, est_, grid_)
GRMarker(std::string name, std::string component, Estimates const& est,
std::shared_ptr<Grid> const& grid)
: Marker<Traits>(std::move(name), std::move(component), est, grid)
{}
/// Implementation of Marker::markElement().
virtual void markElement(AdaptInfo& adaptInfo, const Element& elem)
virtual void markElement(AdaptInfo& adaptInfo, Element const& elem)
{
if (this->refineAllowed)
if (this->refineAllowed_)
this->mark(elem, 1);
}
};
......@@ -210,26 +206,23 @@ namespace AMDiS {
* \brief
* Maximum strategy.
*/
template <class Traits>
class MSMarker : public Marker<Traits>
class MSMarker
: public Marker<Traits>
{
using GlobalBasis = typename Traits::GlobalBasis;
using GridView = typename GlobalBasis::GridView;
using Grid = typename GridView::Grid;
using IndexSet = typename GridView::IndexSet;
using Element = typename GridView::template Codim<0>::Entity;
using EstType = std::vector<double>;
using Super = Marker<Traits>;
using Grid = typename Super::Grid;
using Estimates = typename Super::Estimates;
public:
/// Constructor.
MSMarker(std::string name_, std::string component_, const EstType& est_, std::shared_ptr<Grid> grid_)
: Marker<Traits>(name_, component_, est_, grid_),
MSGamma(0.5),
MSGammaC(0.1)
MSMarker(std::string name, std::string component, Estimates const& est,
std::shared_ptr<Grid> const& grid)
: Super{std::move(name), std::move(component), est, grid}
{
Parameters::get(this->name + "->MSGamma", MSGamma);
Parameters::get(this->name + "->MSGammaC", MSGammaC);
Parameters::get(this->name_ + "->MSGamma", msGamma_);
Parameters::get(this->name_ + "->MSGammaC", msGammaC_);
}
/// Implementation of Marker::initMarking().
......@@ -237,10 +230,10 @@ namespace AMDiS {
protected:
/// Marking parameter.
double MSGamma;
double msGamma_ = 0.5;
/// Marking parameter.
double MSGammaC;
double msGammaC_ = 0.1;
};
......@@ -250,26 +243,23 @@ namespace AMDiS {
* \brief
* Equidistribution strategy
*/
template <class Traits>
class ESMarker : public Marker<Traits>
class ESMarker
: public Marker<Traits>
{
using GlobalBasis = typename Traits::GlobalBasis;
using GridView = typename GlobalBasis::GridView;
using Grid = typename GridView::Grid;
using IndexSet = typename GridView::IndexSet;
using Element = typename GridView::template Codim<0>::Entity;
using EstType = std::vector<double>;
using Super = Marker<Traits>;
using Grid = typename Super::Grid;
using Estimates = typename Super::Estimates;
public:
/// Constructor.
ESMarker(std::string name_, std::string component_, const EstType& est_, std::shared_ptr<Grid> grid_)
: Marker<Traits>(name_, component_, est_, grid_),
ESTheta(0.9),
ESThetaC(0.2)
ESMarker(std::string name, std::string component, Estimates const& est,
std::shared_ptr<Grid> const& grid)
: Super{std::move(name), std::move(component), est, grid}
{
Parameters::get(this->name + "->ESTheta", ESTheta);
Parameters::get(this->name + "->ESThetaC", ESThetaC);
Parameters::get(this->name_ + "->ESTheta", esTheta_);
Parameters::get(this->name_ + "->ESThetaC", esThetaC_);
}
/// Implementation of Marker::initMarking().
......@@ -277,10 +267,10 @@ namespace AMDiS {
protected:
/// Marking parameter.
double ESTheta;
double esTheta_ = 0.9;
/// Marking parameter.
double ESThetaC;
double esThetaC_ = 0.2;
};
......@@ -290,29 +280,25 @@ namespace AMDiS {
* \brief
* Guaranteed error reduction strategy
*/
template <class Traits>
class GERSMarker : public Marker<Traits>
class GERSMarker
: public Marker<Traits>
{
using GlobalBasis = typename Traits::GlobalBasis;
using GridView = typename GlobalBasis::GridView;
using Grid = typename GridView::Grid;
using IndexSet = typename GridView::IndexSet;
using Element = typename GridView::template Codim<0>::Entity;
using EstType = std::vector<double>;
using Super = Marker<Traits>;
using Grid = typename Super::Grid;
using Element = typename Super::Element;
using Estimates = typename Super::Estimates;
public:
/// Constructor.
GERSMarker(std::string name_, std::string component_, const EstType& est_, std::shared_ptr<Grid> grid_)
: Marker<Traits>(name_, component_, est_, grid_),
oldErrSum(0.0),
GERSThetaStar(0.6),
GERSNu(0.1),
GERSThetaC(0.1)
GERSMarker(std::string name, std::string component, Estimates const& est,
std::shared_ptr<Grid> const& grid)
: Super{std::move(name), std::move(component), est, grid}
{
Parameters::get(this->name + "->GERSThetaStar", GERSThetaStar);
Parameters::get(this->name + "->GERSNu", GERSNu);
Parameters::get(this->name + "->GERSThetaC", GERSThetaC);
Parameters::get(this->name_ + "->GERSThetaStar", gersThetaStar_);
Parameters::get(this->name_ + "->GERSNu", gersNu_);
Parameters::get(this->name_ + "->GERSThetaC", gersThetaC_);
}
/// Implementation of Marker::markGrid().
......@@ -320,28 +306,28 @@ namespace AMDiS {
protected:
/// Refinement marking function.
void markElementForRefinement(AdaptInfo& adaptInfo, const Element& elem);
void markElementForRefinement(AdaptInfo& adaptInfo, Element const& elem);
/// Coarsening marking function.
void markElementForCoarsening(AdaptInfo& adaptInfo, const Element& elem);
void markElementForCoarsening(AdaptInfo& adaptInfo, Element const& elem);
protected:
/// Marking parameter.
double GERSSum;
double gersSum_ = 0.0;
/// Marking parameter.
double oldErrSum;
double oldErrSum_ = 0.0;
/// Marking parameter.
double GERSThetaStar;
double gersThetaStar_ = 0.6;
/// Marking parameter.
double GERSNu;
double gersNu_ = 0.1;
/// Marking parameter.
double GERSThetaC;
double gersThetaC_ = 0.1;
};
}
#include "Marker.inc.hpp"
\ No newline at end of file
#include "Marker.inc.hpp"
// TODO: Cleanup of copied comments
#include <amdis/common/Math.hpp>
#include <cmath>
namespace AMDiS {
using std::pow;
template <class Traits>
std::shared_ptr<Marker<Traits> > Marker<Traits>::createMarker(std::string name, std::string component, const EstType& est, std::shared_ptr<Grid> grid)
{
int strategy = 0;
Parameters::get(name + "->strategy", strategy);
std::shared_ptr<Marker<Traits> > marker;
switch (strategy) {
case 0: // no refinement/coarsening
break;
case 1:
msg("Creating global refinement (GR) marker");
marker = std::make_shared<GRMarker<Traits> >(name, component, est, grid);
break;
case 2:
msg("Creating maximum strategy (MS) marker");
marker = std::make_shared<MSMarker<Traits> >(name, component, est, grid);
break;
case 3:
msg("Creating equidistribution strategy (ES) marker");
marker = std::make_shared<ESMarker<Traits> >(name, component, est, grid);
break;
case 4:
msg("Creating guaranteed error reduction strategy (GERS) marker");
marker = std::make_shared<GERSMarker<Traits> >(name, component, est, grid);
break;
default:
error_exit("invalid strategy");
}
return marker;
template <class Traits>
std::shared_ptr<Marker<Traits> > Marker<Traits>::createMarker(std::string name, std::string component, Estimates const& est, std::shared_ptr<Grid> const& grid)
{
int strategy = 0;
Parameters::get(name + "->strategy", strategy);
std::shared_ptr<Marker<Traits> > marker;
switch (strategy) {
case 0: // no refinement/coarsening
break;
case 1:
msg("Creating global refinement (GR) marker");
marker = std::make_shared<GRMarker<Traits> >(name, component, est, grid);
break;
case 2:
msg("Creating maximum strategy (MS) marker");
marker = std::make_shared<MSMarker<Traits> >(name, component, est, grid);
break;
case 3:
msg("Creating equidistribution strategy (ES) marker");
marker = std::make_shared<ESMarker<Traits> >(name, component, est, grid);
break;
case 4:
msg("Creating guaranteed error reduction strategy (GERS) marker");
marker = std::make_shared<GERSMarker<Traits> >(name, component, est, grid);
break;
default:
error_exit("invalid strategy");
}
return marker;
}
template <class Traits>
void Marker<Traits>::initMarking(AdaptInfo& adaptInfo)
{
elMarkRefine = 0;
elMarkCoarsen = 0;
estSum = pow(adaptInfo.getEstSum(component), p);
estMax = adaptInfo.getEstMax(component);
refineAllowed = adaptInfo.isRefinementAllowed(component);
coarsenAllowed = adaptInfo.isCoarseningAllowed(component);
}
template <class Traits>
void Marker<Traits>::initMarking(AdaptInfo& adaptInfo)
{
elMarkRefine_ = 0;
elMarkCoarsen_ = 0;
estSum_ = std::pow(adaptInfo.getEstSum(component_), p_);
estMax_ = adaptInfo.getEstMax(component_);
refineAllowed_ = adaptInfo.isRefinementAllowed(component_);
coarsenAllowed_ = adaptInfo.isCoarseningAllowed(component_);
}
template <class Traits>
void Marker<Traits>::finishMarking(AdaptInfo& adaptInfo)
{
msg(elMarkRefine, " elements marked for refinement");
msg(elMarkCoarsen, " elements marked for coarsening");
}
template <class Traits>
void Marker<Traits>::finishMarking(AdaptInfo& adaptInfo)
{
msg(elMarkRefine_, " elements marked for refinement");
msg(elMarkCoarsen_, " elements marked for coarsening");
}
template <class Traits>
void Marker<Traits>::markElement(AdaptInfo& adaptInfo, const Element& elem)
{
const auto& index = grid->leafIndexSet().index(elem);
double lError = est[index];
if (lError > markRLimit && refineAllowed
&& (maxRefineLevel == -1 || elem.level() < maxRefineLevel)) {
this->mark(elem, 1);
} else {
if (lError <= markCLimit && coarsenAllowed
&& (minRefineLevel == -1 || elem.level() > minRefineLevel))
if (/*!elem->getElementData()->getElementData(COARSENABLE) ||*/
lError /*+ elem->getCoarseningEstimation(index)*/ <= markCLimit)
this->mark(elem, -1);
}
template <class Traits>
void Marker<Traits>::markElement(AdaptInfo& adaptInfo, const Element& elem)
{
const auto& index = grid_->leafIndexSet().index(elem);
double lError = est_[index];
if (lError > markRLimit_ && refineAllowed_
&& (maxRefineLevel_ == -1 || elem.level() < maxRefineLevel_)) {
this->mark(elem, 1);
} else if (lError <= markCLimit_ && coarsenAllowed_
&& (minRefineLevel_ == -1 || elem.level() > minRefineLevel_)) {
this->mark(elem, -1);
}
}
template <class Traits>
Flag Marker<Traits>::markGrid(AdaptInfo& adaptInfo)
{
if (!(grid))
error_exit("No grid!");
template <class Traits>
Flag Marker<Traits>::markGrid(AdaptInfo& adaptInfo)
{
if (!grid_)
error_exit("No grid!");
initMarking(adaptInfo);
initMarking(adaptInfo);
if (!coarsenAllowed &&
!refineAllowed) {
return 0;
}
for (const auto& elem : Dune::elements(grid->leafGridView())) {
markElement(adaptInfo, elem);
}
if (!coarsenAllowed_ && !refineAllowed_)
return 0;
finishMarking(adaptInfo);
for (const auto& elem : Dune::elements(grid_->leafGridView()))
markElement(adaptInfo, elem);
Flag markFlag;
if (elMarkRefine)
markFlag = 1;
if (elMarkCoarsen)
markFlag |= 2;