27 #ifndef VTK_TENSOR_FUNCTION_HH
28 #define VTK_TENSOR_FUNCTION_HH
32 #include <dune/grid/io/file/vtk/function.hh>
33 #include <dune/common/fvector.hh>
34 #include <dune/common/version.hh>
36 #include <opm/common/Exceptions.hpp>
37 #include <opm/common/ErrorMacros.hpp>
48 template <
class Gr
idView,
class Mapper>
51 enum { dim = GridView::dimension };
52 typedef typename GridView::ctype ctype;
53 typedef typename GridView::template Codim<0>::Entity Element;
55 typedef BaseOutputWriter::TensorBuffer TensorBuffer;
59 const GridView& gridView,
61 const TensorBuffer& buf,
63 unsigned matrixColumnIdx)
69 , matrixColumnIdx_(matrixColumnIdx)
70 { assert(
int(buf_.size()) == mapper_.size()); }
72 virtual std::string name()
const
75 virtual int ncomps()
const
76 {
return static_cast<int>(buf_[0].M()); }
78 virtual double evaluate(
int mycomp,
80 const Dune::FieldVector<ctype, dim>& xi)
const
85 idx =
static_cast<size_t>(mapper_.index(e));
87 else if (codim_ == dim) {
92 Dune::GeometryType gt = e.type();
93 int n =
static_cast<int>(e.subEntities(dim));
94 for (
int i = 0; i < n; ++i) {
95 Dune::FieldVector<ctype, dim> local =
96 Dune::ReferenceElements<ctype, dim>::general(gt).position(i, dim);
99 if (local.infinity_norm() < min) {
100 min = local.infinity_norm();
106 idx =
static_cast<size_t>(mapper_.subIndex(e, imin, codim_));
109 OPM_THROW(std::logic_error,
110 "Only element and vertex based tensor fields are supported so far.");
112 unsigned i =
static_cast<unsigned>(mycomp);
113 unsigned j =
static_cast<unsigned>(matrixColumnIdx_);
115 return static_cast<double>(
static_cast<float>(buf_[idx][i][j]));
119 const std::string name_;
120 const GridView gridView_;
121 const Mapper& mapper_;
122 const TensorBuffer& buf_;
124 unsigned matrixColumnIdx_;
The base class for all output writers.
Provides a tensor-valued function using Dune::FieldMatrix objects as elements.
Definition: vtktensorfunction.hh:49