[petsc-users] GAMG Indefinite
Sanjay Govindjee
s_g at berkeley.edu
Sun May 22 16:38:48 CDT 2016
Mark,
Can you give me the full option line that you want me to use? I
currently have:
-ksp_type cg -ksp_monitor -ksp_chebyshev_esteig_random -log_view
-pc_type gamg -pc_gamg_type agg -pc_gamg_agg_nsmooths 1 -options_left
-sanjay
On 5/22/16 2:29 PM, Mark Adams wrote:
> Humm, maybe we have version mixup:
>
> src/ksp/ksp/impls/cheby/cheby.c: ierr =
> PetscOptionsBool("-ksp_chebyshev_esteig_random","Use random right hand
> side for
> estimate","KSPChebyshevEstEigSetUseRandom",cheb->userandom,&cheb->userandom
>
> Also, you should use CG. These other options are the defaults but CG
> is not:
>
> -mg_levels_esteig_ksp_type cg
> -mg_levels_esteig_ksp_max_it 10
> -mg_levels_ksp_chebyshev_esteig 0,0.05,0,1.05
>
> Anyway. you can also run with -info, which will be very noisy, but
> just grep for GAMG and send me that.
>
> Mark
>
>
>
> On Sat, May 21, 2016 at 6:03 PM, Sanjay Govindjee <s_g at berkeley.edu
> <mailto:s_g at berkeley.edu>> wrote:
>
> Mark,
> I added the option you mentioned but it seems not to use it;
> -options_left reports:
>
> #PETSc Option Table entries:
> -ksp_chebyshev_esteig_random
> -ksp_monitor
> -ksp_type cg
> -log_view
> -options_left
> -pc_gamg_agg_nsmooths 1
> -pc_gamg_type agg
> -pc_type gamg
> #End of PETSc Option Table entries
> There is one unused database option. It is:
> Option left: name:-ksp_chebyshev_esteig_random (no value)
>
>
>
> On 5/21/16 12:36 PM, Mark Adams wrote:
>> Barry, this is probably the Chebyshev problem.
>>
>> Sanjay, this is fixed but has not yet been moved to the master
>> branch. You can fix this now with with
>> -ksp_chebyshev_esteig_random. This should recover v3.5 semantics.
>>
>> Mark
>>
>> On Thu, May 19, 2016 at 2:42 PM, Barry Smith <bsmith at mcs.anl.gov
>> <mailto:bsmith at mcs.anl.gov>> wrote:
>>
>>
>> We see this occasionally, there is nothing in the
>> definition of GAMG that guarantees a positive definite
>> preconditioner even if the operator was positive definite so
>> we don't think this is a bug in the code. We've found using a
>> slightly stronger smoother, like one more smoothing step
>> seems to remove the problem.
>>
>> Barry
>>
>> > On May 19, 2016, at 1:07 PM, Sanjay Govindjee
>> <s_g at berkeley.edu <mailto:s_g at berkeley.edu>> wrote:
>> >
>> > I am trying to solve a very ordinary nonlinear elasticity
>> problem
>> > using -ksp_type cg -pc_type gamg in PETSc 3.7.0, which
>> worked fine
>> > in PETSc 3.5.3.
>> >
>> > The problem I am seeing is on my first Newton iteration,
>> the Ax=b
>> > solve returns with and Indefinite Preconditioner error
>> (KSPGetConvergedReason == -8):
>> > (log_view.txt output also attached)
>> >
>> > 0 KSP Residual norm 8.411630828687e-02
>> > 1 KSP Residual norm 2.852209578900e-02
>> > NO CONVERGENCE REASON: Indefinite Preconditioner
>> > NO CONVERGENCE REASON: Indefinite Preconditioner
>> >
>> > On the next and subsequent Newton iterations, I see
>> perfectly normal
>> > behavior and the problem converges quadratically. The
>> results look fine.
>> >
>> > I tried the same problem with -pc_type jacobi as well as
>> super-lu, and mumps
>> > and they all work without complaint.
>> >
>> > My run line for GAMG is:
>> > -ksp_type cg -ksp_monitor -log_view -pc_type gamg
>> -pc_gamg_type agg -pc_gamg_agg_nsmooths 1 -options_left
>> >
>> > The code flow looks like:
>> >
>> > ! If no matrix allocation yet
>> > if(Kmat.eq.0) then
>> > call MatCreate(PETSC_COMM_WORLD,Kmat,ierr)
>> > call
>> MatSetSizes(Kmat,numpeq,numpeq,PETSC_DETERMINE,PETSC_DETERMINE,ierr)
>> > call MatSetBlockSize(Kmat,nsbk,ierr)
>> > call MatSetFromOptions(Kmat, ierr)
>> > call MatSetType(Kmat,MATAIJ,ierr)
>> > call
>> MatMPIAIJSetPreallocation(Kmat,PETSC_NULL_INTEGER,mr(np(246)),PETSC_NULL_INTEGER,mr(np(247)),ierr)
>> > call
>> MatSeqAIJSetPreallocation(Kmat,PETSC_NULL_INTEGER,mr(np(246)),ierr)
>> > endif
>> >
>> > call MatZeroEntries(Kmat,ierr)
>> >
>> > ! Code to set values in matrix
>> >
>> > call MatAssemblyBegin(Kmat, MAT_FINAL_ASSEMBLY, ierr)
>> > call MatAssemblyEnd(Kmat, MAT_FINAL_ASSEMBLY, ierr)
>> > call
>> MatSetOption(Kmat,MAT_NEW_NONZERO_LOCATIONS,PETSC_TRUE,ierr)
>> >
>> > ! If no rhs allocation yet
>> > if(rhs.eq.0) then
>> > call VecCreate (PETSC_COMM_WORLD, rhs, ierr)
>> > call VecSetSizes (rhs, numpeq, PETSC_DECIDE, ierr)
>> > call VecSetFromOptions(rhs, ierr)
>> > endif
>> >
>> > ! Code to set values in RHS
>> >
>> > call VecAssemblyBegin(rhs, ierr)
>> > call VecAssemblyEnd(rhs, ierr)
>> >
>> > if(kspsol_exists) then
>> > call KSPDestroy(kspsol,ierr)
>> > endif
>> >
>> > call KSPCreate(PETSC_COMM_WORLD, kspsol ,ierr)
>> > call KSPSetOperators(kspsol, Kmat, Kmat, ierr)
>> > call KSPSetFromOptions(kspsol,ierr)
>> > call KSPGetPC(kspsol, pc , ierr)
>> >
>> > call PCSetCoordinates(pc,ndm,numpn,hr(np(43)),ierr)
>> >
>> > call KSPSolve(kspsol, rhs, sol, ierr)
>> > call KSPGetConvergedReason(kspsol,reason,ierr)
>> >
>> > ! update solution, go back to the top
>> >
>> > reason is coming back as -8 on my first Ax=b solve and 2 or
>> 3 after that
>> > (with gamg). With the other solvers it is coming back as 2
>> or 3 for
>> > iterative options and 4 if I use one of the direct solvers.
>> >
>> > Any ideas on what is causing the Indefinite PC on the first
>> iteration with GAMG?
>> >
>> > Thanks in advance,
>> > -sanjay
>> >
>> > --
>> > -----------------------------------------------
>> > Sanjay Govindjee, PhD, PE
>> > Professor of Civil Engineering
>> >
>> > 779 Davis Hall
>> > University of California
>> > Berkeley, CA 94720-1710
>> >
>> > Voice: +1 510 642 6060 <tel:%2B1%20510%20642%206060>
>> > FAX: +1 510 643 5264 <tel:%2B1%20510%20643%205264>
>> >
>> > s_g at berkeley.edu <mailto:s_g at berkeley.edu>
>> > http://www.ce.berkeley.edu/~sanjay
>> <http://www.ce.berkeley.edu/%7Esanjay>
>> >
>> > -----------------------------------------------
>> >
>> > Books:
>> >
>> > Engineering Mechanics of Deformable
>> > Solids: A Presentation with Exercises
>> >
>> >
>> http://www.oup.com/us/catalog/general/subject/Physics/MaterialsScience/?view=usa&ci=9780199651641
>> > http://ukcatalogue.oup.com/product/9780199651641.do
>> > http://amzn.com/0199651647
>> >
>> >
>> > Engineering Mechanics 3 (Dynamics) 2nd Edition
>> >
>> > http://www.springer.com/978-3-642-53711-0
>> > http://amzn.com/3642537111
>> >
>> >
>> > Engineering Mechanics 3, Supplementary Problems: Dynamics
>> >
>> > http://www.amzn.com/B00SOXN8JU
>> >
>> >
>> > -----------------------------------------------
>> >
>> > <log_view.txt>
>>
>>
>
> --
> -----------------------------------------------
> Sanjay Govindjee, PhD, PE
> Professor of Civil Engineering
>
> 779 Davis Hall
> University of California
> Berkeley, CA 94720-1710
>
> Voice:+1 510 642 6060 <tel:%2B1%20510%20642%206060>
> FAX:+1 510 643 5264 <tel:%2B1%20510%20643%205264>
> s_g at berkeley.edu <mailto:s_g at berkeley.edu>
> http://www.ce.berkeley.edu/~sanjay
> <http://www.ce.berkeley.edu/%7Esanjay>
> -----------------------------------------------
>
> Books:
>
> Engineering Mechanics of Deformable
> Solids: A Presentation with Exercises
> http://www.oup.com/us/catalog/general/subject/Physics/MaterialsScience/?view=usa&ci=9780199651641
> http://ukcatalogue.oup.com/product/9780199651641.do
> http://amzn.com/0199651647
>
> Engineering Mechanics 3 (Dynamics) 2nd Edition
> http://www.springer.com/978-3-642-53711-0
> http://amzn.com/3642537111
>
> Engineering Mechanics 3, Supplementary Problems: Dynamics
> http://www.amzn.com/B00SOXN8JU
>
> -----------------------------------------------
>
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.mcs.anl.gov/pipermail/petsc-users/attachments/20160522/11c3b6e4/attachment-0001.html>
More information about the petsc-users
mailing list