I was not assuming that the first iteration in the N-R method for the Poisson equation was fully quasi-neutral becuase db_dphi was calculated as if we had a distribution of electrons different from ions. It should be fixed now