Skip to content
Snippets Groups Projects
Commit bd0b7ed4 authored by Oliver Sander's avatar Oliver Sander
Browse files

Remove class VintageBasisGridFunction

The class `VintageBasisGridFunction` has been removed,
because it relied on the old `dune-fufem` interface
for function space bases.

Use `Dune::Functions::DiscreteGlobalBasisFunction` instead.
parent 0d4253a6
No related tags found
1 merge request!145Remove class VintageBasisGridFunction
Pipeline #59649 failed
......@@ -40,6 +40,10 @@
- The class `NamedFunctionMap` based on `Dune::VirtualFunction` has been removed.
Use `std::map<:std::string, std::function<Range(Domain)>>` as a replacement.
- The class `VintageBasisGridFunction` has been removed, because it relied on the
old `dune-fufem` interface for function space bases.
Use `Dune::Functions::DiscreteGlobalBasisFunction` instead.
# 2.9 Release
......
......@@ -8,7 +8,6 @@ install(FILES
gridfunctionhelper.hh
portablegreymap.hh
portablepixelmap.hh
vintagebasisgridfunction.hh
virtualdifferentiablefunction.hh
virtualgridfunction.hh
DESTINATION ${CMAKE_INSTALL_INCLUDEDIR}/dune/fufem/functions)
// -*- tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 2 -*-
// vi: set et ts=4 sw=2 sts=2:
#ifndef DUNE_FUFEM_FUNCTIONS_VINTAGEBASISGRIDFUNCTION_HH
#define DUNE_FUFEM_FUNCTIONS_VINTAGEBASISGRIDFUNCTION_HH
#include <type_traits>
#include <optional>
#include <memory>
#include <dune/common/shared_ptr.hh>
#include <dune/common/exceptions.hh>
#include <dune/matrix-vector/axpy.hh>
#include <dune/functions/gridfunctions/gridviewentityset.hh>
#include <dune/functions/gridfunctions/gridfunction.hh>
namespace Dune::Fufem {
namespace Impl {
template<class Basis, class Coefficients, class DerivativeRangeType>
struct LocalDerivativeEvaluator;
template<class Basis, class Coefficients, class RangeType>
struct LocalEvaluator
{
template<class Element, class LocalDomainType>
static void evaluateLocal(const Basis& basis,
const Coefficients& coefficients,
const Element& e,
const LocalDomainType& x,
RangeType& y)
{
typedef typename Basis::LinearCombination LinearCombination;
typedef typename Basis::LocalFiniteElement::Traits::LocalBasisType LB;
typedef typename LB::Traits::RangeType LBRangeType;
const LB& localBasis = basis.getLocalFiniteElement(e).localBasis();
std::vector<LBRangeType> values;
localBasis.evaluateFunction(x, values);
y = 0;
for (size_t i=0; i<localBasis.size(); ++i)
{
int index = basis.index(e, i);
if (basis.isConstrained(index))
{
const LinearCombination& constraints = basis.constraints(index);
for (size_t w=0; w<constraints.size(); ++w)
Dune::MatrixVector::addProduct(y, values[i] * constraints[w].factor, coefficients[constraints[w].index]);
}
else
Dune::MatrixVector::addProduct(y, values[i], coefficients[index]);
}
}
};
template<class Basis, class Coefficients, class DerivativeRangeType>
struct LocalDerivativeEvaluator
{
template<class Element, class LocalDomainType>
static void evaluateLocal(const Basis& basis,
const Coefficients& c,
const Element& e,
const LocalDomainType& x,
DerivativeRangeType& d)
{
DUNE_THROW(Dune::NotImplemented, "LocalDerivativeEvaluator not implemented for this DerivativeRangeType");
}
};
template<class Basis, class Coefficients, class RangeFieldType, int n, int m>
struct LocalDerivativeEvaluator<Basis, Coefficients, typename Dune::FieldMatrix<RangeFieldType,n,m> >
{
using DerivativeRangeType = typename Dune::FieldMatrix<RangeFieldType,n,m>;
//! Compute y+= a*A*B
template<class Weight, class Jacobian>
static void gaxpy(const Weight& a,
const typename Coefficients::value_type& A,
const Jacobian& B,
DerivativeRangeType& y)
{
for (int i=0; i<DerivativeRangeType::rows; ++i)
for (int j=0; j<DerivativeRangeType::cols; ++j)
y[i][j] += a * A[i] * B[0][j];
}
//! Compute y+= A*B
template<class Jacobian>
static void gaxpy(const typename Coefficients::value_type& A,
const Jacobian& B,
DerivativeRangeType& y)
{
for (int i=0; i<DerivativeRangeType::rows; ++i)
for (int j=0; j<DerivativeRangeType::cols; ++j)
y[i][j] += A[i] * B[0][j];
}
template<class Element, class LocalDomainType>
static void evaluateLocal(const Basis& basis,
const Coefficients& coefficients,
const Element& e,
const LocalDomainType& x,
DerivativeRangeType& d)
{
typedef typename Basis::LinearCombination LinearCombination;
typedef typename Basis::LocalFiniteElement::Traits::LocalBasisType LB;
typedef typename LB::Traits::JacobianType LBJacobianType;
const LB& localBasis = basis.getLocalFiniteElement(e).localBasis();
std::vector<LBJacobianType> localJacobians;
localBasis.evaluateJacobian(x, localJacobians);
DerivativeRangeType localD(0.0);
for (size_t i=0; i<localBasis.size(); ++i)
{
int index = basis.index(e, i);
if (basis.isConstrained(index))
{
const LinearCombination& constraints = basis.constraints(index);
for (size_t w=0; w<constraints.size(); ++w)
gaxpy(constraints[w].factor, coefficients[constraints[w].index], localJacobians[i], localD);
}
else
gaxpy(coefficients[index], localJacobians[i], localD);
}
typename Element::Geometry::JacobianInverseTransposed J = e.geometry().jacobianInverseTransposed(x);
for (int j=0; j<DerivativeRangeType::rows; ++j)
J.mv(localD[j], d[j]);
}
};
/**
* \brief A dune-functions grid function implemented using a dune-fufem basis
*
* This allows to use a dune-fufem basis as dune-functions
* grid function. Moreover you can customize the local function
* evaluation.
*/
template<typename B, typename V, typename LEV, typename R = typename V::value_type>
class VintageBasisGridFunction
{
public:
using Basis = B;
using Vector = V;
using LocalEvaluator = LEV;
using Coefficient = std::decay_t<decltype(std::declval<Vector>()[0])>;
using GridView = typename Basis::GridView;
using EntitySet = Dune::Functions::GridViewEntitySet<GridView, 0>;
using Domain = typename EntitySet::GlobalCoordinate;
using Range = R;
using LocalDomain = typename EntitySet::LocalCoordinate;
using Element = typename EntitySet::Element;
using Traits = Dune::Functions::Imp::GridFunctionTraits<Range(Domain), EntitySet, Dune::Functions::DefaultDerivativeTraits, 16>;
class LocalFunction
{
public:
using GlobalFunction = VintageBasisGridFunction;
using Domain = GlobalFunction::LocalDomain;
using Range = GlobalFunction::Range;
using Element = GlobalFunction::Element;
LocalFunction(const VintageBasisGridFunction& globalFunction)
: globalFunction_(&globalFunction)
{}
LocalFunction(const LocalFunction& other) = default;
LocalFunction(LocalFunction&& other) = default;
LocalFunction& operator=(const LocalFunction& other) = default;
LocalFunction& operator=(LocalFunction&& other) = default;
/**
* \brief Bind LocalFunction to grid element.
*
* You must call this method before operator()
* and after changes to the coefficient vector.
*/
void bind(const Element& element) {
element_ = element;
}
void unbind() {
element_.reset();
}
bool bound() const {
return static_cast<bool>(element_);
}
Range operator()(const Domain& x) const {
auto y = Range(0);
globalFunction_->localEvaluator()(
globalFunction_->basis(),
globalFunction_->dofs(),
*element_,
x,
y);
return y;
}
const Element& localContext() const {
return *element_;
}
friend typename Traits::LocalFunctionTraits::DerivativeInterface derivative(const LocalFunction& t) {
DUNE_THROW(NotImplemented,"not implemented");
}
private:
const VintageBasisGridFunction* globalFunction_;
std::optional<Element> element_;
};
template<class B_T, class V_T, class LEV_T>
VintageBasisGridFunction(B_T&& basis, V_T&& v, LEV_T&& localEvaluator)
: entitySet_(basis.getGridView())
, basis_(wrap_or_move(std::forward<B_T>(basis)))
, coefficients_(wrap_or_move(std::forward<V_T>(v)))
, localEvaluator_(std::forward<LEV_T>(localEvaluator))
{}
template<class B_T, class V_T>
VintageBasisGridFunction(B_T&& basis, V_T&& v)
: entitySet_(basis.getGridView())
, basis_(wrap_or_move(std::forward<B_T>(basis)))
, coefficients_(wrap_or_move(std::forward<V_T>(v)))
, localEvaluator_()
{}
const Basis& basis() const {
return *basis_;
}
const V& dofs() const {
return *coefficients_;
}
const LEV& localEvaluator() const {
return localEvaluator_;
}
Range operator() (const Domain& x) const {
DUNE_THROW(NotImplemented,"not implemented");
}
friend typename Traits::DerivativeInterface derivative(const VintageBasisGridFunction& t) {
DUNE_THROW(NotImplemented,"not implemented");
}
friend LocalFunction localFunction(const VintageBasisGridFunction& t) {
return LocalFunction(t);
}
const EntitySet& entitySet() const {
return entitySet_;
}
private:
EntitySet entitySet_;
std::shared_ptr<const Basis> basis_;
std::shared_ptr<const V> coefficients_;
LocalEvaluator localEvaluator_;
};
/**
* \brief Construction of local functions from a temporary VintageBasisGridFunction (forbidden)
*
* Since a VintageBasisGridFunction::LocalFunction stores a reference
* to the global VintageBasisGridFunction its life time is bound to
* the latter. Hence construction from a temporary VintageBasisGridFunction
* would lead to a dangling reference and is thus forbidden/deleted.
*/
template<typename... TT>
void localFunction(VintageBasisGridFunction<TT...>&& t) = delete;
} // end namespace Impl
/**
* \brief Create a VintageBasisGridFunction
*
* This allows to use a dune-fufem basis as dune-functions
* grid function. Moreover you can customize the local function
* evaluation.
*
* \tparam R Range type of function
* \tparam B Type of dune-fufem basis
* \tparam V Type of coefficient vector
* \tparam LEV Type of local evaluator
*
* \param basis Dune-fufem basis to interpret the coefficient vector
* \param v Coefficient vector of FE function wrt the given basis
* \param localEvaluator A function doing the local evaluation.
*
* The actual local function implementation is forwarded to the
* callback localEvaluator whose signature
* localEvaluator(basis, coefficients, element, localX, y)
* mimics the LocalEvaluator::evaluateLocal() function.
* By implementing a custom callback you can e.g.
* represent partial derivatives.
*/
template<typename R, typename B, typename V, typename LEV>
auto makeVintageBasisGridFunction(B&& basis, V&& v, LEV&& localEvaluator)
{
using Basis = std::decay_t<B>;
using Vector = std::decay_t<V>;
using LocalEvaluator = std::decay_t<LEV>;
return Impl::VintageBasisGridFunction<Basis, Vector, LocalEvaluator, R>(std::forward<B>(basis), std::forward<V>(v), std::forward<LEV>(localEvaluator));
}
/**
* \brief Create a VintageBasisGridFunction
*
* This allows to use a dune-fufem basis as dune-functions
* grid function.
*
* \tparam R Range type of function
* \tparam B Type of dune-fufem basis
* \tparam V Type of coefficient vector
*
* \param basis Dune-fufem basis to interpret the coefficient vector
* \param v Coefficient vector of FE function wrt the given basis
*/
template<typename R, typename B, typename V>
auto makeVintageBasisGridFunction(B&& basis, V&& v)
{
using Basis = std::decay_t<B>;
using Vector = std::decay_t<V>;
auto localEvaluator = [](const auto& b, const auto& v, const auto&e, const auto& x, auto& y) {
Dune::Fufem::Impl::LocalEvaluator<Basis, Vector, R>::evaluateLocal(b, v, e, x, y);
};
using LocalEvaluator = decltype(localEvaluator);
return Impl::VintageBasisGridFunction<Basis, Vector, LocalEvaluator, R>(std::forward<B>(basis), std::forward<V>(v), std::move(localEvaluator));
}
} // end namespace Dune::Fufem
#endif // DUNE_FUFEM_FUNCTIONS_VINTAGEBASISGRIDFUNCTION_HH
......@@ -31,7 +31,6 @@ dune_add_test(SOURCES symmetricmatrixtest.cc)
dune_add_test(SOURCES tensortest.cc)
dune_add_test(SOURCES test-polyhedral-minimisation.cc)
dune_add_test(SOURCES transferoperatorassemblertest.cc)
dune_add_test(SOURCES vintagebasisgridfunctiontest.cc)
# PYTHONLIBS_FOUND is just placed for backward compatibility with 2.7 Core tests
# and can be removed once tests against 2.7 Core are disabled
......
// -*- tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 2 -*-
// vi: set et ts=4 sw=2 sts=2:
#include <config.h>
#include <iostream>
#include <dune/common/exceptions.hh>
#include <dune/common/parallel/mpihelper.hh>
#include <dune/grid/yaspgrid.hh>
#include <dune/functions/functionspacebases/lagrangebasis.hh>
#include <dune/functions/functionspacebases/interpolate.hh>
#include <dune/functions/gridfunctions/discreteglobalbasisfunction.hh>
#include <dune/functions/gridfunctions/test/gridfunctiontest.hh>
#include <dune/fufem/functions/vintagebasisgridfunction.hh>
#include <dune/fufem/functionspacebases/dunefunctionsbasis.hh>
using namespace Dune;
using namespace Dune::Functions;
using namespace Dune::Functions::Test;
int main (int argc, char* argv[]) try
{
using namespace Dune;
using namespace Dune::Functions;
using namespace Dune::Functions::BasisFactory;
Dune::MPIHelper::instance(argc, argv);
bool passed = true;
// Generate grid for testing
const int dim = 2;
using Grid = Dune::YaspGrid<dim>;
auto l = Dune::FieldVector<double,dim>(1);
auto elements = std::array<int,dim>{{10, 10}};
auto grid = Grid(l,elements);
grid.globalRefine(2);
using GridView = Grid::LeafGridView;
using Domain = GridView::template Codim<0>::Geometry::GlobalCoordinate;
auto gridView = grid.leafGridView();
{
auto duneFunctionsBasis = makeBasis(gridView, lagrange<3>());
auto duneFufemBasis = DuneFunctionsBasis<decltype(duneFunctionsBasis)>(duneFunctionsBasis);
using Range = double;
using Vector = std::vector<Range>;
auto f = [](const Domain& x){
return x[0];
};
double exactIntegral = 0.5;
// Construct some coefficient vector by interpolation
auto x = Vector{};
Dune::Functions::interpolate(duneFunctionsBasis, x, f);
// Create and use VintageBasisGridFunction
auto gf = Dune::Fufem::makeVintageBasisGridFunction<Range>(duneFufemBasis, x);
auto y = Vector{};
Dune::Functions::interpolate(duneFunctionsBasis, y, gf);
passed = passed and (x==y);
passed = passed and Dune::Functions::Test::checkGridViewFunction(gridView, gf, exactIntegral);
}
return passed ? 0: 1;
} catch ( Dune::Exception &e )
{
std::cerr << "Dune reported error: " << e << std::endl;
return 1;
}
catch(...)
{
std::cerr << "Unknown exception thrown!" << std::endl;
return 1;
}
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment