[petsc-users] use Superlu as ilu preconditioner
Hong Zhang
hzhang at mcs.anl.gov
Tue Aug 2 15:32:40 CDT 2011
Ping,
I'm sorry, the code is too difficult to read and is not executable by itself :-(
I do not understand what you want to do.
Do you want solve
A*x=rhs with same A and multiple rhs vectors?
>From your previous report, it seems that original matrix A is
overwritten by matrix factorization when superlu's ilu is used.
Please correct me if I'm wrong.
Hong
On Tue, Aug 2, 2011 at 1:50 PM, Ping Rong <ping.rong at tuhh.de> wrote:
> Hong:
> sorry for the delay. I have put the code pieces together, search for the
> keyword ATTENTION, I have put some comments there.
> the overwrite happens when the system matrix is passed to the solver by
>
> linear_solver->solve (*matrix, *matrix, *solution, *rhs, tol,
> maxits);
>
> In the file you should also find the a init() function, which is called when
> the linear_solver object is instantiated.
> By now I created an identical precond_matrix, and called the solve function
> by
>
> linear_solver->solve (*matrix, *precond_matrix, *solution, *rhs, tol,
> maxits);
>
> The overwritten only happens when using superlu as ilu preconditioner. petsc
> default ilu preconditioner doesn't have this issue.
> Please tell me, if you need any more information.
>
> Ping
>
> Am 29.07.2011 16:42, schrieb Hong Zhang:
>>
>> Ping :
>> Please send me your code. I do not understand why the system matrix is
>> overwritten
>> by factored matrix when superlu is used. We do not support in-place
>> factorization
>> for external packages. Something is incorrect here.
>>
>> Hong
>>>
>>> I have found the problem. I don't know whether this is a bug.
>>>
>>> ierr = KSPSetOperators(_ksp, matrix->mat(), precond->mat(),
>>> this->same_preconditioner ?
>>> SAME_PRECONDITIONER : DIFFERENT_NONZERO_PATTERN);
>>>
>>>
>>> I used my system matrix as the precondition matrix, where matrix->mat()
>>> and
>>> precond->mat() are in deed the same matrix. I have left the
>>> KSPSetOperators(ksp) out before the KSPSolve(..) function, so
>>> KSPSetOperators() is called only once while the petsc_linear_solver
>>> object
>>> is initialized, which eliminated the petsc error ("object is in wrong
>>> state").
>>> Also when I used petsc default ilu preconditioner, everything is fine,
>>> also
>>> the solutions are correct. However, as long as I apply the superlu as my
>>> solver package, it overwrites the system matrix. I compared the system
>>> matrix before and after solve(). With petsc ilu, it remains identical,
>>> but
>>> with superlu ilu, the system matrix is changed to the preconditioner
>>> matrix
>>> itself. Now I simply duplicate the system matrix, and use this copy as
>>> precondition matrix. Everything works perfectly.
>>>
>>> The code I used is pretty standard, just as in ex2. It took me quite a
>>> while
>>> to realize the problem, because when the system matrix is replaced by the
>>> preconditioner matrix, you still get a "correct" solution, which is
>>> however
>>> incorrectly scaled.
>>>
>>> Ping
>>>
>>>
>>> Am 28.07.2011 18:21, schrieb Hong Zhang:
>>>>
>>>> Ping :
>>>> I tested similar calling procedure using
>>>> petsc-dev/src/ksp/ksp/examples/tutorials/ex5.c
>>>> successfully.
>>>>
>>>> Try following:
>>>> 1) switch to petsc-dev. See
>>>> http://www.mcs.anl.gov/petsc/petsc-as/developers/index.html on how to
>>>> get it.
>>>> 2) send us a short and stand-alone code that reproduce the error. We
>>>> can investigate it.
>>>>
>>>> Hong
>>>>
>>>>> Dear Hong,
>>>>>
>>>>> thank you very much for your time. Can you tell me when should
>>>>> KSPSetFromOptions(ksp) be called?
>>>>> well, libMesh has a class called "petsc_linear_solver", in its original
>>>>> code
>>>>> KSPSetFromOptions(ksp) is called while the object is created. I
>>>>> couldn't
>>>>> get
>>>>> any output, when I run the program with the option
>>>>>
>>>>> -pc_type ilu -pc_factor_mat_solver_package superlu
>>>>> -mat_superlu_ilu_droptol 0.1 -ksp_view
>>>>>
>>>>> The actual solve() function of the "petsc_linear_solver" contains the
>>>>> similar code as in the ex2.c
>>>>> ...
>>>>> ierr = KSPSetOperators(_ksp, matrix->mat(), precond->mat(),
>>>>> this->same_preconditioner ?
>>>>> SAME_PRECONDITIONER : DIFFERENT_NONZERO_PATTERN);
>>>>> CHKERRABORT(libMesh::COMM_WORLD,ierr);
>>>>>
>>>>> // Set the tolerances for the iterative solver. Use the user-supplied
>>>>> // tolerance for the relative residual& leave the others at default
>>>>> values.
>>>>> ierr = KSPSetTolerances (_ksp, tol, PETSC_DEFAULT, PETSC_DEFAULT,
>>>>> max_its);
>>>>> CHKERRABORT(libMesh::COMM_WORLD,ierr);
>>>>>
>>>>> ierr = KSPSetFromOptions (_ksp);<-------------- I added these two
>>>>> lines,
>>>>> CHKERRABORT(libMesh::COMM_WORLD,ierr);<-------------- then I got
>>>>> proper
>>>>> output.
>>>>>
>>>>> // Solve the linear system
>>>>> ierr = KSPSolve (_ksp, rhs->vec(), solution->vec());
>>>>> CHKERRABORT(libMesh::COMM_WORLD,ierr);
>>>>> ....
>>>>>
>>>>> The problem is that I have to solve a multiple frequency system, and
>>>>> for
>>>>> each frequency I have several RHS. The same ksp context is used to
>>>>> solve
>>>>> the
>>>>> systems repeatedly. I always call KSPSetOperators() with the option
>>>>> DIFFERENT_NONZERO_PATTERN for the first RHS and SAME_PRECONDITIONER for
>>>>> the
>>>>> rest, since the system matrix remains the same. For the first
>>>>> frequency,
>>>>> everything is fine, both DIFFERENT_NONZERO_PATTERN and
>>>>> SAME_PRECONDITIONER
>>>>> options. I got following output from the -ksp_view.
>>>>> ...
>>>>> PC Object:
>>>>> type: ilu
>>>>> ILU: out-of-place factorization
>>>>> 0 levels of fill
>>>>> tolerance for zero pivot 1e-12
>>>>> using diagonal shift to prevent zero pivot
>>>>> matrix ordering: natural
>>>>> factor fill ratio given 0, needed 0
>>>>> Factored matrix follows:
>>>>> Matrix Object:
>>>>> type=seqaij, rows=2883, cols=2883
>>>>> package used to perform factorization:
>>>>> superlu<----------------
>>>>> superlu is properly set
>>>>> ....
>>>>>
>>>>> but when the second frequency comes up, KSPSetOperators() is called
>>>>> again
>>>>> with DIFFERENT_NONZERO_PATTERN option for the first RHS. the command
>>>>> line
>>>>> option doesn't seem to work anymore.
>>>>> ...
>>>>> PC Object:
>>>>> type: ilu
>>>>> ILU: out-of-place factorization
>>>>> 0 levels of fill
>>>>> tolerance for zero pivot 1e-12
>>>>> using diagonal shift to prevent zero pivot
>>>>> matrix ordering: natural
>>>>> factor fill ratio given 1, needed 1
>>>>> Factored matrix follows:
>>>>> Matrix Object:
>>>>> type=seqaij, rows=2883, cols=2883
>>>>> package used to perform factorization: petsc<------------ not
>>>>> superlu anymore
>>>>> ...
>>>>>
>>>>> and then for the second RHS, when it calls KSPSetOperators(ksp) in the
>>>>> solve() routine, petsc throws the following error message.
>>>>>
>>>>> [0]PETSC ERROR: --------------------- Error Message
>>>>> ------------------------------------
>>>>> [0]PETSC ERROR: Object is in wrong state!
>>>>> [0]PETSC ERROR: Cannot change solver matrix package after PC has been
>>>>> setup
>>>>> or used!
>>>>> [0]PETSC ERROR:
>>>>>
>>>>> ------------------------------------------------------------------------
>>>>> [0]PETSC ERROR: Petsc Release Version 3.1.0, Patch 8, Thu Mar 17
>>>>> 13:37:48
>>>>> CDT 2011
>>>>> [0]PETSC ERROR: See docs/changes/index.html for recent updates.
>>>>> [0]PETSC ERROR: See docs/faq.html for hints about trouble shooting.
>>>>> [0]PETSC ERROR: See docs/index.html for manual pages.
>>>>> [0]PETSC ERROR:
>>>>>
>>>>> ------------------------------------------------------------------------
>>>>> [0]PETSC ERROR: ./PLTMainTest.opt on a linux-gnu named abe-vm by mubpr
>>>>> Wed
>>>>> Jul 27 19:41:22 2011
>>>>> [0]PETSC ERROR: Libraries linked from
>>>>> /home/mubpr/libs/petsc/linux-gnu-acml-shared-opt/lib
>>>>> [0]PETSC ERROR: Configure run at Tue Jul 26 22:09:39 2011
>>>>> [0]PETSC ERROR: Configure options
>>>>>
>>>>>
>>>>> --with-blas-lapack-lib=/home/mubpr/libs/acml/gfortran64_mp/lib/libacml_mp.so
>>>>> --with-blas-lapack-include=/home/mubpr/libs/acml/gfortran64_mp/include
>>>>> --with-scalar-type=complex --with-clanguage=c++
>>>>> --with-mpi-dir=/home/mubpr/libs/mpich2 --download-superlu=yes
>>>>> --with-superlu=1 --download-parmetis=yes --with-parmetis=1
>>>>> --download-superlu_dist=yes --with-superlu_dist=1 --download-mumps=yes
>>>>> --with-mumps=1 --download-blacs=yes --with-blacs=1
>>>>> --download-scalapack=yes
>>>>> --with-scalapack=1 --download-spooles=yes --with-spooles=1
>>>>> --download-umfpack=yes --with-umfpack=1 --with-debugging=no
>>>>> --with-shared=1
>>>>> [0]PETSC ERROR:
>>>>>
>>>>> ------------------------------------------------------------------------
>>>>> [0]PETSC ERROR: PCFactorSetMatSolverPackage_Factor() line 183 in
>>>>> src/ksp/pc/impls/factor/factimpl.c
>>>>> [0]PETSC ERROR: PCFactorSetMatSolverPackage() line 276 in
>>>>> src/ksp/pc/impls/factor/factor.c
>>>>> [0]PETSC ERROR: PCSetFromOptions_Factor() line 275 in
>>>>> src/ksp/pc/impls/factor/factimpl.c
>>>>> [0]PETSC ERROR: PCSetFromOptions_ILU() line 112 in
>>>>> src/ksp/pc/impls/factor/ilu/ilu.c
>>>>> [0]PETSC ERROR: PCSetFromOptions() line 185 in
>>>>> src/ksp/pc/interface/pcset.c
>>>>> [0]PETSC ERROR: KSPSetFromOptions() line 323 in
>>>>> src/ksp/ksp/interface/itcl.c
>>>>> [0]PETSC ERROR: User provided function() line 696 in
>>>>> "unknowndirectory/"src/solvers/petsc_linear_solver.C
>>>>>
>>>>> Any idea what the problem may be? I would be very grateful, if you can
>>>>> give
>>>>> me any hint. thanks alot
>>>>>
>>>>> best regards
>>>>> Ping
>>>>>
>>>>>
>>>>>
>>>>> Am 27.07.2011 04:13, schrieb Hong Zhang:
>>>>>>
>>>>>> Did you call KSPSetFromOptions(ksp)?
>>>>>> Run your code with '-options_table' to dump list of options inputted
>>>>>> or '-options_left' to dump list of unused options.
>>>>>>
>>>>>> I tested with petsc-3.1/src/ksp/ksp/examples/tutorials/ex2.c:
>>>>>> ./ex2
>>>>>> ...
>>>>>> SuperLU run parameters:
>>>>>> ...
>>>>>> ILU_DropTol: 0.1
>>>>>> ILU_FillTol: 0.01
>>>>>> ILU_FillFactor: 10
>>>>>> ILU_DropRule: 9
>>>>>> ILU_Norm: 2
>>>>>> ILU_MILU: 2
>>>>>>
>>>>>> Hong
>>>>>>
>>>>>> On Tue, Jul 26, 2011 at 3:53 PM, Ping Rong<ping.rong at tuhh.de>
>>>>>> wrote:
>>>>>>>
>>>>>>> Dear developers,
>>>>>>>
>>>>>>> I have compiled petsc-3.1-p8 for a while. Now I would like to use
>>>>>>> superlu
>>>>>>> as
>>>>>>> an ilu preconditioner, since it offers the drop tolerance option. I
>>>>>>> have
>>>>>>> read in a thread
>>>>>>>
>>>>>>>
>>>>>>>
>>>>>>> (https://lists.mcs.anl.gov/mailman/htdig/petsc-users/2010-December/007439.html)
>>>>>>> that one can run the code with option
>>>>>>>
>>>>>>> -pc_type ilu -pc_factor_mat_solver_package superlu
>>>>>>> -mat_superlu_ilu_droptol
>>>>>>> <>
>>>>>>>
>>>>>>> to setup the ilu preconditioner. I also use the option "-help |grep
>>>>>>> superlu"
>>>>>>> to check the settings. However, no matter how I change the value of
>>>>>>> -mat_superlu_ilu_droptol, I always got the same result
>>>>>>> ...
>>>>>>> -mat_superlu_lwork<0>: size of work array in bytes used by
>>>>>>> factorization
>>>>>>> (None)
>>>>>>> -mat_superlu_ilu_droptol<0.0001>: ILU_DropTol (None)
>>>>>>> -mat_superlu_ilu_filltol<0.01>: ILU_FillTol (None)
>>>>>>> ...
>>>>>>> I dont know whether I understand it correctly. But it seems to me the
>>>>>>> value
>>>>>>> of the droptol has never been changed. In that maillist thread, it
>>>>>>> was
>>>>>>> also
>>>>>>> mentioned that the dev-version is recommended, because of some bugs
>>>>>>> in
>>>>>>> the
>>>>>>> superlu interface. I have then compiled the current dev-version. but
>>>>>>> the
>>>>>>> problem is my program is based on another library (libMesh) which
>>>>>>> uses
>>>>>>> petsc
>>>>>>> as a solver package. Since some of the interfaces have been changed
>>>>>>> in
>>>>>>> the
>>>>>>> dev-version, I would not be able to compile the libMesh with petsc
>>>>>>> anymore.
>>>>>>> Is there anyway I can correct the superlu interface in the 3.1-p8
>>>>>>> version
>>>>>>> directly? Any help will be appreciated!!
>>>>>>>
>>>>>>> Best regards
>>>>>>>
>>>>>>> --
>>>>>>> Ping Rong, M.Sc.
>>>>>>> Hamburg University of Technology
>>>>>>> Institut of modelling and computation
>>>>>>> Denickestraße 17 (Room 3031)
>>>>>>> 21073 Hamburg
>>>>>>>
>>>>>>> Tel.: ++49 - (0)40 42878 2749
>>>>>>> Fax: ++49 - (0)40 42878 43533
>>>>>>> Email: ping.rong at tuhh.de
>>>>>>>
>>>>>>>
>>>>> --
>>>>> Ping Rong, M.Sc.
>>>>> Hamburg University of Technology
>>>>> Institut of modelling and computation
>>>>> Denickestraße 17 (Room 3031)
>>>>> 21073 Hamburg
>>>>>
>>>>> Tel.: ++49 - (0)40 42878 2749
>>>>> Fax: ++49 - (0)40 42878 43533
>>>>> Email: ping.rong at tuhh.de
>>>>>
>>>>>
>>>
>>> --
>>> Ping Rong, M.Sc.
>>> Hamburg University of Technology
>>> Institut of modelling and computation
>>> Denickestraße 17 (Room 3031)
>>> 21073 Hamburg
>>>
>>> Tel.: ++49 - (0)40 42878 2749
>>> Fax: ++49 - (0)40 42878 43533
>>> Email: ping.rong at tuhh.de
>>>
>>>
>
>
> --
> Ping Rong, M.Sc.
> Hamburg University of Technology
> Institut of modelling and computation
> Denickestraße 17 (Room 3031)
> 21073 Hamburg
>
> Tel.: ++49 - (0)40 42878 2749
> Fax: ++49 - (0)40 42878 43533
> Email: ping.rong at tuhh.de
>
>
More information about the petsc-users
mailing list