Rollback ival-exact + fix tunnel vision changes

Author: obround

ops/arith.rkt

@@ -161,8 +161,8 @@
 (define (ival-hypot! out x y)
   (define err? (or (ival-err? x) (ival-err? y)))
   (define err (or (ival-err x) (ival-err y)))
-  (define x* (ival-fabs x))
-  (define y* (ival-fabs y))
+  (define x* (ival-exact-fabs x))
+  (define y* (ival-exact-fabs y))
   (ival (eplinear! (ival-lo-val out) mpfr-hypot! (ival-lo x*) (ival-lo y*) 'down)
         (eplinear! (ival-hi-val out) mpfr-hypot! (ival-hi x*) (ival-hi y*) 'up)
         err?

ops/core.rkt

@@ -17,6 +17,8 @@
          epfn
          split-ival
          ival-max-prec
+         ival-exact-neg
+         ival-exact-fabs
          bf-return-exact?
          ival-lo-fixed?
          ival-hi-fixed?
@@ -191,18 +193,14 @@
 (define (endpoint-min2 e1 e2 rnd)
   (match-define (endpoint x x!) e1)
   (match-define (endpoint y y!) e2)
-  (define out
-    (parameterize ([bf-precision (bf-precision)])
-      (bf 0)))
+  (define out (bf 0))
   (mpfr-min! out x y rnd)
   (endpoint out (or (and (bf=? out x) x!) (and (bf=? out y) y!))))
 
 (define (endpoint-max2 e1 e2 rnd)
   (match-define (endpoint x x!) e1)
   (match-define (endpoint y y!) e2)
-  (define out
-    (parameterize ([bf-precision (bf-precision)])
-      (bf 0)))
+  (define out (bf 0))
   (mpfr-max! out x y rnd)
   (endpoint out (or (and (bf=? out x) x!) (and (bf=? out y) y!))))
 
@@ -349,6 +347,14 @@
 (define (ival-max-prec x)
   (max (bigfloat-precision (ival-lo-val x)) (bigfloat-precision (ival-hi-val x))))
 
+(define (ival-exact-fabs x)
+  (parameterize ([bf-precision (ival-max-prec x)])
+    (ival-fabs x)))
+
+(define (ival-exact-neg x)
+  (parameterize ([bf-precision (ival-max-prec x)])
+    (ival-neg x)))
+
 ;; Since MPFR has a cap on exponents, no value can be more than twice MAX_VAL
 (define exp-overflow-threshold (bfadd (bflog (bfprev +inf.bf)) 1.bf))
 (define exp2-overflow-threshold (bfadd (bflog2 (bfprev +inf.bf)) 1.bf))
@@ -369,7 +375,7 @@
 (define* ival-log2 (compose (monotonic-mpfr mpfr-log2!) (clamp-strict 0.bf +inf.bf)))
 (define* ival-log10 (compose (monotonic-mpfr mpfr-log10!) (clamp-strict 0.bf +inf.bf)))
 (define* ival-log1p (compose (monotonic-mpfr mpfr-log1p!) (clamp-strict -1.bf +inf.bf)))
-[define* ival-logb (compose ival-floor ival-log2 ival-fabs)]
+[define* ival-logb (compose ival-floor ival-log2 ival-exact-fabs)]
 
 (define* ival-sqrt (compose (monotonic-mpfr mpfr-sqrt!) (clamp 0.bf +inf.bf)))
 (define* ival-cbrt (monotonic-mpfr mpfr-cbrt!))
@@ -427,7 +433,7 @@
          (compose (overflows-at (monotonic-mpfr mpfr-cosh!)
                                 (bfneg acosh-overflow-threshold)
                                 acosh-overflow-threshold)
-                  ival-fabs))
+                  ival-exact-fabs))
 (define*
  ival-sinh
  (overflows-at (monotonic-mpfr mpfr-sinh!) (bfneg sinh-overflow-threshold) sinh-overflow-threshold))

ops/fmod.rkt

@@ -63,14 +63,14 @@
         (and (bflte? (ival-lo-val y) 0.bf) (bfgte? (ival-hi-val y) 0.bf))))
   (define err
     (or (ival-err x) (ival-err y) (and (bf=? (ival-lo-val y) 0.bf) (bf=? (ival-hi-val y) 0.bf))))
-  (define y* (ival-fabs y))
+  (define y* (ival-exact-fabs y))
   (cond
-    [(= (mpfr-sign (ival-hi-val x)) -1) (ival-neg (ival-fmod-pos (ival-neg x) y* err? err))]
+    [(= (mpfr-sign (ival-hi-val x)) -1) (ival-neg (ival-fmod-pos (ival-exact-neg x) y* err? err))]
     [(= (mpfr-sign (ival-lo-val x)) 1) (ival-fmod-pos x y* err? err)]
     [else
      (define-values (neg pos) (split-ival x 0.bf))
      (ival-union (ival-fmod-pos pos y* err? err)
-                 (ival-neg (ival-fmod-pos (ival-neg neg) y* err? err)))]))
+                 (ival-neg (ival-fmod-pos (ival-exact-neg neg) y* err? err)))]))
 
 (define (ival-remainder-pos x y err? err)
   ;; Assumes both `x` and `y` are entirely positive
@@ -106,11 +106,12 @@
         (and (bflte? (ival-lo-val y) 0.bf) (bfgte? (ival-hi-val y) 0.bf))))
   (define err
     (or (ival-err x) (ival-err y) (and (bf=? (ival-lo-val y) 0.bf) (bf=? (ival-hi-val y) 0.bf))))
-  (define y* (ival-fabs y))
+  (define y* (ival-exact-fabs y))
   (cond
-    [(= (mpfr-sign (ival-hi-val x)) -1) (ival-neg (ival-remainder-pos (ival-neg x) y* err? err))]
+    [(= (mpfr-sign (ival-hi-val x)) -1)
+     (ival-neg (ival-remainder-pos (ival-exact-neg x) y* err? err))]
     [(= (mpfr-sign (ival-lo-val x)) 1) (ival-remainder-pos x y* err? err)]
     [else
      (define-values (neg pos) (split-ival x 0.bf))
      (ival-union (ival-remainder-pos pos y* err? err)
-                 (ival-neg (ival-remainder-pos (ival-neg neg) y* err? err)))]))
+                 (ival-neg (ival-remainder-pos (ival-exact-neg neg) y* err? err)))]))

ops/pow.rkt

@@ -129,17 +129,17 @@
           ; If y is an integer point interval, there's no error,
           ; because it's always valid to raise to an integer power.
           (if (bfodd? (ival-lo-val y))
-              (ival-neg (ival-pow-pos (ival-fabs x) y)) ; Use fabs in case of [x, 0]
-              (ival-pow-pos (ival-fabs x) y))
+              (ival-neg (ival-pow-pos (ival-exact-fabs x) y)) ; Use fabs in case of [x, 0]
+              (ival-pow-pos (ival-exact-fabs x) y))
           ; If y is non-integer point interval, it must be an even
           ; fraction (because all bigfloats are) so we always error
           ival-illegal)
       ; Moreover, if we have (-x)^y, that's basically x^y U -(x^y).
-      (let ([pospow (ival-pow-pos (ival-fabs x) y)])
+      (let ([pospow (ival-pow-pos (ival-exact-fabs x) y)])
         (ival-then (ival-assert ival-maybe) (ival-union (ival-neg pospow) pospow)))))
 
 (define (ival-pow2 x)
-  ((monotonic->ival (lambda (x) (bfmul x x))) (ival-fabs x)))
+  ((monotonic->ival (lambda (x) (bfmul x x))) (ival-exact-fabs x)))
 
 (define (ival-pow x y)
   (cond