Add fadd, fsub, fmul, fdiv instructions

* libguile/vm-engine.c (fadd, fsub, fmul, fdiv): New instructions.

* module/language/cps/effects-analysis.scm:
* module/language/cps/types.scm: Wire up support for new instructions.

* module/system/vm/assembler.scm: Export emit-fadd and friends.

4 files changed, 133 insertions, 43 deletions

diff --git a/libguile/vm-engine.c b/libguile/vm-engine.c
index d7320059a..d33878d20 100644
--- a/libguile/vm-engine.c
+++ b/libguile/vm-engine.c
@@ -3258,10 +3258,58 @@ VM_NAME (scm_i_thread *thread, struct scm_vm *vp,
       NEXT (1);
     }
 
-  VM_DEFINE_OP (138, unused_138, NULL, NOP)
-  VM_DEFINE_OP (139, unused_139, NULL, NOP)
-  VM_DEFINE_OP (140, unused_140, NULL, NOP)
-  VM_DEFINE_OP (141, unused_141, NULL, NOP)
+  /* fadd dst:8 a:8 b:8
+   *
+   * Add A to B, and place the result in DST.  The operands and the
+   * result are unboxed double-precision floating-point numbers.
+   */
+  VM_DEFINE_OP (138, fadd, "fadd", OP1 (X8_S8_S8_S8) | OP_DST)
+    {
+      scm_t_uint8 dst, a, b;
+      UNPACK_8_8_8 (op, dst, a, b);
+      SP_SET_F64 (dst, SP_REF_F64 (a) + SP_REF_F64 (b));
+      NEXT (1);
+    }
+
+  /* fsub dst:8 a:8 b:8
+   *
+   * Subtract B from A, and place the result in DST.  The operands and
+   * the result are unboxed double-precision floating-point numbers.
+   */
+  VM_DEFINE_OP (139, fsub, "fsub", OP1 (X8_S8_S8_S8) | OP_DST)
+    {
+      scm_t_uint8 dst, a, b;
+      UNPACK_8_8_8 (op, dst, a, b);
+      SP_SET_F64 (dst, SP_REF_F64 (a) - SP_REF_F64 (b));
+      NEXT (1);
+    }
+
+  /* fmul dst:8 a:8 b:8
+   *
+   * Multiply A and B, and place the result in DST.  The operands and
+   * the result are unboxed double-precision floating-point numbers.
+   */
+  VM_DEFINE_OP (140, fmul, "fmul", OP1 (X8_S8_S8_S8) | OP_DST)
+    {
+      scm_t_uint8 dst, a, b;
+      UNPACK_8_8_8 (op, dst, a, b);
+      SP_SET_F64 (dst, SP_REF_F64 (a) * SP_REF_F64 (b));
+      NEXT (1);
+    }
+
+  /* fdiv dst:8 a:8 b:8
+   *
+   * Divide A by B, and place the result in DST.  The operands and the
+   * result are unboxed double-precision floating-point numbers.
+   */
+  VM_DEFINE_OP (141, fdiv, "fdiv", OP1 (X8_S8_S8_S8) | OP_DST)
+    {
+      scm_t_uint8 dst, a, b;
+      UNPACK_8_8_8 (op, dst, a, b);
+      SP_SET_F64 (dst, SP_REF_F64 (a) / SP_REF_F64 (b));
+      NEXT (1);
+    }
+
   VM_DEFINE_OP (142, unused_142, NULL, NOP)
   VM_DEFINE_OP (143, unused_143, NULL, NOP)
   VM_DEFINE_OP (144, unused_144, NULL, NOP)
diff --git a/module/language/cps/effects-analysis.scm b/module/language/cps/effects-analysis.scm
index 3542a1e74..ae7a1a614 100644
--- a/module/language/cps/effects-analysis.scm
+++ b/module/language/cps/effects-analysis.scm
@@ -414,6 +414,10 @@ is or might be a read or a write to the same location as A."
   ((mul . _)                       &type-check)
   ((sub . _)                       &type-check)
   ((div . _)                       &type-check)
+  ((fadd . _))
+  ((fsub . _))
+  ((fmul . _))
+  ((fdiv . _))
   ((sub1 . _)                      &type-check)
   ((add1 . _)                      &type-check)
   ((quo . _)                       &type-check)
diff --git a/module/language/cps/types.scm b/module/language/cps/types.scm
index 8a2cc86d3..dac29f7d2 100644
--- a/module/language/cps/types.scm
+++ b/module/language/cps/types.scm
@@ -842,18 +842,48 @@ minimum, and maximum."
                  min* max*))))))
 
 (define-simple-type-checker (add &number &number))
+(define-type-checker (fadd a b) #t)
 (define-type-inferrer (add a b result)
   (define-binary-result! a b result #t
                          (+ (&min a) (&min b))
                          (+ (&max a) (&max b))))
+(define-type-inferrer (fadd a b result)
+  (define! result &f64
+    (+ (&min a) (&min b))
+    (+ (&max a) (&max b))))
 
 (define-simple-type-checker (sub &number &number))
+(define-type-checker (fsub a b) #t)
 (define-type-inferrer (sub a b result)
   (define-binary-result! a b result #t
                          (- (&min a) (&max b))
                          (- (&max a) (&min b))))
+(define-type-inferrer (fsub a b result)
+  (define! result &f64
+    (- (&min a) (&max b))
+    (- (&max a) (&min b))))
 
 (define-simple-type-checker (mul &number &number))
+(define-type-checker (fmul a b) #t)
+(define (mul-result-range same? nan-impossible? min-a max-a min-b max-b)
+  (define (nan* a b)
+    (if (and (or (and (inf? a) (zero? b))
+                 (and (zero? a) (inf? b)))
+             nan-impossible?)
+        0 
+        (* a b)))
+  (let ((-- (nan* min-a min-b))
+        (-+ (nan* min-a max-b))
+        (++ (nan* max-a max-b))
+        (+- (nan* max-a min-b)))
+    (let ((has-nan? (or (nan? --) (nan? -+) (nan? ++) (nan? +-))))
+      (values (cond
+               (same? 0)
+               (has-nan? -inf.0)
+               (else (min -- -+ ++ +-)))
+              (if has-nan?
+                  +inf.0
+                  (max -- -+ ++ +-))))))
 (define-type-inferrer (mul a b result)
   (let ((min-a (&min a)) (max-a (&max a))
         (min-b (&min b)) (max-b (&max b))
@@ -863,25 +893,20 @@ minimum, and maximum."
         ;; range inference time is 0 and not +nan.0.
         (nan-impossible? (not (logtest (logior (&type a) (&type b))
                                        (logior &flonum &complex)))))
-    (define (nan* a b)
-      (if (and (or (and (inf? a) (zero? b))
-                   (and (zero? a) (inf? b)))
-               nan-impossible?)
-          0 
-          (* a b)))
-    (let ((-- (nan* min-a min-b))
-          (-+ (nan* min-a max-b))
-          (++ (nan* max-a max-b))
-          (+- (nan* max-a min-b)))
-      (let ((has-nan? (or (nan? --) (nan? -+) (nan? ++) (nan? +-))))
-        (define-binary-result! a b result #t
-                               (cond
-                                ((eqv? a b) 0)
-                                (has-nan? -inf.0)
-                                (else (min -- -+ ++ +-)))
-                               (if has-nan?
-                                   +inf.0
-                                   (max -- -+ ++ +-)))))))
+    (call-with-values (lambda ()
+                        (mul-result-range (eqv? a b) nan-impossible?
+                                          min-a max-a min-b max-b))
+      (lambda (min max)
+        (define-binary-result! a b result #t min max)))))
+(define-type-inferrer (fmul a b result)
+  (let ((min-a (&min a)) (max-a (&max a))
+        (min-b (&min b)) (max-b (&max b))
+        (nan-impossible? #f))
+    (call-with-values (lambda ()
+                        (mul-result-range (eqv? a b) nan-impossible?
+                                          min-a max-a min-b max-b))
+      (lambda (min max)
+        (define! result &f64 min max)))))
 
 (define-type-checker (div a b)
   (and (check-type a &number -inf.0 +inf.0)
@@ -889,31 +914,40 @@ minimum, and maximum."
        ;; We only know that there will not be an exception if b is not
        ;; zero.
        (not (<= (&min b) 0 (&max b)))))
+(define-type-checker (fdiv a b) #t)
+(define (div-result-range min-a max-a min-b max-b)
+  (if (<= min-b 0 max-b)
+      ;; If the range of the divisor crosses 0, the result spans
+      ;; the whole range.
+      (values -inf.0 +inf.0)
+      ;; Otherwise min-b and max-b have the same sign, and cannot both
+      ;; be infinity.
+      (let ((--- (if (inf? min-b) 0 (floor/ min-a min-b)))
+            (-+- (if (inf? max-b) 0 (floor/ min-a max-b)))
+            (++- (if (inf? max-b) 0 (floor/ max-a max-b)))
+            (+-- (if (inf? min-b) 0 (floor/ max-a min-b)))
+            (--+ (if (inf? min-b) 0 (ceiling/ min-a min-b)))
+            (-++ (if (inf? max-b) 0 (ceiling/ min-a max-b)))
+            (+++ (if (inf? max-b) 0 (ceiling/ max-a max-b)))
+            (+-+ (if (inf? min-b) 0 (ceiling/ max-a min-b))))
+        (values (min (min --- -+- ++- +--)
+                     (min --+ -++ +++ +-+))
+                (max (max --- -+- ++- +--)
+                     (max --+ -++ +++ +-+))))))
 (define-type-inferrer (div a b result)
   (let ((min-a (&min a)) (max-a (&max a))
         (min-b (&min b)) (max-b (&max b)))
-    (call-with-values
-        (lambda ()
-          (if (<= min-b 0 max-b)
-              ;; If the range of the divisor crosses 0, the result spans
-              ;; the whole range.
-              (values -inf.0 +inf.0)
-              ;; Otherwise min-b and max-b have the same sign, and cannot both
-              ;; be infinity.
-              (let ((--- (if (inf? min-b) 0 (floor/ min-a min-b)))
-                    (-+- (if (inf? max-b) 0 (floor/ min-a max-b)))
-                    (++- (if (inf? max-b) 0 (floor/ max-a max-b)))
-                    (+-- (if (inf? min-b) 0 (floor/ max-a min-b)))
-                    (--+ (if (inf? min-b) 0 (ceiling/ min-a min-b)))
-                    (-++ (if (inf? max-b) 0 (ceiling/ min-a max-b)))
-                    (+++ (if (inf? max-b) 0 (ceiling/ max-a max-b)))
-                    (+-+ (if (inf? min-b) 0 (ceiling/ max-a min-b))))
-                (values (min (min --- -+- ++- +--)
-                             (min --+ -++ +++ +-+))
-                        (max (max --- -+- ++- +--)
-                             (max --+ -++ +++ +-+))))))
+    (call-with-values (lambda ()
+                        (div-result-range min-a max-a min-b max-b))
       (lambda (min max)
         (define-binary-result! a b result #f min max)))))
+(define-type-inferrer (fdiv a b result)
+  (let ((min-a (&min a)) (max-a (&max a))
+        (min-b (&min b)) (max-b (&max b)))
+    (call-with-values (lambda ()
+                        (div-result-range min-a max-a min-b max-b))
+      (lambda (min max)
+        (define! result &f64 min max)))))
 
 (define-simple-type-checker (add1 &number))
 (define-type-inferrer (add1 a result)
diff --git a/module/system/vm/assembler.scm b/module/system/vm/assembler.scm
index 9cb04bbed..ae54d1303 100644
--- a/module/system/vm/assembler.scm
+++ b/module/system/vm/assembler.scm
@@ -136,6 +136,10 @@
             (emit-rem* . emit-rem)
             (emit-mod* . emit-mod)
             (emit-ash* . emit-ash)
+            (emit-fadd* . emit-fadd)
+            (emit-fsub* . emit-fsub)
+            (emit-fmul* . emit-fmul)
+            (emit-fdiv* . emit-fdiv)
             (emit-logand* . emit-logand)
             (emit-logior* . emit-logior)
             (emit-logxor* . emit-logxor)