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 > >
