Commit 520d1838 authored by Praetorius, Simon's avatar Praetorius, Simon

merge intel bug workaround

parents 1492541a 18ecdb30
......@@ -93,7 +93,6 @@ void Assembler<Traits>::assemble(
}
template <class Traits>
template <class ElementContainer, class Container, class Operators, class Geometry, class LocalView>
void Assembler<Traits>::assembleElementOperators(
......@@ -135,13 +134,13 @@ void Assembler<Traits>::assembleElementOperators(
template <class Traits>
template <class ElementContainer, class Container, class Operators, class Geometry, class RowLocalViews, class ColLocalView>
template <class ElementContainer, class Container, class Operators, class Geometry, class RowLocalView, class ColLocalView>
void Assembler<Traits>::assembleElementOperators(
ElementContainer& elementContainer,
Container& container,
Operators& operators,
Geometry const& geometry,
RowLocalViews const& rowLocalView, ColLocalView const& colLocalView) const
RowLocalView const& rowLocalView, ColLocalView const& colLocalView) const
{
auto const& element = getElement(rowLocalView, colLocalView);
auto const& gridView = getGridView(rowLocalView, colLocalView);
......
......@@ -265,7 +265,7 @@ namespace AMDiS
{
using RawExpr = Underlying_t<Expr>;
static_assert(Concepts::HasGridFunctionOrder<RawExpr>,
"Polynomial degree of expression can not be deduced. You need to provide an explicit value for polynomial order or a quadrature rule in `makeOperator()`.");
"Polynomial degree of expression can not be deduced. You need to provide an explicit value for polynomial order or a quadrature rule in 'makeOperator()'.");
return ExpressionPreOperator<Tag, Expr, tag::deduce>{t, expr};
}
......
......@@ -137,6 +137,11 @@ void ProblemStat<Traits>::addMatrixOperator(
RowTreePath row, ColTreePath col,
double* factor, double* estFactor)
{
static_assert( Concept::PreTreePath<RowTreePath>,
"row must be a valid treepath, or an integer/index-constant");
static_assert( Concept::PreTreePath<ColTreePath>,
"col must be a valid treepath, or an integer/index-constant");
auto i = child(globalBasis->localView().tree(), makeTreePath(row));
auto j = child(globalBasis->localView().tree(), makeTreePath(col));
......@@ -156,6 +161,11 @@ void ProblemStat<Traits>::addMatrixOperator(
RowTreePath row, ColTreePath col,
double* factor, double* estFactor)
{
static_assert( Concept::PreTreePath<RowTreePath>,
"row must be a valid treepath, or an integer/index-constant");
static_assert( Concept::PreTreePath<ColTreePath>,
"col must be a valid treepath, or an integer/index-constant");
auto i = child(globalBasis->localView().tree(), makeTreePath(row));
auto j = child(globalBasis->localView().tree(), makeTreePath(col));
using Intersection = typename GridView::Intersection;
......@@ -175,6 +185,9 @@ void ProblemStat<Traits>::addVectorOperator(
TreePath path,
double* factor, double* estFactor)
{
static_assert( Concept::PreTreePath<TreePath>,
"path must be a valid treepath, or an integer/index-constant");
auto i = child(globalBasis->localView().tree(), makeTreePath(path));
auto op = makeGridOperator(preOp, globalBasis->gridView());
......@@ -193,6 +206,9 @@ void ProblemStat<Traits>::addVectorOperator(
TreePath path,
double* factor, double* estFactor)
{
static_assert( Concept::PreTreePath<TreePath>,
"path must be a valid treepath, or an integer/index-constant");
auto i = child(globalBasis->localView().tree(), makeTreePath(path));
using Intersection = typename GridView::Intersection;
......
......@@ -116,7 +116,7 @@ namespace AMDiS
}
/// \brief Return the LocalFunction of the AnalyticGridFunction.
friend auto localFunction(AnalyticGridFunction const& gf)
friend LocalFunction localFunction(AnalyticGridFunction const& gf)
{
return LocalFunction{gf.fct_};
}
......
......@@ -47,19 +47,23 @@ namespace AMDiS
return 0;
}
/// \relates ConstantLocalFunction
friend auto derivative(ConstantLocalFunction const& lf)
{
using RawSignature = typename Dune::Functions::SignatureTraits<R(D)>::RawSignature;
using DerivativeRange = typename Dune::Functions::DefaultDerivativeTraits<RawSignature>::Range;
DerivativeRange diff(0);
return ConstantLocalFunction<DerivativeRange(D),LocalContext,DerivativeRange>{diff};
}
template <class R_, class D_, class L_, class T_>
friend auto derivative(ConstantLocalFunction<R_(D_), L_, T_> const& lf);
private:
T value_;
};
/// \relates ConstantLocalFunction
template <class R, class D, class LocalContext, class T>
auto derivative(ConstantLocalFunction<R(D), LocalContext, T> const& lf)
{
using RawSignature = typename Dune::Functions::SignatureTraits<R(D)>::RawSignature;
using DerivativeRange = typename Dune::Functions::DefaultDerivativeTraits<RawSignature>::Range;
DerivativeRange diff(0);
return ConstantLocalFunction<DerivativeRange(D),LocalContext,DerivativeRange>{diff};
}
/**
* \addtogroup GridFunctions
......@@ -101,7 +105,7 @@ namespace AMDiS
}
/// Return the LocalFunction of the \ref ConstantGridFunction. \relates ConstantGridFunction
friend auto localFunction(ConstantGridFunction const& gf)
friend LocalFunction localFunction(ConstantGridFunction const& gf)
{
return LocalFunction{gf.value_};
}
......
......@@ -8,6 +8,68 @@
namespace AMDiS
{
namespace Concepts
{
/** \addtogroup Concepts
* @{
**/
namespace Definition
{
template <class TP>
struct IsPreTreePath
: any_of_t<std::is_same<TP, int>::value, std::is_same<TP,std::size_t>::value>
{};
template <int I>
struct IsPreTreePath<int_t<I>>
: std::true_type
{};
template <std::size_t I>
struct IsPreTreePath<index_t<I>>
: std::true_type
{};
template <std::size_t... I>
struct IsPreTreePath<Indices<I...>>
: std::true_type
{};
} // end namespace Definition
template <class TP>
constexpr bool PreTreePath = Dune::TypeTree::IsTreePath<TP>::value || Definition::IsPreTreePath<TP>::value;
/** @} **/
} // end namespace Concepts
#ifdef DOXYGEN
/// \brief Converts a (pre)TreePath into a HybridTreePath
/**
* Converts an integer, an integralconstant or a Dune TreePath into an
* \ref Dune::TypeTree::HybridTreePath that is used in GlobalBasis traversal.
*
* **Requirements:**
* - `PreTreePath` one of
* + integer type (`int, std::size_t`),
* + integral constant (`std::integral_constant<[int|std::size_t], i>`)
* + any Dune TreePath (`TreePath<std::size_t...>, DynamicTreePath, HybridTreePath<class... T>`)
*
* **Example:**
* ```
* makeTreePath(1),
* makeTreePath(int_<2>),
* makeTreePath(treepath(1, int_<2>))
* ```
**/
template <class PreTreePath>
auto makeTreePath(PreTreePath tp);
#else // DOXYGEN
auto makeTreePath(int i) { return Dune::TypeTree::hybridTreePath(std::size_t(i)); }
auto makeTreePath(std::size_t i) { return Dune::TypeTree::hybridTreePath(i); }
......@@ -28,61 +90,46 @@ namespace AMDiS
template <class TP>
auto makeTreePath(TP const&) {
static_assert( Dune::TypeTree::IsTreePath<TP>::value,
static_assert( Concepts::PreTreePath<TP>,
"Argument must be a valid treepath, or an integer/index-constant");
return Dune::TypeTree::hybridTreePath();
}
#endif // DOXYGEN
namespace Impl
{
template <std::size_t I, std::size_t N>
struct PrintTreePath
{
template <class... T>
static void apply(std::ostream& os, Dune::TypeTree::HybridTreePath<T...> const& tp)
{
os << "," << treePathIndex(tp, index_<I>);
PrintTreePath<I+1,N>::apply(os, tp);
}
};
template <std::size_t N>
struct PrintTreePath<N,N>
{
template <class... T>
static void apply(std::ostream&, Dune::TypeTree::HybridTreePath<T...> const&) {}
};
template <std::size_t N>
struct PrintTreePath<0,N>
{
template <class... T>
static void apply(std::ostream& os, Dune::TypeTree::HybridTreePath<T...> const& tp)
{
os << treePathIndex(tp, _0);
PrintTreePath<1,N>::apply(os, tp);
}
};
template <>
struct PrintTreePath<0,0>
template <class TreePath, std::size_t... I>
void printTreePath(std::ostream& os, TreePath const& tp, std::index_sequence<I...>)
{
template <class... T>
static void apply(std::ostream&, Dune::TypeTree::HybridTreePath<T...> const&) {}
};
(void)std::initializer_list<int>{
((void)(os << treePathIndex(tp, index_<I>) << ","), 0)...
};
}
}
template <class... T>
std::string to_string(Dune::TypeTree::HybridTreePath<T...> const& tp)
template <class T0, class... T>
std::string to_string(Dune::TypeTree::HybridTreePath<T0,T...> const& tp)
{
std::stringstream ss;
Impl::PrintTreePath<0,sizeof...(T)>::apply(ss, tp);
Impl::printTreePath(ss, tp, std::make_index_sequence<sizeof...(T)>{});
ss << back(tp);
return ss.str();
}
// alias for hybridTreePath(...)
std::string to_string(Dune::TypeTree::HybridTreePath<> const& tp)
{
return "";
}
/// \brief Generate a TreePath from a sequence of integers and integral-constants
/**
* Converts a sequence of arguments to a \ref Dune::TypeTree::HybridTreePath.
* The arguments can be one of
* - integer type (`int, std::size_t`)
* - integral constant (`std::integral_constant<std::size_t,i>, index_t<i>`)
**/
template <class... T>
constexpr Dune::TypeTree::HybridTreePath<T...> treepath(T const&... t)
{
......
......@@ -50,53 +50,42 @@ int main(int argc, char** argv)
auto _p = 1_c;
// <1/tau * u, v>
Op opTime;
opTime.addZOT(density);
auto opTime = makeOperator(tag::testvec_trialvec{}, density);
prob.addMatrixOperator(opTime, _v, _v);
// <1/tau * u^old, v>
Op opTimeOld;
opTimeOld.addZOT(valueOf(prob.getSolution(_v), density));
auto opTimeOld = makeOperator(tag::testvec{}, density * prob.getSolution(_v));
prob.addVectorOperator(opTimeOld, _v);
// sum_i <grad(u_i),grad(v_i)>
Op opL;
opL.addSOT(viscosity);
prob.addMatrixOperator(opL, _v, _v);
for (std::size_t i = 0; i < DOW; ++i) {
auto opL = makeOperator(tag::gradtest_gradtrial{}, viscosity);
prob.addMatrixOperator(opL, treepath(_v,i), treepath(_v,i));
}
// <p, div(v)>
Op opP;
opP.addFOT(1.0, GRD_PSI);
auto opP = makeOperator(tag::divtestvec_trial{}, 1.0);
prob.addMatrixOperator(opP, _v, _p);
// <div(u), q>
Op opDiv;
opDiv.addFOT(1.0, GRD_PHI);
auto opDiv = makeOperator(tag::test_divtrialvec{}, 1.0);
prob.addMatrixOperator(opDiv, _p, _v);
for (std::size_t i = 0; i < DOW; ++i) {
// <(u * nabla)u_i^old, v_i>
for (std::size_t j = 0; j < DOW; ++j) {
Op opNonlin1;
opNonlin1.addZOT(derivativeOf(prob.getSolution(treepath(_v,i)), j, density));
prob.addMatrixOperator(opNonlin1, treepath(_v,i), treepath(_v,j));
}
// <(u * nabla)u_i^old, v_i>
auto opNonlin1 = makeOperator(tag::testvec_trialvec{}, density * trans(gradientAtQP(prob.getSolution(_v))));
prob.addMatrixOperator(opNonlin1, _v, _v);
for (std::size_t i = 0; i < DOW; ++i) {
// <(u^old * nabla)u_i, v_i>
Op opNonlin2;
opNonlin2.addFOT(valueOf(prob.getSolution(_v)), GRD_PHI);
auto opNonlin2 = makeOperator(tag::test_gradtrial{}, density * prob.getSolution(_v));
prob.addMatrixOperator(opNonlin2, treepath(_v,i), treepath(_v,i));
// <(u^old * grad(u_i^old)), v_i>
Op opNonlin3;
opNonlin3.addZOT( func([density](Dune::FieldVector<double,DOW> const& u, Dune::FieldVector<double,DOW> const& grad_ui)
{
return density * (u * grad_ui);
}, valueOf(prob.getSolution(_v)), gradientOf(prob.getSolution(treepath(_v,i)))) );
prob.addVectorOperator(opNonlin, treepath(_v,i));
}
// <(u^old * grad(u_i^old)), v_i>
auto opNonlin3 = makeOperator(tag::testvec{}, trans(gradientAtQP(prob.getSolution(_v))) * prob.getSolution(_v));
prob.addVectorOperator(opNonlin3, _v);
// define boundary regions
auto left = [](auto const& x) { return x[0] < 1.e-8; };
auto not_left = [](auto const& x) { return x[0] > 1.0 - 1.e-8 || x[1] < 1.e-8 || x[1] > 1.0 - 1.e-8; };
......@@ -119,6 +108,9 @@ int main(int argc, char** argv)
// set point constraint for pressure
prob.addDirichletBC([](auto const& x) { return x[0] < 1.e-8 && x[1] < 1.e-8; }, _p, _p, 0.0);
// set initial conditions
prob.getSolution(_v).interpolate(parabolic_y);
prob.getSolution(_p).interpolate(0.0);
AdaptInstationary adapt("adapt", prob, adaptInfo, probInstat, adaptInfo);
adapt.adapt();
......
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