Skip to content
Snippets Groups Projects
Commit e9b08a28 authored by oliver.sander_at_tu-dresden.de's avatar oliver.sander_at_tu-dresden.de
Browse files

Do not use deprecated VirtualFunction in tests

parent ff3a165c
Branches
No related tags found
1 merge request!46Remove deprecated stuff from tests
Pipeline #31415 passed
......@@ -7,7 +7,6 @@
#include <dune/common/fvector.hh>
#include <dune/common/fmatrix.hh>
#include <dune/common/function.hh>
#include <dune/istl/matrixindexset.hh>
......@@ -30,23 +29,6 @@
#include <dune/solvers/norms/energynorm.hh>
#include <dune/solvers/norms/twonorm.hh>
template <class DomainType, class RangeType>
class ConstantFunction :
public Dune::VirtualFunction<DomainType, RangeType>
{
public:
ConstantFunction(double c):
c_(c)
{}
void evaluate(const DomainType& x, RangeType& y) const
{
y = c_;
}
private:
double c_;
};
template<class GridView, class Matrix>
void constructPQ1Pattern(const GridView& gridView, Matrix& matrix)
......@@ -231,7 +213,6 @@ void assemblePQ1RHS(const GridView& gridView, Vector& r, const Function& f)
typedef typename FiniteElementCache::FiniteElementType FiniteElement;
typedef typename FiniteElement::Traits::LocalBasisType::Traits::RangeType RangeType;
typedef typename Function::RangeType FunctionRangeType;
const auto& indexSet = gridView.indexSet();
FiniteElementCache cache;
......@@ -265,8 +246,7 @@ void assemblePQ1RHS(const GridView& gridView, Vector& r, const Function& f)
fe.localBasis().evaluateFunction(quadPos, values);
// evaluate function
FunctionRangeType fAtPos;
f.evaluate(geometry.global(quadPos), fAtPos);
auto fAtPos = f(geometry.global(quadPos));
// add vector entries
double z = pt.weight() * integrationElement;
......
......@@ -123,7 +123,9 @@ bool checkWithGrid(const GridType& grid, const std::string fileName="")
Vector rhs(A.N());
rhs = 0;
ConstantFunction<DomainType, RangeType> f(50);
auto f = [](const DomainType& x) -> RangeType
{return 50.0; };
assemblePQ1RHS(gridView, rhs, f);
Vector u(A.N());
......
......@@ -114,7 +114,9 @@ bool checkWithGrid(const GridType& grid, const std::string fileName="")
Vector rhs(A.N());
rhs = 0;
ConstantFunction<DomainType, RangeType> f(50);
auto f = [](const DomainType& x) -> RangeType
{return 50.0; };
assemblePQ1RHS(gridView, rhs, f);
Vector u(A.N());
......
......@@ -84,7 +84,9 @@ bool checkWithGrid(const GridType& grid, const std::string fileName, int maxIter
Vector rhs(A.N());
rhs = 0;
ConstantFunction<DomainType, RangeType> f(50);
auto f = [](const DomainType& x) -> RangeType
{return 50.0; };
assemblePQ1RHS(gridView, rhs, f);
Vector u(A.N());
......
......@@ -77,7 +77,9 @@ bool checkWithGrid(const GridType& grid, const std::string fileName="")
Vector rhs(A.N());
rhs = 0;
ConstantFunction<DomainType, RangeType> f(50);
auto f = [](const DomainType& x) -> RangeType
{return 50.0; };
assemblePQ1RHS(gridView, rhs, f);
Vector u(A.N());
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment