mini-gmp and mpq
Bradley Lucier
lucier at math.purdue.edu
Sun Feb 25 20:43:35 UTC 2018
> On Feb 18, 2018, at 3:51 AM, Niels Möller <nisse at lysator.liu.se> wrote:
>
> But I might be mistaken, e.g, mpq/aors.c look a lot more more complex
> than in your sketch of mini-mpq.c.
I’m not a gmp developer, but it seems from examining the gmp/mpq code that nearly all the complexity in rational arithmetic is to limit the size or even the need to do gcd. This is covered in the first section of Chapter 4.5 of Knuth's Seminumerical Algorithms.
For example, in the Gambit Scheme runtime you have
(define-prim (##ratnum.* x y)
(##declare (mostly-fixnum)
(standard-bindings))
(let ((p (macro-ratnum-numerator x))
(q (macro-ratnum-denominator x))
(r (macro-ratnum-numerator y))
(s (macro-ratnum-denominator y)))
(if (eq? x y)
(macro-ratnum-make (square p) (square q)) ;; already in lowest form
(let* ((gcd-ps (##gcd p s))
(gcd-rq (##gcd r q))
(num (* (quotient p gcd-ps) (quotient r gcd-rq)))
(den (* (quotient q gcd-rq) (quotient s gcd-ps))))
(if (eqv? den 1)
num
(macro-ratnum-make num den))))))
and
(define-prim (##ratnum.+ x y)
(##declare (mostly-fixnum)
(standard-bindings))
(let ((p (macro-ratnum-numerator x))
(q (macro-ratnum-denominator x))
(r (macro-ratnum-numerator y))
(s (macro-ratnum-denominator y)))
(let ((d1 (##gcd q s)))
(if (eqv? d1 1)
(macro-ratnum-make (+ (* p s)
(* r q))
(* q s))
(let* ((s-prime (quotient s d1))
(t (+ (* p s-prime)
(* r (quotient q d1))))
(d2 (##gcd d1 t))
(num (quotient t d2))
(den (* (quotient q d2)
s-prime)))
(if (eqv? den 1)
num
(macro-ratnum-make num den)))))))
Even in aors.c it seems that what complexity is there would be worth the trouble even in a mini-gmp.
Brad
More information about the gmp-devel
mailing list