Permalink
Browse files

1.0.29.44: Complex float improvements

* On all platforms:
 - Slightly more stable complex-complex float (double and single)
   division;
 - New transform for real-complex division;
 - complex-real and real-complex float addition and subtraction
   behave as though the real was first upgraded to a complex, thus
   losing the sign of any imaginary zero.

* On x86-64
 - Complexes floats are represented packed in a single SSE register;
 - VOPs for all four arithmetic operations, complex-complex, but also
   complex-real and real-complex, except for complex-complex and
   real-complex division;
 - VOPs for =, negate and conjugate of complexes (complex-real and
   complex-complex);
 - VOPs for EQL of floats (real and complexes).
 - Full register moves for float values in SSE registers should also
   speed scalar operations up.
  • Loading branch information...
1 parent aecec1d commit a157ed0be79751f85b8243c06102eea95af06aa3 @pkhuong pkhuong committed Jun 25, 2009
@@ -155,6 +155,18 @@
;;
; :cycle-counter
+ ;; Enabled automatically for platforms which implement complex arithmetic
+ ;; VOPs. Such platforms should implement real-complex, complex-real and
+ ;; complex-complex addition and subtractions (for complex-single-float
+ ;; and complex-double-float). They should also also implement complex-real
+ ;; and real-complex multiplication, complex-real division, and
+ ;; sb!vm::swap-complex, which swaps the real and imaginary parts.
+ ;; Finally, they should implement conjugate and complex-real, real-complex
+ ;; and complex-complex CL:= (complex-complex EQL would usually be a good
+ ;; idea).
+ ;;
+ ; :complex-float-vops
+
;; Peter Van Eynde's increase-bulletproofness code for CMU CL
;;
;; Some of the code which was #+high-security before the fork has now
View
@@ -296,7 +296,7 @@ elif [ "$sbcl_arch" = "x86-64" ]; then
printf ' :compare-and-swap-vops :unwind-to-frame-and-call-vop :raw-instance-init-vops' >> $ltf
printf ' :stack-allocatable-closures :stack-allocatable-vectors' >> $ltf
printf ' :stack-allocatable-lists :stack-allocatable-fixed-objects' >> $ltf
- printf ' :alien-callbacks :cycle-counter' >> $ltf
+ printf ' :alien-callbacks :cycle-counter :complex-float-vops' >> $ltf
elif [ "$sbcl_arch" = "mips" ]; then
printf ' :linkage-table' >> $ltf
printf ' :stack-allocatable-closures :stack-allocatable-vectors' >> $ltf
@@ -2460,7 +2460,10 @@ structure representations"
#!+long-float "COMPLEX-LONG-FLOAT-WIDETAG"
#!+long-float "COMPLEX-LONG-REG-SC-NUMBER"
#!+long-float "COMPLEX-LONG-STACK-SC-NUMBER"
+ #!-x86-64 #!-x86-64
"COMPLEX-SINGLE-FLOAT-IMAG-SLOT" "COMPLEX-SINGLE-FLOAT-REAL-SLOT"
+ #!+x86-64
+ "COMPLEX-SINGLE-FLOAT-DATA-SLOT"
"COMPLEX-SINGLE-FLOAT-SIZE" "COMPLEX-SINGLE-FLOAT-WIDETAG"
"COMPLEX-SINGLE-REG-SC-NUMBER" "COMPLEX-SINGLE-STACK-SC-NUMBER"
"COMPLEX-SIZE" "COMPLEX-BASE-STRING-WIDETAG"
View
@@ -107,7 +107,7 @@
#!+sb-dyncount
(collect-dynamic-statistics nil))
(sb!c::defprinter (segment)
- name)
+ type)
(declaim (inline segment-current-index))
(defun segment-current-index (segment)
@@ -1178,27 +1178,45 @@
;;; of complex operation VOPs.
(macrolet ((frob (type)
`(progn
+ (deftransform complex ((r) (,type))
+ '(complex r ,(coerce 0 type)))
+ (deftransform complex ((r i) (,type (and real (not ,type))))
+ '(complex r (truly-the ,type (coerce i ',type))))
+ (deftransform complex ((r i) ((and real (not ,type)) ,type))
+ '(complex (truly-the ,type (coerce r ',type)) i))
;; negation
+ #!-complex-float-vops
(deftransform %negate ((z) ((complex ,type)) *)
'(complex (%negate (realpart z)) (%negate (imagpart z))))
;; complex addition and subtraction
+ #!-complex-float-vops
(deftransform + ((w z) ((complex ,type) (complex ,type)) *)
'(complex (+ (realpart w) (realpart z))
(+ (imagpart w) (imagpart z))))
+ #!-complex-float-vops
(deftransform - ((w z) ((complex ,type) (complex ,type)) *)
'(complex (- (realpart w) (realpart z))
(- (imagpart w) (imagpart z))))
;; Add and subtract a complex and a real.
+ #!-complex-float-vops
(deftransform + ((w z) ((complex ,type) real) *)
- '(complex (+ (realpart w) z) (imagpart w)))
+ `(complex (+ (realpart w) z)
+ (+ (imagpart w) ,(coerce 0 ',type))))
+ #!-complex-float-vops
(deftransform + ((z w) (real (complex ,type)) *)
- '(complex (+ (realpart w) z) (imagpart w)))
+ `(complex (+ (realpart w) z)
+ (+ (imagpart w) ,(coerce 0 ',type))))
;; Add and subtract a real and a complex number.
+ #!-complex-float-vops
(deftransform - ((w z) ((complex ,type) real) *)
- '(complex (- (realpart w) z) (imagpart w)))
+ `(complex (- (realpart w) z)
+ (- (imagpart w) ,(coerce 0 ',type))))
+ #!-complex-float-vops
(deftransform - ((z w) (real (complex ,type)) *)
- '(complex (- z (realpart w)) (- (imagpart w))))
+ `(complex (- z (realpart w))
+ (- ,(coerce 0 ',type) (imagpart w))))
;; Multiply and divide two complex numbers.
+ #!-complex-float-vops
(deftransform * ((x y) ((complex ,type) (complex ,type)) *)
'(let* ((rx (realpart x))
(ix (imagpart x))
@@ -1207,39 +1225,81 @@
(complex (- (* rx ry) (* ix iy))
(+ (* rx iy) (* ix ry)))))
(deftransform / ((x y) ((complex ,type) (complex ,type)) *)
+ #!-complex-float-vops
'(let* ((rx (realpart x))
(ix (imagpart x))
(ry (realpart y))
(iy (imagpart y)))
(if (> (abs ry) (abs iy))
(let* ((r (/ iy ry))
- (dn (* ry (+ 1 (* r r)))))
+ (dn (+ ry (* r iy))))
(complex (/ (+ rx (* ix r)) dn)
(/ (- ix (* rx r)) dn)))
(let* ((r (/ ry iy))
- (dn (* iy (+ 1 (* r r)))))
+ (dn (+ iy (* r ry))))
(complex (/ (+ (* rx r) ix) dn)
- (/ (- (* ix r) rx) dn))))))
+ (/ (- (* ix r) rx) dn)))))
+ #!+complex-float-vops
+ `(let* ((cs (conjugate (sb!vm::swap-complex x)))
+ (ry (realpart y))
+ (iy (imagpart y)))
+ (if (> (abs ry) (abs iy))
+ (let* ((r (/ iy ry))
+ (dn (+ ry (* r iy))))
+ (/ (+ x (* cs r)) dn))
+ (let* ((r (/ ry iy))
+ (dn (+ iy (* r ry))))
+ (/ (+ (* x r) cs) dn)))))
;; Multiply a complex by a real or vice versa.
+ #!-complex-float-vops
(deftransform * ((w z) ((complex ,type) real) *)
'(complex (* (realpart w) z) (* (imagpart w) z)))
+ #!-complex-float-vops
(deftransform * ((z w) (real (complex ,type)) *)
'(complex (* (realpart w) z) (* (imagpart w) z)))
- ;; Divide a complex by a real.
+ ;; Divide a complex by a real or vice versa.
+ #!-complex-float-vops
(deftransform / ((w z) ((complex ,type) real) *)
'(complex (/ (realpart w) z) (/ (imagpart w) z)))
+ (deftransform / ((x y) (,type (complex ,type)) *)
+ #!-complex-float-vops
+ '(let* ((ry (realpart y))
+ (iy (imagpart y)))
+ (if (> (abs ry) (abs iy))
+ (let* ((r (/ iy ry))
+ (dn (+ ry (* r iy))))
+ (complex (/ x dn)
+ (/ (- (* x r)) dn)))
+ (let* ((r (/ ry iy))
+ (dn (+ iy (* r ry))))
+ (complex (/ (* x r) dn)
+ (/ (- x) dn)))))
+ #!+complex-float-vops
+ '(let* ((ry (realpart y))
+ (iy (imagpart y)))
+ (if (> (abs ry) (abs iy))
+ (let* ((r (/ iy ry))
+ (dn (+ ry (* r iy))))
+ (/ (complex x (- (* x r))) dn))
+ (let* ((r (/ ry iy))
+ (dn (+ iy (* r ry))))
+ (/ (complex (* x r) (- x)) dn)))))
;; conjugate of complex number
+ #!-complex-float-vops
(deftransform conjugate ((z) ((complex ,type)) *)
'(complex (realpart z) (- (imagpart z))))
;; CIS
(deftransform cis ((z) ((,type)) *)
'(complex (cos z) (sin z)))
;; comparison
+ #!-complex-float-vops
(deftransform = ((w z) ((complex ,type) (complex ,type)) *)
'(and (= (realpart w) (realpart z))
(= (imagpart w) (imagpart z))))
+ #!-complex-float-vops
(deftransform = ((w z) ((complex ,type) real) *)
'(and (= (realpart w) z) (zerop (imagpart w))))
+ #!-complex-float-vops
(deftransform = ((w z) (real (complex ,type)) *)
'(and (= (realpart z) w) (zerop (imagpart z)))))))
@@ -735,10 +735,17 @@ core and return a descriptor to it."
(let ((des (allocate-unboxed-object *dynamic* sb!vm:n-word-bits
(1- sb!vm:complex-single-float-size)
sb!vm:complex-single-float-widetag)))
- (write-wordindexed des sb!vm:complex-single-float-real-slot
- (make-random-descriptor (single-float-bits (realpart num))))
- (write-wordindexed des sb!vm:complex-single-float-imag-slot
- (make-random-descriptor (single-float-bits (imagpart num))))
+ #!-x86-64
+ (progn
+ (write-wordindexed des sb!vm:complex-single-float-real-slot
+ (make-random-descriptor (single-float-bits (realpart num))))
+ (write-wordindexed des sb!vm:complex-single-float-imag-slot
+ (make-random-descriptor (single-float-bits (imagpart num)))))
+ #!+x86-64
+ (write-wordindexed des sb!vm:complex-single-float-data-slot
+ (make-random-descriptor
+ (logior (ldb (byte 32 0) (single-float-bits (realpart num)))
+ (ash (single-float-bits (imagpart num)) 32))))
des))
(defun complex-double-float-to-core (num)
@@ -342,13 +342,17 @@
(define-primitive-object (complex-single-float
:lowtag other-pointer-lowtag
:widetag complex-single-float-widetag)
+ #!+x86-64
+ (data :c-type "struct { float data[2]; } ")
+ #!-x86-64
(real :c-type "float")
+ #!-x86-64
(imag :c-type "float"))
(define-primitive-object (complex-double-float
:lowtag other-pointer-lowtag
:widetag complex-double-float-widetag)
- #!-x86-64 (filler)
+ (filler)
(real :c-type "double" :length #!-x86-64 2 #!+x86-64 1)
(imag :c-type "double" :length #!-x86-64 2 #!+x86-64 1))
@@ -584,10 +584,11 @@
(values)))
;;;; transforms for EQL of floating point values
-
+#!-x86-64
(deftransform eql ((x y) (single-float single-float))
'(= (single-float-bits x) (single-float-bits y)))
+#!-x86-64
(deftransform eql ((x y) (double-float double-float))
'(and (= (double-float-low-bits x) (double-float-low-bits y))
(= (double-float-high-bits x) (double-float-high-bits y))))
@@ -3499,7 +3499,13 @@
(cond ((or (and (csubtypep x-type (specifier-type 'float))
(csubtypep y-type (specifier-type 'float)))
(and (csubtypep x-type (specifier-type '(complex float)))
- (csubtypep y-type (specifier-type '(complex float)))))
+ (csubtypep y-type (specifier-type '(complex float))))
+ #!+complex-float-vops
+ (and (csubtypep x-type (specifier-type '(or single-float (complex single-float))))
+ (csubtypep y-type (specifier-type '(or single-float (complex single-float)))))
+ #!+complex-float-vops
+ (and (csubtypep x-type (specifier-type '(or double-float (complex double-float))))
+ (csubtypep y-type (specifier-type '(or double-float (complex double-float))))))
;; They are both floats. Leave as = so that -0.0 is
;; handled correctly.
(give-up-ir1-transform))
Oops, something went wrong.

0 comments on commit a157ed0

Please sign in to comment.