| blob | 04636b388548230ae7b03d644ede1b50746e93de |
1 ;------------------------------------------------------------------------------
3 (define-type num
4 lo16
5 shr16
6 )
8 (define min-fix -2)
9 (define max-fix 2)
11 (define fixnums (make-vector (+ (- max-fix min-fix) 1)))
13 (define fixnum
14 (lambda (n)
15 (vector-ref fixnums (- n min-fix))))
17 (let loop ((n min-fix))
18 (if (<= n max-fix)
19 (begin
20 (vector-set! fixnums (- n min-fix) (make-num (modulo n #x10000) #f))
21 (loop (+ n 1)))))
23 (let loop ((n min-fix))
24 (if (<= n max-fix)
25 (begin
26 (num-shr16-set!
27 (fixnum n)
28 (if (< n 0) (fixnum -1) (fixnum 0)))
29 (loop (+ n 1)))))
31 (define minus1 (fixnum -1))
32 (define zero (fixnum 0))
33 (define one (fixnum 1))
35 (define make-num
36 (let ((mn make-num))
37 (lambda (lo16 shr16)
38 (if (or (< lo16 0) (> lo16 #xffff))
39 (error "lo16 out of range"))
40 (mn lo16 shr16))))
42 (define check-norm
43 (lambda (name result)
44 (cond ((and (eq? (num-shr16 result) zero)
45 (<= (num-lo16 result) max-fix))
46 (if (not (eq? result (fixnum (num-lo16 result))))
47 (pp (list name "nonneg number is not normalized" (dec result) (num-lo16 result)))))
48 ((and (eq? (num-shr16 result) minus1)
49 (>= (num-lo16 result) (+ #x10000 min-fix)))
50 (if (not (eq? result (fixnum (- (num-lo16 result) #x10000))))
51 (pp (list name "neg number is not normalized" (dec result) (- (num-lo16 result) #x10000))))))))
53 (define check-normalized
54 (lambda (name f)
55 (lambda args
56 (let ((result (apply f args)))
57 (check-norm name result)
58 result))))
60 (define invert
61 (lambda (x)
62 (if (eq? x zero) minus1 zero)))
64 (define sign?
65 (lambda (x)
66 (check-norm 'sign?.x x)
67 (or (eq? x zero) (eq? x minus1))))
69 (define neg?
70 (lambda (x)
71 (check-norm 'neg?.x x)
72 (cond ((eq? x zero) #f)
73 ((eq? x minus1) #t)
74 (else (neg? (num-shr16 x))))))
76 (define lt?
77 (lambda (x y)
78 (check-norm 'lt?.x x)
79 (check-norm 'lt?.y y)
80 (< (cmp x y) 0)))
82 (define le?
83 (lambda (x y)
84 (check-norm 'le?.x x)
85 (check-norm 'le?.y y)
86 (<= (cmp x y) 0)))
88 (define cmp
89 (lambda (x y)
90 (check-norm 'cmp.x x)
91 (check-norm 'cmp.y y)
92 (let loop ((x x) (y y) (c 0))
93 (cond ((sign? x)
94 (cond ((eq? x y) c)
95 ((neg? y) 1)
96 (else -1)))
97 ((sign? y)
98 (cond ;((eq? x y) c)
99 ((neg? x) -1)
100 (else 1)))
101 (else
102 (loop (num-shr16 x)
103 (num-shr16 y)
104 (let ((x-lo16 (num-lo16 x))
105 (y-lo16 (num-lo16 y)))
106 (cond ((= x-lo16 y-lo16) c)
107 ((< x-lo16 y-lo16) -1)
108 (else 1)))))))))
110 (define int-len-nonneg
111 (lambda (x)
112 (check-norm 'int-len-nonneg.x x)
113 (let loop1 ((x x) (n 0))
114 (let ((shr16 (num-shr16 x)))
115 (if (eq? shr16 zero)
116 (let loop2 ((lo16 (num-lo16 x)) (n n))
117 (if (= lo16 0)
118 n
119 (loop2 (arithmetic-shift lo16 -1)
120 (+ n 1))))
121 (loop1 shr16 (+ 16 n)))))))
123 (define shift-left
124 (check-normalized 'shift-left
125 (lambda (x n)
126 (check-norm 'shift-left.x x)
127 (if (< n 0) (error "xxxxxxxxxxx"))
128 (if (eq? x zero)
129 zero
130 (let loop1 ((x x) (n n))
131 (if (= (bitwise-and n 15) 0)
132 (let loop2 ((x x) (n n))
133 (if (= n 0)
134 x
135 (loop2 (make-num #x0000 x)
136 (- n 16))))
137 (loop1 (shl x)
138 (- n 1))))))))
140 (define shl
141 (check-normalized 'shl
142 (lambda (x)
143 (check-norm 'shl.x x)
144 (let loop ((x x) (minus-c zero) (r #f))
145 (cond ((eq? x minus-c)
146 (norm r x))
147 (else
148 (let ((z (+ (* 2 (num-lo16 x)) (if (eq? minus-c zero) 0 1))))
149 (loop (num-shr16 x)
150 (if (>= z #x10000) minus1 zero)
151 (make-num (modulo z #x10000) r)))))))))
153 (define shr
154 (check-normalized 'shr
155 (lambda (x)
156 (check-norm 'shr.x x)
157 (let loop ((x x) (r #f))
158 (cond ((eq? x zero)
159 (norm r x))
160 ((eq? x minus1)
161 (norm r x))
162 (else
163 (let ((z (+ (arithmetic-shift (num-lo16 x) -1)
164 (if (even? (num-lo16 (num-shr16 x))) #x0000 #x8000))))
165 (loop (num-shr16 x)
166 (make-num z r)))))))))
168 (define add
169 (check-normalized 'add
170 (lambda (x y)
171 (check-norm 'add.x x)
172 (check-norm 'add.y y)
173 (let loop ((x x) (y y) (minus-c zero) (r #f))
174 (cond ((eq? x minus-c)
175 (norm r y))
176 ((eq? y minus-c)
177 (norm r x))
178 (else
179 (let ((z (+ (num-lo16 x) (num-lo16 y) (if (eq? minus-c zero) 0 1))))
180 (loop (num-shr16 x)
181 (num-shr16 y)
182 (if (>= z #x10000) minus1 zero)
183 (make-num (modulo z #x10000) r)))))))))
185 (define sub
186 (check-normalized 'sub
187 (lambda (x y)
188 (check-norm 'sub.x x)
189 (check-norm 'sub.y y)
190 (let loop ((x x) (y y) (minus-c minus1) (r #f))
191 (cond ((and (eq? x minus-c) (sign? y))
192 (norm r (invert y)))
193 ((eq? y (invert minus-c))
194 (norm r x))
195 (else
196 (let ((z (+ (num-lo16 x) (- #xffff (num-lo16 y)) (if (eq? minus-c zero) 0 1))))
197 (loop (num-shr16 x)
198 (num-shr16 y)
199 (if (>= z #x10000) minus1 zero)
200 (make-num (modulo z #x10000) r)))))))))
202 (define scale
203 (check-normalized 'scale
204 (lambda (n16 x)
205 (check-norm 'scale.x x)
206 (cond ((or (= n16 0)
207 (eq? x zero))
208 zero)
209 ((eq? n16 1)
210 x)
211 (else
212 (let loop ((x x) (carry 0) (r #f))
213 (cond ((eq? x zero)
214 (norm r (build-nonneg carry)))
215 ((eq? x minus1)
216 (norm r (build-neg (+ (modulo (- n16) #x10000) carry))))
217 (else
218 (let ((m (+ (* n16 (num-lo16 x)) carry)))
219 (loop (num-shr16 x)
220 (arithmetic-shift m -16)
221 (make-num (modulo m #x10000) r)))))))))))
223 (define neg
224 (check-normalized 'neg
225 (lambda (x)
226 (check-norm 'neg.x x)
227 (sub zero x))))
229 (define build-nonneg
230 (check-normalized 'build-nonneg
231 (lambda (lo16)
232 (if (<= lo16 max-fix) (fixnum lo16) (make-num lo16 zero)))))
234 (define build-neg
235 (check-normalized 'build-neg
236 (lambda (lo16)
237 (if (>= lo16 (+ #x10000 min-fix)) (fixnum (- lo16 #x10000)) (make-num lo16 minus1)))))
239 (define mul
240 (check-normalized 'mul
241 (lambda (x y)
242 (check-norm 'mul.x x)
243 (check-norm 'mul.y y)
244 (if (neg? x)
245 (neg (mulnonneg (neg x) y))
246 (mulnonneg x y)))))
248 (define mulnonneg
249 (lambda (x y) ; x is nonnegative
250 (check-norm 'mulnonneg.x x)
251 (check-norm 'mulnonneg.y y)
252 (let loop ((x x) (s zero) (r #f))
253 (if (eq? x zero)
254 (norm r s)
255 (let ((a (add s (scale (num-lo16 x) y))))
256 (loop (num-shr16 x)
257 (num-shr16 a)
258 (make-num (num-lo16 a) r)))))))
260 (define rem
261 (check-normalized 'rem
262 (lambda (x y)
263 (check-norm 'rem.x x)
264 (check-norm 'rem.y y)
265 (if (neg? x)
266 (if (neg? y)
267 (neg (cdr (divnonneg (neg x) (neg y))))
268 (neg (cdr (divnonneg (neg x) y))))
269 (if (neg? y)
270 (cdr (divnonneg x (neg y)))
271 (cdr (divnonneg x y)))))))
273 (define div
274 (check-normalized 'div
275 (lambda (x y)
276 (check-norm 'div.x x)
277 (check-norm 'div.y y)
278 (if (neg? x)
279 (if (neg? y)
280 (car (divnonneg (neg x) (neg y)))
281 (neg (car (divnonneg (neg x) y))))
282 (if (neg? y)
283 (neg (car (divnonneg x (neg y))))
284 (car (divnonneg x y)))))))
286 (define divnonneg2
287 (lambda (x y)
288 (let ((x (decode x)) (y (decode y)))
289 (let ((lx (integer-length x))
290 (ly (integer-length y)))
291 (let ((b (- lx ly)))
292 (let loop ((i b) (a x) (b (arithmetic-shift y b)) (r 0))
293 ; (pp (list i a b r))
294 (if (< i 0)
295 (cons (encode r) (encode a))
296 (begin
297 (if (<= b a)
298 (loop (- i 1)
299 (- a b)
300 (arithmetic-shift b -1)
301 (+ 1 (* 2 r)))
302 (loop (- i 1)
303 a
304 (arithmetic-shift b -1)
305 (* 2 r)))))))))))
307 (define divnonneg
308 (lambda (x y)
309 (check-norm 'divnonneg.x x)
310 (check-norm 'divnonneg.y y)
311 (let ((lx (int-len-nonneg x))
312 (ly (int-len-nonneg y)))
313 (if (not (and (= lx (integer-length (decode x)))
314 (= ly (integer-length (decode y)))))
315 (error "!!!!"))
316 (let ((b (- lx ly)))
317 (if (< b 0)
318 (cons zero x)
319 (let loop ((i b) (a x) (b (shift-left y b)) (r zero))
320 ; (pp (list i (decode a) (decode b) (decode r)))
321 (if (< i 0)
322 (cons r a)
323 (if (le? b a)
324 (loop (- i 1)
325 (sub a b)
326 (shr b)
327 (add one (shl r)))
328 (loop (- i 1)
329 a
330 (shr b)
331 (shl r))))))))))
333 (define norm2
334 (check-normalized 'norm2
335 (lambda (shr16 r)
336 (if r
337 (let ((lo16 (num-lo16 r))
338 (newr (num-shr16 r)))
339 (norm newr
340 (if (cond ((eq? shr16 zero) (= lo16 #x0000))
341 ((eq? shr16 minus1) (= lo16 #xffff))
342 (else #f))
343 shr16
344 (begin
345 (num-shr16-set! r shr16)
346 r))))
347 shr16))))
349 (define norm
350 (check-normalized 'norm
351 (lambda (r shr16)
352 (check-norm 'norm.shr16 shr16)
353 (if r
354 (let ((lo16 (num-lo16 r))
355 (newr (num-shr16 r)))
356 (cond ((eq? shr16 zero)
357 (cond ((<= lo16 max-fix)
358 (norm newr (fixnum lo16)))
359 (else
360 (num-shr16-set! r shr16)
361 (norm newr r))))
362 ((eq? shr16 minus1)
363 (cond ((>= lo16 (+ #x10000 min-fix))
364 (norm newr (fixnum (- lo16 #x10000))))
365 (else
366 (num-shr16-set! r shr16)
367 (norm newr r))))
368 (else
369 (num-shr16-set! r shr16)
370 (norm newr r))))
371 shr16))))
373 (define #%list->string
374 (lambda (lst)
375 (list->string (map integer->char (map decode lst)))))
377 (define #%cons cons)
379 (define #%< lt?)
380 (define #%+ add)
381 (define #%neg neg)
382 (define #%quotient div)
383 (define #%remainder rem)
385 (define #%number->string
386 (lambda (n)
387 (#%list->string
388 (if (#%< n zero)
389 (#%cons (make-num 45 zero) (#%number->string-aux (#%neg n) '()))
390 (#%number->string-aux n '())))))
392 (define #%number->string-aux
393 (lambda (n lst)
394 (let ((rest (#%cons (#%+ (make-num 48 zero) (#%remainder n (make-num 10 zero))) lst)))
395 (if (#%< n (make-num 10 zero))
396 rest
397 (#%number->string-aux (#%quotient n (make-num 10 zero)) rest)))))
399 ;------------------------------------------------------------------------------
401 (define fn
402 (lambda (i)
403 (list (- (expt 2 (+ i -1)) 3)
404 (- (expt 2 (+ i -1)) 2)
405 (- (expt 2 (+ i -1)) 1)
406 (expt 2 (+ i -1))
407 (+ (expt 2 (+ i -1)) 1)
408 (+ (expt 2 (+ i -1)) 2)
409 (+ (expt 2 (+ i -1)) 3)
411 (- (expt 2 (+ i 0)) 3)
412 (- (expt 2 (+ i 0)) 2)
413 (- (expt 2 (+ i 0)) 1)
414 (expt 2 (+ i 0))
415 (+ (expt 2 (+ i 0)) 1)
416 (+ (expt 2 (+ i 0)) 2)
417 (+ (expt 2 (+ i 0)) 3)
419 (- (expt 2 (+ i 1)) 3)
420 (- (expt 2 (+ i 1)) 2)
421 (- (expt 2 (+ i 1)) 1)
422 (expt 2 (+ i 1))
423 (+ (expt 2 (+ i 1)) 1)
424 (+ (expt 2 (+ i 1)) 2)
425 (+ (expt 2 (+ i 1)) 3)
426 )))
428 (define nums1
429 (append (fn (* 16 1))
430 (fn (* 16 2))
431 (fn (* 16 3))))
433 (define nums
434 (append (map - (reverse nums1))
435 (list -10 -9 -8 -7 -6 -5 -4 -3 -2 -1 0 1 2 3 4 5 6 7 8 9 10)
436 nums1))
438 (define dec
439 (lambda (n)
440 (cond ((eq? n zero) '(0))
441 ((eq? n minus1) '(-1))
442 (else (cons (num-lo16 n) (dec (num-shr16 n)))))))
444 (define decode
445 (lambda (n)
446 (cond ((eq? n zero) 0)
447 ((eq? n minus1) -1)
448 (else (+ (* #x10000 (decode (num-shr16 n))) (num-lo16 n))))))
450 (define encode
451 (lambda (n)
452 (cond ((and (>= n min-fix)
453 (<= n max-fix))
454 (fixnum n))
455 (else
456 (make-num (modulo n #x10000)
457 (encode (arithmetic-shift n -16)))))))
459 ;------------------------------------------------------------------------------
461 (pp 'cmp)
462 (for-each
463 (lambda (x)
464 (for-each
465 (lambda (y)
466 (let ((a (encode x))
467 (b (encode y))
468 (c (cond ((< x y) -1) ((= x y) 0) (else 1))))
469 (if (not (equal? c (cmp a b)))
470 (error "incorrect" x y (dec a) (dec b) c (cmp a b)))))
471 nums))
472 nums)
474 (pp 'int-len-nonneg)
475 (for-each
476 (lambda (x)
477 (if (>= x 0)
478 (let ((a (encode x))
479 (c (integer-length x)))
480 (if (not (equal? c (int-len-nonneg a)))
481 (error "incorrect" x (dec a) c (int-len-nonneg a))))))
482 nums)
484 (pp 'shl)
485 (for-each
486 (lambda (x)
487 (let ((a (encode x))
488 (c (encode (arithmetic-shift x 1))))
489 (if (not (equal? (dec c) (dec (shl a))))
490 (error "incorrect" x (dec a) (dec c) (dec (shl a))))))
491 nums)
493 (pp 'shr)
494 (for-each
495 (lambda (x)
496 (let ((a (encode x))
497 (c (encode (arithmetic-shift x -1))))
498 (if (not (equal? (dec c) (dec (shr a))))
499 (error "incorrect" x (dec a) (dec c) (dec (shr a))))))
500 nums)
502 (pp 'add)
503 (for-each
504 (lambda (x)
505 (for-each
506 (lambda (y)
507 (let ((a (encode x))
508 (b (encode y))
509 (c (encode (+ x y))))
510 (if (not (equal? (dec c) (dec (add a b))))
511 (error "incorrect" x y (dec a) (dec b) (dec c) (dec (add a b))))))
512 nums))
513 nums)
515 (pp 'sub)
516 (for-each
517 (lambda (x)
518 (for-each
519 (lambda (y)
520 (let ((a (encode x))
521 (b (encode y))
522 (c (encode (- x y))))
523 (if (not (equal? (dec c) (dec (sub a b))))
524 (error "incorrect" x y (dec a) (dec b) (dec c) (dec (sub a b))))))
525 nums))
526 nums)
528 (pp 'scale)
529 (for-each
530 (lambda (n16)
531 (for-each
532 (lambda (x)
533 (let ((a (encode x))
534 (c (encode (* n16 x))))
535 (if (not (equal? (dec c) (dec (scale n16 a))))
536 (error "incorrect" x (dec a) (dec c) (dec (scale n16 a))))))
537 nums))
538 '(0 1 2 3 4 5 6 7 8 9 10 32767 32768 65535))
540 (pp 'mul)
541 (for-each
542 (lambda (x)
543 (for-each
544 (lambda (y)
545 (let ((a (encode x))
546 (b (encode y))
547 (c (encode (* x y))))
548 (if (not (equal? (dec c) (dec (mul a b))))
549 (error "incorrect" x y (dec a) (dec b) (dec c) (dec (mul a b))))))
550 nums))
551 nums)
553 (pp 'div)
554 (for-each
555 (lambda (x)
556 (for-each
557 (lambda (y)
558 (if (not (= y 0))
559 (let ((a (encode x))
560 (b (encode y))
561 (c (encode (quotient x y))))
562 (if (not (equal? (dec c) (dec (div a b))))
563 (error "incorrect" x y (dec a) (dec b) (dec c) (dec (div a b)))))))
564 nums))
565 nums)
567 (pp 'rem)
568 (for-each
569 (lambda (x)
570 (for-each
571 (lambda (y)
572 (if (not (= y 0))
573 (let ((a (encode x))
574 (b (encode y))
575 (c (encode (remainder x y))))
576 (if (not (equal? (dec c) (dec (rem a b))))
577 (error "incorrect" x y (dec a) (dec b) (dec c) (dec (rem a b)))))))
578 nums))
579 nums)
581 (pp 'fact)
582 (define fact
583 (lambda (n)
584 (if (eq? n zero)
585 one
586 (mul n (fact (sub n one))))))
588 (define n (time (fact (encode 50))))
589 (pp (time (#%number->string n)))