diff --git a/collects/math/private/flonum/expansion/expansion-base.rkt b/collects/math/private/flonum/expansion/expansion-base.rkt
index b597b9785f..44af901f75 100644
--- a/collects/math/private/flonum/expansion/expansion-base.rkt
+++ b/collects/math/private/flonum/expansion/expansion-base.rkt
@@ -10,16 +10,20 @@ Discrete & Computational Geometry 18(3):305–363, October 1997
 Other parts shamelessly stolen from crlibm (which is LGPL)
 |#
 
-(require "../flonum-functions.rkt"
+(require racket/math
+         "../flonum-functions.rkt"
          "../flonum-bits.rkt"
          "../flonum-error.rkt"
          "../flonum-constants.rkt"
          "../utils.rkt")
 
-(provide fl2? fl2 fl2->real
-         fl2ulp fl2ulp-error
-         +max-fl2-subnormal.0 -max-fl2-subnormal.0
-         fl2abs fl2+ fl2- fl2*split-fl fl2* fl2sqr fl2/
+(provide fl2? fl2zero? fl2rational? fl2positive? fl2negative? fl2infinite? fl2nan?
+         fl2 fl2->real
+         fl2ulp fl2ulp-error fl2step fl2next fl2prev
+         +max.hi +max.lo -max.hi -max.lo
+         +max-subnormal.hi -max-subnormal.hi
+         fl2abs fl2+ fl2- fl2= fl2> fl2< fl2>= fl2<=
+         fl2*split-fl fl2* fl2sqr fl2/
          fl2sqrt flsqrt/error)
 
 (: floverlapping? (Flonum Flonum -> Boolean))
@@ -38,30 +42,58 @@ Other parts shamelessly stolen from crlibm (which is LGPL)
   (cond [(flrational? x2)
          (cond [(flrational? x1)
                 (cond [((flabs x2) . < . (flabs x1))  #f]
-                      [else  (not (floverlapping? x2 x1))])]
+                      [else (not (floverlapping? x2 x1))])]
                [else  #f])]
-        [else  (fl= x1 0.0)]))
+        [else
+         (fl= x1 0.0)]))
+
+(define-syntax-rule (define-simple-fl2-predicate fl2pred? flpred?)
+  (begin
+    (: fl2pred? (Flonum Flonum -> Boolean))
+    (define (fl2pred? x2 x1)
+      (flpred? (fl+ x2 x1)))))
+
+(define-simple-fl2-predicate fl2zero? (λ (x) (fl= x 0.0)))
+(define-simple-fl2-predicate fl2positive? (λ (x) (fl> x 0.0)))
+(define-simple-fl2-predicate fl2negative? (λ (x) (fl< x 0.0)))
+(define-simple-fl2-predicate fl2rational? flrational?)
+(define-simple-fl2-predicate fl2nan? flnan?)
+(define-simple-fl2-predicate fl2infinite? flinfinite?)
 
 ;; ===================================================================================================
 ;; Conversion
 
-(: fl2 (Real -> (Values Flonum Flonum)))
-(define (fl2 x)
-  (cond [(flonum? x)  (values x 0.0)]
-        [(single-flonum? x)  (values (fl x) 0.0)]
-        [else  (let ([x2  (fl x)])
-                 (cond [(flrational? x2)
-                        (let ([x2  (fl+ x2 (fl (- x (inexact->exact x2))))])
-                          (cond [(flrational? x2)
-                                 (values x2 (fl (- x (inexact->exact x2))))]
-                                [else  (values x2 0.0)]))]
-                       [else  (values x2 0.0)]))]))
+(: fl2 (case-> (Real -> (Values Flonum Flonum))
+               (Flonum Flonum -> (Values Flonum Flonum))))
+(define fl2
+  (case-lambda
+    [(x)
+     (cond [(flonum? x)  (values x 0.0)]
+           [(single-flonum? x)  (values (fl x) 0.0)]
+           [else
+            (define x2 (fl x))
+            (if (flinfinite? x2)
+                (values x2 0.0)
+                (let* ([x  (- x (inexact->exact x2))]
+                       [x1  (fl x)]
+                       [x  (- x (inexact->exact x1))])
+                  (let-values ([(x2 x1)  (fl+/error x2 x1)])
+                    (values x2 (fl+ x1 (fl x))))))])]
+    [(x2 x1)
+     (if (and (fl= x2 0.0) (fl= x1 0.0))
+         (values x2 0.0)
+         (fl+/error x2 x1))]))
+
+(: fl2eqv? (case-> (Flonum Flonum Flonum -> Boolean)
+                   (Flonum Flonum Flonum Flonum -> Boolean)))
+(define (fl2eqv? x2 x1 y2 [y1 0.0])
+  (and (eqv? x2 y2) (fl= x1 y1)))
 
 (: fl2->real (Flonum Flonum -> Real))
 (define (fl2->real x2 x1)
-  (cond [(and (flrational? x2) (flrational? x1))
-         (+ (inexact->exact x2) (inexact->exact x1))]
-        [else  (fl+ x1 x2)]))
+  (if (flrational? x2)
+      (+ (inexact->exact x2) (inexact->exact x1))
+      x2))
 
 (: fl4->fl2 (Flonum Flonum Flonum Flonum -> (Values Flonum Flonum)))
 (define (fl4->fl2 e4 e3 e2 e1)
@@ -73,7 +105,7 @@ Other parts shamelessly stolen from crlibm (which is LGPL)
 (: fl2ulp (Flonum Flonum -> Flonum))
 (define (fl2ulp x2 x1)
   (cond [(fl= x2 0.0)  0.0]
-        [else  (flmax +min.0 (fl* (flulp x2) (flexpt 2.0 -52.0)))]))
+        [else  (flmax +min.0 (fl* (flulp x2) epsilon.0))]))
 
 (: fl2ulp-error (Flonum Flonum Real -> Flonum))
 (define (fl2ulp-error x2 x1 r)
@@ -89,20 +121,47 @@ Other parts shamelessly stolen from crlibm (which is LGPL)
                        (inexact->exact (flmax +min.0 (fl2ulp r2 r1))))))]
         [else  +inf.0]))
 
-(define +max-fl2-subnormal.0 1.0020841800044863e-292)  ; (flprev (flexpt 2.0 -970.0))
-(define -max-fl2-subnormal.0 (- +max-fl2-subnormal.0))
+(define-values (+max.hi +max.lo)
+  (values +max.0 (flprev (* 0.5 (flulp +max.0)))))
+
+(define-values (-max.hi -max.lo)
+  (values (- +max.hi) (- +max.lo)))
+
+(: fl2step (Flonum Flonum Integer -> (Values Flonum Flonum)))
+(define (fl2step x2 x1 n)
+  (let-values ([(x2 x1)  (fast-fl+/error x2 x1)])
+    (cond [(flnan? x2)  (values +nan.0 0.0)]
+          [(fl= x2 +inf.0)  (fl+/error +max.hi (flstep +max.lo (+ n 1)))]
+          [(fl= x2 -inf.0)  (fl+/error -max.hi (flstep -max.lo (- n 1)))]
+          [else  (fl+/error x2 (flstep x1 n))])))
+
+(: fl2next (Flonum Flonum -> (Values Flonum Flonum)))
+(define (fl2next x2 x1) (fl2step x2 x1 1))
+
+(: fl2prev (Flonum Flonum -> (Values Flonum Flonum)))
+(define (fl2prev x2 x1) (fl2step x2 x1 -1))
+
+(define +min-normal.hi (fl/ (flnext +max-subnormal.0) epsilon.0))
+
+(define-values (+max-subnormal.hi +max-subnormal.lo)
+  (fl2prev +min-normal.hi 0.0))
+
+(define-values (-max-subnormal.hi -max-subnormal.lo)
+  (values (- +max-subnormal.hi) (- +max-subnormal.lo)))
 
 ;; ===================================================================================================
 ;; Absolute value
 
 (: fl2abs (case-> (Flonum -> (Values Flonum Flonum))
                   (Flonum Flonum -> (Values Flonum Flonum))))
-(define (fl2abs x2 [x1 0.0])
-  (define x (fl+ x2 x1))
-  (cond [(not (flrational? x))  (values (flabs x2) 0.0)]
-        [(fl= x 0.0)  (values 0.0 0.0)]
-        [(fl> x 0.0)  (values x2 x1)]
-        [else  (values (- x2) (- x1))]))
+(define fl2abs
+  (case-lambda
+    [(x)  (values (flabs x) 0.0)]
+    [(x2 x1)
+     (cond [(flnan? x2)  (values +nan.0 0.0)]
+           [(fl= x2 0.0)  (values 0.0 0.0)]
+           [(fl> x2 0.0)  (values x2 x1)]
+           [else  (values (- x2) (- x1))])]))
 
 ;; ===================================================================================================
 ;; Addition and subtraction
@@ -125,6 +184,22 @@ Other parts shamelessly stolen from crlibm (which is LGPL)
 (define (fl2- x2 x1 y2 [y1 0.0])
   (fl2+ x2 x1 (- y2) (- y1)))
 
+;; ===================================================================================================
+;; Comparison
+
+(define-syntax-rule (define-fl2-comparison name flcomp)
+  (begin
+    (: name (Flonum Flonum Flonum Flonum -> Boolean))
+    (define (name x2 x1 y2 y1)
+      (let-values ([(z2 z1)  (fl2- x2 x1 y2 y1)])
+        ((fl+ z2 z1) . flcomp . 0.0)))))
+
+(define-fl2-comparison fl2= fl=)
+(define-fl2-comparison fl2> fl>)
+(define-fl2-comparison fl2< fl<)
+(define-fl2-comparison fl2>= fl>=)
+(define-fl2-comparison fl2<= fl<=)
+
 ;; ===================================================================================================
 ;; Multiplication and division
 
@@ -164,10 +239,13 @@ Other parts shamelessly stolen from crlibm (which is LGPL)
                 (Flonum Flonum Flonum Flonum -> (Values Flonum Flonum))))
 (define (fl2* x2 x1 y2 [y1 0.0])
   (define z (fl* x2 y2))
-  (cond [(and (flrational? z) (not (fl= z 0.0)) (not (flsubnormal? z)))
+  (cond [(fl= z 0.0)  (values z 0.0)]
+        [(flsubnormal? z)  (values z 0.0)]
+        [(and (flrational? x2) (flrational? y2) (z . fl>= . -inf.0) (z . fl<= . +inf.0))
          (define dx (near-pow2 x2))
          (define dy (near-pow2 y2))
          (define d (fl* dx dy))
+         (define d? (and (d . fl> . 0.0) (d . fl< . +inf.0)))
          (let* ([x2  (fl/ x2 dx)]
                 [x1  (fl/ x1 dx)]
                 [y2  (fl/ y2 dy)]
@@ -187,8 +265,8 @@ Other parts shamelessly stolen from crlibm (which is LGPL)
                           (fl* x1 y2))]
                 [z2  (fl+ m2 m1)]
                 [z1  (fl+ (fl- m2 z2) m1)]
-                [z2  (fl* d z2)])
-           (values z2 (if (flrational? z2) (fl* d z1) 0.0)))]
+                [z2  (if d? (fl* z2 d) (fl* (fl* z2 dx) dy))])
+           (values z2 (if (flrational? z2) (if d? (fl* z1 d) (fl* (fl* z1 dx) dy)) 0.0)))]
         [else
          (values z 0.0)]))
 
@@ -200,11 +278,14 @@ Other parts shamelessly stolen from crlibm (which is LGPL)
     [(x)  (flsqr/error x)]
     [(x2 x1)
      (define z (fl* x2 x2))
-     (cond [(and (flrational? z) (not (fl= z 0.0)) (not (flsubnormal? z)))
-            (define d (near-pow2 x2))
-            (define d^2 (fl* d d))
-            (let* ([x2  (fl/ x2 d)]
-                   [x1  (fl/ x1 d)]
+     (cond [(fl= z 0.0)  (values z 0.0)]
+           [(flsubnormal? z)  (values z 0.0)]
+           [(and (flrational? x2) (z . fl>= . -inf.0) (z . fl<= . +inf.0))
+            (define dx (near-pow2 x2))
+            (define d (fl* dx dx))
+            (define d? (and (d . fl> . 0.0) (d . fl< . +inf.0)))
+            (let* ([x2  (fl/ x2 dx)]
+                   [x1  (fl/ x1 dx)]
                    [up  (fl* x2 (fl+ 1.0 (flexpt 2.0 27.0)))]
                    [u1  (fl+ (fl- x2 up) up)]
                    [u2  (fl- x2 u1)]
@@ -215,8 +296,8 @@ Other parts shamelessly stolen from crlibm (which is LGPL)
                              (fl* 2.0 (fl* x2 x1)))]
                    [z2  (fl+ m2 m1)]
                    [z1  (fl+ (fl- m2 z2) m1)]
-                   [z2  (fl* d^2 z2)])
-              (values z2 (if (flrational? z2) (fl* d^2 z1) 0.0)))]
+                   [z2  (if d? (fl* z2 d) (fl* (fl* z2 dx) dx))])
+              (values z2 (if (flrational? z2) (if d? (fl* z1 d) (fl* (fl* z1 dx) dx)) 0.0)))]
            [else
             (values z 0.0)])]))
 
@@ -245,8 +326,8 @@ Other parts shamelessly stolen from crlibm (which is LGPL)
 (define (fl2sqrt x2 [x1 0.0])
   (cond [(and (flrational? x2) (not (fl= x2 0.0)))
          (define-values (d^2 d)
-           (cond [(x2 . fl<= . +max-fl2-subnormal.0)  (values (flexpt 2.0 -104.0)
-                                                              (flexpt 2.0 -52.0))]
+           (cond [(x2 . fl<= . +max-subnormal.hi)  (values (flexpt 2.0 -104.0)
+                                                           (flexpt 2.0 -52.0))]
                  [(x2 . fl> . 1e300)  (values (flexpt 2.0 104.0)
                                               (flexpt 2.0 52.0))]
                  [else  (values 1.0 1.0)]))
diff --git a/collects/math/private/flonum/expansion/expansion-exp-reduction.rkt b/collects/math/private/flonum/expansion/expansion-exp-reduction.rkt
index 7a5dc7c388..6ba36a5356 100644
--- a/collects/math/private/flonum/expansion/expansion-exp-reduction.rkt
+++ b/collects/math/private/flonum/expansion/expansion-exp-reduction.rkt
@@ -88,6 +88,7 @@ negative.
                     (cons y2 y1)))))
 |#
 
+(: positive-ds (Vectorof (Pair Flonum Flonum)))
 (define positive-ds
   #((0.00048511293333186416 . 7.287053259455084e-21)
     (0.0005002765039165264 . 5.174812552711603e-20)
@@ -553,6 +554,7 @@ negative.
     (1.0134858165681662e+305 . -1.1960732643227e+288)
     (+inf.0 . -inf.0)))
 
+(: negative-ds (Vectorof (Pair Flonum Flonum)))
 (define negative-ds
   #((-0.0004848777128822606 . -2.106959040416988e-20)
     (-0.0005000263524810411 . 1.7625505823099235e-20)
diff --git a/collects/math/private/flonum/expansion/expansion-exp.rkt b/collects/math/private/flonum/expansion/expansion-exp.rkt
index a2df696b2b..45aa069377 100644
--- a/collects/math/private/flonum/expansion/expansion-exp.rkt
+++ b/collects/math/private/flonum/expansion/expansion-exp.rkt
@@ -72,14 +72,13 @@ Note to self: exp(x+y)-1 = (exp(x)-1) * (exp(y)-1) + (exp(x)-1) + (exp(y)-1)
 
 ;; See "expansion-exp-reduction.rkt" for details on the argument reduction here
 
-(: flexp/error (Flonum -> (Values Flonum Flonum)))
-(define (flexp/error x)
-  (let: loop : (Values Flonum Flonum) ([x : Flonum  x] [n : Nonnegative-Fixnum  0])
+(: flexp/positive/error (Flonum -> (Values Flonum Flonum)))
+(define (flexp/positive/error x)
+  (let: loop ([x x] [n : Nonnegative-Fixnum  0])
     (cond
       [(or ((flabs x) . fl< . exp-min) (n . fx> . 15))  ; even n > 5 should never happen
        (let-values ([(y2 y1)  (flexpm1-small/error x)])
          (fl2+ y2 y1 1.0))]
-      [(x . fl< . (- (fllog +max.0)))  (values (flexp x) 0.0)]
       [(x . fl> . (fllog +max.0))  (values +inf.0 0.0)]
       [(rational? x)
        (let*-values ([(x d2 d1)  (flexpm1-reduction x)]
@@ -87,14 +86,29 @@ Note to self: exp(x+y)-1 = (exp(x)-1) * (exp(y)-1) + (exp(x)-1) + (exp(y)-1)
                      [(y2 y1)  (loop x (fx+ n 1))])
          (fl2* d2 d1 y2 y1))]
       [else
-       (values +nan.0 +nan.0)])))
+       (values +nan.0 0.0)])))
+
+(: flexp/error (Flonum -> (Values Flonum Flonum)))
+(define (flexp/error x)
+  (cond [(x . fl> . 0.0)  (flexp/positive/error x)]
+        [(x . fl= . 0.0)  (values 1.0 0.0)]
+        [(x . fl> . (- (fllog +max.0)))
+         (let*-values ([(y2 y1)  (flexp/positive/error (- x))])
+           (fl2/ 1.0 0.0 y2 y1))]
+        [(x . fl>= . -inf.0)  (values (flexp x) 0.0)]
+        [else  (values +nan.0 0.0)]))
 
 (: fl2exp (Flonum Flonum -> (Values Flonum Flonum)))
 (define (fl2exp x2 x1)
-  (let*-values ([(a2 a1)  (flexp/error x2)]
-                [(b2 b1)  (flexpm1-tiny/error x1)]
-                [(b2 b1)  (fl2+ b2 b1 1.0)])
-    (fl2* a2 a1 b2 b1)))
+  (cond [(x2 . fl> . -746.0)
+         (let*-values ([(a2 a1)  (flexp/error x2)]
+                       [(b2 b1)  (flexpm1-tiny/error x1)]
+                       [(b2 b1)  (fl2+ b2 b1 1.0)])
+           (fl2* a2 a1 b2 b1))]
+        [(x2 . fl>= . -inf.0)
+         (values 0.0 0.0)]
+        [else
+         (values +nan.0 0.0)]))
 
 ;; ===================================================================================================
 ;; expm1
@@ -118,12 +132,8 @@ The computation of R+1 definitely loses precision, but it doesn't seem to matter
 (: flexpm1/error (Flonum -> (Values Flonum Flonum)))
 (define (flexpm1/error x)
   (cond
-    [(x . < . -1.0)
-     ;; exp(x) is near zero here, so this is more accurate
-     (let-values ([(y2 y1)  (flexp/error x)])
-       (fl2+ y2 y1 -1.0))]
-    [else
-     (let: loop : (Values Flonum Flonum) ([x : Flonum  x] [n : Nonnegative-Fixnum  0])
+    [(x . fl> . -1.0)
+     (let: loop ([x x] [n : Nonnegative-Fixnum  0])
        (cond [(or ((flabs x) . fl< . expm1-min) (n . fx> . 15))  ; even n > 5 should never happen
               (flexpm1-small/error x)]
              [(x . fl< . (- (fllog +max.0)))  (values -1.0 0.0)]
@@ -135,14 +145,29 @@ The computation of R+1 definitely loses precision, but it doesn't seem to matter
                             [(y2 y1)  (fl2* d2 d1 w2 w1)])
                 (fl2+ y2 y1 r2 r1))]
              [else
-              (values +nan.0 +nan.0)]))]))
+              (values +nan.0 0.0)]))]
+    [(x . fl> . (fllog (* 0.25 epsilon.0)))
+     ;; exp(x) is near zero here, so this is more accurate
+     (let-values ([(y2 y1)  (flexp/error x)])
+       (fl2+ y2 y1 -1.0))]
+    [else
+     ;; The high-order flonum is -1 here, so return -1 + exp(x) directly
+     (values -1.0 (flexp x))]))
+
+(define-values (expm1-max-hi expm1-max-lo)
+  (values 709.782712893384 2.3691528222554853e-14))
 
 (: fl2expm1 (Flonum Flonum -> (Values Flonum Flonum)))
 (define (fl2expm1 x2 x1)
+  (define x (fl+ x2 x1))
   (cond
-    [((fl+ x1 x2) . < . -1.0)
+    [(x . fl< . -1.0)
      (let-values ([(y2 y1)  (fl2exp x2 x1)])
        (fl2+ y2 y1 -1.0))]
+    [(x2 x1 . fl2>= . expm1-max-hi expm1-max-lo)
+     (values +inf.0 0.0)]
+    [(x . fl= . 0.0)
+     (values x2 0.0)]
     [else
      (let*-values ([(a2 a1)  (flexpm1/error x2)]
                    [(b2 b1)  (flexpm1-tiny/error x1)]
diff --git a/collects/math/private/flonum/expansion/expansion-log.rkt b/collects/math/private/flonum/expansion/expansion-log.rkt
index 19e697f850..e91af5413e 100644
--- a/collects/math/private/flonum/expansion/expansion-log.rkt
+++ b/collects/math/private/flonum/expansion/expansion-log.rkt
@@ -50,7 +50,8 @@ A value of k that reduces any x to [0.5, 1.0] is
 (define (fl2log x2 x1)
   (define x (fl+ x1 x2))
   (cond [(x . fl<= . 0.0)  (cond [(fl= x 0.0)  (values -inf.0 0.0)]
-                                 [else  (values +nan.0 +nan.0)])]
+                                 [else  (values +nan.0 0.0)])]
+        [(x . fl= . +inf.0)  (values +inf.0 0.0)]
         [(or (x . fl< . 0.5) (x . fl> . 2.5))
          ;; Reduce arguments
          (let*-values ([(k x2 x1)  (fl2log-reduction x2 x1)]
@@ -80,10 +81,8 @@ A `k'  that reduces any argument `x' to (-1/2,1/2) is
   k = round(log1p(x)/log(2))
 |#
 
-(: fl2log1p-reduction (Flonum Flonum -> (Values Flonum Flonum Flonum)))
-(define (fl2log1p-reduction x2 x1)
-  (define-values (a2 a1) (fl2+ x2 x1 1.0))
-  (define y (fllog+ a2 a1))
+(: fl2log1p-reduction (Flonum Flonum Flonum Flonum Flonum -> (Values Flonum Flonum Flonum)))
+(define (fl2log1p-reduction x2 x1 a2 a1 y)
   (define k (flround (fl/ y (fllog 2.0))))
   (define 2^k (flexpt 2.0 k))
   (define-values (j2 j1) (fast-fl-/error (/ 1.0 2^k) 1.0))
@@ -96,9 +95,10 @@ A `k'  that reduces any argument `x' to (-1/2,1/2) is
   (define-values (a2 a1) (fl2+ x2 x1 1.0))
   (define y (fllog+ a2 a1))
   (cond
-    [(y . fl< . -0.5)  (fl2log a2 a1)]
+    [(or (y . fl< . -0.5) (a2 . fl> . 2.0))  (fl2log a2 a1)]
+    [(fl= (fl+ x2 x1) 0.0)  (values x2 0.0)]
     [(y . fl> . 0.5)
-     (let*-values ([(k x2 x1)  (fl2log1p-reduction x2 x1)]
+     (let*-values ([(k x2 x1)  (fl2log1p-reduction x2 x1 a2 a1 y)]
                    [(y2 y1)  (fl2log1p x2 x1)]
                    [(z2 z1)  (fl2* log2-hi log2-lo k)])
        (fl2+ y2 y1 z2 z1))]
diff --git a/collects/math/private/flonum/flonum-error.rkt b/collects/math/private/flonum/flonum-error.rkt
index 45b4cc810a..7933074e11 100644
--- a/collects/math/private/flonum/flonum-error.rkt
+++ b/collects/math/private/flonum/flonum-error.rkt
@@ -146,12 +146,15 @@
         (values x2 (if (and (flrational? x2) (not (flsubnormal? x2)))
                        (let*-values ([(da db)  (values (near-pow2 a) (near-pow2 b))]
                                      [(d)  (fl* da db)]
+                                     [(d?)  (and (d . fl> . 0.0) (d . fl< . +inf.0))]
                                      [(a2 a1)  (flsplit (fl/ a da))]
-                                     [(b2 b1)  (flsplit (fl/ b db))])
-                         (fl* d (- (fl- (fl- (fl- (fl- (fl/ x2 d) (fl* a2 b2))
-                                                  (fl* a1 b2))
-                                             (fl* a2 b1))
-                                        (fl* a1 b1)))))
+                                     [(b2 b1)  (flsplit (fl/ b db))]
+                                     [(x2)  (if d? (fl/ x2 d) (fl/ (fl/ x2 da) db))]
+                                     [(x1)  (- (fl- (fl- (fl- (fl- x2 (fl* a2 b2))
+                                                              (fl* a1 b2))
+                                                         (fl* a2 b1))
+                                                    (fl* a1 b1)))])
+                         (if d? (fl* x1 d) (fl* (fl* x1 da) db)))
                        0.0))))
     
     (: flsqr/error (Flonum -> (Values Flonum Flonum)))
@@ -160,10 +163,13 @@
         (values x2 (if (and (flrational? x2) (not (flsubnormal? x2)))
                        (let*-values ([(d)  (near-pow2 a)]
                                      [(d^2)  (fl* d d)]
-                                     [(a2 a1)  (flsplit (fl/ a d))])
-                         (fl* d^2 (- (fl- (fl- (fl- (fl/ x2 d^2) (fl* a2 a2))
-                                               (fl* 2.0 (fl* a1 a2)))
-                                          (fl* a1 a1)))))
+                                     [(d^2?)  (and (d^2 . fl> . 0.0) (d^2 . fl< . +inf.0))]
+                                     [(a2 a1)  (flsplit (fl/ a d))]
+                                     [(x2)  (if d^2? (fl/ x2 d^2) (fl/ (fl/ x2 d) d))]
+                                     [(x1)  (- (fl- (fl- (fl- x2 (fl* a2 a2))
+                                                         (fl* 2.0 (fl* a1 a2)))
+                                                    (fl* a1 a1)))])
+                         (if d^2? (fl* x1 d^2) (fl* (fl* x1 d) d)))
                        0.0))))
     
     (: fl//error (Flonum Flonum -> (Values Flonum Flonum)))
diff --git a/collects/math/private/flonum/flonum-functions.rkt b/collects/math/private/flonum/flonum-functions.rkt
index 0fafa850e4..73357d663a 100644
--- a/collects/math/private/flonum/flonum-functions.rkt
+++ b/collects/math/private/flonum/flonum-functions.rkt
@@ -9,7 +9,7 @@
 
 (provide (all-from-out racket/flonum)
          fl
-         flsubnormal? flrational? flnan? flinteger?
+         flsubnormal? flrational? flinfinite? flnan? flinteger?
          flnext* flprev*
          flulp-error
          fleven? flodd? flsgn flhypot fllog/base
@@ -40,6 +40,10 @@
     (λ: ([x : Flonum])
       (and (x . fl> . -inf.0) (x . fl< . +inf.0))))
   
+  (define flinfinite?
+    (λ: ([x : Flonum])
+      (or (x . fl= . -inf.0) (x . fl= . +inf.0))))
+  
   (define flnan?
     (λ: ([x : Flonum])
       (not (and (x . fl>= . -inf.0) (x . fl<= . +inf.0)))))
diff --git a/collects/math/private/flonum/utils.rkt b/collects/math/private/flonum/utils.rkt
index f7e374af1d..480df81fe7 100644
--- a/collects/math/private/flonum/utils.rkt
+++ b/collects/math/private/flonum/utils.rkt
@@ -9,7 +9,7 @@
   (provide (all-defined-out))
   
   (define-syntax-rule (near-pow2 a)
-    (flexpt 2.0 (flmax -511.0 (flmin 511.0 (flround (fl/ (fllog (flabs a)) (fllog 2.0)))))))
+    (flexpt 2.0 (flmax -1023.0 (flmin 1023.0 (flround (fl/ (fllog (flabs a)) (fllog 2.0)))))))
   
   (define-syntax-rule (near-pow2/div a b)
     ;; Clamping both values makes this work properly when a or b is infinite or zero
diff --git a/collects/math/private/utils/flonum-tests.rkt b/collects/math/private/utils/flonum-tests.rkt
index 3f41f3e355..25d056d5db 100644
--- a/collects/math/private/utils/flonum-tests.rkt
+++ b/collects/math/private/utils/flonum-tests.rkt
@@ -6,6 +6,8 @@
          typed/rackunit
          "../base/base-random.rkt"
          "../flonum/expansion/expansion-base.rkt"
+         "../flonum/expansion/expansion-exp.rkt"
+         "../flonum/expansion/expansion-log.rkt"
          "../flonum/flonum-functions.rkt"
          "../flonum/flonum-constants.rkt"
          "../flonum/flonum-bits.rkt"
@@ -26,8 +28,15 @@
 ;; Allowable error for different kinds of functions, in ulps
 (define flonum-fun-ulps 0.5)
 (define flonum/error-fun-ulps 0.5)
+(define flexp/error-fun-ulps 2.0)
+(define fl2-conversion-ulps 0.5)
 (define unary-fl2-fun-ulps 1.0)
 (define binary-fl2-fun-ulps 8.0)
+(define fl2exp-fun-ulps 2.5)
+(define fl2log-fun-ulps 1.5)
+
+(: current-max-ulp-error (Parameterof Nonnegative-Flonum))
+(define current-max-ulp-error (make-parameter 0.0))
 
 ;; ===================================================================================================
 ;; Helpers
@@ -94,13 +103,13 @@
    ;; Test squaring
    (- (flsqrt +min.0)) (flsqrt +min.0)
    (- (flsqrt +max-subnormal.0)) (flsqrt +max-subnormal.0)
-   (- (flsqrt +max-fl2-subnormal.0)) (flsqrt +max-fl2-subnormal.0)
+   (- (flsqrt +max-subnormal.hi)) (flsqrt +max-subnormal.hi)
    (- (flsqrt +max.0)) (flsqrt +max.0)
    ;; Test exp limits
-   (fllog +min.0) (fllog +max-subnormal.0) (fllog +max-fl2-subnormal.0) (fllog +max.0)
+   (fllog +min.0) (fllog +max-subnormal.0) (fllog +max-subnormal.hi) (fllog +max.0)
    ;; Standard special values
-   -inf.0 -max.0 -1.0 -max-fl2-subnormal.0 -max-subnormal.0 -min.0 -0.0
-   +inf.0 +max.0 +1.0 +max-fl2-subnormal.0 +max-subnormal.0 +min.0 +0.0
+   -inf.0 -max.0 -1.0 -max-subnormal.hi -max-subnormal.0 -min.0 -0.0
+   +inf.0 +max.0 +1.0 +max-subnormal.hi +max-subnormal.0 +min.0 +0.0
    +nan.0))
 
 (define standard-rs
@@ -123,8 +132,8 @@
 (define min-subnormal-ord (flonum->ordinal -max-subnormal.0))
 (define max-subnormal-ord (+ 1 (flonum->ordinal +max-subnormal.0)))
 
-(define min-fl2-subnormal-ord (flonum->ordinal -max-fl2-subnormal.0))
-(define max-fl2-subnormal-ord (+ 1 (flonum->ordinal +max-fl2-subnormal.0)))
+(define min-fl2-subnormal-ord (flonum->ordinal -max-subnormal.hi))
+(define max-fl2-subnormal-ord (+ 1 (flonum->ordinal +max-subnormal.hi)))
 
 (: sample-flonum (case-> (Integer -> (Listof Flonum))
                          (Integer Flonum Flonum -> (Listof Flonum))))
@@ -183,7 +192,7 @@
                                                           [i  (in-naturals 1)])
      (maybe-print-progress name i m)
      (list (list name x) (unary-flonum-fun-error f g x)))
-   flonum-fun-ulps))
+   (current-max-ulp-error)))
 
 (: binary-flonum-fun-error
    ((Flonum Flonum -> Flonum) (Bigfloat Bigfloat -> Bigfloat) Flonum Flonum -> Any))
@@ -205,7 +214,7 @@
                                                                  [i  (in-naturals 1)])
      (maybe-print-progress name i m)
      (list (list name x y) (binary-flonum-fun-error f g x y)))
-   flonum-fun-ulps))
+   (current-max-ulp-error)))
 
 ;; ===================================================================================================
 ;; fl2 conversion
@@ -230,7 +239,7 @@
                                                       [i  (in-naturals 1)])
      (maybe-print-progress 'fl2 i m)
      (list (list 'fl2 x) (fl2-conversion-error x)))
-   flonum/error-fun-ulps))
+   (current-max-ulp-error)))
 
 ;; ===================================================================================================
 ;; Flonum arithmetic with error
@@ -262,7 +271,7 @@
                                                           [i  (in-naturals 1)])
      (maybe-print-progress name i m)
      (list (list name x) (unary-flonum/error-fun-error f g x)))
-   flonum/error-fun-ulps))
+   (current-max-ulp-error)))
 
 (: test-binary-flonum/error-fun
    (Symbol (Flonum Flonum -> (Values Flonum Flonum)) (Bigfloat Bigfloat -> Bigfloat) Integer
@@ -278,7 +287,7 @@
                                                                  [i  (in-naturals 1)])
      (maybe-print-progress name i m)
      (list (list name x y) (binary-flonum/error-fun-error f g x y)))
-   flonum/error-fun-ulps))
+   (current-max-ulp-error)))
 
 ;; ===================================================================================================
 ;; Flonum expansions
@@ -302,7 +311,7 @@
      (maybe-print-progress name i m)
      (define-values (x2 x1) (fl2 x))
      (list (list name x2 x1) (unary-fl2-fun-error f g x2 x1)))
-   unary-fl2-fun-ulps))
+   (current-max-ulp-error)))
 
 (: binary-fl2-fun-error ((Flonum Flonum Flonum Flonum -> (Values Flonum Flonum))
                          (Bigfloat Bigfloat -> Bigfloat)
@@ -330,83 +339,110 @@
      (define-values (x2 x1) (fl2 x))
      (define-values (y2 y1) (fl2 y))
      (list (list name x2 x1 y2 y1) (binary-fl2-fun-error f g x2 x1 y2 y1)))
-   binary-fl2-fun-ulps))
+   (current-max-ulp-error)))
 
 ;; ===================================================================================================
 
 (: test-fpu-arith (Natural -> Any))
 (define (test-fpu-arith n)
-  (check-equal? (test-unary-flonum-fun 'flabs flabs bfabs n -inf.0 +inf.0)
-                '())
-  (check-equal? (test-binary-flonum-fun 'fl+ fl+ bf+ n)
-                '())
-  (check-equal? (test-binary-flonum-fun 'fl- fl- bf- n)
-                '())
-  (check-equal? (test-binary-flonum-fun 'fl* fl* bf* n)
-                '())
-  (check-equal? (test-binary-flonum-fun 'fl/ fl/ bf/ n)
-                '()))
+  (parameterize ([current-max-ulp-error  flonum-fun-ulps])
+    (check-equal? (test-unary-flonum-fun 'flabs flabs bfabs n -inf.0 +inf.0)
+                  '())
+    (check-equal? (test-binary-flonum-fun 'fl+ fl+ bf+ n)
+                  '())
+    (check-equal? (test-binary-flonum-fun 'fl- fl- bf- n)
+                  '())
+    (check-equal? (test-binary-flonum-fun 'fl* fl* bf* n)
+                  '())
+    (check-equal? (test-binary-flonum-fun 'fl/ fl/ bf/ n)
+                  '())))
 
 (: test-fpu-trig (Natural -> Any))
 (define (test-fpu-trig n)
-  (check-equal? (test-unary-flonum-fun 'flsin flsin bfsin n -inf.0 +inf.0)
-                '())
-  (check-equal? (test-unary-flonum-fun 'flcos flcos bfcos n -inf.0 +inf.0)
-                '())
-  (check-equal? (test-unary-flonum-fun 'fltan fltan bftan n -inf.0 +inf.0)
-                '())
-  (check-equal? (test-unary-flonum-fun 'flasin flasin bfasin n -1.0 1.0)
-                '())
-  (check-equal? (test-unary-flonum-fun 'flacos flacos bfacos n -1.0 1.0)
-                '())
-  (check-equal? (test-unary-flonum-fun 'flatan flatan bfatan n -inf.0 +inf.0)
-                '()))
+  (parameterize ([current-max-ulp-error  flonum-fun-ulps])
+    (check-equal? (test-unary-flonum-fun 'flsin flsin bfsin n -inf.0 +inf.0)
+                  '())
+    (check-equal? (test-unary-flonum-fun 'flcos flcos bfcos n -inf.0 +inf.0)
+                  '())
+    (check-equal? (test-unary-flonum-fun 'fltan fltan bftan n -inf.0 +inf.0)
+                  '())
+    (check-equal? (test-unary-flonum-fun 'flasin flasin bfasin n -1.0 1.0)
+                  '())
+    (check-equal? (test-unary-flonum-fun 'flacos flacos bfacos n -1.0 1.0)
+                  '())
+    (check-equal? (test-unary-flonum-fun 'flatan flatan bfatan n -inf.0 +inf.0)
+                  '())))
 
 (: test-fpu-non-trig (Natural -> Any))
 (define (test-fpu-non-trig n)
-  (check-equal? (test-unary-flonum-fun 'flsqrt flsqrt bfsqrt n 0.0 +inf.0)
-                '())
-  (check-equal? (test-unary-flonum-fun 'fllog fllog bflog n 0.0 +inf.0)
-                '())
-  (check-equal? (test-unary-flonum-fun 'flexp flexp bfexp n -746.0 710.0)
-                '())
-  (check-equal? (test-binary-flonum-fun 'flexpt flexpt bfexpt n)
-                '()))
+  (parameterize ([current-max-ulp-error  flonum-fun-ulps])
+    (check-equal? (test-unary-flonum-fun 'flsqrt flsqrt bfsqrt n 0.0 +inf.0)
+                  '())
+    (check-equal? (test-unary-flonum-fun 'fllog fllog bflog n 0.0 +inf.0)
+                  '())
+    (check-equal? (test-unary-flonum-fun 'flexp flexp bfexp n -746.0 710.0)
+                  '())
+    (check-equal? (test-binary-flonum-fun 'flexpt flexpt bfexpt n)
+                  '())))
 
 (: test-fpu-arith/error (Natural -> Any))
 (define (test-fpu-arith/error n)
-  (check-equal? (test-binary-flonum/error-fun 'fl+/error fl+/error bf+ n)
-                '())
-  (check-equal? (test-binary-flonum/error-fun 'fl-/error fl-/error bf- n)
-                '())
-  (check-equal? (test-binary-flonum/error-fun 'fl*/error fl*/error bf* n)
-                '())
-  (check-equal? (test-unary-flonum/error-fun 'flsqr/error flsqr/error bfsqr n)
-                '())
-  (check-equal? (test-binary-flonum/error-fun 'fl//error fl//error bf/ n)
-                '()))
+  (parameterize ([current-max-ulp-error  flonum/error-fun-ulps])
+    (check-equal? (test-binary-flonum/error-fun 'fl+/error fl+/error bf+ n)
+                  '())
+    (check-equal? (test-binary-flonum/error-fun 'fl-/error fl-/error bf- n)
+                  '())
+    (check-equal? (test-binary-flonum/error-fun 'fl*/error fl*/error bf* n)
+                  '())
+    (check-equal? (test-unary-flonum/error-fun 'flsqr/error flsqr/error bfsqr n)
+                  '())
+    (check-equal? (test-binary-flonum/error-fun 'fl//error fl//error bf/ n)
+                  '())))
+
+(: test-fpu-non-trig/error (Natural -> Any))
+(define (test-fpu-non-trig/error n)
+  (parameterize ([current-max-ulp-error  flexp/error-fun-ulps])
+    (check-equal? (test-unary-flonum/error-fun 'flexp/error flexp/error bfexp n)
+                  '())))
 
 (: test-fpu-arith/fl2 (Natural -> Any))
 (define (test-fpu-arith/fl2 n)
-  (check-equal? (test-fl2-conversion n)
-                '())
-  (check-equal? (test-unary-fl2-fun 'fl2abs fl2abs bfabs n)
-                '())
-  (check-equal? (test-binary-fl2-fun 'fl2+ fl2+ bf+ n)
-                '())
-  (check-equal? (test-binary-fl2-fun 'fl2- fl2- bf- n)
-                '())
-  (check-equal? (test-binary-fl2-fun 'fl2* fl2* bf* n)
-                '())
-  (check-equal? (test-unary-fl2-fun 'fl2sqr fl2sqr bfsqr n)
-                '())
-  (check-equal? (test-binary-fl2-fun 'fl2/ fl2/ bf/ n)
-                '()))
+  (parameterize ([current-max-ulp-error  fl2-conversion-ulps])
+    (check-equal? (test-fl2-conversion n)
+                  '()))
+  (parameterize ([current-max-ulp-error  unary-fl2-fun-ulps])
+    (check-equal? (test-unary-fl2-fun 'fl2abs fl2abs bfabs n)
+                  '())
+    (check-equal? (test-unary-fl2-fun 'fl2sqr fl2sqr bfsqr n)
+                  '()))
+  (parameterize ([current-max-ulp-error  binary-fl2-fun-ulps])
+    (check-equal? (test-binary-fl2-fun 'fl2+ fl2+ bf+ n)
+                  '())
+    (check-equal? (test-binary-fl2-fun 'fl2- fl2- bf- n)
+                  '())
+    (check-equal? (test-binary-fl2-fun 'fl2* fl2* bf* n)
+                  '())
+    (check-equal? (test-binary-fl2-fun 'fl2/ fl2/ bf/ n)
+                  '())))
 
 (: test-fpu-non-trig/fl2 (Natural -> Any))
 (define (test-fpu-non-trig/fl2 n)
-  (check-equal? (test-unary-fl2-fun 'fl2sqrt fl2sqrt bfsqrt n)
-                '()))
+  (parameterize ([current-max-ulp-error  unary-fl2-fun-ulps])
+    (check-equal? (test-unary-fl2-fun 'fl2sqrt fl2sqrt bfsqrt n)
+                  '()))
+  (parameterize ([current-max-ulp-error  fl2exp-fun-ulps])
+    (check-equal? (test-unary-fl2-fun 'fl2exp fl2exp bfexp n)
+                  '()))
+  (parameterize ([current-max-ulp-error  fl2exp-fun-ulps])
+    (check-equal? (test-unary-fl2-fun 'fl2expm1 fl2expm1 bfexpm1 n)
+                  '()))
+  (parameterize ([current-max-ulp-error  fl2log-fun-ulps])
+    (check-equal? (test-unary-fl2-fun 'fl2log fl2log bflog n)
+                  '()))
+  (parameterize ([current-max-ulp-error  fl2log-fun-ulps])
+    (check-equal? (test-unary-fl2-fun 'fl2log1p fl2log1p bflog1p n)
+                  '()))
+  )
 
 (: test-fpu (Natural -> Any))
 (define (test-fpu n)
@@ -414,5 +450,43 @@
   (test-fpu-trig n)
   (test-fpu-non-trig n)
   (test-fpu-arith/error n)
+  (test-fpu-non-trig/error n)
   (test-fpu-arith/fl2 n)
   (test-fpu-non-trig/fl2 n))
+
+(for*: ([x2  (list -inf.0 -max.0 -1.0 -min.0 -0.0 0.0 +min.0 +1.0 +max.0 +inf.0 +nan.0)]
+        [x1  (list -inf.0 -max.0 -1.0 -min.0 -0.0 0.0 +min.0 +1.0 +max.0 +inf.0 +nan.0)])
+  (define n
+    (count (λ: ([b : Boolean]) b)
+           (map (λ: ([f : (Flonum Flonum -> Boolean)])
+                  (f x2 x1))
+                (list fl2rational? fl2infinite? fl2nan?))))
+  (unless (= n 1) (printf "x2 = ~v  x1 = ~v~n" x2 x1)))
+
+#|
+Tests to add
+
+(for*: ([x2  (list -inf.0 -max.0 -1.0 -min.0 -0.0 0.0 +min.0 +1.0 +max.0 +inf.0 +nan.0)]
+        [x1  (list -inf.0 -max.0 -1.0 -min.0 -0.0 0.0 +min.0 +1.0 +max.0 +inf.0 +nan.0)])
+  (define n
+    (count (λ: ([b : Boolean]) b)
+           (map (λ: ([f : (Flonum Flonum -> Boolean)])
+                  (f x2 x1))
+                (list fl2rational? fl2infinite? fl2nan?))))
+  (unless (= n 1) (printf "x2 = ~v  x1 = ~v~n" x2 x1)))
+
+fl2=
+fl2>
+fl2<
+fl2>=
+fl2<=
+
+(fl2step x2 x1 n/2) twice = (fl2step x2 x1 n)
+
+|#
+
+(check-true (let-values ([(y2 y1)  (fl+/error +max.hi +max.lo)])
+              (fl2= y2 y1 +max.hi +max.lo)))
+
+(check-true (let*-values ([(y2 y1)  (fl2next +max.hi +max.lo)])
+              (fl2infinite? y2 y1)))