Skip to content
Snippets Groups Projects
Commit 6e9f47b6 authored by Max Kahnt's avatar Max Kahnt
Browse files

Remove necessity to specify fallback QuadRuleKey in SubgridL2FunctionalAssembler.

Provide two separate constructors instead.
parent 6c5c5563
Branches
No related tags found
No related merge requests found
...@@ -50,22 +50,41 @@ class SubgridL2FunctionalAssembler : ...@@ -50,22 +50,41 @@ class SubgridL2FunctionalAssembler :
* It can assemble functionals on the subgrid given by grid * It can assemble functionals on the subgrid given by grid
* functions on the underlying hostgrid exactly. * functions on the underlying hostgrid exactly.
* *
* The QuadratureRuleKey given here does only specify what is * \param f the basisgrid function representing the functional
* needed to integrate f if f is not a BasisGridFunction on
* the hostgrid. Otherwise a quadrature rule is determined
* automatically.
*
* \param f the (hostgrid) function representing the functional
* \param grid the subgrid (!) * \param grid the subgrid (!)
* \param fQuadKey A fallback QuadratureRuleKey that specifies how to integrate f if its not a BasisGridFunction on the hostgrid
*/ */
SubgridL2FunctionalAssembler(const HostGridFunction& f, const GridType& grid, const QuadratureRuleKey& fQuadKey) : template<typename B, typename C, std::enable_if<std::is_base_of<HostGridFunction, BasisGridFunction<B,C>>::value_type>>
functionQuadKey_(fQuadKey), SubgridL2FunctionalAssembler(const BasisGridFunction<B, C>& f, const GridType& grid) :
functionQuadKey_(0),
grid_(grid), grid_(grid),
f_(f), f_(f),
hostBasisGridFunctionInfo_(dynamic_cast<const HostBasisGridFunctionInfo*>(&f)) hostBasisGridFunctionInfo_(&f)
{} {}
/** \brief constructor
*
* Creates a local functional assembler for an L2-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. It can be used to override the
* automatically determined QuadratureRule if f is a BasisGridFunction.
*
* \param f the hostgrid function representing the functional
* \param grid the subgrid (!)
* \param fQuadKey A QuadratureRuleKey that specifies how to integrate f
*/
SubgridL2FunctionalAssembler(const HostGridFunction& f, const GridType& grid, const QuadratureRuleKey& fQuadKey) :
functionQuadKey_(&fQuadKey),
grid_(grid),
f_(f),
hostBasisGridFunctionInfo_(0)
{
if(dynamic_cast<const HostBasisGridFunctionInfo*>(&f))
Dune::dwarn << "Overriding HostBasis Quadrature rule." << std::endl;
}
/** \copydoc L2FunctionalAssembler::assemble /** \copydoc L2FunctionalAssembler::assemble
*/ */
void assemble(const Element& element, LocalVector& localVector, const TrialLocalFE& tFE) const void assemble(const Element& element, LocalVector& localVector, const TrialLocalFE& tFE) const
...@@ -90,7 +109,7 @@ class SubgridL2FunctionalAssembler : ...@@ -90,7 +109,7 @@ class SubgridL2FunctionalAssembler :
if (hostBasisGridFunctionInfo_) if (hostBasisGridFunctionInfo_)
quadKey = quadKey.product(hostBasisGridFunctionInfo_->quadratureRuleKey(*baseHostElement)); quadKey = quadKey.product(hostBasisGridFunctionInfo_->quadratureRuleKey(*baseHostElement));
else else
quadKey = quadKey.product(functionQuadKey_); quadKey = quadKey.product(*functionQuadKey_);
const Dune::template QuadratureRule<double, dim>& quad = QuadratureRuleCache<double, dim>::rule(quadKey); const Dune::template QuadratureRule<double, dim>& quad = QuadratureRuleCache<double, dim>::rule(quadKey);
// loop over quadrature points // loop over quadrature points
...@@ -136,7 +155,7 @@ class SubgridL2FunctionalAssembler : ...@@ -136,7 +155,7 @@ class SubgridL2FunctionalAssembler :
if (hostBasisGridFunctionInfo_) if (hostBasisGridFunctionInfo_)
quadKey = quadKey.product(hostBasisGridFunctionInfo_->quadratureRuleKey(hostElement)); quadKey = quadKey.product(hostBasisGridFunctionInfo_->quadratureRuleKey(hostElement));
else else
quadKey = quadKey.product(functionQuadKey_); quadKey = quadKey.product(*functionQuadKey_);
const Dune::template QuadratureRule<double, dim>& quad = QuadratureRuleCache<double, dim>::rule(quadKey); const Dune::template QuadratureRule<double, dim>& quad = QuadratureRuleCache<double, dim>::rule(quadKey);
...@@ -171,7 +190,7 @@ class SubgridL2FunctionalAssembler : ...@@ -171,7 +190,7 @@ class SubgridL2FunctionalAssembler :
} }
private: private:
const QuadratureRuleKey functionQuadKey_; const QuadratureRuleKey* functionQuadKey_;
const GridType& grid_; const GridType& grid_;
const HostGridFunction& f_; const HostGridFunction& f_;
const HostBasisGridFunctionInfo* hostBasisGridFunctionInfo_; const HostBasisGridFunctionInfo* hostBasisGridFunctionInfo_;
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment