From b2d36ca67eb81ffec9e7d57ec51b46eda2df250c Mon Sep 17 00:00:00 2001
From: Phillip Weinberg <weinbe58@gmail.com>
Date: Tue, 2 Jan 2024 13:32:47 -0500
Subject: [PATCH] adding error code for when step size is `NaN`

---
 src/dp5/solver.jl |  9 ++++-----
 src/dp8/solver.jl | 11 +++++------
 src/types.jl      |  1 +
 test/debug.jl     | 19 +++++++++++++------
 4 files changed, 23 insertions(+), 17 deletions(-)

diff --git a/src/dp5/solver.jl b/src/dp5/solver.jl
index 3403161..6628fee 100644
--- a/src/dp5/solver.jl
+++ b/src/dp5/solver.jl
@@ -98,11 +98,10 @@ function dopri5(
 
     ###### Basic Integration Step
     for _ in 1:solver.options.maximum_allowed_steps
-        # if nstep > solver.options.maximum_allowed_steps
-        #     # GOTO 78
-        #     # println(" MORE THAN NMAX = ", solver.options.maximum_allowed_steps, " STEPS ARE NEEDED")
-        #     return h, Report(solver.vars.x, DormandPrince.INPUT_CHECKS_SUCCESSFUL, DormandPrince.LARGER_NMAX_NEEDED , 0, 0, 0, 0)
-        # end
+        if isnan(h)
+            idid = DormandPrince.STEP_SIZE_BECOMES_NAN
+            break
+        end
         
         if (0.10 * abs(h)) <= abs(solver.vars.x)*solver.options.uround 
             # GOTO 77
diff --git a/src/dp8/solver.jl b/src/dp8/solver.jl
index 57d53eb..9377eea 100644
--- a/src/dp8/solver.jl
+++ b/src/dp8/solver.jl
@@ -97,12 +97,11 @@ function dop853(
 
     ###### Basic Integration Step
     for _ in 1:solver.options.maximum_allowed_steps
-        # if nstep > solver.options.maximum_allowed_steps
-        #     # GOTO 78
-        #     # println(" MORE THAN NMAX = ", solver.options.maximum_allowed_steps, " STEPS ARE NEEDED")
-        #     return h, Report(solver.vars.x, DormandPrince.INPUT_CHECKS_SUCCESSFUL, DormandPrince.LARGER_NMAX_NEEDED , 0, 0, 0, 0)
-        # end
-        
+        if isnan(h)
+            idid = DormandPrince.STEP_SIZE_BECOMES_NAN
+            break
+        end
+
         if (0.10 * abs(h)) <= abs(solver.vars.x)*solver.options.uround 
             # GOTO 77
             # println("STEP SIZE TOO SMALL, H = ", h)
diff --git a/src/types.jl b/src/types.jl
index 5e2149f..1ad85e2 100644
--- a/src/types.jl
+++ b/src/types.jl
@@ -6,6 +6,7 @@ abstract type AbstractDPSolver{T <: Real, StateType <: AbstractVector, F} end
     INPUT_NOT_CONSISTENT = -1 # use for check failures in the beginning of  core_integrator call
     LARGER_NMAX_NEEDED = -2
     STEP_SIZE_BECOMES_TOO_SMALL = -3
+    STEP_SIZE_BECOMES_NAN = -4
 end
 
 @enum Checks begin
diff --git a/test/debug.jl b/test/debug.jl
index cab8887..2a61684 100644
--- a/test/debug.jl
+++ b/test/debug.jl
@@ -1,5 +1,5 @@
 using LinearAlgebra
-using DormandPrince: DP5Solver, DP8Solver, integrate
+using DormandPrince: DP8Solver, DP5Solver, integrate
 
 function evolution_operator(t::Float64)
     ϕ = 2.2 * sin(π * t)^2
@@ -33,13 +33,20 @@ function run()
         ComplexF64[1.0, 0.0]
     )
     
-    report = integrate(solver, 1.0)
+    report = integrate(solver, 2π)
     
-    println(report.num_rejected_steps)
-    println(norm(solver.y))
+    if !(solver.y ≈ solution(2π))
+        println("failed, $(solver.y) != $(solution(2π)), $(report.idid)")
+    end
+    @assert solver.y ≈ solution(2π)
+    # println(report.num_rejected_steps)
+    # println(norm(solver.y))
     solver.vars.h
-    @assert solver.y ≈ solution(1.0)
+    solver.y
 
 end
 
-run()
\ No newline at end of file
+for i in 1:1000
+    run()
+    println("-------------------------------------------------------------------------------")
+end
\ No newline at end of file