[petsc-users] Using PCREDISTRIBUTE together with PCFIELDSPLIT
Matthew Knepley
knepley at gmail.com
Mon Jan 23 11:46:06 CST 2023
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);
>
> 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.
>
Let me reply to the substance of the question. I will restate it to make
sure I understand. You would like to eliminate the diagonal rows/cols from
your operator before the solve, which PCREDISTRIBUTE can do. So you could
have
-ksp_type preonly -pc_type redistribute
-redistribute_ksp_type bcgs -redistribute_pc_type fieldsplit
-redistribute_pc_fieldsplit_detect_saddle_point
-redistribute_pc_fieldsplit_type schur
-redistribute_pc_fieldsplit_schur_fact_type lower
-redistribute_pc_fieldsplit_schur_precondition selfp
Using options solves the ordering problems you are having in setup and this
is what we recommend. Satisfying the ordering constraints by hand, rather
than programmatically, is very hard.
Also note the selfp is not a good preconditioner for Stokes. The Schur
complement looks like
div (lap)^{-1} grad ~ id
so you want the mass matrix as the preconditioning matrix, not the
diagonally weighted pressure Laplacian, which is what you get with selfp.
What you can do is make a preconditioning matrix for the whole problem that
looks like
/ A B \
\ C M /
where M is the mass matrix, or leave out B and C and use -pc_amat to
cleverly get them from the system matrix. You can see this in SNES ex69.
Thanks,
Matt
> 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/84471b37/attachment.html>
More information about the petsc-users
mailing list