4 #ifndef DUNE_GRID_IO_FILE_VTK_FUNCTION_HH 5 #define DUNE_GRID_IO_FILE_VTK_FUNCTION_HH 9 #include <dune/common/exceptions.hh> 10 #include <dune/common/fvector.hh> 12 #include <dune/geometry/type.hh> 13 #include <dune/geometry/referenceelements.hh> 14 #include <dune/geometry/multilineargeometry.hh> 37 template<
class Gr
idView >
43 typedef typename GridView::template Codim< 0 >::Entity
Entity;
47 virtual int ncomps ()
const = 0;
57 virtual double evaluate (
int comp,
const Entity& e,
58 const Dune::FieldVector<ctype,dim>& xi)
const = 0;
61 virtual std::string
name ()
const = 0;
87 template<
typename GV,
typename V>
121 const Dune::FieldVector<ctype,dim>&)
const 123 return v[mapper.
index(e)*ncomps_+mycomp_];
127 virtual std::string
name ()
const 150 int ncomps=1,
int mycomp=0 )
157 if (v.size()!=(
unsigned int)(mapper.
size()*ncomps_))
158 DUNE_THROW(IOError,
"P0VTKFunction: size mismatch");
185 template<
typename GV,
typename V>
218 virtual double evaluate (
int comp,
const Entity& e,
219 const Dune::FieldVector<ctype,dim>& xi)
const 222 const unsigned int nVertices = e.subEntities(dim);
224 std::vector<FieldVector<ctype,1> > cornerValues(nVertices);
225 for (
unsigned i=0; i<nVertices; ++i)
226 cornerValues[i] = v[mapper.
subIndex(e,i,dim)*ncomps_+mycomp_];
229 const MultiLinearGeometry<ctype,dim,1> interpolation(e.type(), cornerValues);
230 return interpolation.global(xi);
234 virtual std::string
name ()
const 257 int ncomps=1,
int mycomp=0 )
264 if (v.size()!=(
unsigned int)(mapper.
size()*ncomps_))
265 DUNE_THROW(IOError,
"P1VTKFunction: size mismatch");
276 #endif // DUNE_GRID_IO_FILE_VTK_FUNCTION_HH Base::Entity Entity
Definition: function.hh:207
virtual int ncomps() const
return number of components
Definition: function.hh:114
Base::Entity Entity
Definition: function.hh:109
GridView::template Codim< 0 >::Entity Entity
Definition: function.hh:43
virtual double evaluate(int comp, const Entity &e, const Dune::FieldVector< ctype, dim > &xi) const =0
evaluate single component comp in the entity e at local coordinates xi
virtual std::string name() const =0
get name
P0VTKFunction(const GV &gv, const V &v_, const std::string &s_, int ncomps=1, int mycomp=0)
construct from a vector and a name
Definition: function.hh:149
GridView::ctype ctype
Definition: function.hh:41
Definition: function.hh:42
virtual double evaluate(int, const Entity &e, const Dune::FieldVector< ctype, dim > &) const
evaluate
Definition: function.hh:120
virtual int ncomps() const
return number of components
Definition: function.hh:212
P1VTKFunction(const GV &gv, const V &v_, const std::string &s_, int ncomps=1, int mycomp=0)
construct from a vector and a name
Definition: function.hh:256
virtual ~VTKFunction()
virtual destructor
Definition: function.hh:64
Take a vector and interpret it as point data for the VTKWriter.
Definition: function.hh:186
Base::ctype ctype
Definition: function.hh:208
virtual int ncomps() const =0
A base class for grid functions with any return type and dimension.
Definition: function.hh:38
virtual std::string name() const
get name
Definition: function.hh:127
Dimensionality of the reference element of the entity.
Definition: common/entity.hh:112
Grid::ctype ctype
type used for coordinates in grid
Definition: common/gridview.hh:127
virtual ~P0VTKFunction()
destructor
Definition: function.hh:162
virtual ~P1VTKFunction()
destructor
Definition: function.hh:269
Index index(const EntityType &e) const
Map entity to array index.
Definition: mcmgmapper.hh:162
virtual std::string name() const
get name
Definition: function.hh:234
Include standard header files.
Definition: agrid.hh:59
Mapper for multiple codim and multiple geometry types.
Take a vector and interpret it as cell data for the VTKWriter.
Definition: function.hh:88
Base::ctype ctype
Definition: function.hh:110
The dimension of the grid.
Definition: common/gridview.hh:130
int size() const
Return total number of entities in the entity set managed by the mapper.
Definition: mcmgmapper.hh:205
Index subIndex(const typename GV::template Codim< 0 >::Entity &e, int i, unsigned int codim) const
Map subentity of codim 0 entity to array index.
Definition: mcmgmapper.hh:190
virtual double evaluate(int comp, const Entity &e, const Dune::FieldVector< ctype, dim > &xi) const
evaluate
Definition: function.hh:218