On Tue, Jun 7, 2011 at 10:44 AM, Danesh Daroui <danesh.daroui at ltu.se> wrote:
> > So, it is not implemented? That's too bad! ILU with drop threshold is > very useful so I am wondering why it is not implemented yet? Do you have > any plan to implement it? > It is not very scalable, so it is only really useful in a serial code. We tend to focus on parallel algorithms which are scalable. We don't have any plans to implement it, but we do take code submissions. > Setting levels did work OK for me. I have another question. If I want to > calculate the preconditioner once, and use it later to solve other > equations, since I know that the preconditioner would be OK because I > solve several equations in several iterations and I know that my > coefficient matrix will not change that much so calculating the > preconditioner at once would be OK. What routines should I call once and > how can I apply the pre-calculated preconditioner each time I want to > solve the equation? > You would use KSPSetOperators() with SAME_PRECONDITIONER. Matt > > Thanks, > > D. > > > > On Tue, 2011-06-07 at 10:28 -0500, Barry Smith wrote: > > Sorry for the confusion, it is our fault. We don't have a drop tolerance > ILU hence PCFactorSetDropTolerance() isn't doing anything. > > > > You can run with -pc_factor_levels M or call PCFactorSetLevels(pc,M) > in the code where M is the size of the matrix. > > > > Barry > > > > On Jun 7, 2011, at 10:20 AM, Danesh Daroui wrote: > > > > > > > > Hi, > > > > > > I want to use incomplete LU preconditioner. My code works OK when I use > > > no preconditioner and direct solver for test: > > > > > > > > > ierr=KSPCreate(PETSC_COMM_WORLD, &ksp); > > > ierr=KSPSetOperators(ksp, Mp, Mp, DIFFERENT_NONZERO_PATTERN); > > > ierr=KSPSetTolerances(ksp, 1.e-2/Msize, 1.e-50, PETSC_DEFAULT, > > > PETSC_DEFAULT); > > > > > > // GMRES iterative solver is used > > > ierr=KSPSetType(ksp, KSPGMRES); > > > > > > // no PETSc preconditioner is used > > > ierr=KSPGetPC(ksp, &prec); > > > ierr=PCSetType(prec, PCLU); > > > > > > // set up the solver according to the options > > > ierr=KSPSetFromOptions(ksp); > > > > > > ierr=KSPSetUp(ksp); > > > > > > // solve the equation using an iterative solver > > > ierr=KSPSolve(ksp, bp, xp); > > > > > > > > > Again for some tests, at the first step, I want to use ILU but create > > > the exact LU and then change the parameters to tune the preconditioner > > > and I almost know what values to set according to my tests in MATLAB. I > > > use the code below to create the get exact LU, using ILU > preconditioner: > > > > > > > > > > > > > > > > > > ierr=KSPCreate(PETSC_COMM_WORLD, &ksp); > > > ierr=KSPSetOperators(ksp, Mp, Mp, DIFFERENT_NONZERO_PATTERN); > > > ierr=KSPSetTolerances(ksp, 1.e-2/Msize, 1.e-50, PETSC_DEFAULT, > > > PETSC_DEFAULT); > > > > > > // GMRES iterative solver is used > > > ierr=KSPSetType(ksp, KSPGMRES); > > > > > > // no PETSc preconditioner is used > > > ierr=KSPGetPC(ksp, &prec); > > > ierr=PCSetType(prec, PCILU); > > > > > > PetscReal dt=1e-10; > > > PetscReal dtcol=0.1; > > > PetscInt maxrowcount=M; > > > PCFactorSetDropTolerance(prec, dt, dtcol, maxrowcount); > > > > > > // set up the solver according to the options > > > ierr=KSPSetFromOptions(ksp); > > > > > > ierr=KSPSetUp(ksp); > > > > > > // solve the equation using an iterative solver > > > ierr=KSPSolve(ksp, bp, xp); > > > > > > > > > and my coefficient matrix is a square "MxM" matrix and I pass "M" as > > > "maxrowcount". Using this code, the solver never converges while it > > > solved the equation using previous code (direct solver). Can anybody > let > > > me know how can I get exact LU using ILU and what parameters should I > > > change in the above code? > > > > > > Regards, > > > > > > D. > > > > > > > > > > > > > -- 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/20110607/ff246b8c/attachment.htm>
