diff --git a/src/equilibrium.f90 b/src/equilibrium.f90
index ad8cb61..43dedbb 100644
--- a/src/equilibrium.f90
+++ b/src/equilibrium.f90
@@ -1155,7 +1155,12 @@ contains
           f(1) = frhotor(x(1)) - rho_t
         else
           ! return f'(x), computed numerically
-          df(1,1) = (frhotor(x(1) + e) - frhotor(x(1) - e)) / (2*e)
+          if (x(1) - e > 0) then
+            df(1,1) = (frhotor(x(1) + e) - frhotor(x(1) - e)) / (2*e)
+          else
+            ! single-sided when close to ρ=0
+            df(1,1) = (frhotor(x(1) + e) - frhotor(x(1))) / e
+          end if
         end if
       end subroutine