[petsc-users] Using PCREDISTRIBUTE together with PCFIELDSPLIT

[Matthew Knepley]

On Mon, Jan 23, 2023 at 9:18 AM Jonas Lundgren via petsc-users <
petsc-users at mcs.anl.gov> wrote:

> Hi!
>
>
>
> (Sorry for a long message, I have tried to cover the essentials only. I am
> happy to provide further details and logs if necessary, but I have tried to
> keep it as short as possible.)
>
>
>
> I am trying to solve a Stokes flow problem with PCFIELDSPLIT as
> preconditioner (and KSPBCGS as solver). I have successfully set up the
> preconditioner and solved several examples by applying the following
> options when setting up the solver:
>
>
>
> // Preliminaries
>
> KSP ksp;
>
> PC pc;
>
> KSPCreate(PETSC_COMM_WORLD, ksp);
>
> KSPSetType(ksp, KSPBCGS);
>
> KSPSetFromOptions(ksp);  // here, “-pc_type redistribute” is read from
> command line
>
> KSPSetOperators(ksp, K, K);
>
> KSPGetPC(ksp, &pc);
>
>
>
> // Preconditioner
>
> PCSetType(pc, PCFIELDSPLIT);
>
> PCFieldSplitSetType(pc, PC_COMPOSITE_SCHUR);
>
> PCFieldSplitSetSchurFactType(pc, PC_FIELDSPLIT_SCHUR_FACT_LOWER);
>
> PCFieldSplitSetSchurPre(pc, PC_FIELDSPLIT_SCHUR_PRE_SELFP, K);
>

This seems wrong. K is your system matrix, but the matrix here should be
the size of the Schur complement, and
furthermore it is overwritten with the Schur complement preconditioner. I
think you should be passing NULL.

  Thanks,

    Matt


> PCFieldSplitSetIS(pc,"0",isu);
>
>
>
> where “isu” in the last row is an IS containing all flow indices (no
> pressure indices), and “K” is my system matrix, which is created using
> DMCreateMatrix(dm_state, &K); and “dm_state” is a DMStag object. (However,
> K is not assembled yet.)
>
>
>
> I want to try to use PCREDISTRIBUTE on top of this, since I believe that I
> have a lot of degrees-of-freedom (DOF) locked, and if I can reduce the size
> of my linear problem, the solution time can decrease as well. My idea was
> to use PCREDISTRIBUTE as the main preconditioner (with KSPPREONLY as the
> solver), and to move the KSPBCGS and PCFIELDSPLIT down one level, to act as
> a sub-KSP and sub-PC, by introducing the following between the two code
> blocks:
>
> PC ipc;
>
> PCRedistributeGetKSP(pc, &iksp);
>
> KSPGetPC(iksp, &ipc);
>
> and letting “ipc” replace “pc” in the second code block (“Preconditioner”).
>
>
>
> The two last rows of in the “Preconditioner” block refers to “K” and
> “isu”, which are defined in terms of the entire system, not the reduced. I
> believe that I have to replace these two variables with a Mat (“A”) that
> has locked DOFs removed, and an IS (“isfree”) that contains only free flow
> DOF indices (no locked flow indices and no pressure indices).
>
>
>
> And here comes my problem: the Mat “A” and the IS “isfree” will only be
> available when the original KSP (“ksp”) – or possibly the PC “pc” – are set
> up (using KSPSetUp(ksp); or PCSetUp(pc); ). I have no way of knowing how
> “A” or “isfree” looks before it is automatically extracted from “K” and
> “isu” during set up – they will be different depending on the specifics of
> the example I am running.
>
>
>
> I think I can solve the matrix problem by inserting the following (instead
> of the second to last row above):
>
> PCFieldSplitSetSchurPre(ipc,PC_FIELDSPLIT_SCHUR_PRE_SELFP, ipc->mat);
>
> Even though “ipc->mat” is not created yet, I believe that this might lead
> to the use of the extracted sub-flow-matrix as (a part of the)
> preconditioner for the flow block in PCFIELDSPLIT. I do not know for sure,
> since I have other issues that have prevented me from testing this
> approach. Does this sound reasonable?
>
>
>
> (One might come up with a better alternative than the system matrix to use
> for this sub-PC, but that work is yet to be done…)
>
>
>
> Furthermore, the second issue I have is regarding the IS “isfree”. After
> the “Preconditioner” code block I am trying to set the Fieldsplit sub-PCs
> by calling PCFieldSplitSchurGetSubKSP(); It appears as if I need to set up
> the outer KSP “ksp” or PC “pc” before I call PCFieldSplitSchurGetSubKSP();
> or otherwise I get an error message looking like:
>
> [0]PETSC ERROR: Object is in wrong state
>
> [0]PETSC ERROR: Must call KSPSetUp() or PCSetUp() before calling
> PCFieldSplitSchurGetSubKSP()
>
>
>
> However, when I try to call KSPSetUp() as suggested above (which is
> possible thanks to assembling the system matrix before the call to set up
> the KSP), I just get another error message indicating that I have the wrong
> IS (“isu”) (most likely because the indices “isu” refers to the entire
> system, not the reduced):
>
> [1]PETSC ERROR: Argument out of range
>
> [1]PETSC ERROR: Index 0's value 3864 is larger than maximum given 2334
>
> [1]PETSC ERROR: #1 ISComplement() at
> /mnt/c/Users/jonlu89/install/petsc-v3-18-1/src/vec/is/is/utils/iscoloring.c:804
>
> [1]PETSC ERROR: #2 PCFieldSplitSetDefaults() at
> /mnt/c/Users/jonlu89/install/petsc-v3-18-1/src/ksp/pc/impls/fieldsplit/fieldsplit.c:544
>
> [1]PETSC ERROR: #3 PCSetUp_FieldSplit() at
> /mnt/c/Users/jonlu89/install/petsc-v3-18-1/src/ksp/pc/impls/fieldsplit/fieldsplit.c:588
>
> [1]PETSC ERROR: #4 PCSetUp() at
> /mnt/c/Users/jonlu89/install/petsc-v3-18-1/src/ksp/pc/interface/precon.c:994
>
> [1]PETSC ERROR: #5 KSPSetUp() at
> /mnt/c/Users/jonlu89/install/petsc-v3-18-1/src/ksp/ksp/interface/itfunc.c:406
>
> [1]PETSC ERROR: #6 PCSetUp_Redistribute() at
> /mnt/c/Users/jonlu89/install/petsc-v3-18-1/src/ksp/pc/impls/redistribute/redistribute.c:240
>
> [1]PETSC ERROR: #7 PCSetUp() at
> /mnt/c/Users/jonlu89/install/petsc-v3-18-1/src/ksp/pc/interface/precon.c:994
>
> [1]PETSC ERROR: #8 KSPSetUp() at
> /mnt/c/Users/jonlu89/install/petsc-v3-18-1/src/ksp/ksp/interface/itfunc.c:406
>
>
>
> Okay, so now I want to obtain the reduced IS “isfree” that contains only
> free flow DOFs (not locked flow DOFs or pressure DOFs). One option is to
> try to obtain it from the same IS that extracts the sub-matrix when setting
> up the PCREDISTRIBUTE preconditioner, see red->is in MatCreateSubMatrix()
> in
> https://petsc.org/main/src/ksp/pc/impls/redistribute/redistribute.c.html
>
>
>
> And annoyingly, I need to have called PCSetUp(pc) or KSPSetUp(ksp) before
> using red->is (thus, before calling PCFieldSplitSetIS() ), which is
> something that is impossible it seems:
>
> [0]PETSC ERROR: Petsc has generated inconsistent data
>
> [0]PETSC ERROR: Unhandled case, must have at least two fields, not 1
>
> [0]PETSC ERROR: #1 PCFieldSplitSetDefaults() at
> /mnt/c/Users/jonlu89/install/petsc-v3-18-1/src/ksp/pc/impls/fieldsplit/fieldsplit.c:550
>
> [0]PETSC ERROR: #2 PCSetUp_FieldSplit() at
> /mnt/c/Users/jonlu89/install/petsc-v3-18-1/src/ksp/pc/impls/fieldsplit/fieldsplit.c:588
>
> [0]PETSC ERROR: #3 PCSetUp() at
> /mnt/c/Users/jonlu89/install/petsc-v3-18-1/src/ksp/pc/interface/precon.c:994
>
> [0]PETSC ERROR: #4 KSPSetUp() at
> /mnt/c/Users/jonlu89/install/petsc-v3-18-1/src/ksp/ksp/interface/itfunc.c:406
>
> [0]PETSC ERROR: #5 PCSetUp_Redistribute() at
> /mnt/c/Users/jonlu89/install/petsc-v3-18-1/src/ksp/pc/impls/redistribute/redistribute.c:240
>
> [0]PETSC ERROR: #6 PCSetUp() at
> /mnt/c/Users/jonlu89/install/petsc-v3-18-1/src/ksp/pc/interface/precon.c:994
>
> [0]PETSC ERROR: #7 KSPSetUp() at
> /mnt/c/Users/jonlu89/install/petsc-v3-18-1/src/ksp/ksp/interface/itfunc.c:406
>
>
>
> So, to summarize, I seem to be stuck in a Catch 22; I cannot set up the
> KSP (“ksp”) before I have the correct IS (“isfree”) set to my Fieldsplit
> precondtitioner, but I cannot obtain the correct IS preconditioner
> (“isfree”) before I set up the (outer) KSP (“ksp”).
>
>
>
> Is there any possible ideas on how to resolve this issue? Do I have to
> hard code the IS “isfree” in order to achieve what I want? (I don’t know
> how…) Or maybe this approach is wrong all together. Can I achieve my goal
> of using PCFIELDSPLIT and PCREDISTRIBUTE together in another way?
>
>
>
> Is this even a sane way of looking at this problem?
>
>
>
> Best regards,
>
> Jonas Lundgren
>


-- 
What most experimenters take for granted before they begin their
experiments is infinitely more interesting than any results to which their
experiments lead.
-- Norbert Wiener

https://www.cse.buffalo.edu/~knepley/ <http://www.cse.buffalo.edu/~knepley/>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.mcs.anl.gov/pipermail/petsc-users/attachments/20230123/147417f9/attachment.html>