[petsc-users] PetscDSSetJacobianPreconditioner causing DIVERGED_LINE_SEARCH for multi-field problem
Sander Arens
Sander.Arens at UGent.be
Sun Apr 3 09:26:40 CDT 2016
No problem. I think I solved the problem with the block size, so I wanted
to know if I should also add this to my pull request or if I should make a
separate one?
Thanks,
Sander
On 3 April 2016 at 15:17, Matthew Knepley <knepley at gmail.com> wrote:
> On Sun, Apr 3, 2016 at 6:55 AM, Sander Arens <Sander.Arens at ugent.be>
> wrote:
>
>> Any thoughts on this, Matt?
>>
>
> Shoot, this fell out of my working set, and I am behind on incorporating
> your pull request. I will
> look at both this week.
>
> Thanks,
>
> Matt
>
>
>> Thanks,
>> Sander
>>
>> On 21 March 2016 at 09:49, Sander Arens <Sander.Arens at ugent.be> wrote:
>>
>>> Matt,
>>>
>>> I think the problem is in DMCreateSubDM_Section_Private. The block size
>>> can be counted there (very similar as in DMCreateMatrix_Plex) and then
>>> passed down to the IS. I can create a pull request for this if you'd like
>>> (or should I just add this to my other pull request for the Neumann bc's)?
>>>
>>> On 3 March 2016 at 18:02, Sander Arens <Sander.Arens at ugent.be> wrote:
>>>
>>>> Ok, thanks for clearing that up!
>>>>
>>>> On 3 March 2016 at 18:00, Matthew Knepley <knepley at gmail.com> wrote:
>>>>
>>>>> On Thu, Mar 3, 2016 at 10:48 AM, Sander Arens <Sander.Arens at ugent.be>
>>>>> wrote:
>>>>>
>>>>>> So is the block size here refering to the dimension of the
>>>>>> near-nullspace or has it something to do with the coarsening?
>>>>>>
>>>>>> I would have thought to also see the block size in this part:
>>>>>>
>>>>>
>>>>> Dang, that is a problem. It should have the correct block size when FS
>>>>> pulls it out of the overall matrix. This should be
>>>>> passed down using the block size of the IS. Something has been broken
>>>>> here. I will put it on my list of things to look at.
>>>>>
>>>>> Thanks,
>>>>>
>>>>> Matt
>>>>>
>>>>>
>>>>>> Mat Object: 1 MPI
>>>>>> processes
>>>>>> type: seqaij
>>>>>> rows=1701, cols=1701
>>>>>> total: nonzeros=109359, allocated nonzeros=109359
>>>>>> total number of mallocs used during MatSetValues
>>>>>> calls =0
>>>>>> has attached near null space
>>>>>> using I-node routines: found 567 nodes, limit
>>>>>> used is 5
>>>>>>
>>>>>> Thanks,
>>>>>> Sander
>>>>>>
>>>>>>
>>>>>> On 3 March 2016 at 17:29, Matthew Knepley <knepley at gmail.com> wrote:
>>>>>>
>>>>>>> You can see the block size in the output
>>>>>>>
>>>>>>> Mat Object: 1 MPI processes
>>>>>>> type: seqaij
>>>>>>> rows=12, cols=12, bs=6
>>>>>>> total: nonzeros=144, allocated nonzeros=144
>>>>>>> total number of mallocs used during MatSetValues
>>>>>>> calls =0
>>>>>>> using I-node routines: found 3 nodes, limit
>>>>>>> used is 5
>>>>>>> linear system matrix = precond matrix:
>>>>>>> Mat Object: 1 MPI processes
>>>>>>> type: seqaij
>>>>>>> rows=12, cols=12, bs=6
>>>>>>> total: nonzeros=144, allocated nonzeros=144
>>>>>>> total number of mallocs used during MatSetValues
>>>>>>> calls =0
>>>>>>> using I-node routines: found 3 nodes, limit used
>>>>>>> is 5
>>>>>>>
>>>>>>> Thanks,
>>>>>>>
>>>>>>> Matt
>>>>>>>
>>>>>>> On Thu, Mar 3, 2016 at 10:28 AM, Sander Arens <Sander.Arens at ugent.be
>>>>>>> > wrote:
>>>>>>>
>>>>>>>> Sure, here it is.
>>>>>>>>
>>>>>>>> Thanks,
>>>>>>>> Sander
>>>>>>>>
>>>>>>>> On 3 March 2016 at 15:33, Matthew Knepley <knepley at gmail.com>
>>>>>>>> wrote:
>>>>>>>>
>>>>>>>>> On Thu, Mar 3, 2016 at 8:05 AM, Sander Arens <
>>>>>>>>> Sander.Arens at ugent.be> wrote:
>>>>>>>>>
>>>>>>>>>> Yes, the matrix is created by one of the assembly functions in
>>>>>>>>>> Plex, so that answers my question. I was confused by this because I
>>>>>>>>>> couldn't see this information when using -snes_view.
>>>>>>>>>>
>>>>>>>>>
>>>>>>>>> Hmm, then something is wrong because the block size should be
>>>>>>>>> printed for the matrix at the end of the solver in -snes_view, Can you
>>>>>>>>> send that to me?
>>>>>>>>>
>>>>>>>>> Thanks,
>>>>>>>>>
>>>>>>>>> Matt
>>>>>>>>>
>>>>>>>>>
>>>>>>>>>> Thanks for the helpful replies,
>>>>>>>>>> Sander
>>>>>>>>>>
>>>>>>>>>> On 3 March 2016 at 14:52, Matthew Knepley <knepley at gmail.com>
>>>>>>>>>> wrote:
>>>>>>>>>>
>>>>>>>>>>> On Thu, Mar 3, 2016 at 7:49 AM, Sander Arens <
>>>>>>>>>>> Sander.Arens at ugent.be> wrote:
>>>>>>>>>>>
>>>>>>>>>>>> And how can I do this? Because when I look at all the options
>>>>>>>>>>>> with -help I can strangely enough only find
>>>>>>>>>>>> -fieldsplit_pressure_mat_block_size and not
>>>>>>>>>>>> -fieldsplit_displacement_mat_block_size.
>>>>>>>>>>>>
>>>>>>>>>>>
>>>>>>>>>>> Maybe I am missing something here. The matrix from which you
>>>>>>>>>>> calculate the preconditioner using GAMG must be created somewhere. Why
>>>>>>>>>>> is the block size not specified at creation time? If the matrix
>>>>>>>>>>> is created by one of the assembly functions in Plex, and the PetscSection
>>>>>>>>>>> has
>>>>>>>>>>> a number of field components, then the block size will already
>>>>>>>>>>> be set.
>>>>>>>>>>>
>>>>>>>>>>> Thanks,
>>>>>>>>>>>
>>>>>>>>>>> Matt
>>>>>>>>>>>
>>>>>>>>>>>
>>>>>>>>>>>> Thanks,
>>>>>>>>>>>> Sander
>>>>>>>>>>>>
>>>>>>>>>>>> On 3 March 2016 at 14:21, Matthew Knepley <knepley at gmail.com>
>>>>>>>>>>>> wrote:
>>>>>>>>>>>>
>>>>>>>>>>>>> On Thu, Mar 3, 2016 at 6:20 AM, Sander Arens <
>>>>>>>>>>>>> Sander.Arens at ugent.be> wrote:
>>>>>>>>>>>>>
>>>>>>>>>>>>>> Ok, I forgot to call SNESSetJacobian(snes, J, P, NULL, NULL)
>>>>>>>>>>>>>> with J != P, which caused to write the mass matrix into the (otherwise
>>>>>>>>>>>>>> zero) (1,1) block of the Jacobian and which was the reason for the
>>>>>>>>>>>>>> linesearch to fail.
>>>>>>>>>>>>>> However, after fixing that and trying to solve it with
>>>>>>>>>>>>>> FieldSplit with LU factorization for the (0,0) block it failed because
>>>>>>>>>>>>>> there were zero pivots for all rows.
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> Anyway, I found out that attaching the mass matrix to the
>>>>>>>>>>>>>> Lagrange multiplier field also worked.
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> Another related question for my elasticity problem: after
>>>>>>>>>>>>>> creating the rigid body modes with DMPlexCreateRigidBody and attaching it
>>>>>>>>>>>>>> to the displacement field, does the matrix block size of the (0,0) block
>>>>>>>>>>>>>> still have to be set for good performance with gamg? If so, how can I do
>>>>>>>>>>>>>> this?
>>>>>>>>>>>>>>
>>>>>>>>>>>>>
>>>>>>>>>>>>> Yes, it should be enough to set the block size of the
>>>>>>>>>>>>> preconditioner matrix.
>>>>>>>>>>>>>
>>>>>>>>>>>>> Matt
>>>>>>>>>>>>>
>>>>>>>>>>>>>
>>>>>>>>>>>>>> Thanks,
>>>>>>>>>>>>>> Sander
>>>>>>>>>>>>>>
>>>>>>>>>>>>>> On 2 March 2016 at 12:25, Matthew Knepley <knepley at gmail.com>
>>>>>>>>>>>>>> wrote:
>>>>>>>>>>>>>>
>>>>>>>>>>>>>>> On Wed, Mar 2, 2016 at 5:13 AM, Sander Arens <
>>>>>>>>>>>>>>> Sander.Arens at ugent.be> wrote:
>>>>>>>>>>>>>>>
>>>>>>>>>>>>>>>> Hi,
>>>>>>>>>>>>>>>>
>>>>>>>>>>>>>>>> I'm trying to set a mass matrix preconditioner for the
>>>>>>>>>>>>>>>> Schur complement of an incompressible finite elasticity problem. I tried
>>>>>>>>>>>>>>>> using the command PetscDSSetJacobianPreconditioner(prob, 1, 1,
>>>>>>>>>>>>>>>> g0_pre_mass_pp, NULL, NULL, NULL) (field 1 is the Lagrange multiplier
>>>>>>>>>>>>>>>> field).
>>>>>>>>>>>>>>>> However, this causes a DIVERGED_LINE_SEARCH due to to Nan
>>>>>>>>>>>>>>>> or Inf in the function evaluation after Newton iteration 1. (Btw, I'm using
>>>>>>>>>>>>>>>> the next branch).
>>>>>>>>>>>>>>>>
>>>>>>>>>>>>>>>> Is this because I didn't use
>>>>>>>>>>>>>>>> PetscDSSetJacobianPreconditioner for the other blocks (which uses the
>>>>>>>>>>>>>>>> Jacobian itself for preconditioning)? If so, how can I tell Petsc to use
>>>>>>>>>>>>>>>> the Jacobian for those blocks?
>>>>>>>>>>>>>>>>
>>>>>>>>>>>>>>>
>>>>>>>>>>>>>>> 1) I put that code in very recently, and do not even have
>>>>>>>>>>>>>>> sufficient test, so it may be buggy
>>>>>>>>>>>>>>>
>>>>>>>>>>>>>>> 2) If you are using FieldSplit, you can control which blocks
>>>>>>>>>>>>>>> come from A and which come from the preconditioner P
>>>>>>>>>>>>>>>
>>>>>>>>>>>>>>>
>>>>>>>>>>>>>>> http://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/PC/PCFieldSplitSetDiagUseAmat.html#PCFieldSplitSetDiagUseAmat
>>>>>>>>>>>>>>>
>>>>>>>>>>>>>>> http://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/PC/PCFieldSplitSetOffDiagUseAmat.html#PCFieldSplitSetOffDiagUseAmat
>>>>>>>>>>>>>>>
>>>>>>>>>>>>>>>
>>>>>>>>>>>>>>>> I guess when using PetscDSSetJacobianPreconditioner the
>>>>>>>>>>>>>>>> preconditioner is recomputed at every Newton step, so for a constant mass
>>>>>>>>>>>>>>>> matrix this might not be ideal. How can I avoid recomputing this at every
>>>>>>>>>>>>>>>> Newton iteration?
>>>>>>>>>>>>>>>>
>>>>>>>>>>>>>>>
>>>>>>>>>>>>>>> Maybe we need another flag like
>>>>>>>>>>>>>>>
>>>>>>>>>>>>>>>
>>>>>>>>>>>>>>> http://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/SNES/SNESSetLagPreconditioner.html
>>>>>>>>>>>>>>>
>>>>>>>>>>>>>>> or we need to expand
>>>>>>>>>>>>>>>
>>>>>>>>>>>>>>>
>>>>>>>>>>>>>>> http://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/SNES/SNESSetLagJacobian.html
>>>>>>>>>>>>>>>
>>>>>>>>>>>>>>> to separately cover the preconditioner matrix. However, both
>>>>>>>>>>>>>>> matrices are computed by one call so this would
>>>>>>>>>>>>>>> involve interface changes to user code, which we do not like
>>>>>>>>>>>>>>> to do. Right now it seems like a small optimization.
>>>>>>>>>>>>>>> I would want to wait and see whether it would really be
>>>>>>>>>>>>>>> maningful.
>>>>>>>>>>>>>>>
>>>>>>>>>>>>>>> Thanks,
>>>>>>>>>>>>>>>
>>>>>>>>>>>>>>> Matt
>>>>>>>>>>>>>>>
>>>>>>>>>>>>>>>
>>>>>>>>>>>>>>>> Thanks,
>>>>>>>>>>>>>>>> Sander
>>>>>>>>>>>>>>>>
>>>>>>>>>>>>>>>
>>>>>>>>>>>>>>>
>>>>>>>>>>>>>>>
>>>>>>>>>>>>>>> --
>>>>>>>>>>>>>>> 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
>>>>>>>>>>>>>>>
>>>>>>>>>>>>>>
>>>>>>>>>>>>>>
>>>>>>>>>>>>>
>>>>>>>>>>>>>
>>>>>>>>>>>>> --
>>>>>>>>>>>>> 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
>>>>>>>>>>>>>
>>>>>>>>>>>>
>>>>>>>>>>>>
>>>>>>>>>>>
>>>>>>>>>>>
>>>>>>>>>>> --
>>>>>>>>>>> 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
>>>>>>>>>>>
>>>>>>>>>>
>>>>>>>>>>
>>>>>>>>>
>>>>>>>>>
>>>>>>>>> --
>>>>>>>>> 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
>>>>>>>>>
>>>>>>>>
>>>>>>>>
>>>>>>>
>>>>>>>
>>>>>>> --
>>>>>>> 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
>>>>>>>
>>>>>>
>>>>>>
>>>>>
>>>>>
>>>>> --
>>>>> 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
>>>>>
>>>>
>>>>
>>>
>>
>
>
> --
> 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
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.mcs.anl.gov/pipermail/petsc-users/attachments/20160403/c005714e/attachment-0001.html>
More information about the petsc-users
mailing list