...
 
......@@ -112,4 +112,5 @@ namespace AMDiS {
(*vector)[dofIndices[i]] = value;
}
}
}
} // end namespace AMDiS
......@@ -29,10 +29,10 @@
#include "BoundaryCondition.h"
#include "AbstractFunction.h"
#include "FixVec.h"
#include "Expressions.h"
namespace AMDiS
{
// some tags used for tag-dispatching
struct _value_by_abstractfunction {};
struct _value_by_dofvector {};
......@@ -93,102 +93,133 @@ namespace AMDiS
} // end namespace detail
/**
* \ingroup Assembler
*
* \brief
* Sub class of BoundaryCondition. Implements Dirichlet boundary conditions.
* A DOFVectors is set to a given value at a Dirichlet dof and in a DOFMatrix
* the row corresponding to a Dirichlet dof is replaced by a row containing
* only a 1.0 in the diagonal.
*/
template <class ValueTag>
class DirichletBC {};
// specialization for AbstractFunctions as value container
template <>
class DirichletBC<_value_by_abstractfunction> : public detail::DirichletBC
{
typedef detail::DirichletBC super;
public:
/// Constructor.
DirichletBC(BoundaryType type,
AbstractFunction<double, WorldVector<double> > *fct,
const FiniteElemSpace *rowFeSpace,
const FiniteElemSpace *colFeSpace = NULL,
bool apply = true)
: super(type, rowFeSpace, colFeSpace, apply),
container(*fct)
{ }
/// Implementation of BoundaryCondition::fillBoundaryCondition().
virtual void fillBoundaryCondition(DOFVectorBase<double>* vector,
ElInfo* elInfo,
const DegreeOfFreedom* dofIndices,
const BoundaryType* localBound,
int nBasFcts) override;
protected:
AbstractFunction<double, WorldVector<double> > &container;
};
// specialization for DOFVectors as value container
template <>
class DirichletBC<_value_by_dofvector> : public detail::DirichletBC
{
typedef detail::DirichletBC super;
public:
/// Constructor.
DirichletBC(BoundaryType type,
DOFVectorBase<double> *vec,
bool apply = true);
/**
* \ingroup Assembler
*
* \brief
* Sub class of BoundaryCondition. Implements Dirichlet boundary conditions.
* A DOFVectors is set to a given value at a Dirichlet dof and in a DOFMatrix
* the row corresponding to a Dirichlet dof is replaced by a row containing
* only a 1.0 in the diagonal.
*/
template <class ValueTag>
class DirichletBC {};
// specialization for AbstractFunctions as value container
template <>
class DirichletBC<_value_by_abstractfunction> : public detail::DirichletBC
{
typedef detail::DirichletBC super;
public:
/// Constructor.
DirichletBC(BoundaryType type,
AbstractFunction<double, WorldVector<double> > *fct,
const FiniteElemSpace *rowFeSpace,
const FiniteElemSpace *colFeSpace = NULL,
bool apply = true)
: super(type, rowFeSpace, colFeSpace, apply),
container(*fct)
{ }
/// Implementation of BoundaryCondition::fillBoundaryCondition().
virtual void fillBoundaryCondition(DOFVectorBase<double>* vector,
ElInfo* elInfo,
const DegreeOfFreedom* dofIndices,
const BoundaryType* localBound,
int nBasFcts) override;
protected:
AbstractFunction<double, WorldVector<double> > &container;
};
// specialization for DOFVectors as value container
template <>
class DirichletBC<_value_by_dofvector> : public detail::DirichletBC
{
typedef detail::DirichletBC super;
public:
/// Constructor.
DirichletBC(BoundaryType type,
DOFVectorBase<double> *vec,
bool apply = true);
/// Implementation of BoundaryCondition::fillBoundaryCondition().
virtual void fillBoundaryCondition(DOFVectorBase<double>* vector,
ElInfo* elInfo,
const DegreeOfFreedom* dofIndices,
const BoundaryType* localBound,
int nBasFcts) override;
/// Implementation of BoundaryCondition::fillBoundaryCondition().
virtual void fillBoundaryCondition(DOFVectorBase<double>* vector,
ElInfo* elInfo,
const DegreeOfFreedom* dofIndices,
const BoundaryType* localBound,
int nBasFcts) override;
protected:
DOFVectorBase<double> *container;
};
protected:
DOFVectorBase<double> *container;
};
#if HAS_CXX11
// specialization for std::function or lambdas as value container
template <>
class DirichletBC<_value_by_function> : public detail::DirichletBC
{
typedef detail::DirichletBC super;
public:
/// Constructor.
DirichletBC(BoundaryType type,
std::function<double(WorldVector<double>)> fct,
const FiniteElemSpace *rowFeSpace,
const FiniteElemSpace *colFeSpace = NULL,
bool apply = true)
: super(type, rowFeSpace, colFeSpace, apply),
container(fct)
{ }
/// Implementation of BoundaryCondition::fillBoundaryCondition().
virtual void fillBoundaryCondition(DOFVectorBase<double>* vector,
ElInfo* elInfo,
const DegreeOfFreedom* dofIndices,
const BoundaryType* localBound,
int nBasFcts) override;
protected:
std::function<double(WorldVector<double>)> container;
};
// specialization for std::function or lambdas as value container
template <>
class DirichletBC<_value_by_function> : public detail::DirichletBC
{
typedef detail::DirichletBC super;
public:
/// Constructor.
DirichletBC(BoundaryType type,
std::function<double(WorldVector<double>)> fct,
const FiniteElemSpace *rowFeSpace,
const FiniteElemSpace *colFeSpace = NULL,
bool apply = true)
: super(type, rowFeSpace, colFeSpace, apply),
container(fct)
{ }
/// Implementation of BoundaryCondition::fillBoundaryCondition().
virtual void fillBoundaryCondition(DOFVectorBase<double>* vector,
ElInfo* elInfo,
const DegreeOfFreedom* dofIndices,
const BoundaryType* localBound,
int nBasFcts) override;
protected:
std::function<double(WorldVector<double>)> container;
};
#endif
}
// specialization for Expressions as value container
template <class Term>
class DirichletExpressionBC : public detail::DirichletBC
{
typedef detail::DirichletBC super;
public:
/// Constructor.
DirichletExpressionBC(BoundaryType type,
Term const& term,
const FiniteElemSpace *rowFeSpace,
const FiniteElemSpace *colFeSpace = NULL,
bool apply = true)
: super(type, rowFeSpace, colFeSpace, apply),
term(term),
ot(term)
{ }
/// Implementation of BoundaryCondition::fillBoundaryCondition().
virtual void fillBoundaryCondition(DOFVectorBase<double>* vector,
ElInfo* elInfo,
const DegreeOfFreedom* dofIndices,
const BoundaryType* localBound,
int nBasFcts) override;
protected:
Term term;
GenericOperatorTerm<Term> ot;
};
} // end namespace AMDiS
#include "DirichletBC.hh"
#endif
/******************************************************************************
*
* AMDiS - Adaptive multidimensional simulations
*
* Copyright (C) 2013 Dresden University of Technology. All Rights Reserved.
* Web: https://fusionforge.zih.tu-dresden.de/projects/amdis
*
* Authors:
* Simon Vey, Thomas Witkowski, Andreas Naumann, Simon Praetorius, et al.
*
* This file is provided AS IS with NO WARRANTY OF ANY KIND, INCLUDING THE
* WARRANTY OF DESIGN, MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE.
*
*
* This file is part of AMDiS
*
* See also license.opensource.txt in the distribution.
*
******************************************************************************/
namespace AMDiS
{
template <class Term>
void DirichletExpressionBC<Term>::fillBoundaryCondition(
DOFVectorBase<double>* vector,
ElInfo* elInfo,
const DegreeOfFreedom* dofIndices,
const BoundaryType* localBound,
int nBasFcts)
{
bool isBoundaryElement = false;
for (int i = 0; i < nBasFcts; i++) {
if (localBound[i] == boundaryType) {
isBoundaryElement = true;
break;
}
}
if (!isBoundaryElement)
return;
const BasisFunction *basFcts = rowFeSpace->getBasisFcts();
term.initElement(&ot, elInfo, NULL, NULL, basFcts);
for (int i = 0; i < nBasFcts; i++) {
if (localBound[i] == boundaryType) {
typename Term::value_type value = term(i);
vector->setDirichletDofValue(dofIndices[i], value);
(*vector)[dofIndices[i]] = value;
}
}
}
} // end namespace AMDiS
......@@ -22,6 +22,7 @@
#define AMDIS_FUNCTORS_H
#include "AbstractFunction.h"
#include "Traits.h"
#include <boost/math/special_functions/cbrt.hpp>
#include <boost/math/special_functions/pow.hpp>
......
......@@ -34,11 +34,12 @@
#include "Boundary.h"
#include "MatrixVector.h"
#include "StandardProblemIteration.h"
#include "io/ElementFileWriter.h"
#include "ComponentTraverseInfo.h"
#include "AbstractFunction.h"
#include "solver/SolverMatrix.h"
#include "SystemVector.h"
#include "expressions/expr_traits.hpp"
#include "io/ElementFileWriter.h"
#include "solver/SolverMatrix.h"
namespace AMDiS {
......@@ -204,6 +205,11 @@ namespace AMDiS {
std::function<double(WorldVector<double>)> b);
#endif
template <class Term>
typename boost::enable_if<typename traits::is_expr<Term>::type>::type
addDirichletBC(BoundaryType type, int row, int col,
Term const& term);
/// Adds a Dirichlet boundary condition, where the rhs is given by a DOF
/// vector.
virtual void addDirichletBC(BoundaryType type, int row, int col,
......@@ -749,4 +755,6 @@ namespace AMDiS {
#endif
}
#include "ProblemStat.hh"
#endif
/******************************************************************************
*
* AMDiS - Adaptive multidimensional simulations
*
* Copyright (C) 2013 Dresden University of Technology. All Rights Reserved.
* Web: https://fusionforge.zih.tu-dresden.de/projects/amdis
*
* Authors:
* Simon Vey, Thomas Witkowski, Andreas Naumann, Simon Praetorius, et al.
*
* This file is provided AS IS with NO WARRANTY OF ANY KIND, INCLUDING THE
* WARRANTY OF DESIGN, MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE.
*
*
* This file is part of AMDiS
*
* See also license.opensource.txt in the distribution.
*
******************************************************************************/
#include "DirichletBC.h"
namespace AMDiS
{
template <class Term>
typename boost::enable_if<typename traits::is_expr<Term>::type>::type
ProblemStatSeq::addDirichletBC(BoundaryType type, int row, int col,
Term const& term)
{
FUNCNAME("ProblemStat::addDirichletBC()");
TEST_EXIT(row >= 0 && row < nComponents)("Wrong row number: %d\n", row);
TEST_EXIT(col >= 0 && col < nComponents)("Wrong col number: %d\n", col);
boundaryConditionSet = true;
DirichletExpressionBC<Term> *dirichletApply = new DirichletExpressionBC<Term>(type, term, componentSpaces[row], componentSpaces[col], true);
DirichletExpressionBC<Term> *dirichletNotApply = new DirichletExpressionBC<Term>(type, term, componentSpaces[row], componentSpaces[col], false);
for (int i = 0; i < nComponents; i++)
if (systemMatrix && (*systemMatrix)[row][i]) {
if (i == col)
(*systemMatrix)[row][i]->getBoundaryManager()->addBoundaryCondition(dirichletApply);
else
(*systemMatrix)[row][i]->getBoundaryManager()->addBoundaryCondition(dirichletNotApply);
}
if (rhs)
rhs->getDOFVector(row)->getBoundaryManager()->addBoundaryCondition(dirichletApply);
if (solution)
solution->getDOFVector(col)->getBoundaryManager()->addBoundaryCondition(dirichletApply);
}
} // end namespace AMDiS
......@@ -26,7 +26,9 @@
#define AMDIS_COORDS_EXPRESSION_HPP
#include "AMDiS_fwd.h"
#include "ElInfo.h"
#include "LazyOperatorTerm.h"
#include "Mesh.h"
namespace AMDiS
{
......
......@@ -27,6 +27,8 @@
#include "LazyOperatorTerm.h"
#include <boost/lexical_cast.hpp>
namespace AMDiS
{
namespace expressions
......
......@@ -190,10 +190,10 @@ namespace AMDiS {
newTimestep = adaptInfo->getTimestep();
} else {
if (firstTimestep || succRejection) {
fac = pow((timeTol / errorEst), 1.0 / order);
fac = std::pow((timeTol / errorEst), 1.0 / order);
} else {
fac = adaptInfo->getTimestep() / tauAcc *
pow((timeTol * estAcc / (errorEst * errorEst)), 1.0 / order);
std::pow((timeTol * estAcc / (errorEst * errorEst)), 1.0 / order);
}
fac = std::max(0.1, std::min(fac, 3.0));
newTimestep = fac * adaptInfo->getTimestep();
......@@ -210,13 +210,13 @@ namespace AMDiS {
succRejection = true;
double reducedP =
log(errorEst / estRej) / log(adaptInfo->getTimestep() / tauRej);
std::log(errorEst / estRej) / std::log(adaptInfo->getTimestep() / tauRej);
if (reducedP < order && reducedP > 0.0)
order = reducedP;
}
newTimestep = pow((timeTol / errorEst), 1.0 / order) * adaptInfo->getTimestep();
newTimestep = std::pow((timeTol / errorEst), 1.0 / order) * adaptInfo->getTimestep();
newTimestep *= 0.95;
tauRej = adaptInfo->getTimestep();
......
#include "AMDiS.h"
#include "boost/date_time/posix_time/posix_time.hpp"
using namespace AMDiS;
using namespace boost::posix_time;
// ===========================================================================
// ===== function definitions ================================================
// ===========================================================================
/// Dirichlet boundary function
class G : public AbstractFunction<double, WorldVector<double> >
{
public:
/// Implementation of AbstractFunction::operator().
double operator()(const WorldVector<double>& x) const
{
return exp(-10.0 * (x * x));
}
};
/// RHS function
class F : public AbstractFunction<double, WorldVector<double> >
{
......@@ -60,18 +45,18 @@ int main(int argc, char* argv[])
// ===== create matrix operator =====
Operator matrixOperator(ellipt.getFeSpace());
matrixOperator.addTerm(new Simple_SOT);
addSOT(matrixOperator, 1.0);
ellipt.addMatrixOperator(matrixOperator, 0, 0);
int degree = ellipt.getFeSpace()->getBasisFcts()->getDegree();
// ===== create rhs operator =====
Operator rhsOperator(ellipt.getFeSpace());
rhsOperator.addTerm(new CoordsAtQP_ZOT(new F(degree)));
addZOT(rhsOperator, eval(new F(degree)));
ellipt.addVectorOperator(rhsOperator, 0);
// ===== add boundary conditions =====
ellipt.addDirichletBC(1, 0, 0, new G);
ellipt.addDirichletBC(1, 0, 0, exp(-10.0 * (X() * X())) );
// ===== start adaption loop =====
......