3.1-git
DUNE for Multi-{Phase, Component, Scale, Physics, ...} flow and transport in porous media
vtknestedfunction.hh
Go to the documentation of this file.
1// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
2// vi: set et ts=4 sw=4 sts=4:
3/*****************************************************************************
4 * See the file COPYING for full copying permissions. *
5 * *
6 * This program is free software: you can redistribute it and/or modify *
7 * it under the terms of the GNU General Public License as published by *
8 * the Free Software Foundation, either version 3 of the License, or *
9 * (at your option) any later version. *
10 * *
11 * This program is distributed in the hope that it will be useful, *
12 * but WITHOUT ANY WARRANTY; without even the implied warranty of *
13 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the *
14 * GNU General Public License for more details. *
15 * *
16 * You should have received a copy of the GNU General Public License *
17 * along with this program. If not, see <http://www.gnu.org/licenses/>. *
18 *****************************************************************************/
25#ifndef VTK_NESTED_FUNCTION_HH
26#define VTK_NESTED_FUNCTION_HH
27
28#include <string>
29
30#include <dune/common/fvector.hh>
31#include <dune/istl/bvector.hh>
32#include <dune/grid/io/file/vtk/function.hh>
33
34namespace Dumux {
35
41template <class GridView, class Mapper, class Buffer>
42class VtkNestedFunction : public Dune::VTKFunction<GridView>
43{
44 enum { dim = GridView::dimension };
45 using ctype = typename GridView::ctype;
46 using Element = typename GridView::template Codim<0>::Entity;
47 using ReferenceElements = Dune::ReferenceElements<ctype, dim>;
48
49public:
51 const GridView &gridView,
52 const Mapper &mapper,
53 const Buffer &buf,
54 int codim,
55 int numComp)
56 : name_(name)
57 , gridView_(gridView)
58 , mapper_(mapper)
59 , buf_(buf)
60 , codim_(codim)
61 , numComp_(numComp)
62 {
63 assert(std::size_t(buf_.size()) == std::size_t(mapper_.size()));
64 }
65
66 virtual std::string name () const
67 { return name_; }
68
69 virtual int ncomps() const
70 { return numComp_; }
71
72 virtual double evaluate(int mycomp,
73 const Element &element,
74 const Dune::FieldVector< ctype, dim > &xi) const
75 {
76 int idx;
77 if (codim_ == 0) {
78 // cells. map element to the index
79 idx = mapper_.index(element);
80 }
81 else if (codim_ == dim) {
82 // find vertex which is closest to xi in local
83 // coordinates. This code is based on Dune::P1VTKFunction
84 double min=1e100;
85 int imin=-1;
86 Dune::GeometryType geomType = element.type();
87 int n = element.subEntities(dim);
88
89 for (int i=0; i < n; ++i)
90 {
91 Dune::FieldVector<ctype,dim> local =
92 ReferenceElements::general(geomType).position(i,dim);
93 local -= xi;
94 if (local.infinity_norm()<min)
95 {
96 min = local.infinity_norm();
97 imin = i;
98 }
99 }
100
101 // map vertex to an index
102 idx = mapper_.subIndex(element, imin, codim_);
103 }
104 else
105 DUNE_THROW(Dune::InvalidStateException,
106 "Only element and vertex based vector "
107 " fields are supported so far.");
108
109 double val = buf_[idx][mycomp];
110 using std::abs;
111 if (abs(val) < std::numeric_limits<float>::min())
112 val = 0;
113
114 return val;
115 }
116
117private:
118 const std::string name_;
119 const GridView gridView_;
120 const Mapper &mapper_;
121 const Buffer &buf_;
122 int codim_;
123 int numComp_;
124};
125
126} // end namespace Dumux
127
128#endif
make the local view function available whenever we use the grid geometry
Definition: adapt.hh:29
Provides a vector-valued function using Dune::FieldVectors as elements.
Definition: vtknestedfunction.hh:43
virtual int ncomps() const
Definition: vtknestedfunction.hh:69
virtual double evaluate(int mycomp, const Element &element, const Dune::FieldVector< ctype, dim > &xi) const
Definition: vtknestedfunction.hh:72
virtual std::string name() const
Definition: vtknestedfunction.hh:66
VtkNestedFunction(std::string name, const GridView &gridView, const Mapper &mapper, const Buffer &buf, int codim, int numComp)
Definition: vtknestedfunction.hh:50