Skip to content
Snippets Groups Projects
Commit 5ff69653 authored by Youett, Jonathan's avatar Youett, Jonathan Committed by akbib@FU-BERLIN.DE
Browse files

don't let the class inherit from gridfunction because there is no derivative...

don't let the class inherit from gridfunction because there is no derivative type implemented for functions with values in manifolds

[[Imported from SVN: r7985]]
parent 9138aa46
No related branches found
No related tags found
No related merge requests found
...@@ -8,16 +8,14 @@ ...@@ -8,16 +8,14 @@
#include <dune/gfe/localgeodesicfefunction.hh> #include <dune/gfe/localgeodesicfefunction.hh>
#include <dune/fufem/functions/virtualgridfunction.hh>
/** \brief Global geodesic finite element function. /** \brief Global geodesic finite element function.
* *
* \tparam B - The global basis type. * \tparam B - The global basis type.
* \tparam TargetSpace - The manifold that this functions takes its values in. * \tparam TargetSpace - The manifold that this functions takes its values in.
* \tparam CoefficientType - The coefficient vector type.
*/ */
template<class B, class TargetSpace, class CoefficientType> template<class B, class TargetSpace>
GlobalGeodesicFEFunction : public VirtualGridFunction<typename Basis::GridView::Grid, TargetSpace> { class GlobalGeodesicFEFunction
{
public: public:
typedef B Basis; typedef B Basis;
...@@ -38,15 +36,16 @@ public: ...@@ -38,15 +36,16 @@ public:
//! Create global function by a global basis and the corresponding coefficient vector //! Create global function by a global basis and the corresponding coefficient vector
GlobalGeodesicFEFunction(const Basis& basis, const CoefficientType& coefficients) : GlobalGeodesicFEFunction(const Basis& basis, const std::vector<TargetSpace>& coefficients) :
basis_(basis), basis_(basis),
coefficients_(coefficients) coefficients_(coefficients)
{} {}
/** \brief Evaluate the function at local coordinates. */ /** \brief Evaluate the function at local coordinates. */
void evaluateLocal(const Element& element, const Dune::FieldVector<ctype,gridDim>& local, TargetSpace& out) const void evaluateLocal(const Element& element, const Dune::FieldVector<ctype,gridDim>& local, TargetSpace& out) const
{ {
int numOfBasisFct = basis_.getLocalFiniteElement(element).size(); int numOfBaseFct = basis_.getLocalFiniteElement(element).localBasis().size();
// Extract local coefficients // Extract local coefficients
std::vector<TargetSpace> localCoeff(numOfBaseFct); std::vector<TargetSpace> localCoeff(numOfBaseFct);
...@@ -63,7 +62,7 @@ public: ...@@ -63,7 +62,7 @@ public:
void evaluateDerivativeLocal(const Element& element, const Dune::FieldVector<ctype,gridDim>& local, void evaluateDerivativeLocal(const Element& element, const Dune::FieldVector<ctype,gridDim>& local,
Dune::FieldMatrix<ctype, embeddedDim, gridDim>& out) const Dune::FieldMatrix<ctype, embeddedDim, gridDim>& out) const
{ {
int numOfBasisFct = basis_.getLocalFiniteElement(element).size(); int numOfBaseFct = basis_.getLocalFiniteElement(element).localBasis().size();
// Extract local coefficients // Extract local coefficients
std::vector<TargetSpace> localCoeff(numOfBaseFct); std::vector<TargetSpace> localCoeff(numOfBaseFct);
...@@ -88,6 +87,6 @@ private: ...@@ -88,6 +87,6 @@ private:
//! The global basis //! The global basis
const Basis& basis_; const Basis& basis_;
//! The coefficient vector //! The coefficient vector
const CoefficientType& coefficients_; const std::vector<TargetSpace>& coefficients_;
}; };
#endif #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