[petsc-users] Question concerning ilu and bcgs

Matthew Knepley knepley at gmail.com
Wed Feb 18 10:54:54 CST 2015


On Wed, Feb 18, 2015 at 10:47 AM, Sun, Hui <hus003 at ucsd.edu> wrote:

>  The matrix is from a 3D fluid problem, with complicated irregular
> boundary conditions. I've tried using direct solvers such as UMFPACK,
> SuperLU_dist and MUMPS. It seems that SuperLU_dist does not solve for my
> linear system; UMFPACK solves the system but would run into memory issue
> even with small size matrices and it cannot parallelize; MUMPS does solve
> the system but it also fails when the size is big and it takes much time.
> That's why I'm seeking an iterative method.
>
>  I guess the direct method is faster than an iterative method for a small
> A, but that may not be true for bigger A.
>

If this is a Stokes flow, you should use PCFIELDSPLIT and multigrid. If it
is advection dominated, I know of nothing better
than sparse direct or perhaps Block-Jacobi with sparse direct blocks. Since
MUMPS solved your system, I would consider
using BJacobi/ASM and MUMPS or UMFPACK as the block solver.

  Thanks,

     Matt


>
>  ------------------------------
> *From:* Matthew Knepley [knepley at gmail.com]
> *Sent:* Wednesday, February 18, 2015 8:33 AM
> *To:* Sun, Hui
> *Cc:* hong at aspiritech.org; petsc-users at mcs.anl.gov
> *Subject:* Re: [petsc-users] Question concerning ilu and bcgs
>
>    On Wed, Feb 18, 2015 at 10:31 AM, Sun, Hui <hus003 at ucsd.edu> wrote:
>
>>  So far I just try around, I haven't looked into literature yet.
>>
>>  However, both MATLAB's ilu+gmres and ilu+bcgs work. Is it possible that
>> some parameter or options need to be tuned in using PETSc's ilu or hypre's
>> ilu? Besides, is there a way to view how good the performance of the pc is
>> and output the matrices L and U, so that I can do some test in MATLAB?
>>
>
>  1) Its not clear exactly what Matlab is doing
>
>  2) PETSc uses ILU(0) by default (you can set it to use ILU(k))
>
>  3) I don't know what Hypre's ILU can do
>
>  I would really discourage from using ILU. I cannot imagine it is faster
> than sparse direct factorization
> for your system, such as from SuperLU or MUMPS.
>
>    Thanks,
>
>       Matt
>
>
>>   Hui
>>
>>
>>  ------------------------------
>> *From:* Matthew Knepley [knepley at gmail.com]
>> *Sent:* Wednesday, February 18, 2015 8:09 AM
>> *To:* Sun, Hui
>> *Cc:* hong at aspiritech.org; petsc-users at mcs.anl.gov
>> *Subject:* Re: [petsc-users] Question concerning ilu and bcgs
>>
>>    On Wed, Feb 18, 2015 at 10:02 AM, Sun, Hui <hus003 at ucsd.edu> wrote:
>>
>>>  Yes I've tried other solvers, gmres/ilu does not work, neither does
>>> bcgs/ilu. Here are the options:
>>>
>>> -pc_type ilu -pc_factor_nonzeros_along_diagonal -pc_factor_levels 0
>>> -pc_factor_reuse_ordering -ksp_ty\
>>>
>>> pe bcgs -ksp_rtol 1e-6 -ksp_max_it 10 -ksp_monitor_short -ksp_view
>>>
>>
>>  Note here that ILU(0) is an unreliable and generally crappy
>> preconditioner. Have you looked in the
>> literature for the kinds of preconditioners that are effective for your
>> problem?
>>
>>    Thanks,
>>
>>       Matt
>>
>>
>>>   Here is the output:
>>>
>>>   0 KSP Residual norm 211292
>>>
>>>   1 KSP Residual norm 13990.2
>>>
>>>   2 KSP Residual norm 9870.08
>>>
>>>   3 KSP Residual norm 9173.9
>>>
>>>   4 KSP Residual norm 9121.94
>>>
>>>   5 KSP Residual norm 7386.1
>>>
>>>   6 KSP Residual norm 6222.55
>>>
>>>   7 KSP Residual norm 7192.94
>>>
>>>   8 KSP Residual norm 33964
>>>
>>>   9 KSP Residual norm 33960.4
>>>
>>>  10 KSP Residual norm 1068.54
>>>
>>> KSP Object: 1 MPI processes
>>>
>>>   type: bcgs
>>>
>>>   maximum iterations=10, initial guess is zero
>>>
>>>   tolerances:  relative=1e-06, absolute=1e-50, divergence=10000
>>>
>>>   left preconditioning
>>>
>>>   using PRECONDITIONED norm type for convergence test
>>>
>>> PC Object: 1 MPI processes
>>>
>>>   type: ilu
>>>
>>>     ILU: out-of-place factorization
>>>
>>>     ILU: Reusing reordering from past factorization
>>>
>>>     0 levels of fill
>>>
>>>     tolerance for zero pivot 2.22045e-14
>>>
>>>     using diagonal shift on blocks to prevent zero pivot [INBLOCKS]
>>>
>>>     matrix ordering: natural
>>>
>>>     factor fill ratio given 1, needed 1
>>>
>>>       Factored matrix follows:
>>>
>>>         Mat Object:         1 MPI processes
>>>
>>>           type: seqaij
>>>
>>>           rows=62500, cols=62500
>>>
>>>           package used to perform factorization: petsc
>>>
>>>           total: nonzeros=473355, allocated nonzeros=473355
>>>
>>>           total number of mallocs used during MatSetValues calls =0
>>>
>>>             not using I-node routines
>>>
>>>   linear system matrix = precond matrix:
>>>
>>>   Mat Object:   1 MPI processes
>>>
>>>     type: seqaij
>>>
>>>     rows=62500, cols=62500
>>>
>>>     total: nonzeros=473355, allocated nonzeros=7.8125e+06
>>>
>>>     total number of mallocs used during MatSetValues calls =0
>>>
>>>       not using I-node routines
>>>
>>> Time cost: 0.307149,  0.268402,  0.0990018
>>>
>>>
>>>
>>>
>>>  ------------------------------
>>> *From:* hong at aspiritech.org [hong at aspiritech.org]
>>> *Sent:* Wednesday, February 18, 2015 7:49 AM
>>> *To:* Sun, Hui
>>> *Cc:* Matthew Knepley; petsc-users at mcs.anl.gov
>>> *Subject:* Re: [petsc-users] Question concerning ilu and bcgs
>>>
>>>    Have you tried other solvers, e.g., PETSc default gmres/ilu,
>>> bcgs/ilu etc.
>>> The matrix is small. If it is ill-conditioned, then pc_type lu would
>>> work the best.
>>>
>>>  Hong
>>>
>>> On Wed, Feb 18, 2015 at 9:34 AM, Sun, Hui <hus003 at ucsd.edu> wrote:
>>>
>>>>  With options:
>>>>
>>>>  -pc_type hypre -pc_hypre_type pilut -pc_hypre_pilut_maxiter 1000
>>>> -pc_hypre_pilut_tol 1e-3 -ksp_type bcgs -ksp_rtol 1e-10 -ksp_max_it 10
>>>> -ksp_monitor_short -ksp_converged_reason -ksp_view
>>>>
>>>>  Here is the full output:
>>>>
>>>>    0 KSP Residual norm 1404.62
>>>>
>>>>   1 KSP Residual norm 88.9068
>>>>
>>>>   2 KSP Residual norm 64.73
>>>>
>>>>   3 KSP Residual norm 71.0224
>>>>
>>>>   4 KSP Residual norm 69.5044
>>>>
>>>>   5 KSP Residual norm 455.458
>>>>
>>>>   6 KSP Residual norm 174.876
>>>>
>>>>   7 KSP Residual norm 183.031
>>>>
>>>>   8 KSP Residual norm 650.675
>>>>
>>>>   9 KSP Residual norm 79.2441
>>>>
>>>>  10 KSP Residual norm 84.1985
>>>>
>>>> Linear solve did not converge due to DIVERGED_ITS iterations 10
>>>>
>>>> KSP Object: 1 MPI processes
>>>>
>>>>   type: bcgs
>>>>
>>>>   maximum iterations=10, initial guess is zero
>>>>
>>>>   tolerances:  relative=1e-10, absolute=1e-50, divergence=10000
>>>>
>>>>   left preconditioning
>>>>
>>>>   using PRECONDITIONED norm type for convergence test
>>>>
>>>> PC Object: 1 MPI processes
>>>>
>>>>   type: hypre
>>>>
>>>>     HYPRE Pilut preconditioning
>>>>
>>>>     HYPRE Pilut: maximum number of iterations 1000
>>>>
>>>>     HYPRE Pilut: drop tolerance 0.001
>>>>
>>>>     HYPRE Pilut: default factor row size
>>>>
>>>>   linear system matrix = precond matrix:
>>>>
>>>>   Mat Object:   1 MPI processes
>>>>
>>>>     type: seqaij
>>>>
>>>>     rows=62500, cols=62500
>>>>
>>>>     total: nonzeros=473355, allocated nonzeros=7.8125e+06
>>>>
>>>>     total number of mallocs used during MatSetValues calls =0
>>>>
>>>>       not using I-node routines
>>>>
>>>> Time cost: 0.756198,  0.662984,  0.105672
>>>>
>>>>
>>>>
>>>>
>>>>  ------------------------------
>>>> *From:* Matthew Knepley [knepley at gmail.com]
>>>> *Sent:* Wednesday, February 18, 2015 3:30 AM
>>>> *To:* Sun, Hui
>>>> *Cc:* petsc-users at mcs.anl.gov
>>>> *Subject:* Re: [petsc-users] Question concerning ilu and bcgs
>>>>
>>>>     On Wed, Feb 18, 2015 at 12:33 AM, Sun, Hui <hus003 at ucsd.edu> wrote:
>>>>
>>>>>  I have a matrix system Ax = b, A is of type MatSeqAIJ or MatMPIAIJ,
>>>>> depending on the number of cores.
>>>>>
>>>>>  I try to solve this problem by pc_type ilu and ksp_type bcgs, it
>>>>> does not converge. The options I specify are:
>>>>>
>>>>> -pc_type hypre -pc_hypre_type pilut -pc_hypre_pilut_maxiter 1000
>>>>> -pc_hypre_pilut_tol 1e-3 -ksp_type b\
>>>>>
>>>>> cgs -ksp_rtol 1e-10 -ksp_max_it 1000 -ksp_monitor_short
>>>>> -ksp_converged_reason
>>>>>
>>>>
>>>>  1) Run with -ksp_view, so we can see exactly what was used
>>>>
>>>>  2) ILUT is unfortunately not a well-defined algorithm, and I believe
>>>> the parallel version makes different decisions
>>>>     than the serial version.
>>>>
>>>>    Thanks,
>>>>
>>>>      Matt
>>>>
>>>>
>>>>>   The first a few lines of the output are:
>>>>>
>>>>>   0 KSP Residual norm 1404.62
>>>>>
>>>>>   1 KSP Residual norm 88.9068
>>>>>
>>>>>   2 KSP Residual norm 64.73
>>>>>
>>>>>   3 KSP Residual norm 71.0224
>>>>>
>>>>>   4 KSP Residual norm 69.5044
>>>>>
>>>>>   5 KSP Residual norm 455.458
>>>>>
>>>>>   6 KSP Residual norm 174.876
>>>>>
>>>>>   7 KSP Residual norm 183.031
>>>>>
>>>>>   8 KSP Residual norm 650.675
>>>>>
>>>>>   9 KSP Residual norm 79.2441
>>>>>
>>>>>  10 KSP Residual norm 84.1985
>>>>>
>>>>>
>>>>>  This clearly indicates non-convergence. However, I output the sparse
>>>>> matrix A and vector b to MATLAB, and run the following command:
>>>>>
>>>>> [L,U] = ilu(A,struct('type','ilutp','droptol',1e-3));
>>>>>
>>>>> [ux1,fl1,rr1,it1,rv1] = bicgstab(A,b,1e-10,1000,L,U);
>>>>>
>>>>>
>>>>>  And it converges in MATLAB, with flag fl1=0, relative residue
>>>>> rr1=8.2725e-11, and iteration it1=89.5. I'm wondering how can I figure out
>>>>> what's wrong.
>>>>>
>>>>>
>>>>>  Best,
>>>>>
>>>>> Hui
>>>>>
>>>>
>>>>
>>>>
>>>>  --
>>>> 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/20150218/0d4b5a8f/attachment-0001.html>


More information about the petsc-users mailing list