vtkscalarfunction.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  This file is part of the Open Porous Media project (OPM).
5 
6  OPM 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 2 of the License, or
9  (at your option) any later version.
10 
11  OPM 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 OPM. If not, see <http://www.gnu.org/licenses/>.
18 
19  Consult the COPYING file in the top-level source directory of this
20  module for the precise wording of the license and the list of
21  copyright holders.
22 */
27 #ifndef VTK_SCALAR_FUNCTION_HH
28 #define VTK_SCALAR_FUNCTION_HH
29 
31 
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>
36 
37 #include <opm/common/Unused.hpp>
38 
39 #include <opm/common/Exceptions.hpp>
40 #include <opm/common/ErrorMacros.hpp>
41 
42 #include <string>
43 #include <limits>
44 #include <vector>
45 
46 namespace Ewoms {
47 
52 template <class GridView, class Mapper>
53 class VtkScalarFunction : public Dune::VTKFunction<GridView>
54 {
55  enum { dim = GridView::dimension };
56  typedef typename GridView::ctype ctype;
57  typedef typename GridView::template Codim<0>::Entity Element;
58 
59  typedef BaseOutputWriter::ScalarBuffer ScalarBuffer;
60 
61 public:
62  VtkScalarFunction(std::string name,
63  const GridView& gridView,
64  const Mapper& mapper,
65  const ScalarBuffer& buf,
66  unsigned codim)
67  : name_(name)
68  , gridView_(gridView)
69  , mapper_(mapper)
70  , buf_(buf)
71  , codim_(codim)
72  { assert(int(buf_.size()) == mapper_.size()); }
73 
74  virtual std::string name() const
75  { return name_; }
76 
77  virtual int ncomps() const
78  { return 1; }
79 
80  virtual double evaluate(int mycomp OPM_UNUSED,
81  const Element& e,
82  const Dune::FieldVector<ctype, dim>& xi) const
83  {
84  unsigned idx;
85  if (codim_ == 0) {
86  // cells. map element to the index
87  idx = static_cast<unsigned>(mapper_.index(e));
88  }
89  else if (codim_ == dim) {
90  // find vertex which is closest to xi in local
91  // coordinates. This code is based on Dune::P1VTKFunction
92  double min = 1e100;
93  int imin = -1;
94  Dune::GeometryType gt = e.type();
95  int n = static_cast<int>(e.subEntities(dim));
96  for (int i = 0; i < n; ++i) {
97  Dune::FieldVector<ctype, dim> local =
98  Dune::ReferenceElements<ctype, dim>::general(gt).position(i, dim);
99 
100  local -= xi;
101  if (local.infinity_norm() < min) {
102  min = local.infinity_norm();
103  imin = static_cast<int>(i);
104  }
105  }
106 
107  // map vertex to an index
108  idx = static_cast<unsigned>(mapper_.subIndex(e, imin, codim_));
109  }
110  else
111  OPM_THROW(std::logic_error, "Only element and vertex based vector "
112  " fields are supported so far.");
113 
114  return static_cast<double>(static_cast<float>(buf_[idx]));
115  }
116 
117 private:
118  const std::string name_;
119  const GridView gridView_;
120  const Mapper& mapper_;
121  const ScalarBuffer& buf_;
122  unsigned codim_;
123 };
124 
125 } // namespace Ewoms
126 
127 #endif
Definition: baseauxiliarymodule.hh:37
Provides a vector-valued function using Dune::FieldVectors as elements.
Definition: vtkscalarfunction.hh:53
The base class for all output writers.