Fix 0/0 = NaN in the rise/set parabola fit

skyfield/almanac.py

@@ -336,7 +336,11 @@ def _transit_ha(latitude, declination, altitude_radians):
 
 def _q(a, b, c, sign):
     discriminant = np.maximum(b*b - 4*a*c, 0.0)  # avoid tiny negative results
-    return - 2*c / (b + sign * sqrt(discriminant))
+    with np.errstate(invalid='ignore', divide='ignore'):
+        x = - 2*c / (b + sign * sqrt(discriminant))
+    # If c is zero, then x=0 is exactly a root, but the quotient above
+    # can be 0/0 = NaN when the denominator is also zero (see #1114).
+    return np.where(c == 0.0, 0.0, x)
 
 def _intersection(y0, y1, v0, v1):
     # Return x at which a curve reaches y=0, given its position and

skyfield/tests/test_almanac.py

@@ -136,6 +136,15 @@ def test_dark_twilight_day():
 
 # Logic.
 
+def test_intersection_when_curve_starts_exactly_at_zero():
+    # If the rise/set Newton iteration happens to land, to the last bit,
+    # exactly on the horizon, then the final parabola fit was dividing
+    # 0/0 = NaN (#1114).  These operands were captured live from a Mars
+    # setting at latitude 52.0N on 2029-11-04.
+    v = -1.5990198773309593e-16
+    x = almanac._intersection(0.0, 0.0, v, v)
+    assert x == 0.0
+
 def test_close_start_and_end():
     ts = api.load.timescale()
     t0 = ts.utc(2018, 9, 23, 1)