diff --git a/dune/vtkgridfunction.hh b/dune/vtkgridfunction.hh new file mode 100644 index 0000000000000000000000000000000000000000..f127fd1d50c964287c44383a6e08f9e9d9636e4e --- /dev/null +++ b/dune/vtkgridfunction.hh @@ -0,0 +1,78 @@ +#ifndef VTK_BASIS_GRID_FUNCTION_HH +#define VTK_BASIS_GRID_FUNCTION_HH + +#include <dune/grid/io/file/vtk/function.hh> + +/** \brief A VTK basis grid function. + * + * This function "evaluates" by evaluating the global basis and + *interpolating the function values. + * \tparam Basis The global basis. + * \tparam CoefficientType The vector type for the coefficients. +*/ +template <class Basis, class CoefficientType> +class VTKBasisGridFunction + : public VTKFunction<typename Basis::GridView::Grid> { + + typedef VTKFunction<typename Basis::GridView::Grid> Base; + typedef typename Basis::LocalFiniteElement::Traits::LocalBasisType::Traits:: + RangeType RangeType; + +public: + typedef typename Base::Entity Entity; + typedef typename Base::ctype ctype; + using Base::dim; + + /** \brief Construct from given global basis, coefficient vector and name. + * + * \param basis The global basis. + * \param v A corresponding vector of coefficients. + * \param s A name of the function. + */ + VTKBasisGridFunction(const Basis& basis, const CoefficientType& v, + const std::string& s) + : basis_(basis), coeffs_(v), s(s_) { + if (v.size() != basis_.size()) + DUNE_THROW( + IOError, + "VTKGridFunction: Coefficient vector is not matching the basis"); + } + + /** \brief Get the number of components the function has. */ + virtual int ncomps() const { return CoefficientType::block_type::dimension; } + + /** \brief Locally evaluate a component of the function. + * + * \param comp The component to evaluate. + * \param e The element the local coordinates are taken from. + * \param xi The local coordinates where to evaluate the function. + */ + virtual double evaluate(int comp, const Entity& e, + const Dune::FieldVector<ctype, dim>& xi) const { + // evaluate the local shapefunctions + const typename Basis::LocalFiniteElement& localFE = + basis_.getLocalFiniteElement(e); + std::vector<RangeType> values(localFE.localBasis().size()); + localFE.localBasis().evaluateFunction(xi, values); + + // get the comp'th component + double rt = 0; + for (int i = 0; i < localFE.localBasis().size(); i++) + rt += values[i] * coeffs_[basis_.index(e, i)][comp]; + + return rt; + } + + /** \brief Stupid function to return the name. */ + virtual std::string name() const { return s_; } + + /** \brief Destructor. */ + virtual ~VTKGridFunction() {} + +private: + const Basis& basis_; + const CoefficientType& coeffs_; + std::string s_; +}; + +#endif