vertexdata.hh 2.08 KB
Newer Older
Peter Bastian's avatar
Peter Bastian committed
1 2
// -*- tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 2 -*-
// vi: set et ts=4 sw=2 sts=2:
3 4 5
#ifndef __DUNE_GRID_HOWTO_VERTEXDATA_HH__
#define __DUNE_GRID_HOWTO_VERTEXDATA_HH__

6 7
#include <dune/grid/common/mcmgmapper.hh>
#include <dune/grid/io/file/vtk/vtkwriter.hh>
Peter Bastian's avatar
Peter Bastian committed
8
#if HAVE_GRAPE
9
#include <dune/grid/io/visual/grapedatadisplay.hh>
Peter Bastian's avatar
Peter Bastian committed
10 11 12 13 14 15
#endif

// demonstrate attaching data to elements
template<class G, class F>
void vertexdata (const G& grid, const F& f)
{
16
  // get dimension from Grid
Peter Bastian's avatar
Peter Bastian committed
17
  const int dim = G::dimension;
18
  typedef typename G::LeafGridView GridView;
19
  // determine type of LeafIterator for codimension = dimension
20 21 22
  typedef typename GridView::template Codim<dim>::Iterator VertexLeafIterator;

  // get grid view on the leaf part
23
  GridView gridView = grid.leafGridView();
Peter Bastian's avatar
Peter Bastian committed
24 25

  // make a mapper for codim 0 entities in the leaf grid
26 27
  Dune::LeafMultipleCodimMultipleGeomTypeMapper<G>
  mapper(grid, Dune::mcmgVertexLayout());
Peter Bastian's avatar
Peter Bastian committed
28 29 30 31

  // allocate a vector for the data
  std::vector<double> c(mapper.size());

32
  // iterate through all entities of codim 0 at the leaves
33 34
  for (VertexLeafIterator it = gridView.template begin<dim>();
       it!=gridView.template end<dim>(); ++it)
Peter Bastian's avatar
Peter Bastian committed
35 36
  {
    // evaluate functor and store value
37
    c[mapper.index(*it)] = f(it->geometry().corner(0));
Peter Bastian's avatar
Peter Bastian committed
38 39 40
  }

  // generate a VTK file
41
  Dune::VTKWriter<typename G::LeafGridView> vtkwriter(grid.leafGridView());
42
  vtkwriter.addVertexData(c,"data");
43
  vtkwriter.write( "vertexdata", Dune::VTK::appendedraw );
Peter Bastian's avatar
Peter Bastian committed
44 45 46

  // online visualization with Grape
#if HAVE_GRAPE
47 48 49 50 51 52 53 54
  {
    const int polynomialOrder = 1; // we piecewise linear data
    const int dimRange = 1; // we have scalar data here
    // create instance of data display
    Dune::GrapeDataDisplay<G> grape(grid);
    // display data
    grape.displayVector("concentration", // name of data that appears in grape
                        c,  // data vector
55
                        gridView.indexSet(), // used index set
56 57 58
                        polynomialOrder, // polynomial order of data
                        dimRange); // dimRange of data
  }
Peter Bastian's avatar
Peter Bastian committed
59 60
#endif
}
61
#endif // __DUNE_GRID_HOWTO_VERTEXDATA_HH__