[petsc-users] Near nullspace lost in fieldsplit

[Jeremy Theler (External)]

> > Because of a combination of settings, our code passes through this line:
> >
> > https://gitlab.com/petsc/petsc/-/blob/main/src/ksp/pc/impls/fieldsplit/fieldsplit.c?ref_type=heads#L692
> >
> > i.e. the matrices associated with each of the sub-KSPs of a fieldsplit are destroyed and then re-created later.
> > The thing is that one of these destroyed matrices had a near nullspace attached, which is lost because the new matrix does > not have it anymore.
> >
> > Is this a bug or are we missing something?
>
> I just want to get a clear picture. You create a PCFIELDSPLIT, set it up, then pull out the matrices and attach a nullspace before the solve.

We need to solve an SNES. We use dmplex so we have the jacobian allocated before starting the solve.
At setup time we

 1. define the PC of the KSP of the SNES to be fieldsplit
 2. define the fields with ISes
 3. call PCSetup() to create the sub-KSPs
 4. retrieve the matrix attached to the sub-KSP that needs the near nullspace and attach it to that matrix

> At a later time, you start another solve with this PC, and it has the DIFFERENT_NONZERO_PATTERN flag, so it recreates these matrices and loses your attached nullspace.

At a later time, in the jacobian evaluation we populate the global matrix (i.e. not the matrices attached to each sub-KSPs) and then we set dirichlet bcs with MatZeroRowsColumns() on that same global matrix.
For some reason, in serial the near nullspace is not lost but in parallel the call to MatZeroRowsColumns() does change the non-zero structure (even though the manual says it does not) and then the code goes through that line 692 in fieldsplit.c and the near nullspace is lost.

> First, does the matrix really change?

Well, the matrix during setup is not filled in, just allocated.
The thing is that if we set MAT_KEEP_NONZERO_PATTERN to true with MatSetOption() before setting the dirichlet BCs, then the near nullspace is not lost (because the code does not go through line 692 of fieldsplit.c).

So there are (at least) two issues:

 1. Code going through line 692 looses the near nullspace of the matrices attached to the sub-KSPs
 2. The call to MatZeroRowsColumns() changes then non-zero structure for MPIAIJ but not for SEQAIJ


> Second, I had the same problem. I added callbacks to DM which allow a nullspace to be automatically attached if you extract a certain subfield. Are you using a DM?


Yes. Can you give us an example?

Regards
--
jeremy



-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.mcs.anl.gov/pipermail/petsc-users/attachments/20240209/6cffd1b5/attachment.html>