This repository has been archived on 2022-08-10. You can view files and clone it, but cannot push or open issues or pull requests.
chez-openbsd/mats/cfl.ms

378 lines
13 KiB
Scheme
Raw Normal View History

2022-07-29 15:12:07 +02:00
;;; cfl.ms
;;; Copyright 1984-2017 Cisco Systems, Inc.
;;;
;;; Licensed under the Apache License, Version 2.0 (the "License");
;;; you may not use this file except in compliance with the License.
;;; You may obtain a copy of the License at
;;;
;;; http://www.apache.org/licenses/LICENSE-2.0
;;;
;;; Unless required by applicable law or agreed to in writing, software
;;; distributed under the License is distributed on an "AS IS" BASIS,
;;; WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
;;; See the License for the specific language governing permissions and
;;; limitations under the License.
(define *fuzz* 1e-14)
(define ~=
(lambda (x y)
(or (= x y)
(and (fl~= (inexact (real-part x))
(inexact (real-part y)))
(fl~= (inexact (imag-part x))
(inexact (imag-part y)))))))
(define fl~=
(lambda (x y)
(cond
[(and (fl>= (flabs x) 2.0) (fl>= (flabs y) 2.0))
(fl~= (fl/ x 2.0) (fl/ y 2.0))]
[(and (fl< 0.0 (flabs x) 1.0) (fl< 0.0 (flabs y) 1.0))
(fl~= (fl* x 2.0) (fl* y 2.0))]
[else (let ([d (flabs (fl- x y))])
(or (fl<= d *fuzz*)
(begin (printf "fl~~=: ~s~%" d) #f)))])))
(define cfl~=
(lambda (x y)
(and (fl~= (cfl-real-part x) (cfl-real-part y))
(fl~= (cfl-imag-part x) (cfl-imag-part y)))))
(define zero 0.0)
(define a 1.1)
(define b +1.1i)
(define c 1.1+1.1i)
(define aa 1.21)
(define ab +1.21i)
(define ac 1.21+1.21i)
(define bb -1.21)
(define bc -1.21+1.21i)
(define cc +2.42i)
(mat cflonum?
(not (cflonum? 3))
(not (cflonum? 18/2))
(not (cflonum? 1+0i))
(not (cflonum? 23084982309482034820348023423048230482304))
(not (cflonum? 203480234802384/23049821))
(not (cflonum? -3/4))
(not (cflonum? -1))
(not (cflonum? 0))
(not (cflonum? -12))
(cflonum? 3.5)
(cflonum? 1.8e-10)
(cflonum? -3e5)
(cflonum? -1231.2344)
(cflonum? 3+5.0i)
(cflonum? 1.8e10@10)
(cflonum? -3e5+1.0i)
(cflonum? -1.0i)
(cflonum? +1.0i)
(not (cflonum? 'a))
(not (cflonum? "hi"))
(not (cflonum? (cons 3 4)))
(cflonum? a)
(cflonum? b)
(cflonum? c)
)
(mat fl-make-rectangular
(error? (fl-make-rectangular 3 'a))
(error? (fl-make-rectangular 'b 4))
(error? (fl-make-rectangular 3 -4))
(eqv? (fl-make-rectangular 3.0 -4.0) 3.0-4.0i)
(eqv? (fl-make-rectangular a a) c)
)
(mat cfl-real-part
(error? (cfl-real-part 'a))
(error? (cfl-real-part 3/2))
(eqv? (cfl-real-part 3.2) 3.2)
(eqv? (cfl-real-part -1.0+2.0i) -1.0)
(eqv? (cfl-real-part a) a)
(eqv? (cfl-real-part c) a)
(eqv? (cfl-real-part b) zero)
)
(mat cfl-imag-part
(error? (cfl-imag-part 'a))
(error? (cfl-imag-part -3))
(eqv? (cfl-imag-part 3.2) zero)
(eqv? (cfl-imag-part -1.0+2.0i) 2.0)
(eqv? (cfl-imag-part a) zero)
(eqv? (cfl-imag-part c) a)
(eqv? (cfl-imag-part b) a)
)
(mat cfl-conjugate
(error? (cfl-conjugate 'a))
(eqv? (cfl-conjugate 3.2) 3.2)
(eqv? (cfl-conjugate 3.2+2.0i) 3.2-2.0i)
(eqv? (cfl-conjugate a) a)
(eqv? (cfl-conjugate c) (+ a (- b)))
(eqv? (cfl-conjugate b) -1.1i)
)
(mat conjugate
(error? (conjugate 'a))
(eqv? (conjugate 3.2) 3.2)
(eqv? (conjugate 3.2+2.0i) 3.2-2.0i)
)
(mat cfl-magnitude-squared
(error? (cfl-magnitude-squared 'a))
(eqv? (cfl-magnitude-squared 3.2) (fl* 3.2 3.2))
(eqv? (cfl-magnitude-squared 3.5-2.0i) 16.25)
(fl~= (cfl-magnitude-squared 3.5@2.0) 12.25)
)
(mat magnitude-squared
(error? (magnitude-squared 'a))
(eqv? (magnitude-squared 3.5) 12.25)
(eqv? (magnitude-squared 3.5-2.0i) 16.25)
(fl~= (magnitude-squared 3.5@2.0) 12.25)
)
(mat cfl+
(error? (cfl+ 'a))
(error? (cfl+ 'a 3))
(error? (cfl+ 'a 3 4))
(eqv? (cfl+) zero)
(eqv? (cfl+ a) a)
(eqv? (cfl+ b) b)
(eqv? (cfl+ c) c)
(eqv? (cfl+ a b) c)
(cfl~= (cfl+ a b c) (cfl+ a (cfl+ b c)))
(cfl~= (cfl+ a b c a b c) (cfl+ (cfl+ a b c) (cfl+ a b c)))
(cfl~= (cfl+ 1+2.0i 3.0) 4.0+2.0i)
(cfl~= (cfl+ 1.0+2.2i -3.7+5.3i) -2.7+7.5i)
(cfl~= (cfl+ 1.0+2.2i -3.7) -2.7+2.2i)
(cfl~= (cfl+ 1.0 -3.7+5.3i) -2.7+5.3i)
(cfl~= (cfl+ 1.0+2.2i +5.3i) 1.0+7.5i)
(cfl~= (cfl+ +2.2i -3.7+5.3i) -3.7+7.5i)
(cfl~= (cfl+ 26.0 2.0) 28.0)
(test-cp0-expansion eqv? '(cfl+) zero)
(test-cp0-expansion eqv? `(cfl+ ,a) a)
(test-cp0-expansion eqv? `(cfl+ ,b) b)
(test-cp0-expansion eqv? `(cfl+ ,c) c)
(test-cp0-expansion eqv? `(cfl+ ,a ,b) c)
(test-cp0-expansion cfl~= `(cfl+ ,a ,b ,c) (cfl+ a (cfl+ b c)))
(test-cp0-expansion cfl~= `(cfl+ ,a ,b ,c ,a ,b ,c) (cfl+ (cfl+ a b c) (cfl+ a b c)))
(test-cp0-expansion cfl~= '(cfl+ 1+2.0i 3.0) 4.0+2.0i)
(test-cp0-expansion cfl~= '(cfl+ 1.0+2.2i -3.7+5.3i) -2.7+7.5i)
(test-cp0-expansion cfl~= '(cfl+ 1.0+2.2i -3.7) -2.7+2.2i)
(test-cp0-expansion cfl~= '(cfl+ 1.0 -3.7+5.3i) -2.7+5.3i)
(test-cp0-expansion cfl~= '(cfl+ 1.0+2.2i +5.3i) 1.0+7.5i)
(test-cp0-expansion cfl~= '(cfl+ +2.2i -3.7+5.3i) -3.7+7.5i)
(test-cp0-expansion cfl~= '(cfl+ 26.0 2.0) 28.0)
)
(mat cfl-
(error? (cfl- 'a))
(error? (cfl- 'a 3))
(error? (cfl- 'a 3 4))
(error? (cfl-))
(eqv? (cfl- a) -1.1)
(eqv? (cfl- b) -0.0-1.1i)
(eqv? (cfl- c) -1.1-1.1i)
(eqv? (cfl- a a) zero)
(cfl~= (cfl- b b) zero)
(cfl~= (cfl- c c) zero)
(eqv? (cfl- c a) b)
(cfl~= (cfl- c b) a)
(cfl~= (cfl- a b c) (cfl- (cfl- a b) c))
(cfl~= (cfl- a b c a b c) (cfl- a (cfl+ b c a b c)))
(cfl~= (cfl- 1+2.0i 3.0) -2.0+2.0i)
(cfl~= (cfl- 1.0+2.2i -3.7+5.3i) 4.7-3.1i)
(cfl~= (cfl- 1.0+2.2i -3.7) 4.7+2.2i)
(cfl~= (cfl- 1.0 -3.7+5.3i) 4.7-5.3i)
(cfl~= (cfl- 1.0+2.2i +5.3i) 1.0-3.1i)
(cfl~= (cfl- +2.2i -3.7+5.3i) 3.7-3.1i)
(cfl~= (cfl- 26.0 2.0) 24.0)
(andmap
(lambda (a)
(andmap
(lambda (b)
(andmap
(lambda (c) (eqv? (cfl- a b c) (cfl- (cfl- a b) c)))
'(0.0 -0.0)))
'(0.0 -0.0)))
'(0.0 -0.0))
(let ()
(define-syntax ff
(syntax-rules ()
[(_ k1 k2) (lambda (x) (eqv? (cfl- k1 x k2) (cfl- (cfl- k1 x) k2)))]))
(andmap
(lambda (p) (and (p +0.0) (p -0.0)))
(list (ff +0.0 +0.0) (ff +0.0 -0.0) (ff -0.0 +0.0) (ff -0.0 -0.0))))
(error? (cfl- 3.0 5.4 'a))
(error? (cfl- 'a 3.0 5.4))
(error? (cfl- 3.0 'a 5.4))
(eqv? (cfl- 5.0 4.0 3.0 2.0) -4.0)
(eqv? (cfl- 5.0 4.0 3.0 2.0 1.0 0.0 -1.0 -2.0) -2.0)
(cfl~= (cfl- 1e30 1e30 7.0) -7.0)
(test-cp0-expansion eqv? `(cfl- ,a) -1.1)
(test-cp0-expansion eqv? `(cfl- ,b) -0.0-1.1i)
(test-cp0-expansion eqv? `(cfl- ,c) -1.1-1.1i)
(test-cp0-expansion eqv? `(cfl- ,a ,a) zero)
(test-cp0-expansion cfl~= `(cfl- ,b ,b) zero)
(test-cp0-expansion cfl~= `(cfl- ,c ,c) zero)
(test-cp0-expansion eqv? `(cfl- ,c ,a) b)
(test-cp0-expansion cfl~= `(cfl- ,c ,b) a)
(test-cp0-expansion cfl~= `(cfl- ,a ,b ,c) (cfl- (cfl- a b) c))
(test-cp0-expansion cfl~= `(cfl- ,a ,b ,c ,a ,b ,c) (cfl- a (cfl+ b c a b c)))
(test-cp0-expansion cfl~= '(cfl- 1+2.0i 3.0) -2.0+2.0i)
(test-cp0-expansion cfl~= '(cfl- 1.0+2.2i -3.7+5.3i) 4.7-3.1i)
(test-cp0-expansion cfl~= '(cfl- 1.0+2.2i -3.7) 4.7+2.2i)
(test-cp0-expansion cfl~= '(cfl- 1.0 -3.7+5.3i) 4.7-5.3i)
(test-cp0-expansion cfl~= '(cfl- 1.0+2.2i +5.3i) 1.0-3.1i)
(test-cp0-expansion cfl~= '(cfl- +2.2i -3.7+5.3i) 3.7-3.1i)
(test-cp0-expansion cfl~= '(cfl- 26.0 2.0) 24.0)
(test-cp0-expansion eqv? '(cfl- 5.0 4.0 3.0 2.0) -4.0)
(test-cp0-expansion eqv? '(cfl- 5.0 4.0 3.0 2.0 1.0 0.0 -1.0 -2.0) -2.0)
(test-cp0-expansion cfl~= '(cfl- 1e30 1e30 7.0) -7.0)
)
(mat cfl*
(error? (cfl* 'a))
(error? (cfl* 'a 3))
(error? (cfl* 'a 3 4))
(eqv? (cfl*) 1.0)
(eqv? (cfl* a) a)
(eqv? (cfl* b) b)
(eqv? (cfl* c) c)
(eqv? (cfl* zero a) zero)
(cfl~= (cfl* zero b) zero)
(cfl~= (cfl* zero c) zero)
(cfl~= (cfl* a a) aa)
(cfl~= (cfl* a b) ab)
(cfl~= (cfl* a c) ac)
(cfl~= (cfl* b b) bb)
(cfl~= (cfl* b c) bc)
(cfl~= (cfl* c c) cc)
(cfl~= (cfl* a b c) (cfl* a (cfl* b c)))
(cfl~= (cfl* a b c a b c) (cfl* (cfl* a b c) (cfl* a b c)))
(cfl~= (cfl* 1+2.0i 3.0) 3.0+6.0i)
(cfl~= (cfl* 1.0+2.0i 3.0+4.0i) -5.0+10.0i)
(cfl~= (cfl* 1.0+2.0i 3.0) 3.0+6.0i)
(cfl~= (cfl* -2.0 3.0+4.0i) -6.0-8.0i)
(cfl~= (cfl* 1.0+2.0i +4.0i) -8.0+4.0i)
(cfl~= (cfl* +2.0i 3.0+4.0i) -8.0+6.0i)
(cfl~= (cfl* 26.0 2.0) 52.0)
(test-cp0-expansion eqv? '(cfl*) 1.0)
(test-cp0-expansion eqv? `(cfl* ,a) a)
(test-cp0-expansion eqv? `(cfl* ,b) b)
(test-cp0-expansion eqv? `(cfl* ,c) c)
(test-cp0-expansion eqv? `(cfl* ,zero ,a) zero)
(test-cp0-expansion cfl~= `(cfl* ,zero ,b) zero)
(test-cp0-expansion cfl~= `(cfl* ,zero ,c) zero)
(test-cp0-expansion cfl~= `(cfl* ,a ,a) aa)
(test-cp0-expansion cfl~= `(cfl* ,a ,b) ab)
(test-cp0-expansion cfl~= `(cfl* ,a ,c) ac)
(test-cp0-expansion cfl~= `(cfl* ,b ,b) bb)
(test-cp0-expansion cfl~= `(cfl* ,b ,c) bc)
(test-cp0-expansion cfl~= `(cfl* ,c ,c) cc)
(test-cp0-expansion cfl~= `(cfl* ,a ,b ,c) (cfl* a (cfl* b c)))
(test-cp0-expansion cfl~= `(cfl* ,a ,b ,c ,a ,b ,c) (cfl* (cfl* a b c) (cfl* a b c)))
(test-cp0-expansion cfl~= '(cfl* 1+2.0i 3.0) 3.0+6.0i)
(test-cp0-expansion cfl~= '(cfl* 1.0+2.0i 3.0+4.0i) -5.0+10.0i)
(test-cp0-expansion cfl~= '(cfl* 1.0+2.0i 3.0) 3.0+6.0i)
(test-cp0-expansion cfl~= '(cfl* -2.0 3.0+4.0i) -6.0-8.0i)
(test-cp0-expansion cfl~= '(cfl* 1.0+2.0i +4.0i) -8.0+4.0i)
(test-cp0-expansion cfl~= '(cfl* +2.0i 3.0+4.0i) -8.0+6.0i)
(test-cp0-expansion cfl~= '(cfl* 26.0 2.0) 52.0)
)
(mat cfl/
(error? (cfl/ 'a))
(error? (cfl/ 'a 3))
(error? (cfl/ 'a 3 4))
(error? (cfl/))
(fl~= (cfl/ a) (fl/ a))
(eqv? (cfl/ zero a) zero)
(cfl~= (cfl/ zero b) zero)
(cfl~= (cfl/ zero c) zero)
(cfl~= (cfl/ a a) 1.0)
(cfl~= (cfl/ b b) 1.0)
(cfl~= (cfl/ c c) 1.0)
(cfl~= (cfl/ aa a) a)
(cfl~= (cfl/ ab b) a)
(cfl~= (cfl/ ab a) b)
(cfl~= (cfl/ ac c) a)
(cfl~= (cfl/ ac a) c)
(cfl~= (cfl/ bc c) b)
(cfl~= (cfl/ bc b) c)
(cfl~= (cfl/ cc c) c)
(cfl~= (cfl/ a b c) (cfl/ (cfl/ a b) c))
(cfl~= (cfl/ a b c a b c) (cfl/ a (cfl* b c a b c)))
(cfl~= (cfl/ 3+6.0i 3.0) 1.0+2.0i)
(cfl~= (cfl/ -5.0+10.0i 1.0+2.0i) 3.0+4.0i)
(cfl~= (cfl/ -6.0-8.0i -2.0) 3.0+4.0i)
(cfl~= (cfl/ 26.0 3.0-2.0i) 6.0+4.0i)
(cfl~= (cfl/ -8.0+6.0i +2.0i) 3.0+4.0i)
(cfl~= (cfl/ +26.0i 3.0+2.0i) 4.0+6.0i)
(cfl~= (cfl/ 26.0 2.0) 13.0)
(andmap
(lambda (a)
(andmap
(lambda (b)
(andmap
(lambda (c) (eqv? (cfl/ a b c) (cfl/ (cfl/ a b) c)))
'(1e300 1e250)))
'(1e300 1e250)))
'(1e300 1e250))
(error? (cfl/ 3.0 5.4 'a))
(error? (cfl/ 'a 3.0 5.4))
(error? (cfl/ 3.0 'a 5.4))
(eqv? (cfl/ 16.0 2.0 -2.0 2.0) -2.0)
(eqv? (cfl/ 16.0 2.0 -2.0 2.0 4.0 1.0 -1.0) 0.5)
(test-cp0-expansion eqv? `(cfl/ ,zero ,a) zero)
(test-cp0-expansion eqv? '(cfl/ 16.0 2.0 -2.0 2.0) -2.0)
(test-cp0-expansion eqv? '(cfl/ 16.0 2.0 -2.0 2.0 4.0 1.0 -1.0) 0.5)
(test-cp0-expansion cfl~= `(cfl/ ,zero ,b) zero)
(test-cp0-expansion cfl~= `(cfl/ ,zero ,c) zero)
(test-cp0-expansion cfl~= `(cfl/ ,a ,a) 1.0)
(test-cp0-expansion cfl~= `(cfl/ ,b ,b) 1.0)
(test-cp0-expansion cfl~= `(cfl/ ,c ,c) 1.0)
(test-cp0-expansion cfl~= `(cfl/ ,aa ,a) a)
(test-cp0-expansion cfl~= `(cfl/ ,ab ,b) a)
(test-cp0-expansion cfl~= `(cfl/ ,ab ,a) b)
(test-cp0-expansion cfl~= `(cfl/ ,ac ,c) a)
(test-cp0-expansion cfl~= `(cfl/ ,ac ,a) c)
(test-cp0-expansion cfl~= `(cfl/ ,bc ,c) b)
(test-cp0-expansion cfl~= `(cfl/ ,bc ,b) c)
(test-cp0-expansion cfl~= `(cfl/ ,cc ,c) c)
(test-cp0-expansion cfl~= `(cfl/ ,a ,b ,c) (cfl/ (cfl/ a b) c))
(test-cp0-expansion cfl~= `(cfl/ ,a ,b ,c ,a ,b ,c) (cfl/ a (cfl* b c a b c)))
(test-cp0-expansion cfl~= '(cfl/ 3+6.0i 3.0) 1.0+2.0i)
(test-cp0-expansion cfl~= '(cfl/ -5.0+10.0i 1.0+2.0i) 3.0+4.0i)
(test-cp0-expansion cfl~= '(cfl/ -6.0-8.0i -2.0) 3.0+4.0i)
(test-cp0-expansion cfl~= '(cfl/ 26.0 3.0-2.0i) 6.0+4.0i)
(test-cp0-expansion cfl~= '(cfl/ -8.0+6.0i +2.0i) 3.0+4.0i)
(test-cp0-expansion cfl~= '(cfl/ +26.0i 3.0+2.0i) 4.0+6.0i)
(test-cp0-expansion cfl~= '(cfl/ 26.0 2.0) 13.0)
)
(mat cfl=
(error? (cfl= 'a))
(error? (cfl= 'a 3))
(error? (cfl= 'a 3 4))
(error? (cfl=))
(cfl= a a)
(cfl= b b)
(cfl= c c)
(cfl= (- c c) zero)
(cfl= (+ a b) c)
(not (cfl= a b))
(cfl= 1.1+1.1i c)
(cfl= c 1.1+1.1i c)
(not (cfl= c 1.1+1.1i c a))
(not (cfl= 3+6.0i 3.0))
(not (cfl= 3+6.0i +6.0i))
(cfl= 1.0+2.0i 1.0+2.0i)
(cfl= 5.4 5.4)
)