[petsc-dev] MatNest and FieldSplit

Mon Apr 15 07:17:48 CDT 2019

On Mon, Apr 15, 2019 at 5:03 AM Pierre Jolivet via petsc-dev <
petsc-dev at mcs.anl.gov> wrote:

> OK, my solver is now working properly with a call to
> PetscObjectCompose((PetscObject)isS, "pmat", (PetscObject) myS);
>
> I have two follow-up questions:
> 1) am I supposed to call this, or is it the sign of something done wrong
> in my sequence of SNESSolve/KSPSetUp/KSPSetFromOptions/KSPSetOperators…?
>

Yes, unfortunately. We want to let the user pass in preconditioning
matrices for arbitrarily nested splits without pulling out
the objects explicitly, so we need a way to locate a particular split. I
guess you could imagine using the prefix, which would
be an XML way of doing things. We are doing it by attaching things to the
IS, since the Solver/Mat objects are created on the
fly. In fact, the IS is created on the fly by a DM sometimes, so you can
attach to the DM (I do this for Plex).

> 2) I have currently hardwired the split name I’m using when
> calling PCFieldSplitGetIS to get “isS” (from above) for debugging purposes.
> Could I create a PR with a new function like PCFieldSplitGetISByIndex(PC
> pc,const PetscInt n,IS *is) that will return the nth IS? Right now, I would
> need to get the KSP prefix followed up by some string comparison to get the
> actual IS prefix, whereas I know the position of the KSP in the
> PC_FieldSplitLink.
>

That sounds fine to me.

   Matt

> Thanks,
> Pierre
>
> On 14 Apr 2019, at 10:54 PM, Pierre Jolivet via petsc-dev <
> petsc-dev at mcs.anl.gov> wrote:
>
> I think I figured out what my problem is _exactly_.
> The Mat inside the MATNEST on which I’m using a PCFIELDSPLIT is
> unassembled before the first KSPSolve, except for the last field.
> Matt, you nailed it, when I call KSPSetFromOptions on the global
> PCFIELDSPLIT, then KSPSetUp explicitly on the inner PCFIELDSPLIT, the
> matrices associated to all fields are created here:
> https://www.mcs.anl.gov/petsc/petsc-current/src/ksp/pc/impls/fieldsplit/fieldsplit.c.html#line656
> Then, whatever I do with the matrix from the last field, currently trying
> MatCopy(myAssembledS, generatedS), before the first KSPSolve is reset by
> https://www.mcs.anl.gov/petsc/petsc-current/src/ksp/pc/impls/fieldsplit/fieldsplit.c.html#line688 and
> my solver fails…
>
> So pretty simple question, how do I set a “pmat” for my last assembled
> field so that
> https://www.mcs.anl.gov/petsc/petsc-current/src/ksp/pc/impls/fieldsplit/fieldsplit.c.html#line646
>  and
> https://www.mcs.anl.gov/petsc/petsc-current/src/ksp/pc/impls/fieldsplit/fieldsplit.c.html#line686 will
> return a non null pmat.
>
> This may sound really trivial but I’m lost in limbo right now. When
> everything is not wrapped inside an outer PCFIELDSPLIT, everything just
> work.
>
> Thanks,
> Pierre
>
> On 25 Mar 2019, at 6:57 PM, Pierre Jolivet via petsc-dev <
> petsc-dev at mcs.anl.gov> wrote:
>
> Thanks, this makes (slightly) more sense to me know.
> For some reason my application is still not acting properly but I must be
> screwing somewhere else the nested FieldSplit…
>
> Thank you,
> Pierre
>
> On 24 Mar 2019, at 11:42 PM, Dave May via petsc-dev <petsc-dev at mcs.anl.gov>
> wrote:
>
> Matt is right.
>
> When you defined the operator S, you basically invalidate the operator N
> (in the sense that they are no longer consistent). Hence when you use KSP
> nest to solve your problem your A matrix looks like
>   A = diag[1, 2, 4, 0, 8]
> but the B matrix you have defined looks like
>   B = diag[1, 2, 4, 0.00001]
>
> The only way to obtain the correct answer with your code is thus to use
> the option
> -ksp_type preonly
>
> Thanks
> Dave
>
>
>
> On Sun, 24 Mar 2019 at 22:09, Mark Adams via petsc-dev <
> petsc-dev at mcs.anl.gov> wrote:
>
>> I think he is saying that this line seems to have no effect (and the
>> comment is hence wrong):
>>
>> KSPSetOperators(subksp[nsplits - 1], S, S);
>>
>> // J2 = [[4, 0] ; [0, 0.00001]]
>>
>>
>> J2 is a 2x2 but this block has been changed into two single equation
>> fields. Does this KSPSetOperators supposed to copy this 1x1 S matrix into
>> the (1,1) block of the "J2", or do some sort of correct mixing internally,
>> to get what he wants?
>>
>>
>> BTW, this line does not seem necessary to me so maybe I'm missing something.
>>
>>
>> KSPSetOperators(sub, J2, J2);
>>
>>
>>
>> On Sun, Mar 24, 2019 at 4:33 PM Matthew Knepley via petsc-dev <
>> petsc-dev at mcs.anl.gov> wrote:
>>
>>> On Sun, Mar 24, 2019 at 10:21 AM Pierre Jolivet <
>>> pierre.jolivet at enseeiht.fr> wrote:
>>>
>>>> It’s a 4x4 matrix.
>>>> The first 2x2 diagonal matrix is a field.
>>>> The second 2x2 diagonal matrix is another field.
>>>> In the second field, the first diagonal coefficient is a subfield.
>>>> In the second field, the second diagonal coefficient is another
>>>> subfield.
>>>> I’m changing the operators from the second subfield (last diagonal
>>>> coefficient of the matrix).
>>>> When I solve a system with the complete matrix (2 fields), I get a
>>>> different “partial solution" than when I solve the “partial system” on just
>>>> the second field (with the two subfields in which I modified the operators
>>>> from the second one).
>>>>
>>>
>>> I may understand waht you are doing.
>>> Fieldsplit calls MatGetSubMatrix() which can copy values, depending on
>>> the implementation,
>>> so changing values in the original matrix may or may not change it in
>>> the PC.
>>>
>>>    Matt
>>>
>>> I don’t know if this makes more or less sense… sorry :\
>>>> Thanks,
>>>> Pierre
>>>>
>>>> On 24 Mar 2019, at 8:42 PM, Matthew Knepley <knepley at gmail.com> wrote:
>>>>
>>>> On Sat, Mar 23, 2019 at 9:12 PM Pierre Jolivet via petsc-dev <
>>>> petsc-dev at mcs.anl.gov> wrote:
>>>>
>>>>> I’m trying to figure out why both solutions are not consistent in the
>>>>> following example.
>>>>> Is what I’m doing complete nonsense?
>>>>>
>>>>
>>>> The code does not make clear what you are asking. I can see its a
>>>> nested fieldsplit.
>>>>
>>>>   Thanks,
>>>>
>>>>      Matt
>>>>
>>>>
>>>>> Thanks in advance for your help,
>>>>> Pierre
>>>>>
>>>>>
>>>>
>>>> --
>>>> 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/>
>>>>
>>>>
>>>>
>>>
>>> --
>>> 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/>
>>>
>>
>
>
>

-- 
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-dev/attachments/20190415/ede04808/attachment.html>