2 (:use
"COMMON-LISP" "SB-ALIEN" "SB-BIGNUM")
3 ;; we need a few very internal symbols
4 (:import-from
"SB-BIGNUM"
5 "%BIGNUM-0-OR-PLUSP" "%NORMALIZE-BIGNUM"
6 "NEGATE-BIGNUM-IN-PLACE")
8 ;; bignum integer operations
13 #:mpz-mul-2exp
; shift left
16 #:mpz-fdiv-2exp
; shift right
23 #:mpz-probably-prime-p
26 ;; the following functions are GMP >= 5.1 only
30 ;; number theoretic functions
34 ;; random number generation
45 ;; (un)installer functions
46 ; these insert/remove the runtime patch in SBCL's runtime
49 ; these also load/unload the shared library and setup/clear
50 ; hooks to handle core saves
60 (defvar *gmp-disabled
* nil
)
62 (defconstant +bignum-raw-area-offset
+
63 (- (* sb-vm
:bignum-digits-offset sb-vm
:n-word-bytes
)
64 sb-vm
:other-pointer-lowtag
))
66 (declaim (inline bignum-data-sap
))
67 (defun bignum-data-sap (x)
68 (declare (type bignum x
))
69 (sb-sys:sap
+ (sb-sys:int-sap
(sb-kernel:get-lisp-obj-address x
))
70 +bignum-raw-area-offset
+))
72 (defun try-load-shared-object (pathname)
74 (load-shared-object pathname
:dont-save t
)
80 (or (some #'try-load-shared-object
81 #-
(or win32 darwin
) '("libgmp.so" "libgmp.so.10" "libgmp.so.3")
82 #+darwin
'("libgmp.dylib" "libgmp.10.dylib" "libgmp.3.dylib")
83 #+win32
'("libgmp.dll" "libgmp-10.dll" "libgmp-3.dll"))
84 (warn "GMP not loaded.")))
86 (defvar *gmp-features
* nil
)
87 (defvar *gmp-version
* nil
)
89 ;; We load only the library right now to avoid undefined alien
94 ;;; types and initialization
96 (define-alien-type gmp-limb
97 #-
(and win32 x86-64
) unsigned-long
98 #+(and win32 x86-64
) unsigned-long-long
)
101 #-
(and win32 x86-64
) '(unsigned-byte #.sb-vm
:n-word-bits
)
102 #+(and win32 x86-64
) '(unsigned-byte 32))
105 #-
(and win32 x86-64
) '(signed-byte #.sb-vm
:n-word-bits
)
106 #+(and win32 x86-64
) '(signed-byte 32))
108 (define-alien-type nil
112 (mp_d (* gmp-limb
))))
114 ;; Section 3.6 "Memory Management" of the GMP manual states: "mpz_t
115 ;; and mpq_t variables never reduce their allocated space. Normally
116 ;; this is the best policy, since it avoids frequent
117 ;; reallocation. Applications that need to return memory to the heap
118 ;; at some particular point can use mpz_realloc2, or clear variables
119 ;; no longer needed."
121 ;; We can therefore allocate a bignum of sufficiant size and use the
122 ;; space for GMP computations without the need for memory transfer
123 ;; from C to Lisp space.
124 (declaim (inline z-to-bignum z-to-bignum-neg
))
126 (defun z-to-bignum (b count
)
127 "Convert GMP integer in the buffer of a pre-allocated bignum."
128 (declare (optimize (speed 3) (space 3) (safety 0))
130 (type bignum-length count
))
133 (the unsigned-byte
(%normalize-bignum b count
))))
135 (defun z-to-bignum-neg (b count
)
136 "Convert to twos complement int the buffer of a pre-allocated
138 (declare (optimize (speed 3) (space 3) (safety 0))
140 (type bignum-length count
))
141 (negate-bignum-in-place b
)
142 (the (integer * 0) (%normalize-bignum b count
)))
144 ;;; conversion functions that also copy from GMP to SBCL bignum space
145 (declaim (inline gmp-z-to-bignum
))
147 (defun gmp-z-to-bignum (z b count
)
148 "Convert and copy a positive GMP integer into the buffer of a
149 pre-allocated bignum. The allocated bignum-length must be (1+ COUNT)."
150 (declare (optimize (speed 3) (space 3) (safety 0))
151 (type (alien (* gmp-limb
)) z
)
153 (type bignum-length count
))
154 (dotimes (i count
(%normalize-bignum b
(1+ count
)))
155 (%bignum-set b i
(deref z i
))))
157 (declaim (inline blength bassert
)
158 (ftype (function (integer) (values bignum-length
&optional
)) blength
)
159 (ftype (function (integer) (values bignum
&optional
)) bassert
))
162 (declare (optimize (speed 3) (space 3) (safety 0)))
165 (t (%bignum-length a
))))
168 (declare (optimize (speed 3) (space 3) (safety 0)))
170 (fixnum (make-small-bignum a
))
174 (define-alien-type nil
176 (mp_num (struct gmpint
))
177 (mp_den (struct gmpint
))))
179 ;;; memory initialization functions to support non-alloced results
180 ;;; since an upper bound cannot always correctly predetermined
181 ;;; (e.g. the memory required for the fib function exceed the number
182 ;;; of limbs that are be determined through the infamous Phi-relation
183 ;;; resulting in a memory access error.
185 ;; use these for non-prealloced bignum values, but only when
186 ;; ultimately necessary since copying back into bignum space a the end
187 ;; of the operation is about three times slower than the shared buffer
189 (declaim (inline __gmpz_init __gmpz_clear
))
190 (define-alien-routine __gmpz_init void
191 (x (* (struct gmpint
))))
193 (define-alien-routine __gmpz_clear void
194 (x (* (struct gmpint
))))
197 ;;; integer interface functions
198 (defmacro define-twoarg-mpz-funs
(funs)
199 (loop for i in funs collect
`(define-alien-routine ,i void
200 (r (* (struct gmpint
)))
201 (a (* (struct gmpint
))))
203 finally
(return `(progn
204 (declaim (inline ,@funs
))
207 (defmacro define-threearg-mpz-funs
(funs)
208 (loop for i in funs collect
`(define-alien-routine ,i void
209 (r (* (struct gmpint
)))
210 (a (* (struct gmpint
)))
211 (b (* (struct gmpint
))))
213 finally
(return `(progn
214 (declaim (inline ,@funs
))
217 (defmacro define-fourarg-mpz-funs
(funs)
218 (loop for i in funs collect
`(define-alien-routine ,i void
219 (r (* (struct gmpint
)))
220 (a (* (struct gmpint
)))
221 (b (* (struct gmpint
)))
222 (c (* (struct gmpint
))))
224 finally
(return `(progn
225 (declaim (inline ,@funs
))
229 (define-twoarg-mpz-funs (__gmpz_sqrt
232 (define-threearg-mpz-funs (__gmpz_add
239 (define-fourarg-mpz-funs (__gmpz_cdiv_qr
244 (declaim (inline __gmpz_mul_2exp
247 __gmpz_probab_prime_p
256 (define-alien-routine __gmpz_mul_2exp void
257 (r (* (struct gmpint
)))
258 (b (* (struct gmpint
)))
261 (define-alien-routine __gmpz_fdiv_q_2exp void
262 (r (* (struct gmpint
)))
263 (b (* (struct gmpint
)))
266 (define-alien-routine __gmpz_pow_ui void
267 (r (* (struct gmpint
)))
268 (b (* (struct gmpint
)))
271 (define-alien-routine __gmpz_probab_prime_p int
272 (n (* (struct gmpint
)))
275 (define-alien-routine __gmpz_fac_ui void
276 (r (* (struct gmpint
)))
279 (define-alien-routine __gmpz_2fac_ui void
280 (r (* (struct gmpint
)))
283 (define-alien-routine __gmpz_mfac_uiui void
284 (r (* (struct gmpint
)))
288 (define-alien-routine __gmpz_primorial_ui void
289 (r (* (struct gmpint
)))
292 (define-alien-routine __gmpz_remove unsigned-long
293 (r (* (struct gmpint
)))
294 (x (* (struct gmpint
)))
295 (f (* (struct gmpint
))))
297 (define-alien-routine __gmpz_bin_ui void
298 (r (* (struct gmpint
)))
299 (n (* (struct gmpint
)))
302 (define-alien-routine __gmpz_fib2_ui void
303 (r (* (struct gmpint
)))
304 (a (* (struct gmpint
)))
309 (defmacro define-threearg-mpq-funs
(funs)
310 (loop for i in funs collect
`(define-alien-routine ,i void
311 (r (* (struct gmprat
)))
312 (a (* (struct gmprat
)))
313 (b (* (struct gmprat
))))
315 finally
(return `(progn
316 (declaim (inline ,@funs
))
319 (define-threearg-mpq-funs (__gmpq_add
327 ;;; utility macros for GMP mpz variable and result declaration and
328 ;;; incarnation of associated SBCL bignums
330 (defmacro with-mpz-results
(pairs &body body
)
331 (loop for
(gres size
) in pairs
332 for res
= (gensym "RESULT")
333 collect
`(when (> ,size sb-bignum
::maximum-bignum-length
)
334 (error "Size of result exceeds maxim bignum length")) into checks
335 collect
`(,gres
(struct gmpint
)) into declares
336 collect
`(,res
(%allocate-bignum
,size
))
338 collect
`(setf (slot ,gres
'mp_alloc
) (%bignum-length
,res
)
339 (slot ,gres
'mp_size
) 0
340 (slot ,gres
'mp_d
) (bignum-data-sap ,res
))
342 collect
`(if (minusp (slot ,gres
'mp_size
)) ; check for negative result
343 (z-to-bignum-neg ,res
,size
)
344 (z-to-bignum ,res
,size
))
346 collect res into results
351 (sb-sys:with-pinned-objects
,results
352 (with-alien ,declares
355 (values ,@normlimbs
))))))))
357 (defmacro with-mpz-vars
(pairs &body body
)
358 (loop for
(a ga
) in pairs
359 for length
= (gensym "LENGTH")
360 for plusp
= (gensym "PLUSP")
361 for barg
= (gensym "BARG")
362 for arg
= (gensym "ARG")
363 collect
`(,ga
(struct gmpint
)) into declares
364 collect
`(,barg
(bassert ,a
)) into gmpinits
365 collect
`(,plusp
(%bignum-0-or-plusp
,barg
(%bignum-length
,barg
))) into gmpinits
366 collect
`(,arg
(if ,plusp
,barg
(negate-bignum ,barg nil
))) into gmpinits
367 collect
`(,length
(%bignum-length
,arg
)) into gmpinits
368 collect arg into vars
369 collect
`(setf (slot ,ga
'mp_alloc
) ,length
371 (progn ;; handle twos complements/ulong limbs mismatch
372 (when (zerop (%bignum-ref
,arg
(1- ,length
)))
374 (if ,plusp
,length
(- ,length
)))
375 (slot ,ga
'mp_d
) (bignum-data-sap ,arg
))
378 `(with-alien ,declares
380 (sb-sys:with-pinned-objects
,vars
384 (defmacro with-gmp-mpz-results
(resultvars &body body
)
385 (loop for gres in resultvars
386 for res
= (gensym "RESULT")
387 for size
= (gensym "SIZE")
388 collect size into sizes
389 collect
`(,gres
(struct gmpint
)) into declares
390 collect
`(__gmpz_init (addr ,gres
)) into inits
391 collect
`(,size
(abs (slot ,gres
'mp_size
)))
393 collect
`(when (> ,size
(1- sb-bignum
::maximum-bignum-length
))
394 (error "Size of result exceeds maxim bignum length")) into checks
395 collect
`(,res
(%allocate-bignum
(1+ ,size
)))
397 collect
`(setf ,res
(if (minusp (slot ,gres
'mp_size
)) ; check for negative result
398 (- (gmp-z-to-bignum (slot ,gres
'mp_d
) ,res
,size
))
399 (gmp-z-to-bignum (slot ,gres
'mp_d
) ,res
,size
)))
401 collect
`(__gmpz_clear (addr ,gres
)) into clears
402 collect res into results
404 `(with-alien ,declares
408 (declare (type bignum-length
,@sizes
))
411 ;; copy GMP limbs into result bignum
412 (sb-sys:with-pinned-objects
,results
415 (values ,@results
)))))))
417 ;;; function definition and foreign function relationships
418 (defmacro defgmpfun
(name args
&body body
)
420 (declaim (sb-ext:maybe-inline
,name
))
422 (declare (optimize (speed 3) (space 3) (safety 0))
423 (type integer
,@args
))
427 ;; SBCL/GMP functions
428 (defgmpfun mpz-add
(a b
)
429 (with-mpz-results ((result (1+ (max (blength a
)
431 (with-mpz-vars ((a ga
) (b gb
))
432 (__gmpz_add (addr result
) (addr ga
) (addr gb
)))))
434 (defgmpfun mpz-sub
(a b
)
435 (with-mpz-results ((result (1+ (max (blength a
)
437 (with-mpz-vars ((a ga
) (b gb
))
438 (__gmpz_sub (addr result
) (addr ga
) (addr gb
)))))
440 (defgmpfun mpz-mul
(a b
)
441 (with-mpz-results ((result (+ (blength a
)
443 (with-mpz-vars ((a ga
) (b gb
))
444 (__gmpz_mul (addr result
) (addr ga
) (addr gb
)))))
446 (defgmpfun mpz-mul-2exp
(a b
)
448 (with-mpz-results ((result (+ (1+ (blength a
))
449 (floor b sb-vm
:n-word-bits
))))
450 (with-mpz-vars ((a ga
))
451 (__gmpz_mul_2exp (addr result
) (addr ga
) b
))))
453 (defgmpfun mpz-mod
(a b
)
454 (with-mpz-results ((result (1+ (max (blength a
)
456 (with-mpz-vars ((a ga
) (b gb
))
457 (__gmpz_mod (addr result
) (addr ga
) (addr gb
))
458 (when (and (minusp (slot gb
'mp_size
))
459 (/= 0 (slot result
'mp_size
)))
460 (__gmpz_add (addr result
) (addr result
) (addr gb
))))))
462 (defgmpfun mpz-cdiv
(n d
)
463 (let ((size (1+ (max (blength n
)
465 (with-mpz-results ((quot size
)
467 (with-mpz-vars ((n gn
) (d gd
))
468 (__gmpz_cdiv_qr (addr quot
) (addr rem
) (addr gn
) (addr gd
))))))
470 (defgmpfun mpz-fdiv
(n d
)
471 (let ((size (1+ (max (blength n
)
473 (with-mpz-results ((quot size
)
475 (with-mpz-vars ((n gn
) (d gd
))
476 (__gmpz_fdiv_qr (addr quot
) (addr rem
) (addr gn
) (addr gd
))))))
478 (defgmpfun mpz-fdiv-2exp
(a b
)
480 (with-mpz-results ((result (1+ (- (blength a
)
481 (floor b sb-vm
:n-word-bits
)))))
482 (with-mpz-vars ((a ga
))
483 (__gmpz_fdiv_q_2exp (addr result
) (addr ga
) b
))))
485 (defgmpfun mpz-tdiv
(n d
)
486 (let ((size (max (blength n
)
488 (with-mpz-results ((quot size
)
490 (with-mpz-vars ((n gn
) (d gd
))
491 (__gmpz_tdiv_qr (addr quot
) (addr rem
) (addr gn
) (addr gd
))))))
493 (defgmpfun mpz-pow
(base exp
)
494 (check-type exp
(integer 0 #.most-positive-fixnum
))
495 (with-gmp-mpz-results (rop)
496 (with-mpz-vars ((base gbase
))
497 (__gmpz_pow_ui (addr rop
) (addr gbase
) exp
))))
499 (defgmpfun mpz-powm
(base exp mod
)
500 (with-mpz-results ((rop (1+ (blength mod
))))
501 (with-mpz-vars ((base gbase
) (exp gexp
) (mod gmod
))
502 (__gmpz_powm (addr rop
) (addr gbase
) (addr gexp
) (addr gmod
)))))
504 (defgmpfun mpz-gcd
(a b
)
505 (with-mpz-results ((result (min (blength a
)
507 (with-mpz-vars ((a ga
) (b gb
))
508 (__gmpz_gcd (addr result
) (addr ga
) (addr gb
)))))
510 (defgmpfun mpz-lcm
(a b
)
511 (with-mpz-results ((result (+ (blength a
)
513 (with-mpz-vars ((a ga
) (b gb
))
514 (__gmpz_lcm (addr result
) (addr ga
) (addr gb
)))))
516 (defgmpfun mpz-sqrt
(a)
517 (with-mpz-results ((result (1+ (ceiling (blength a
) 2))))
518 (with-mpz-vars ((a ga
))
519 (__gmpz_sqrt (addr result
) (addr ga
)))))
522 ;;; Functions that use GMP-side allocated integers and copy the result
523 ;;; into a SBCL bignum at the end of the computation when the required
524 ;;; bignum length is known.
525 (defun mpz-probably-prime-p (n &optional
(reps 25))
526 (declare (optimize (speed 3) (space 3) (safety 0))
527 (type integer n reps
))
528 (check-type reps fixnum
)
529 (with-mpz-vars ((n gn
))
530 (__gmpz_probab_prime_p (addr gn
) reps
)))
532 (defgmpfun mpz-nextprime
(a)
533 (with-gmp-mpz-results (prime)
534 (with-mpz-vars ((a ga
))
535 (__gmpz_nextprime (addr prime
) (addr ga
)))))
537 (defgmpfun mpz-fac
(n)
539 (with-gmp-mpz-results (fac)
540 (__gmpz_fac_ui (addr fac
) n
)))
542 (defgmpfun %mpz-2fac
(n)
544 (with-gmp-mpz-results (fac)
545 (__gmpz_2fac_ui (addr fac
) n
)))
547 (defgmpfun %mpz-mfac
(n m
)
550 (with-gmp-mpz-results (fac)
551 (__gmpz_mfac_uiui (addr fac
) n m
)))
553 (defgmpfun %mpz-primorial
(n)
555 (with-gmp-mpz-results (r)
556 (__gmpz_primorial_ui (addr r
) n
)))
558 (defgmpfun mpz-remove-5.1
(n f
)
559 (check-type f unsigned-byte
560 "handled by GMP prior to version 5.1")
562 (res (with-gmp-mpz-results (r)
563 (with-mpz-vars ((n gn
)
565 (setf c
(__gmpz_remove (addr r
) (addr gn
) (addr gf
)))))))
568 (defgmpfun mpz-remove
(n f
)
570 (res (with-gmp-mpz-results (r)
571 (with-mpz-vars ((n gn
)
573 (setf c
(__gmpz_remove (addr r
) (addr gn
) (addr gf
)))))))
576 (defun setup-5.1-stubs
()
577 (macrolet ((stubify (name implementation
&rest arguments
)
578 `(setf (fdefinition ',name
)
579 (if (member :sb-gmp-5.1
*gmp-features
*)
580 (fdefinition ',implementation
)
582 (declare (ignore ,@arguments
))
583 (error "~S is only available in GMP >= 5.1"
585 (stubify mpz-2fac %mpz-2fac n
)
586 (stubify mpz-mfac %mpz-mfac n m
)
587 (stubify mpz-primorial %mpz-primorial n
))
588 (unless (member :sb-gmp-5.1
*gmp-features
*)
589 (setf (fdefinition 'mpz-remove
) #'mpz-remove-5.1
)))
591 (defgmpfun mpz-bin
(n k
)
593 (with-gmp-mpz-results (r)
594 (with-mpz-vars ((n gn
))
595 (__gmpz_bin_ui (addr r
) (addr gn
) k
))))
597 (defgmpfun mpz-fib2
(n)
598 ;; (let ((size (1+ (ceiling (* n (log 1.618034 2)) 64)))))
599 ;; fibonacci number magnitude in bits is assymptotic to n(log_2 phi)
600 ;; This is correct for the result but appears not to be enough for GMP
601 ;; during computation (memory access error), so use GMP-side allocation.
603 (with-gmp-mpz-results (fibn fibn-1
)
604 (__gmpz_fib2_ui (addr fibn
) (addr fibn-1
) n
)))
607 ;;;; Random bignum (mpz) generation
609 ;; we do not actually use the gestalt of the struct but need its size
610 ;; for allocation purposes
611 (define-alien-type nil
613 (mp_seed (struct gmpint
))
617 (declaim (inline __gmp_randinit_mt
618 __gmp_randinit_lc_2exp
624 (define-alien-routine __gmp_randinit_mt void
625 (s (* (struct gmprandstate
))))
627 (define-alien-routine __gmp_randinit_lc_2exp void
628 (s (* (struct gmprandstate
)))
629 (a (* (struct gmpint
)))
633 (define-alien-routine __gmp_randseed void
634 (s (* (struct gmprandstate
)))
635 (sd (* (struct gmpint
))))
637 (define-alien-routine __gmp_randseed_ui void
638 (s (* (struct gmprandstate
)))
641 (define-alien-routine __gmpz_urandomb void
642 (r (* (struct gmpint
)))
643 (s (* (struct gmprandstate
)))
644 (bcnt unsigned-long
))
646 (define-alien-routine __gmpz_urandomm void
647 (r (* (struct gmpint
)))
648 (s (* (struct gmprandstate
)))
649 (n (* (struct gmpint
))))
651 (defstruct (gmp-rstate (:constructor %make-gmp-rstate
))
652 (ref (make-alien (struct gmprandstate
))
653 :type
(alien (* (struct gmprandstate
))) :read-only t
))
655 (declaim (sb-ext:maybe-inline make-gmp-rstate
661 (defun make-gmp-rstate ()
662 "Instantiate a state for the GMP Mersenne-Twister random number generator."
663 (declare (optimize (speed 3) (space 3) (safety 0)))
664 (let* ((state (%make-gmp-rstate
))
665 (ref (gmp-rstate-ref state
)))
666 (__gmp_randinit_mt ref
)
667 (sb-ext:finalize state
(lambda () (free-alien ref
)))
670 (defun make-gmp-rstate-lc (a c m2exp
)
671 "Instantiate a state for the GMP linear congruential random number generator."
672 (declare (optimize (speed 3) (space 3)))
674 (check-type m2exp ui
)
675 (let* ((state (%make-gmp-rstate
))
676 (ref (gmp-rstate-ref state
)))
677 (with-mpz-vars ((a ga
))
678 (__gmp_randinit_lc_2exp ref
(addr ga
) c m2exp
))
679 (sb-ext:finalize state
(lambda () (free-alien ref
)))
682 (defun rand-seed (state seed
)
683 "Initialize a random STATE with SEED."
684 (declare (optimize (speed 3) (space 3) (safety 0)))
685 (check-type state gmp-rstate
)
686 (let ((ref (gmp-rstate-ref state
)))
689 (__gmp_randseed_ui ref seed
))
690 ((typep seed
'(integer 0 *))
691 (with-mpz-vars ((seed gseed
))
692 (__gmp_randseed ref
(addr gseed
))))
694 (error "SEED must be a positive integer")))))
696 (defun random-bitcount (state bitcount
)
697 "Return a random integer in the range 0..(2^bitcount - 1)."
698 (declare (optimize (speed 3) (space 3) (safety 0)))
699 (check-type state gmp-rstate
)
700 (check-type bitcount ui
)
701 (let ((ref (gmp-rstate-ref state
)))
702 (with-mpz-results ((result (+ (ceiling bitcount sb-vm
:n-word-bits
) 2)))
703 (__gmpz_urandomb (addr result
) ref bitcount
))))
705 (defun random-int (state boundary
)
706 "Return a random integer in the range 0..(boundary - 1)."
707 (declare (optimize (speed 3) (space 3)))
708 (check-type state gmp-rstate
)
709 (let ((b (bassert boundary
))
710 (ref (gmp-rstate-ref state
)))
711 (with-mpz-results ((result (1+ (%bignum-length b
))))
712 (with-mpz-vars ((b gb
))
713 (__gmpz_urandomm (addr result
) ref
(addr gb
))))))
716 ;;; Rational functions
717 (declaim (inline %lsize
))
718 (defun %lsize
(minusp n
)
719 (declare (optimize (speed 3) (space 3) (safety 0)))
720 "n must be a (potentially denormalized) bignum"
721 (let ((length (%bignum-length n
)))
722 (when (zerop (%bignum-ref n
(1- length
)))
724 (if minusp
(- length
) length
)))
726 (defmacro with-mpq-var
(pair &body body
)
727 (destructuring-bind (var mpqvar
) pair
728 `(let* ((an (bassert (numerator ,var
)))
729 (ad (bassert (denominator ,var
)))
730 (asign (not (%bignum-0-or-plusp an
(%bignum-length an
)))))
732 (setf an
(negate-bignum an nil
)))
733 (let* ((anlen (%lsize asign an
))
734 (adlen (%lsize NIL ad
)))
735 (with-alien ((,mpqvar
(struct gmprat
)))
736 (sb-sys:with-pinned-objects
(an ad
)
737 (setf (slot (slot ,mpqvar
'mp_num
) 'mp_size
) anlen
738 (slot (slot ,mpqvar
'mp_num
) 'mp_alloc
) (abs anlen
)
739 (slot (slot ,mpqvar
'mp_num
) 'mp_d
)
740 (bignum-data-sap an
))
741 (setf (slot (slot ,mpqvar
'mp_den
) 'mp_size
) adlen
742 (slot (slot ,mpqvar
'mp_den
) 'mp_alloc
) (abs adlen
)
743 (slot (slot ,mpqvar
'mp_den
) 'mp_d
)
744 (bignum-data-sap ad
))
747 (defmacro defmpqfun
(name gmpfun
)
749 (declaim (sb-ext:maybe-inline
,name
))
751 (declare (optimize (speed 3) (space 3) (safety 0)))
752 (let ((size (+ (max (blength (numerator a
))
753 (blength (denominator a
)))
754 (max (blength (numerator b
))
755 (blength (denominator b
)))
757 (with-alien ((r (struct gmprat
)))
758 (let ((num (%allocate-bignum size
))
759 (den (%allocate-bignum size
)))
760 (sb-sys:with-pinned-objects
(num den
)
761 (setf (slot (slot r
'mp_num
) 'mp_size
) 0
762 (slot (slot r
'mp_num
) 'mp_alloc
) size
763 (slot (slot r
'mp_num
) 'mp_d
) (bignum-data-sap num
))
764 (setf (slot (slot r
'mp_den
) 'mp_size
) 0
765 (slot (slot r
'mp_den
) 'mp_alloc
) size
766 (slot (slot r
'mp_den
) 'mp_d
) (bignum-data-sap den
))
767 (let* ((an (bassert (numerator a
)))
768 (ad (bassert (denominator a
)))
769 (asign (not (%bignum-0-or-plusp an
(%bignum-length an
))))
770 (bn (bassert (numerator b
)))
771 (bd (bassert (denominator b
)))
772 (bsign (not (%bignum-0-or-plusp bn
(%bignum-length bn
)))))
774 (setf an
(negate-bignum an nil
)))
776 (setf bn
(negate-bignum bn nil
)))
777 (let* ((anlen (%lsize asign an
))
778 (adlen (%lsize NIL ad
))
779 (bnlen (%lsize bsign bn
))
780 (bdlen (%lsize NIL bd
)))
781 (with-alien ((arga (struct gmprat
))
782 (argb (struct gmprat
)))
783 (sb-sys:with-pinned-objects
(an ad bn bd
)
784 (setf (slot (slot arga
'mp_num
) 'mp_size
) anlen
785 (slot (slot arga
'mp_num
) 'mp_alloc
) (abs anlen
)
786 (slot (slot arga
'mp_num
) 'mp_d
)
787 (bignum-data-sap an
))
788 (setf (slot (slot arga
'mp_den
) 'mp_size
) adlen
789 (slot (slot arga
'mp_den
) 'mp_alloc
) (abs adlen
)
790 (slot (slot arga
'mp_den
) 'mp_d
)
791 (bignum-data-sap ad
))
792 (setf (slot (slot argb
'mp_num
) 'mp_size
) bnlen
793 (slot (slot argb
'mp_num
) 'mp_alloc
) (abs bnlen
)
794 (slot (slot argb
'mp_num
) 'mp_d
)
795 (bignum-data-sap bn
))
796 (setf (slot (slot argb
'mp_den
) 'mp_size
) bdlen
797 (slot (slot argb
'mp_den
) 'mp_alloc
) (abs bdlen
)
798 (slot (slot argb
'mp_den
) 'mp_d
)
799 (bignum-data-sap bd
))
800 (,gmpfun
(addr r
) (addr arga
) (addr argb
)))))
801 (locally (declare (optimize (speed 1)))
802 (sb-kernel::build-ratio
(if (minusp (slot (slot r
'mp_num
) 'mp_size
))
803 (z-to-bignum-neg num size
)
804 (z-to-bignum num size
))
805 (z-to-bignum den size
)))))))))))
807 (defmpqfun mpq-add __gmpq_add
)
808 (defmpqfun mpq-sub __gmpq_sub
)
809 (defmpqfun mpq-mul __gmpq_mul
)
810 (defmpqfun mpq-div __gmpq_div
)
813 ;;;; SBCL interface and integration installation
814 (macrolet ((def (name original
)
815 (let ((special (intern (format nil
"*~A-FUNCTION*" name
))))
817 (declaim (type function
,special
)
819 (defvar ,special
(symbol-function ',original
))
820 (defun ,name
(&rest args
)
821 (apply (load-time-value ,special t
) args
))))))
822 (def orig-mul multiply-bignums
)
823 (def orig-truncate bignum-truncate
)
824 (def orig-gcd bignum-gcd
)
825 (def orig-lcm sb-kernel
:two-arg-lcm
)
826 (def orig-isqrt isqrt
)
827 (def orig-two-arg-
+ sb-kernel
:two-arg-
+)
828 (def orig-two-arg-- sb-kernel
:two-arg--
)
829 (def orig-two-arg-
* sb-kernel
:two-arg-
*)
830 (def orig-two-arg-
/ sb-kernel
:two-arg-
/)
831 (def orig-intexp sb-kernel
::intexp
))
835 (declare (optimize (speed 3) (space 3))
836 (type bignum-type a b
)
838 (if (or (< (min (%bignum-length a
)
845 (defun gmp-truncate (a b
)
846 (declare (optimize (speed 3) (space 3))
847 (type bignum-type a b
)
849 (if (or (< (min (%bignum-length a
)
857 (declare (optimize (speed 3) (space 3))
860 (if (or (and (typep a
'fixnum
)
867 (declare (optimize (speed 3) (space 3))
868 (type unsigned-byte n
)
870 (if (or (typep n
'fixnum
)
876 (defun gmp-two-arg-+ (x y
)
877 (declare (optimize (speed 3) (space 3))
879 (if (and (or (typep x
'ratio
)
883 (not *gmp-disabled
*))
885 (orig-two-arg-+ x y
)))
887 (defun gmp-two-arg-- (x y
)
888 (declare (optimize (speed 3) (space 3))
890 (if (and (or (typep x
'ratio
)
894 (not *gmp-disabled
*))
896 (orig-two-arg-- x y
)))
898 (defun gmp-two-arg-* (x y
)
899 (declare (optimize (speed 3) (space 3))
901 (if (and (or (typep x
'ratio
)
905 (not *gmp-disabled
*))
907 (orig-two-arg-* x y
)))
909 (defun gmp-two-arg-/ (x y
)
910 (declare (optimize (speed 3) (space 3))
912 (if (and (rationalp x
)
915 (not *gmp-disabled
*))
917 (orig-two-arg-/ x y
)))
919 (defun gmp-intexp (base power
)
920 (declare (inline mpz-mul-2exp mpz-pow
))
921 (check-type power
(integer #.
(1+ most-negative-fixnum
) #.most-positive-fixnum
))
923 ((or (and (integerp base
)
925 (< (blength base
) 4))
927 (orig-intexp base power
))
929 (cond ((minusp power
)
930 (/ (the integer
(gmp-intexp base
(- power
)))))
932 (mpz-mul-2exp 1 power
))
934 (sb-kernel::%make-ratio
(gmp-intexp (numerator base
) power
)
935 (gmp-intexp (denominator base
) power
)))
937 (mpz-pow base power
))))))
940 (defun install-gmp-funs ()
941 (sb-ext:without-package-locks
942 (macrolet ((def (destination source
)
943 `(setf (fdefinition ',destination
)
944 (fdefinition ',source
))))
945 (def sb-bignum
:multiply-bignums gmp-mul
)
946 (def sb-bignum
:bignum-truncate gmp-truncate
)
947 (def sb-bignum
:bignum-gcd mpz-gcd
)
948 (def sb-kernel
:two-arg-lcm gmp-lcm
)
949 (def sb-kernel
:two-arg-
+ gmp-two-arg-
+)
950 (def sb-kernel
:two-arg-- gmp-two-arg--
)
951 (def sb-kernel
:two-arg-
* gmp-two-arg-
*)
952 (def sb-kernel
:two-arg-
/ gmp-two-arg-
/)
953 (def isqrt gmp-isqrt
)
954 (def sb-kernel
::intexp gmp-intexp
)))
957 (defun uninstall-gmp-funs ()
958 (sb-ext:without-package-locks
959 (macrolet ((def (destination source
)
960 `(setf (fdefinition ',destination
)
961 ,(intern (format nil
"*~A-FUNCTION*" source
)))))
962 (def multiply-bignums orig-mul
)
963 (def bignum-truncate orig-truncate
)
964 (def bignum-gcd orig-gcd
)
965 (def sb-kernel
:two-arg-lcm orig-lcm
)
966 (def sb-kernel
:two-arg-
+ orig-two-arg-
+)
967 (def sb-kernel
:two-arg-- orig-two-arg--
)
968 (def sb-kernel
:two-arg-
* orig-two-arg-
*)
969 (def sb-kernel
:two-arg-
/ orig-two-arg-
/)
970 (def isqrt orig-isqrt
)
971 (def sb-kernel
::intexp orig-intexp
)))
974 (defun load-gmp (&key
(persistently t
))
975 (setf *gmp-features
* nil
977 *features
* (set-difference *features
* '(:sb-gmp
:sb-gmp-5.0
:sb-gmp-5.1
)))
979 (pushnew 'load-gmp sb-ext
:*init-hooks
*)
980 (pushnew 'uninstall-gmp-funs sb-ext
:*save-hooks
*))
981 (let ((success (%load-gmp
)))
983 (setf *gmp-version
* (extern-alien "__gmp_version" c-string
)))
984 (cond ((null *gmp-version
*))
985 ((string<= *gmp-version
* "5.")
986 (warn "SB-GMP requires at least GMP version 5.0")
989 (pushnew :sb-gmp
*gmp-features
*)
990 (pushnew :sb-gmp-5.0
*gmp-features
*)
991 (when (string>= *gmp-version
* "5.1")
992 (pushnew :sb-gmp-5.1
*gmp-features
*))
993 (setf *features
* (union *features
* *gmp-features
*))))
996 (uninstall-gmp-funs))
1000 (defun unload-gmp ()
1001 (setf sb-ext
:*init-hooks
* (remove 'load-gmp sb-ext
:*init-hooks
*))
1002 (uninstall-gmp-funs)
1003 (setf sb-ext
:*save-hooks
* (remove 'uninstall-gmp-funs sb-ext
:*save-hooks
*))