Skip to content
Snippets Groups Projects
Commit 099b0ff7 authored by Carsten Gräser's avatar Carsten Gräser
Browse files

Switch function interface in `SubgridH1FunctionalAssembler`

* Switch to dune-functions interface for grid functions
* Modernize constructor: You now have to pass a `QuadratureRuleKey`
  describing how to integrade the passed function instead of
  a raw quadrature rule order for the whole integrant.
* Reduce magic:
  * This no longer works like a simple `H1FunctionalAssembler`
    if you pass Subgrid `GridFunction`s. In this use case
    use a plain `H1FunctionalAssembler` instead.

After this patch you can no longer pass a `Dune::VirtualFunction` or `VirtualGridFunction`,
but have to use either the `std::function` or the `Dune::Functions::GridViewFunction` interface.
parent 243a9db5
No related branches found
No related tags found
1 merge request!146Stop using Dune::VirtualFunction interface
......@@ -5,17 +5,12 @@
#include <dune/common/fvector.hh>
#include <dune/common/function.hh>
#include <dune/grid/common/rangegenerators.hh>
#include <dune/istl/bvector.hh>
#include <dune/fufem/functions/virtualgridfunction.hh>
#include <dune/fufem/functions/virtualdifferentiablefunction.hh>
#include <dune/fufem/functions/basisgridfunction.hh>
#include <dune/fufem/quadraturerules/quadraturerulecache.hh>
#include <dune/fufem/functionspacebases/refinedp1nodalbasis.hh>
#include <dune/functions/common/defaultderivativetraits.hh>
#include <dune/fufem/functions/gridfunctionhelper.hh>
#include <dune/fufem/quadraturerules/quadraturerulecache.hh>
#include <dune/fufem/assemblers/localfunctionalassembler.hh>
#include <dune/fufem/assemblers/localassemblers/h1functionalassembler.hh>
......@@ -37,11 +32,9 @@ class SubgridH1FunctionalAssembler :
typedef typename GridType::template Codim<0>::Geometry::GlobalCoordinate GlobalCoordinate;
typedef typename GridType::template Codim<0>::Geometry::LocalCoordinate LocalCoordinate;
typedef typename DerivativeTypefier<GlobalCoordinate,T>::DerivativeType FRangeType;
typedef VirtualGridFunction<GridType, FRangeType> GridFunction;
typedef VirtualGridFunction<HostGrid, FRangeType> HostGridFunction;
typedef BasisGridFunctionInfo<HostGrid> HostBasisGridFunctionInfo;
using FRangeType = typename Dune::Functions::DefaultDerivativeTraits<T(GlobalCoordinate)>::Range;
using LocalHostFunction = Dune::Fufem::Impl::LocalFunctionInterfaceForGrid<HostGrid, FRangeType>;
public:
typedef typename LocalFunctionalAssembler<GridType,TrialLocalFE,T>::Element Element;
......@@ -50,33 +43,31 @@ class SubgridH1FunctionalAssembler :
typedef typename GridType::HostGridType::template Codim<0>::Entity::Geometry HostGeometry;
typedef typename Dune::VirtualFunction<GlobalCoordinate, FRangeType> Function;
/** \brief constructor
*
* Creates a local functional assembler for an H1-functional.
* It can assemble functionals on the subgrid given by grid
* functions on the underlying hostgrid exactly.
*
* The QuadratureRuleKey given here does specify what is
* needed to integrate f.
*
* \param f the (hostgrid) function representing the functional (if assembling e.g. \f$(\nabla u,\nabla v)\f$ you need to provide \f$f=\nabla u\f$)
* \param grid the subgrid (!) type
* \param order the quadrature order (DEFAULT 2)
* \param grid the subgrid (!)
* \param fQuadKey A QuadratureRuleKey that specifies how to integrate f
*/
SubgridH1FunctionalAssembler(const Function& f, const GridType& grid, int order=2) :
H1FunctionalAssembler<GridType,TrialLocalFE,T>(f, order),
f_(f),
template<class FT>
SubgridH1FunctionalAssembler(const FT& f, const GridType& grid, const QuadratureRuleKey& fQuadKey) :
H1FunctionalAssembler<GridType,TrialLocalFE, T>(f, fQuadKey),
grid_(grid),
order_(order),
fAsHostGridFunction_(dynamic_cast<const HostGridFunction*>(&f)),
hostBasisGridFunctionInfo_(dynamic_cast<const HostBasisGridFunctionInfo*>(&f))
localHostF_(Dune::Fufem::Impl::localFunctionForGrid<HostGrid, FRangeType>(f)),
functionQuadKey_(fQuadKey)
{}
/** \copydoc H1FunctionalAssembler::assemble
*/
void assemble(const Element& element, LocalVector& localVector, const TrialLocalFE& tFE) const
{
if (not(fAsHostGridFunction_))
{
H1FunctionalAssembler<GridType,TrialLocalFE,T>::assemble(element, localVector, tFE);
return;
}
typedef typename TrialLocalFE::Traits::LocalBasisType::Traits::JacobianType LocalBasisJacobianType;
// Make sure we got suitable shape functions
......@@ -89,18 +80,19 @@ class SubgridH1FunctionalAssembler :
const auto hostelement = grid_.template getHostEntity<0>(element);
bool trialSpaceIsRefined = IsRefinedLocalFiniteElement<TrialLocalFE>::value(tFE);
// store gradients of shape functions and base functions
std::vector<LocalBasisJacobianType> referenceGradients(tFE.localBasis().size());
std::vector<GlobalCoordinate> gradients(tFE.localBasis().size());
if (hostelement.isLeaf())
{
bool hostGFIsRefined = hostBasisGridFunctionInfo_ and hostBasisGridFunctionInfo_->isRefinedLocalFiniteElement(hostelement);
// get quadrature rule
const auto& quad = QuadratureRuleCache<double, dim>::rule(element.type(), order_, (trialSpaceIsRefined or hostGFIsRefined) );
QuadratureRuleKey quadKey(tFE);
quadKey.setGeometryType(element.type());
quadKey = quadKey.derivative().product(functionQuadKey_);
const auto& quad = QuadratureRuleCache<double, dim>::rule(quadKey);
localHostF_.bind(hostelement);
// loop over quadrature points
for (size_t pt=0; pt < quad.size(); ++pt)
......@@ -122,8 +114,7 @@ class SubgridH1FunctionalAssembler :
invJacobian.mv(referenceGradients[i][0], gradients[i]);
// compute values of function
FRangeType f_pos;
fAsHostGridFunction_->evaluateLocal(hostelement, quadPos, f_pos);
auto f_pos = localHostF_(quadPos);
// and vector entries
for (size_t i=0; i<gradients.size(); ++i)
......@@ -142,10 +133,13 @@ class SubgridH1FunctionalAssembler :
{
const HostGeometry hostGeometry = descElement.geometry();
bool hostGFIsRefined = hostBasisGridFunctionInfo_ and hostBasisGridFunctionInfo_->isRefinedLocalFiniteElement(descElement);
// get quadrature rule
const auto& quad = QuadratureRuleCache<double, dim>::rule(descElement.type(), order_, hostGFIsRefined);
QuadratureRuleKey quadKey(tFE);
quadKey.setGeometryType(descElement.type());
quadKey = quadKey.product(functionQuadKey_);
const auto& quad = QuadratureRuleCache<double, dim>::rule(quadKey);
localHostF_.bind(descElement);
// loop over quadrature points
for (size_t pt=0; pt < quad.size(); ++pt)
......@@ -168,8 +162,7 @@ class SubgridH1FunctionalAssembler :
invJacobian.mv(referenceGradients[i][0], gradients[i]);
// compute values of function
FRangeType f_pos;
fAsHostGridFunction_->evaluateLocal(descElement, quadPos, f_pos);
auto f_pos = localHostF_(quadPos);
// and vector entries
for (size_t i=0; i<gradients.size(); ++i)
......@@ -187,11 +180,9 @@ class SubgridH1FunctionalAssembler :
}
private:
const Function& f_;
const GridType& grid_;
const int order_;
const HostGridFunction* fAsHostGridFunction_;
const HostBasisGridFunctionInfo* hostBasisGridFunctionInfo_;
mutable LocalHostFunction localHostF_;
const QuadratureRuleKey functionQuadKey_;
};
#endif
......
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