Suggestions on PETSc SNES Nonlinear Solver?

[Lex Lee]

Hello all, 

In my FSI problem, there are 4 vector-valued variables + 3 scalar variables. Most of them are strongly coupled with each other. It's hard to decouple the system and design a good preconditoner. Previously, I solved the whole system with the direct solver (UMFPACK)  + Newton's iteration method. 

Several days ago, I added deal.ii's new nonlinear solver (SNES) into my code, and now its default method is jacobian-free newton-krylov. However, it behaves bad and never get converged in my case. Maybe I am using it in a wrong way, as I heard that SNES solver is powerful.

Do you have the experience in using this nonlinear solver? If yes, could you kindly give me some advice for multiple-variable problem? Do I need to create block containers for the matrix, solution and right-hand-side ?

Thanks a lot

Lex  

constraints_hp.distribute(dis_current_solution);
current_solution = dis_current_solution;

const double target_tolerance = 1.e-8;
PETScWrappers::NonlinearSolverData additional_data;
additional_data.absolute_tolerance = target_tolerance;
additional_data.snes_linesearch_type = "basic";
//additional_data.snes_type = "ngmres" ;

NonlinearSolver nonlinear_solver(additional_data,mpi_communicator);
nonlinear_solver.residual = [&](const VectorType &evaluation_point,
VectorType & residual){
compute_hp_residual(evaluation_point, residual);
};

bool user_control = false;
if (user_control)
{
nonlinear_solver.setup_jacobian =
[&](const VectorType &current_u) {
compute_jacobian_and_initialize_preconditioner(current_u);
};
nonlinear_solver.solve_with_jacobian = [&](const VectorType &rhs,
VectorType &dst) {
this->solve(rhs, dst);
};
//nonlinear_solver.set_matrix(jacobian_matrix);
}
else //jacobian-free newton-krylov
{
nonlinear_solver.set_matrix(jacobian_matrix);
nonlinear_solver.jacobian =
[&](const VectorType &current_u, MatrixType &, MatrixType &P) {
Assert(P == jacobian_matrix, ExcInternalError());
compute_hp_jacobian_matrix(current_u);
(void)P;
};
}
nonlinear_solver.monitor =
[&](const VectorType &, unsigned int step, double gnorm) {
pcout << step << " norm=" << gnorm << std::endl;
};

nonlinear_solver.solve(dis_current_solution);
constraints_hp.distribute(dis_current_solution);
current_solution = dis_current_solution;

[Wolfgang Bangerth]

Lex,

> In my FSI problem, there are 4 vector-valued variables + 3 scalar variables.
> Most of them are strongly coupled with each other. It's hard to decouple the
> system and design a good preconditoner. Previously, I solved the whole system
> with the direct solver (UMFPACK)  + Newton's iteration method.
>
> Several days ago, I added deal.ii's new nonlinear solver (SNES) into my code,
> and now its default method is jacobian-free newton-krylov. However, it behaves
> bad and never get converged in my case. Maybe I am using it in a wrong way, as
> I heard that SNES solver is powerful.

My first question is *what* doesn't converge. Is it the linear solver, or the
nonlinear solver?

I am not too familiar with PETSc's SNES, but it should be relatively easy to
use the KINSOL solvers (see step-77) instead where the choice of linear solver
is independent of the nonlinear solver. So you can continue to use the direct
solver for the linear systems (solving with the Jacobian) and still get the
benefit from using an excellent nonlinear solver package. If that works, you
can then independently think about using a better linear solver.

Best
W.

--
------------------------------------------------------------------------
Wolfgang Bangerth email: bang...@colostate.edu
www: http://www.math.colostate.edu/~bangerth/