grid_function_writer.hh 4.81 KB
Newer Older
1 2
#ifndef DUNE_DORIE_GRID_FUNCTION_WRITER_HH
#define DUNE_DORIE_GRID_FUNCTION_WRITER_HH
3

4 5 6 7 8 9
#include <string>
#include <memory>
#include <type_traits>

#include <dune/common/exceptions.hh>
#include <dune/common/parametertree.hh>
Santiago Ospina's avatar
Santiago Ospina committed
10 11

#include <dune/grid/io/file/vtk/vtksequencewriter.hh>
12
#include <dune/grid/io/file/vtk/subsamplingvtkwriter.hh>
Santiago Ospina's avatar
Santiago Ospina committed
13

14 15
#include <dune/pdelab/gridfunctionspace/gridfunctionspaceutilities.hh>

16 17 18
namespace Dune{
namespace Dorie{

19 20 21 22 23 24 25 26 27
/**
 * @brief      Converts PDELab grid function to match a dune-function.
 * @note       This adapter differs from the one provided in PDELab in that this
 *             one has the current interface of dune-grid to write vtk files.
 *             The main difference using this adapter is that it allows to write
 *             more than 3 components.
 *
 * @tparam     GF    PDELab grid function type.
 */
28 29 30 31 32 33 34 35 36 37
template<class GF>
class VTKGridFunctionAdapter 
{
  using GridFunction = typename std::decay<GF>::type;
  using GridView = typename GF::Traits::GridViewType;
  using Entity = typename GridView::template Codim<0>::Entity;
  using Range = typename GF::Traits::RangeType;
  using Domain = typename GF::Traits::DomainType;
public:

38 39 40 41 42
  /**
   * @brief      Constructs the object
   *
   * @param[in]  gf    pointer to a constant grid function
   */
43 44 45 46
  VTKGridFunctionAdapter(std::shared_ptr<const GF> gf)
    : _gf(gf)
  {}

47 48 49 50 51
  /**
   * @brief      Binds an entity to the grid function
   *
   * @param[in]  e     Entity
   */
52 53 54 55 56
  void bind(const Entity& e)
  {
    _entity = &e;
  }

57 58 59
  /**
   * @brief      Unbinds the previously binded entity
   */
60 61 62 63 64
  void unbind()
  {
    _entity = nullptr;
  }

65 66 67 68 69
  /**
   * @brief      Evaluates the grid function at a given coordinate
   *
   * @param[in]  x     Local coordinate with respect to the last binded entity
   */
70 71 72 73 74 75 76 77 78 79 80 81 82
  Range operator()(const Domain& x) const
  {
    Range y;
    _gf->evaluate(*_entity,x,y);
    return y;
  }

private:
  const Entity* _entity;
  std::shared_ptr<const GF> _gf;
};


Santiago Ospina's avatar
Santiago Ospina committed
83
/*-------------------------------------------------------------------------*//**
84
 * @brief      Class for grid fucntion vtk sequence writer. This class works
Santiago Ospina's avatar
Santiago Ospina committed
85 86 87 88 89
 *             exacly the same way as a VTKSequence writer but it receives
 *             directly smart pointers to grid functions
 *
 * @tparam     GV    GridView of the dune grid
 */
90
template<class GV>
Santiago Ospina's avatar
Santiago Ospina committed
91
class GridFunctionVTKSequenceWriter : public Dune::VTKSequenceWriter<GV>
92 93
{
public:
94 95 96
  // export constructors of vtk sequence writer
  using Dune::VTKSequenceWriter<GV>::VTKSequenceWriter;

Santiago Ospina's avatar
Santiago Ospina committed
97 98 99 100 101 102 103 104
  /*-----------------------------------------------------------------------*//**
   * @brief      Adds a vertex data.
   *
   * @param[in]  gf_ptr  A pointer to a grid function of the type GF
   * @param[in]  name    Name of the variable
   *
   * @tparam     GF      Type of the grid function
   */
105
  template<class GF>
106
  void addVertexData(std::shared_ptr<const GF> gf_ptr, std::string name)
107 108
  {
    static_assert(std::is_same<typename GF::Traits::GridViewType,GV>::value, 
Santiago Ospina's avatar
Santiago Ospina committed
109
      "GridFunctionVTKSequenceWriter only works with only one GridView type.");
110

111 112 113 114 115 116
    Dune::Dorie::VTKGridFunctionAdapter<GF> vtk_gf(gf_ptr);

    const int vtk_ncomps = GF::Traits::dimRange;
    const int dim = GF::Traits::dimDomain;
    auto vtk_type = (vtk_ncomps == dim) ? VTK::FieldInfo::Type::vector 
                                        : VTK::FieldInfo::Type::scalar;
117
    auto vtk_info = VTK::FieldInfo(name,vtk_type,vtk_ncomps);
118 119

    this->vtkWriter()->addVertexData(vtk_gf,vtk_info);
120 121
  }

Santiago Ospina's avatar
Santiago Ospina committed
122 123 124 125 126 127 128 129
  /*-----------------------------------------------------------------------*//**
   * @brief      Adds a cell data.
   *
   * @param[in]  gf_ptr  A pointer to a grid function of the type GF
   * @param[in]  name    Name of the variable
   *
   * @tparam     GF      Type of the grid function
   */
130
  template<class GF>
131
  void addCellData(std::shared_ptr<const GF> gf_ptr, std::string name)
132 133
  {
    static_assert(std::is_same<typename GF::Traits::GridViewType,GV>::value, 
Santiago Ospina's avatar
Santiago Ospina committed
134
      "GridFunctionVTKSequenceWriter only works with only one GridView type.");
135
    
136 137 138 139 140 141
    Dune::Dorie::VTKGridFunctionAdapter<GF> vtk_gf(gf_ptr);

    const int vtk_ncomps = GF::Traits::dimRange;
    const int dim = GF::Traits::dimDomain;
    auto vtk_type = (vtk_ncomps == dim) ? VTK::FieldInfo::Type::vector 
                                        : VTK::FieldInfo::Type::scalar;
142
    auto vtk_info = VTK::FieldInfo(name,vtk_type,vtk_ncomps);
143 144

    this->vtkWriter()->addCellData(vtk_gf,vtk_info);
145
  }
146

147
  /*-----------------------------------------------------------------------*//**
148 149 150
   * @brief      Clear all the attached VTKFunctions.
   * @details    This is necessary for each iteration if VTKFunctions are
   *             managed internaly with pointers instead of references
151
   */
152 153 154 155
  void clear()
  {
    this->vtkWriter()->clear();
  }
156 157 158 159 160
};

} // namespace Dorie
} // namespace Dune

161
#endif // DUNE_DORIE_GRID_FUNCTION_WRITER_HH