Setting up blocked data structures

When solving vector-valued PDEs the result at each grid point must be a vector. When doing this the ordering of the global indices can be chosen, e.g. lexicographic or entity blocked. For many solvers, e.g. AMG, the ordering makes a difference in the efficiency of the solver. Further options, such as choosing to aggregate over the blocks and specifying the norm to be used can produce faster results.

First, we have to define one or more scalar grid function spaces to combine

  // Set up scalar grid function space
  typedef Dune::PDELab::GridFunctionSpace<GV, FEM, Constraints, Dune::PDELab::ISTL::VectorBackend<>> SCALAR_GFS;
  SCALAR_GFS U1(grid.leafGridView(),fem); U1.name("U1");
  SCALAR_GFS U2(grid.leafGridView(),fem); U2.name("U2");

There are two ways to define a vector-valued grid function space. If all scalar grid function spaces are the same a power grid function space can be used. In contrast a composite grid function space allows for different scalar grid function spaces. Both can be used with different blocking.

Let's first define a power grid function space with lexiographic ordering

  // Use lexiographical blocked ordering
  typedef Dune::PDELab::LexicographicOrderingTag LexiographicOrderingTag;
 
  // Set up power grid function space
  using VBE = Dune::PDELab::ISTL::VectorBackend<Dune::PDELab::ISTL::Blocking::fixed>; // blocking vector backend
  typedef Dune::PDELab::PowerGridFunctionSpace<SCALAR_GFS,
                                                dim, // block size
                                                VBE, // blocked vector backend
                                                LexiographicOrderingTag> GFS;
  GFS gfs(U1,U2);

And then a composite grid function space with entity blocked ordering

  // Use entity blocked ordering
  typedef Dune::PDELab::EntityBlockedOrderingTag EntityOrderingTag;
 
  // Setting up a composite grid function space with the same scalar grid function space in both components
  using VBE = Dune::PDELab::ISTL::VectorBackend<Dune::PDELab::ISTL::Blocking::fixed>; // blocking vector backend
  typedef Dune::PDELab::CompositeGridFunctionSpace<VBE, // blocked vector backend
                                                    EntityOrderingTag,
                                                    SCALAR_GFS, SCALAR_GFS> GFS;
  GFS gfs(U1,U2);

Full example code: recipe-blocking.cc