[racket] Math library kudos

FYI, log gamma is another fast way to calculate the number of combinations
if you want to deal with really big numbers.

On Tue, Feb 19, 2013 at 7:28 PM, Joe Gilray <jgilray at gmail.com> wrote:

> Racketeers,
>
> Thanks for putting together the fantastic math library.  It will be a
> wonderful resource.  Here are some quick impressions (after playing mostly
> with math/number-theory)
>
> 1) The functions passed all my tests and were very fast.  If you need even
> more speed you can keep a list of primes around and write functions to use
> that, but that should be rarely necessary
>
> 2) I have a couple of functions to donate if you want them:
>
> 2a) Probablistic primality test:
>
> ; function that performs a Miller-Rabin probabalistic primality test k
> times, returns #t if n is probably prime
> ; algorithm from http://rosettacode.org/wiki/Miller-Rabin_primality_test,
> code adapted from Lisp example
> ; (module+ test (check-equal? (is-mr-prime? 1000000000000037 8) #t))
> (define (is-mr-prime? n k)
>   ; function that returns two values r and e such that number = divisor^e
> * r, and r is not divisible by divisor
>   (define (factor-out number divisor)
>     (do ([e 0 (add1 e)] [r number (/ r divisor)])
>       ((not (zero? (remainder r divisor))) (values r e))))
>
>   ; function that performs fast modular exponentiation by repeated squaring
>   (define (expt-mod base exponent modulus)
>     (let expt-mod-iter ([b base] [e exponent] [p 1])
>       (cond
>         [(zero? e) p]
>         [(even? e) (expt-mod-iter (modulo (* b b) modulus) (/ e 2) p)]
>         [else (expt-mod-iter b (sub1 e) (modulo (* b p) modulus))])))
>
>   ; function to return a random, exact number in the passed range
> (inclusive)
>   (define (shifted-rand lower upper)
>     (+ lower (random (add1 (- (modulo upper 4294967088) (modulo lower
> 4294967088))))))
>
>   (cond
>     [(= n 1) #f]
>     [(< n 4) #t]
>     [(even? n) #f]
>     [else
>      (let-values ([(d s) (factor-out (- n 1) 2)]) ; represent n-1 as 2^s-d
>        (let lp ([a (shifted-rand 2 (- n 2))] [cnt k])
>          (if (zero? cnt) #t
>              (let ([x (expt-mod a d n)])
>                (if (or (= x 1) (= x (sub1 n))) (lp (shifted-rand 2 (- n
> 2)) (sub1 cnt))
>                    (let ctestlp ([r 1] [ctest (modulo (* x x) n)])
>                      (cond
>                        [(>= r s) #f]
>                        [(= ctest 1) #f]
>                        [(= ctest (sub1 n)) (lp (shifted-rand 2 (- n 2))
> (sub1 cnt))]
>                        [else (ctestlp (add1 r) (modulo (* ctest ctest)
> n))])))))))]))
>
> 2b) combinations calculator
>
> ; function that returns the number of combinations, not the combinations
> themselves
> ; faster than using n! / (r! (n-r)!)
> (define (combinations n r)
>   (cond
>     [(or (< n 0) (< r 0)) (error "combinations: illegal arguments, n and r
> must be >= 0")]
>     [(> r n) 0]
>     [else
>      (let lp ([mord n] [total 1] [mult #t])
>        (cond
>          [(or (= 0 mord) (= 1 mord)) total]
>          [(and mult (= mord (- n r))) (lp r total #f)]
>          [(and mult (= mord r)) (lp (- n r) total #f)]
>          [mult (lp (sub1 mord) (* total mord) #t)]
>          [else (lp (sub1 mord) (/ total mord) #f)]))]))
>
> Thanks again!
> -Joe
>
> ____________________
>   Racket Users list:
>   http://lists.racket-lang.org/users
>
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.racket-lang.org/users/archive/attachments/20130219/906d577d/attachment.html>