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