Virtual Refinement

Please have a look on the non-virtual Refinement documentation.

General

Refinement can only be used when you know the geometryType of your entities at compile time. You could circumvent this by using a switch(geometryType), but each case would look very much the same. If you have many such switch() statements, or each case contains lots of code, or you simply have many possible geometryTypes, this can be quiet annoying.

VirtualRefinement does all of this switch() statements for you. It defines a common virtual base class per dimension, and derives one class for each geometryType and coerceTo from that class. The derived classes simply wrap the non-virtual classes from Refinement. This makes it possible to treat each geometryType (of a given dimension) the same, and thus eliminates the many repetitions of lots of code.

But the case statements are not totally gone yet. VirtualRefinement does these statements once and for all by wrapping them into the buildRefinement() function.

The user Interface

The VirtualRefinement class

VirtualRefinement is not a set of unrelated specialisations of the same template class. VirtualRefinement is a base class with several virtual methods, which are overloaded by the concrete VirtualRefinement implementation classes. Each implementation class wraps one of the non-virtual Refinement classes.

The user interface is modelled closely after the Refinement interface. The main differences are:

• VirtualRefinement is not a static class, but a singleton. Thus each VirtualRefinement implementation has to be instantiated before use. This is done with the template function buildRefinement (see below).
• Since the methods of VirtualRefinement are virtual (or use virtual methods themself) they have to be called like

  refinementInstace.nElements(Dune::refinementIntervals(n));

  instead of

  RefinementTypedef::nElements(Dune::refinementIntervals(n));

• IndexVector is a std::vector instead of a FieldVector since the number of corners of different geometry types may be different at runtime. The user is responsible to always pass the same coerceTo parameter to buildRefinement() so he always gets the same number of corners.

template<int dimension>
class VirtualRefinement
{
public:
  template<int Codimension>
  struct Codim {
    class SubEntityIterator;
  };
  typedef VertexIterator; // These are aliases for Codim<codim>::SubEntityIterator
  typedef ElementIterator;
 
  typedef IndexVector; // This is a std::vector
  typedef CoordVector; // This is a FieldVector
 
  virtual int nVertices(Dune::RefinementIntervals intervals) const;
  VertexIterator vBegin(Dune::RefinementIntervals intervals) const;
  VertexIterator vEnd(Dune::RefinementIntervals intervals) const;
  virtual int nElements(Dune::RefinementIntervals intervals) const;
  ElementIterator eBegin(Dune::RefinementIntervals intervals) const;
  ElementIterator eEnd(Dune::RefinementIntervals intervals) const;
};

The iterators have the same interface as the Refinement iterators except that IndexVector is a std::vector instead of a FieldVector (see above). Also the restriction that the Iterators are not dereferenceable applies.

template<int dimension>
class VertexIterator
{
public:
  typedef VirtualRefinement<dimension> Refinement;
 
  int index() const;
  Refinement::CoordVector coords() const;
};
 
template<int dimension>
class ElementIterator
{
public:
  typedef VirtualRefinement<dimension> Refinement;
 
  int index() const;
  // Coords of the center of mass of the element
  Refinement::CoordVector coords() const;
  Refinement::IndexVector vertexIndices() const;
};

buildRefinement()

The declaration for buildRefinement is

template<int dimension, class CoordType>
VirtualRefinement<dimension, CoordType> &buildRefinement(GeometryType geometryType, GeometryType coerceTo);

It is expected that you know the dimension and the coordinate type of the elements you want to refine at compile time.

The simple case is that you want to refine, say, quadrilaterals and the subentities should look like quadrilaterals as well. In that case you would call buildRefinement() like

VirtualRefinement<2, CoordType> &refinement = buildRefinement<2, CoordType>(quadrilateral, quadrilateral);

The more complicated case is that your entity is a quadrilateral, but the subentities should look like triangles. In this case call buildRefinement() like

VirtualRefinement<2, CoordType> &refinement = buildRefinement<2, CoordType>(quadrilateral, triangle);

Summary: geometryType is the geometry type of the entity you want to refine, while coerceTo is the geometry type of the subentities.

Implementing a new Refinement type

When you write a Refinement implementation for a new combination of geometryType and coerceTo, you have to tell buildRefinement() about it.

• First, you have to implement the non-virtual part in Refinement, if you have not done so yet.
• Second, visit the end of refinementvirtual.cc, and look for the specialisations of template<int dimension, class CoordType> class RefinementBuilder. There is one specialisation for each dimension, containing the single method build().
• The build() contains two levels of switch statements, the outer for geomentryType and the inner for coerceTo. Each case will either return the correct VirtualRefinement or fall through to the end of the method and throw an error. Insert the cases for your refinement.

Everything else has been done for you automatically.

Namespaces

VirtualRefinement does not use a complicated namespace scheme like Refinement. The complete VirtualRefinement stuff simply lives directly in namespace Dune.

Conceptual layers

VirtualRefinement adds to more layers to the ones already defined in Refinement:

• Layer 3 makes it easy to use several Refinement implementations in the same code, when you only know at run-time, which Refinement implementation you need. It wraps class Refinement and its iterators into a Proxy class, retaining its interface but all deriving from a virtual base class VirtualRefinement<dimension, CoordType>. This is located in refinementvirtual.cc.
• Layer 4 defines function buildRefinement(geometryType, coerceTo), which returns the right refinement for a runtime-determined GeometryType. It is also located in refinementvirtual.cc

Implementation

The interface is defined by the template class VirtualRefinement. It simply defines the CoordVectors and IndexVectors appropriate for this dimension and CoordType, defines which iterators to use, and provides some proxy or pure virtual functions.

For each class Refinement<geometryType, CoordType, coercTo, dim> we provide a class VirtualRefinementImp<geometryType, CoordType, coercTo, dim>, which wraps the matching class Refinement<geometryType, CoordType, coercTo, dim> and derives from the matching base class VirtualRefinement<dimension, CoordType>. Each VirtualRefinementImp is a singleton and has a static instance() method which will return this instance as a reference to the base class VirtualRefinement. All this is done in a single template class.

The iterators

We can't do the same thing with the iterators as we do with class VirtualRefinement. Since they are polymorph we cannot simply pass them into user code. They are not singletons, so we also cannot pass references to them. Passing pointers to iterators would work, but then the programmer has to remember to explicitly delete them. Also, it is uncommon for iterators to be handled by their pointers.

What we do instead is having a wrapper class which conforms to the iterator interface and is the same for all VirtualRefinementIterators of a given dimension. This class contains a pointer to a polymorph backend object implementing the iterator. The various VirtualRefinementImps then derive from the abstract backend class and pass a pointer to a concrete backend object when instantiating an iterator.

buildRefinement()

The template function buildRefinement() has to be specialized for each dimension. It makes no sense to test for geometryType.isPrism() when dimension==2. But this way we run into a limitation of C++: we can't do partial function specialisation.

The workaround is to create a class RefinementBuilder with a lone static method build() and to call that from buildRefinement(). Since RefinementBuilder is a class and not a function we can do partial specialisations.

It is probably possible to automatically generate the switch statements with linked lists of template structs, functions implementing the cases, and a recursive template function that will iterate over the list, but it is probably not worth the effort, as long as buildRefinement() is enough for the job.