diff --git a/CHANGELOG.md b/CHANGELOG.md
index c196ce1..b57826d 100644
--- a/CHANGELOG.md
+++ b/CHANGELOG.md
@@ -12,7 +12,7 @@ and this project adheres to [Semantic Versioning](http://semver.org/).
   need to change it to `BvhTraverseIterator<T, D, Ray, S>`. [#128](https://github.com/svenstaro/bvh/pull/128) (thanks @finnbear)
 - **Breaking change:** Distance-traversal no longer outputs non-intersected shapes, but note that
   `Ray::intersection_slice_for_aabb` now returns `None` instead of `(-1.0, -1.0)` in the case of no 
-  intersection, and `Some((entry, exit))` in the case of intersection. [#133](https://github.com/svenstaro/bvh/pull/133) (thanks @finnbear)
+  intersection, and `Some((entry, exit))` in the case of intersection. [#133](https://github.com/svenstaro/bvh/pull/133) [#142](https://github.com/svenstaro/bvh/pull/142) (thanks @finnbear)
 - Fix panic on empty `DistanceTraverseIterator` [#117](https://github.com/svenstaro/bvh/pull/117) (thanks @finnbear)
 - Fix center() for very large AABBs [#118](https://github.com/svenstaro/bvh/pull/118) (thanks @finnbear)
 - Fix more cases where an empty BVH would panic [#116](https://github.com/svenstaro/bvh/pull/116) (thanks @finnbear)
diff --git a/fuzz/Cargo.toml b/fuzz/Cargo.toml
index 764bcef..406427f 100644
--- a/fuzz/Cargo.toml
+++ b/fuzz/Cargo.toml
@@ -3,7 +3,7 @@ name = "bvh-fuzz"
 version = "0.0.0"
 authors = ["Automatically generated"]
 publish = false
-edition = "2018"
+edition = "2021"
 
 [package.metadata]
 cargo-fuzz = true
diff --git a/fuzz/fuzz_targets/fuzz.rs b/fuzz/fuzz_targets/fuzz.rs
index 492133b..66aec6f 100644
--- a/fuzz/fuzz_targets/fuzz.rs
+++ b/fuzz/fuzz_targets/fuzz.rs
@@ -322,8 +322,19 @@ impl<const D: usize> Workload<D> {
     ///
     /// The contents are explained in the module-level comment up above.
     fn fuzz(mut self) {
-        let mut bvh = Bvh::build(&mut self.shapes);
         let ray = self.ray.ray();
+        let aabb = self.aabb.aabb();
+
+        // Make sure these functions agree and that the latter doesn't return NaN.
+        let intersects = ray.intersects_aabb(&aabb);
+        if let Some((entry, exit)) = ray.intersection_slice_for_aabb(&aabb) {
+            assert!(intersects);
+            assert!(entry <= exit, "{entry} {exit}");
+        } else {
+            assert!(!intersects);
+        }
+
+        let mut bvh = Bvh::build(&mut self.shapes);
 
         if self.shapes.len()
             + self
@@ -367,18 +378,9 @@ impl<const D: usize> Workload<D> {
             bvh.assert_tight();
             let flat_bvh = bvh.flatten();
 
-            let traverse_ray = self.fuzz_traversal(
-                &bvh,
-                &flat_bvh,
-                &self.ray.ray(),
-                assert_ray_traversal_agreement,
-            );
-            self.fuzz_traversal(
-                &bvh,
-                &flat_bvh,
-                &self.aabb.aabb(),
-                assert_aabb_traversal_agreement,
-            );
+            let traverse_ray =
+                self.fuzz_traversal(&bvh, &flat_bvh, &ray, assert_ray_traversal_agreement);
+            self.fuzz_traversal(&bvh, &flat_bvh, &aabb, assert_aabb_traversal_agreement);
             self.fuzz_traversal(
                 &bvh,
                 &flat_bvh,
diff --git a/src/ray/ray_impl.rs b/src/ray/ray_impl.rs
index 9210797..5e79ac6 100644
--- a/src/ray/ray_impl.rs
+++ b/src/ray/ray_impl.rs
@@ -1,6 +1,8 @@
 //! This module defines a Ray structure and intersection algorithms
 //! for axis aligned bounding boxes and triangles.
 
+use std::cmp::Ordering;
+
 use crate::aabb::IntersectsAabb;
 use crate::utils::fast_max;
 use crate::{aabb::Aabb, bounding_hierarchy::BHValue};
@@ -124,15 +126,15 @@ impl<T: BHValue, const D: usize> Ray<T, D> {
 
         let (inf, sup) = lbr.inf_sup(&rtr);
 
-        let tmin = inf.max();
+        let tmin = fast_max(inf.max(), T::zero());
         let tmax = sup.min();
 
-        // no intersection
-        if tmin > tmax || tmax < T::zero() {
+        if matches!(tmin.partial_cmp(&tmax), Some(Ordering::Greater) | None) {
+            // tmin > tmax or either was NaN, meaning no intersection.
             return None;
         }
 
-        Some((fast_max(tmin, T::zero()), tmax))
+        Some((tmin, tmax))
     }
 
     /// Implementation of the