Comment 6 for bug 1088175

Revision history for this message
Anders Logg (logg) wrote : Re: [Bug 1088175] [NEW] LinearOperator does not work in parallel

On Sun, Dec 09, 2012 at 04:02:03PM +0100, Anders Logg wrote:
> On Sun, Dec 09, 2012 at 02:32:19PM -0000, Garth Wells wrote:
> > On Sun, Dec 9, 2012 at 2:01 PM, Anders Logg <email address hidden> wrote:
> > > On Sun, Dec 09, 2012 at 01:27:50PM -0000, Garth Wells wrote:
> > >> On Sun, Dec 9, 2012 at 1:20 PM, Anders Logg <email address hidden> wrote:
> > >> > I can take a look but object to the description. It is indeed written
> > >> > to work in parallel. It has also been tested to work in parallel
> > >> > before.
> > >> >
> > >>
> > >> Tested where?
> > >
> > > I thought the unit test was running fine on the buildbots but I see
> > > now it was never added to the list of tests (in the main test.py).
> > >
> >
> > The unit test is a special case for which it may run on 3 processes
> > because the number of dofs is divisible by 3. In general, the parallel
> > layouts will not match.
>
> It doesn't even work for that case. So what happened is that I thought
> the test was run on the buildbot and therefore worked in parallel. I
> see now why it doesn't work but it should be relatively easy to fix.

I've made an attempt to fix this but it still crashes. I'm
initializing the local range for the mat shell matrix using the same
local range as for the solution vector:

  // Get local range
  std::size_t m_local = M;
  std::size_t n_local = N;
  if (MPI::num_processes() > 1)
  {
    std::pair<std::size_t, std::size_t> local_range = x.local_range();
    m_local = local_range.first;
    n_local = local_range.second;
  }

  // Initialize PETSc matrix
  A.reset(new Mat, PETScMatrixDeleter());
  MatCreateShell(PETSC_COMM_WORLD, m_local, n_local, M, N, (void*) this, A.get());
  MatShellSetOperation(*A, MATOP_MULT, (void (*)()) usermult);

What am I missing? Does the solution vector not have the same parallel
layout as the right-hand side b (which gets multiplied by the matrix
in the Krylov solver)?

--
Anders