dolfin team mailing list archive

[Anders Logg <logg@xxxxxxxxx>]

On Wed, Mar 22, 2006 at 10:17:10AM -0600, Anders Logg wrote:
> On Wed, Mar 22, 2006 at 04:53:16PM +0100, Garth N. Wells wrote:
> > On Wed, 2006-03-22 at 08:47 -0600, Anders Logg wrote:
> > > On Wed, Mar 22, 2006 at 10:13:22AM +0000, Alexander Jarosch wrote:
> > > > I got a bit lost in the current updates and now I just have a short 
> > > > question. I set up my solver like this:
> > > > 
> > > >  // Discretize equation
> > > >  Matrix A;
> > > >  Vector x, b;
> > > >  FEM::assemble(*a, *L, A, b, mesh, bc);
> > > > 
> > > >  // Solve the linear system
> > > >  KrylovSolver solver(Preconditioner::ilu);
> > > >  KSP ksp;
> > > >  ksp = solver.solver();
> > > > 
> > > >  PC pc;
> > > >  KSPGetPC(ksp, &pc);
> > > >  PCFactorSetShiftNonzero(pc, PETSC_DECIDE);
> > > >  //PCFactorSetShiftPd(pc, PETSC_TRUE);
> > > > 
> > > >  solver.set("Krylov relative tolerance", 1e-14);
> > > > 
> > > >  solver.solve(A, x, b);
> > > > 
> > > > but Petsc still complains about zero Pivot. should the line 
> > > > "PCFactorSetShiftNonzero(pc, PETSC_DECIDE);" not take care of that.
> > > 
> > > The problem is that KrylovSolver::solve() reads all parameters from
> > > the (DOLFIN) database, including "Krylov relative tolerance", and sets
> > > them before each solve. The solve() function also sets the
> > > preconditioner. Therefore, setting options manually by PETSc calls
> > > might not always work.
> > > 
> > > The plan is to remove the KrylovSolver::solve() function if no one
> > > objects (Garth?). This will leave two options: either use the simple
> > > (but limited) DOLFIN wrapper or get the PETSc data structures and
> > > solve using PETSc calls:
> > > 
> > > Mat Amat = A.mat();
> > > Vec xvec = x.vec();
> > > 
> > > PETSc calls...
> > > 
> > 
> > Another option is to add PETSc options at the command line. This is
> > really easy,
> > 
> >   ./dolfin -pc_factor_shift_nonzero
> > 
> > Just make sure that you have 
> >   int main(int argc, char* argv[])
> >  	
> > 
> > > That said, I will (re-)add a parameter for "Krylov shift nonzero" so
> > > you can do this natively in DOLFIN.
> > > 
> > > I'll await comments from Garth and others before going ahead.
> > > 
> > > /Anders
> > > 
> > 
> > I think that removing KrylovSolver.solver() is the best approach from a
> > robustness point of view. KrylovSolver then doesn't have to worry that a
> > user might have messed with the PETSc solver directly. I do find the
> > direct access to ksp useful, but I can live without.
> > 
> > Garth  
> 
> ok, I'll make the update.
> 
> /Anders

ok, I have added the new option. Try it and see if it works:

    KrylovSolver solver(Preconditioner::ilu);
    solver.set("Krylov shift nonzero", whatever);
    solver.set("Krylov relative tolerance", 1e-14);
    solver.solve(A, x, b);

/Anders