Skip to content

#658 Cannot get GeometryGrid with DiscreteCoordFunction to work

Metadata

Property Value
Reported by Oliver Sander (oliver.sander@tu-dresden.de)
Reported at Nov 11, 2009 15:23
Type Bug Report
Version 1.2
Operating System Unspecified / All
Last edited by Martin Nolte (nolte@mathematik.uni-freiburg.de)
Last edited at Feb 13, 2011 21:22

Description

I am not sure whether this is a bug or whether I am not understanding the GeometryGrid documentation right. I am trying to deform a grid with a vector given the new position for each vertex. In the documentation it says that the deformation function must derive from either AnalyticalCoordFunction or from DiscreteCoordFunction. Well, my function is not analytic, so I try the second one. No example code is given but I try my best and produce the following deformation function

template <class GridView>
class DeformationFunction
    : public Dune :: DiscreteCoordFunction< double, dim, DeformationFunction<GridView> >
{
    typedef DeformationFunction<GridView> This;
    typedef Dune :: DiscreteCoordFunction< double, dim, This > Base;

  public:

    DeformationFunction(const GridView& gridView,
                        const double* deformedPosition)
        : gridView_(gridView),
          deformedPosition_(deformedPosition)
    {}

    template< class HostEntity >
    void evaluate ( const HostEntity &hostEntity, unsigned int corner,
                    Dune::FieldVector<double,dim> &y ) const
    {

        const typename GridView::IndexSet& indexSet = gridView_.indexSet();

        int idx = indexSet.subIndex(hostEntity, corner,dim);

        for (int i=0; i<dim; i++)
            y[i] = deformedPosition_[idx*dim + i];
    }

private:

    GridView gridView_;

    const double* deformedPosition_;

};

Unfortunately, the code won't compile if I ask for a vertex position directly

std::cout << deformedGrid.leafbegin<dim>()->geometry().corner(0);

because then HostEntity apparently is a Codim entity for which subIndex is not declared. On the other hand, if I only query elements I can compile this, and the corner parameter really takes on different values. That means the corner parameter is not superfluous and I cannot use index instead of subIndex. I tried to read and understand the code, but it is completely undocumented...

Full testcase attached, thank you for your help.

Attachments