diff --git a/NEWS.md b/NEWS.md
index f1f881d00..b7ee6576c 100644
--- a/NEWS.md
+++ b/NEWS.md
@@ -7,15 +7,11 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0
 
 ## [Unreleased]
 
-<<<<<<< fix/nlsolver-convergence-check
 ### Fixed
 
-- `NLSolver` now checks convergence status before returning solution. Previously, when `nlsolve()` failed to converge, the unconverged last iterate was silently returned. Now throws an informative error with iteration count and residual norm. Since PR[#1283](https://github.com/gridap/Gridap.jl/pull/1283). (@Ady0333)
-=======
-## Fixed
-
+- `NewtonRaphsonSolver` now correctly handles already-converged initial states and reports non-convergence properly. Previously crashed with `@unreachable` when starting residual was ~0. Now checks initial convergence and throws informative error instead of assertion failure. Since PR[#1285](https://github.com/gridap/Gridap.jl/pull/1285).
+- `NLSolver` now checks convergence status before returning solution. Previously, when `nlsolve()` failed to converge, the unconverged last iterate was silently returned. Now throws an informative error with iteration count and residual norm. Since PR[#1283](https://github.com/gridap/Gridap.jl/pull/1283).
 - Fix a bug where external types were not correctly parsed with generated functions in Polynomials, see Issue[#1286](https://github.com/gridap/Gridap.jl/issues/1286). Since PR[#1287](https://github.com/gridap/Gridap.jl/pull/1287).
->>>>>>> master
 
 ## [0.20.4] - 2026-04-23
 
diff --git a/src/Algebra/NLSolvers.jl b/src/Algebra/NLSolvers.jl
index cabffef1b..066a0769a 100644
--- a/src/Algebra/NLSolvers.jl
+++ b/src/Algebra/NLSolvers.jl
@@ -64,7 +64,7 @@ function _solve_nr!(x,A,b,dx,ns,nls,op)
     if isconv; return; end
 
     if nliter == nls.max_nliters
-      @unreachable
+      error("NewtonRaphsonSolver failed to converge in $(nls.max_nliters) iterations. Final residual: $(maximum(abs,b))")
     end
 
     # Assemble jacobian (fast in-place version)
@@ -78,7 +78,7 @@ end
 
 function _check_convergence(nls,b)
   m0 = maximum(abs,b)
-  return (false, m0)
+  return (m0 < nls.tol, m0)
 end
 
 function _check_convergence(nls,b,m0)