[petsc-users] use Superlu as ilu preconditioner
Ping Rong
ping.rong at tuhh.de
Tue Aug 2 13:50:51 CDT 2011
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
-------------- next part --------------
An embedded and charset-unspecified text was scrubbed...
Name: libMesh-petsc_linear_solver.cpp
URL: <http://lists.mcs.anl.gov/pipermail/petsc-users/attachments/20110802/a22b2d59/attachment-0001.diff>
More information about the petsc-users
mailing list