[petsc-dev] PC fieldsplit / schur issues / questions

[Jed Brown]

Kai Germaschewski <kai.germaschewski at unh.edu> writes:

> Hi,
>
> we've been looking at using the built-in Schur complement capabilities in
> PETSc, and encountered some issues which I think are bugs, though maybe I
> just can't figure out how to use things properly.
>
> Attached is a somewhat boiled down test program. It's based on
> Crank-Nicolson TS, which uses SNES, which uses KSP, which has a PC that's
> set to type fieldsplit. We'd like to have the KSP solver for the actual
> Schur complement to use multigrid as PC, eventually.
>
> One, kinda minor issue: when ComputeJacobian is called, we have to use
> MatZeroEntries and recompute the Jacobian (of course, that's not really
> needed for the simple problem we're solving right now because it's actually
> linear). The reason for that is that, I think, the TS CN has shifted /
> added a diagonal to the Jacobian matrix, as needed for CN. For this, I'm
> just wondering if that's expected behavior (it's not obvious that one
> wouldn't be passed the same matrix that was calculated the previous time).

There is a conflict between modifying the matrix in-place to save space
and returning a pristine matrix to the user.  You can use
MatStoreValues() and MatRetrieveValues() to manage saving and restoring
the pristine part of your matrix.  If you have suggestions for how to
make this more transparent, I'm interested in hearing them.

> Secondly, running the code with "-ts_view_pre" crashes, because the TS
> isn't fully set up (in particular, the fieldsplit PC doesn't seem to have
> properly created its sub-KSPs yet). I'm pretty sure it shouldn't crash
> because of that. But maybe more importantly, the way things are set up
> doesn't seem quite right. In particular, even after calling TSSetUp(), the
> sub objects in the TS aren't properly set up yet -- that only seems to
> happen in the first call to TSSolve() (or TSStep()). Trying to set up the
> SNES inside the TS explicitly by calling SNESSetUp() doesn't work, either.
> During TSStep, SNESSolve() is called with explicit vectors for b, x, it
> appears only then TSSetUp() will actually work without errors.

PCFieldSplit (and several other preconditioners) need the matrices in
order to create their inner solvers.  With SNES (and TS), those matrices
are not available until the first iteration, unless you use
SNESComputeJacobian (or TSComputeIJacobian).  I hate this, but it's hard
to make a method choose "smart" defaults if it has to configure itself
without any data.  I'd love to have _all_ configuration take place in
XXSetFromOptions(), but if we did that, many methods would frequently be
configuring themselves incorrectly.

It should not crash (thanks for reporting), but it's not going to give
you much useful information.

> I'm not sure when the fieldsplit PC should be set up to have its proper
> sub-KSP, I suppose at least at TSSetUp() time, if not earlier. For us,
> since we want to change the sub-KSPs, we need some way to do that before
> the actual solve, but at that point, PCFieldSplitGetSubKSP() will fail --
> this part of the code is #if 0'd out right now.

I suppose PCFieldSplit could be made to create the MatSchurComplement
up-front, even though it is empty.  Right now, there is an option
processed in PCSetUp_FieldSplit to decide whether to use the same KSP
inside the Schur complement as in the outer factors (the
-fieldsplit_%s_inner_ option).  This option can't be processed without
knowing the field names, which are (automatically) determined in
PCFieldSplitSetDefaults (called by PCSetUp_FieldSplit).  You are
providing index sets so enough information is available, but in the
general case, the matrix is required to extract that information.

> I believe that all of this works properly when using SNES directly, but we
> can't get it to work inside the TS, so any help would be appreciated.

At this point, I think the easiest thing is to take one micro time step,
then pull the solvers apart to reconfigure as you wish.  Or configure
using run-time options, which should always work.
-------------- next part --------------
A non-text attachment was scrubbed...
Name: not available
Type: application/pgp-signature
Size: 835 bytes
Desc: not available
URL: <http://lists.mcs.anl.gov/pipermail/petsc-dev/attachments/20131106/89157ff1/attachment.sig>