27#ifndef VTK_VECTOR_FUNCTION_HH
28#define VTK_VECTOR_FUNCTION_HH
32#include <dune/grid/io/file/vtk/function.hh>
33#include <dune/istl/bvector.hh>
34#include <dune/common/fvector.hh>
35#include <dune/common/version.hh>
48template <
class Gr
idView,
class Mapper>
51 enum { dim = GridView::dimension };
52 using ctype =
typename GridView::ctype;
53 using Element =
typename GridView::template Codim<0>::Entity;
59 const GridView& gridView,
61 const VectorBuffer& buf,
68 { assert(
int(buf_.size()) ==
int(mapper_.size())); }
70 virtual std::string
name()
const
74 {
return static_cast<int>(buf_[0].size()); }
78 const Dune::FieldVector<ctype, dim>& xi)
const
83 idx =
static_cast<unsigned>(mapper_.index(e));
85 else if (codim_ == dim) {
90 Dune::GeometryType gt = e.type();
91 int n =
static_cast<int>(e.subEntities(dim));
92 for (
int i = 0; i < n; ++i) {
93 Dune::FieldVector<ctype, dim> local =
94 Dune::ReferenceElements<ctype, dim>::general(gt).position(i, dim);
97 if (local.infinity_norm() < min) {
98 min = local.infinity_norm();
104 idx =
static_cast<unsigned>(mapper_.subIndex(e, imin, codim_));
107 throw std::logic_error(
"Only element and vertex based vector fields are "
108 "supported so far.");
110 return static_cast<double>(
static_cast<float>(buf_[idx][
static_cast<unsigned>(mycomp)]));
114 const std::string name_;
115 const GridView gridView_;
116 const Mapper& mapper_;
117 const VectorBuffer& buf_;
std::vector< Vector > VectorBuffer
Definition: baseoutputwriter.hh:50
Provides a vector-valued function using Dune::FieldVectors as elements.
Definition: vtkvectorfunction.hh:50
virtual double evaluate(int mycomp, const Element &e, const Dune::FieldVector< ctype, dim > &xi) const
Definition: vtkvectorfunction.hh:76
virtual int ncomps() const
Definition: vtkvectorfunction.hh:73
VtkVectorFunction(std::string name, const GridView &gridView, const Mapper &mapper, const VectorBuffer &buf, unsigned codim)
Definition: vtkvectorfunction.hh:58
virtual std::string name() const
Definition: vtkvectorfunction.hh:70
Definition: blackoilboundaryratevector.hh:37