[Dune] Constraints for parallel overlapping CG with SSORk (ISTLBackend_OVLP_CG_SSORk<GFS, C>)

[Eike Mueller]

Dear Dune list,

I'm trying to use the parallel version of the ISTL backend for the CG  
solver preconditioned with SSOR (overlapping version), i.e.  
Dune::PDELab::ISTLBackend_OVLP_CG_SSORk<GFS,C>  and can't get it to  
work. I'm basically following example 04 in the pdelab howto, i.e. my  
local operator implements the finite volume discretisation. I must get  
the constraints C wrong, the code below shows what I tried. It works  
when I run on one processor, but not on more.

// Extract types from gridview
typedef typename GV::Grid::ctype Coord;
typedef double Real;
const int dim = GV::dimension;

// Construct constant finite element space
typedef Dune::PDELab::P0LocalFiniteElementMap<Coord, Real, dim> FEM;
FEM fem(Dune::GeometryType(Dune::GeometryType::cube,dim));

// Construct grid function space
typedef Dune::PDELab::OverlappingConformingDirichletConstraints CON;
typedef Dune::PDELab::ISTLVectorBackend<1> VBE;
typedef Dune::PDELab::GridFunctionSpace<GV, FEM, CON, VBE> GFS;
GFS gfs(gv,fem);

// Construct local operator for finite volume discretisation
// (see example 4 in pdelab-howto)
typedef HelmholtzLocalOperator<aReference,bReference,fRHS> LOP;
LOP lop(aref,
         bref,
         frhs,
         model_parameters);

// Construct grid operator space
typedef Dune::PDELab::ISTLBCRSMatrixBackend<1,1> MBE;
typedef Dune::PDELab::EmptyTransformation CC;
typedef Dune::PDELab::GridOperatorSpace<GFS, GFS, LOP, CC, CC, MBE> GOS;
GOS gos(gfs,gfs,lop);

// Container for solution
typedef typename GFS::template VectorContainer<Real>::Type U;
U u(gfs, 0.0);

// Constraints for ISTL backend
typedef Dune::PDELab::EmptyTransformation C;
C c;

// ISTL backend: Overlapping CG, preconditioned with SSOR
typedef Dune::PDELab::ISTLBackend_OVLP_CG_SSORk<GFS,C> LS;
LS ls(gfs,c, maxiter, smoothsteps, verbose);
// Solve
Dune::PDELab::StationaryLinearProblemSolver<GOS,LS,U>
   slp(gos,u,ls, model_parameters.solver.resreduction);
slp.apply();

I also tried replacing the EmptyTransformation constraints by an  
implementation of the boundary conditions, i.e. I replaced the  
following lines

typedef Dune::PDELab::EmptyTransformation C;
C c;

by this:
  typedef BCType<GV> B;
  B b(gv,model_parameters);
  typedef typename GFS::template ConstraintsContainer<Real>::Type C;
  C c;
  Dune::PDELab::constraints(b,gfs,c);

where BCType is similar to the example in the pdelab howto. Again, it  
doesn't work on more than one processor.

Also, could someone help me clarify the difference between overlap and  
ghosts elements and why there are overlapping and non-overlapping  
versions of the solvers at all?

Thanks a lot,

Eike