function.hh

Go to the documentation of this file.

// SPDX-FileCopyrightText: Copyright © DUNE Project contributors, see file LICENSE.md in module root
// SPDX-License-Identifier: LicenseRef-GPL-2.0-only-with-DUNE-exception
// -*- tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 2 -*-
// vi: set et ts=4 sw=2 sts=2:

#ifndef DUNE_GRID_IO_FILE_VTK_FUNCTION_HH
#define DUNE_GRID_IO_FILE_VTK_FUNCTION_HH

#include <string>

#include <dune/common/exceptions.hh>
#include <dune/common/fvector.hh>

#include <dune/geometry/type.hh>
#include <dune/geometry/referenceelements.hh>
#include <dune/geometry/multilineargeometry.hh>

#include <dune/grid/common/mcmgmapper.hh>
#include <dune/grid/io/file/vtk/common.hh>

namespace Dune
{

  //
  //  Base VTKFunction
  //

  template< class GridView >
  class VTKFunction
  {
  public:
    typedef typename GridView::ctype ctype;
    constexpr static int dim = GridView::dimension;
    typedef typename GridView::template Codim< 0 >::Entity Entity;

    virtual int ncomps () const = 0;


    virtual double evaluate (int comp, const Entity& e,
                             const Dune::FieldVector<ctype,dim>& xi) const = 0;

    virtual std::string name () const = 0;

    virtual VTK::Precision precision() const
    { return VTK::Precision::float32; }

    virtual ~VTKFunction () {}
  };

  //
  //  P0VTKFunction
  //


  template<typename GV, typename V>
  class P0VTKFunction
    : public VTKFunction< GV >
  {
    typedef VTKFunction< GV > Base;
    typedef MultipleCodimMultipleGeomTypeMapper<GV> Mapper;

    const V& v;
    std::string s;
    int ncomps_;
    int mycomp_;
    VTK::Precision prec_;
    Mapper mapper;

  public:
    typedef typename Base::Entity Entity;
    typedef typename Base::ctype ctype;
    using Base::dim;

    int ncomps () const override
    {
      return 1;
    }

    double evaluate (int, const Entity& e,
                     const Dune::FieldVector<ctype,dim>&) const override
    {
      return v[mapper.index(e)*ncomps_+mycomp_];
    }

    std::string name () const override
    {
      return s;
    }

    VTK::Precision precision() const override
    {
      return prec_;
    }


    P0VTKFunction(const GV &gv, const V &v_, const std::string &s_,
                  int ncomps=1, int mycomp=0, VTK::Precision prec = VTK::Precision::float32)
      : v( v_ ),
        s( s_ ),
        ncomps_(ncomps),
        mycomp_(mycomp),
        prec_(prec),
        mapper( gv, mcmgElementLayout() )
    {
      if (v.size()!=(unsigned int)(mapper.size()*ncomps_))
        DUNE_THROW(IOError, "P0VTKFunction: size mismatch");
    }

    virtual ~P0VTKFunction() {}
  };

  //
  //  P1VTKFunction
  //


  template<typename GV, typename V>
  class P1VTKFunction
    : public VTKFunction< GV >
  {
    typedef VTKFunction< GV > Base;
    typedef MultipleCodimMultipleGeomTypeMapper<GV> Mapper;

    const V& v;
    std::string s;
    int ncomps_;
    int mycomp_;
    VTK::Precision prec_;
    Mapper mapper;

  public:
    typedef typename Base::Entity Entity;
    typedef typename Base::ctype ctype;
    using Base::dim;

    int ncomps () const override
    {
      return 1;
    }

    double evaluate ([[maybe_unused]] int comp, const Entity& e,
                     const Dune::FieldVector<ctype,dim>& xi) const override
    {
      const unsigned int myDim = Entity::mydimension;
      const unsigned int nVertices = e.subEntities(dim);

      std::vector<FieldVector<ctype,1> > cornerValues(nVertices);
      for (unsigned i=0; i<nVertices; ++i)
        cornerValues[i] = v[mapper.subIndex(e,i,myDim)*ncomps_+mycomp_];

      // (Ab)use the MultiLinearGeometry class to do multi-linear interpolation between scalars
      const MultiLinearGeometry<ctype,dim,1> interpolation(e.type(), cornerValues);
      return interpolation.global(xi);
    }

    std::string name () const override
    {
      return s;
    }

    VTK::Precision precision() const override
    {
      return prec_;
    }


    P1VTKFunction(const GV& gv, const V &v_, const std::string &s_,
                  int ncomps=1, int mycomp=0, VTK::Precision prec = VTK::Precision::float32)
      : v( v_ ),
        s( s_ ),
        ncomps_(ncomps),
        mycomp_(mycomp),
        prec_(prec),
        mapper( gv, mcmgVertexLayout() )
    {
      if (v.size()!=(unsigned int)(mapper.size()*ncomps_))
        DUNE_THROW(IOError,"P1VTKFunction: size mismatch");
    }

    virtual ~P1VTKFunction() {}
  };


} // namespace Dune

#endif // DUNE_GRID_IO_FILE_VTK_FUNCTION_HH