I'm trying to numerically solve a simple Laplace equation in 2D, with a nonlinear source term:

$\nabla^2 u = u^2$

with boundary conditions as $u=0$ everywhere except for $y=1$ where $u=u_0$. I'm using scipy's newton_krylov solver (method="lgmres") to minimized the discretized equation using the finite-difference method.

Here's the problem: for relatively small values of $u_0$ such as $u_0=1$, LGMRES converges relatively fast with fewer than 7 iterations. However, when I increase $u_0$ to $10$ or $100$, the the solver converges much slower and needs on the order of O(1000) iterations to converge.

Is this an expected behaviour of Newton-Krylov solvers? If so, what can I do to alleviate the issue?

  • $\begingroup$ Sounds like a conditioning issue. Try preconditioning with a Fast Poisson Solver and that should help a lot $\endgroup$ – whpowell96 Jul 15 '19 at 10:38
  • $\begingroup$ When the problem is linear does it present the same problem? $\endgroup$ – nicoguaro Jul 15 '19 at 14:14
  • $\begingroup$ @nicoguaro, when I omit the nonlinear term, it's much better. But, if I keep increasing the Dirichlet BC, scipy eventually spits out the following error: Jacobian inversion yielded zero vector. This indicates a bug in the Jacobian approximation. $\endgroup$ – Blademaster Jul 15 '19 at 16:38

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Browse other questions tagged or ask your own question.