diff --git a/basis/math/primes/miller-rabin/miller-rabin-docs.factor b/basis/math/primes/miller-rabin/miller-rabin-docs.factor index 2455dafdd5..2d19d51e06 100644 --- a/basis/math/primes/miller-rabin/miller-rabin-docs.factor +++ b/basis/math/primes/miller-rabin/miller-rabin-docs.factor @@ -3,20 +3,6 @@ USING: help.markup help.syntax kernel sequences math ; IN: math.primes.miller-rabin -HELP: find-relative-prime -{ $values - { "n" integer } - { "p" integer } -} -{ $description "Returns a number that is relatively prime to " { $snippet "n" } "." } ; - -HELP: find-relative-prime* -{ $values - { "n" integer } { "guess" integer } - { "p" integer } -} -{ $description "Returns a number that is relatively prime to " { $snippet "n" } ", starting by trying " { $snippet "guess" } "." } ; - HELP: miller-rabin { $values { "n" integer } @@ -33,68 +19,10 @@ HELP: miller-rabin* } { $description "Performs " { $snippet "numtrials" } " trials of the Miller-Rabin probabilistic primality test algorithm and returns true if prime." } ; -HELP: next-prime -{ $values - { "n" integer } - { "p" integer } -} -{ $description "Tests consecutive numbers for primality with " { $link miller-rabin } " and returns the next prime." } ; - -HELP: next-safe-prime -{ $values - { "n" integer } - { "q" integer } -} -{ $description "Tests consecutive numbers and returns the next safe prime. A safe prime is desirable in cryptography applications such as Diffie-Hellman and SRP6." } ; - -HELP: random-bits* -{ $values - { "numbits" integer } - { "n" integer } -} -{ $description "Returns an integer exactly " { $snippet "numbits" } " in length, with the topmost bit set to one." } ; - -HELP: random-prime -{ $values - { "numbits" integer } - { "p" integer } -} -{ $description "Returns a prime number exactly " { $snippet "numbits" } " bits in length, with the topmost bit set to one." } ; - -HELP: random-safe-prime -{ $values - { "numbits" integer } - { "p" integer } -} -{ $description "Returns a safe prime number " { $snippet "numbits" } " bits in length, with the topmost bit set to one." } ; - -HELP: safe-prime? -{ $values - { "q" integer } - { "?" "a boolean" } -} -{ $description "Tests whether the number is a safe prime. A safe prime " { $snippet "p" } " must be prime, as must " { $snippet "(p - 1) / 2" } "." } ; - -HELP: unique-primes -{ $values - { "numbits" integer } { "n" integer } - { "seq" sequence } -} -{ $description "Generates a sequence of " { $snippet "n" } " unique prime numbers with exactly " { $snippet "numbits" } " bits." } ; - ARTICLE: "math.primes.miller-rabin" "Miller-Rabin probabilistic primality test" "The " { $vocab-link "math.primes.miller-rabin" } " vocabulary implements the Miller-Rabin probabilistic primality test and utility words that use it in order to generate random prime numbers." $nl "The Miller-Rabin probabilistic primality test:" { $subsection miller-rabin } -{ $subsection miller-rabin* } -"Generating relative prime numbers:" -{ $subsection find-relative-prime } -{ $subsection find-relative-prime* } -"Generating prime numbers:" -{ $subsection next-prime } -{ $subsection random-prime } -"Generating safe prime numbers:" -{ $subsection next-safe-prime } -{ $subsection random-safe-prime } ; +{ $subsection miller-rabin* } ; ABOUT: "math.primes.miller-rabin" diff --git a/basis/math/primes/miller-rabin/miller-rabin-tests.factor b/basis/math/primes/miller-rabin/miller-rabin-tests.factor index 9c635c8f38..aeae6cac1b 100644 --- a/basis/math/primes/miller-rabin/miller-rabin-tests.factor +++ b/basis/math/primes/miller-rabin/miller-rabin-tests.factor @@ -1,5 +1,6 @@ -USING: math.primes.miller-rabin tools.test kernel sequences -math.primes.miller-rabin.private math ; +USING: kernel math math.primes math.primes.miller-rabin +math.primes.miller-rabin.private math.primes.safe +math.primes.safe.private random sequences tools.test ; IN: math.primes.miller-rabin.tests [ f ] [ 473155932665450549999756893736999469773678960651272093993257221235459777950185377130233556540099119926369437865330559863 miller-rabin ] unit-test diff --git a/basis/math/primes/miller-rabin/miller-rabin.factor b/basis/math/primes/miller-rabin/miller-rabin.factor index 35ee97a897..b0dfc4ed35 100755 --- a/basis/math/primes/miller-rabin/miller-rabin.factor +++ b/basis/math/primes/miller-rabin/miller-rabin.factor @@ -1,18 +1,9 @@ ! Copyright (c) 2008-2009 Doug Coleman. ! See http://factorcode.org/license.txt for BSD license. -USING: combinators kernel locals math math.functions math.ranges -random sequences sets combinators.short-circuit math.bitwise -math math.order ; +USING: combinators combinators.short-circuit kernel locals math +math.functions math.ranges random sequences sets ; IN: math.primes.miller-rabin -: >odd ( n -- int ) 0 set-bit ; foldable - -: >even ( n -- int ) 0 clear-bit ; foldable - -: next-even ( m -- n ) >even 2 + ; - -: next-odd ( m -- n ) dup even? [ 1 + ] [ 2 + ] if ; - } cond ; : miller-rabin ( n -- ? ) 10 miller-rabin* ; - -ERROR: prime-range-error n ; - -: next-prime ( n -- p ) - dup 1 < [ prime-range-error ] when - dup 1 = [ - drop 2 - ] [ - next-odd dup miller-rabin [ next-prime ] unless - ] if ; - -: random-bits* ( numbits -- n ) - 1 - [ random-bits ] keep set-bit ; - -: random-prime ( numbits -- p ) - random-bits* next-prime ; - -ERROR: no-relative-prime n ; - - [ 2 + (find-relative-prime) ] [ nip ] if ; - -PRIVATE> - -: find-relative-prime* ( n guess -- p ) - #! find a prime relative to n with initial guess - >odd (find-relative-prime) ; - -: find-relative-prime ( n -- p ) - dup random find-relative-prime* ; - -ERROR: too-few-primes ; - -: unique-primes ( numbits n -- seq ) - #! generate two primes - swap - dup 5 < [ too-few-primes ] when - 2dup [ random-prime ] curry replicate - dup all-unique? [ 2nip ] [ drop unique-primes ] if ; - -! Safe primes are of the form p = 2q + 1, p,q are prime -! See http://en.wikipedia.org/wiki/Safe_prime - - - -: safe-prime? ( q -- ? ) - { - [ 1 - 2 / dup integer? [ miller-rabin ] [ drop f ] if ] - [ miller-rabin ] - } 1&& ; - -: next-safe-prime ( n -- q ) - next-safe-prime-candidate - dup safe-prime? [ next-safe-prime ] unless ; - -: random-safe-prime ( numbits -- p ) - random-bits* next-safe-prime ; diff --git a/basis/math/primes/primes-docs.factor b/basis/math/primes/primes-docs.factor index c7dbc950e8..fa991e800f 100644 --- a/basis/math/primes/primes-docs.factor +++ b/basis/math/primes/primes-docs.factor @@ -1,10 +1,10 @@ -USING: help.markup help.syntax ; +USING: help.markup help.syntax math sequences ; IN: math.primes { next-prime prime? } related-words HELP: next-prime -{ $values { "n" "an integer not smaller than 2" } { "p" "a prime number" } } +{ $values { "n" integer } { "p" "a prime number" } } { $description "Return the next prime number greater than " { $snippet "n" } "." } ; HELP: prime? @@ -20,3 +20,49 @@ HELP: primes-upto HELP: primes-between { $values { "low" "an integer" } { "high" "an integer" } { "seq" "a sequence" } } { $description "Return a sequence containing all the prime numbers between " { $snippet "low" } " and " { $snippet "high" } "." } ; + +HELP: find-relative-prime +{ $values + { "n" integer } + { "p" integer } +} +{ $description "Returns a number that is relatively prime to " { $snippet "n" } "." } ; + +HELP: find-relative-prime* +{ $values + { "n" integer } { "guess" integer } + { "p" integer } +} +{ $description "Returns a number that is relatively prime to " { $snippet "n" } ", starting by trying " { $snippet "guess" } "." } ; + +HELP: random-prime +{ $values + { "numbits" integer } + { "p" integer } +} +{ $description "Returns a prime number exactly " { $snippet "numbits" } " bits in length, with the topmost bit set to one." } ; + +HELP: unique-primes +{ $values + { "numbits" integer } { "n" integer } + { "seq" sequence } +} +{ $description "Generates a sequence of " { $snippet "n" } " unique prime numbers with exactly " { $snippet "numbits" } " bits." } ; + + +ARTICLE: "math.primes" "Prime numbers" +"The " { $vocab-link "math.primes" } " vocabulary implements words related to prime numbers." $nl +"Testing if a number is prime:" +{ $subsection prime? } +"Generating prime numbers:" +{ $subsection next-prime } +{ $subsection primes-upto } +{ $subsection primes-between } +{ $subsection random-prime } +"Generating relative prime numbers:" +{ $subsection find-relative-prime } +{ $subsection find-relative-prime* } +"Make a sequence of random prime numbers:" +{ $subsection unique-primes } ; + +ABOUT: "math.primes" diff --git a/basis/math/primes/primes-tests.factor b/basis/math/primes/primes-tests.factor index db738399ef..6580f0780e 100644 --- a/basis/math/primes/primes-tests.factor +++ b/basis/math/primes/primes-tests.factor @@ -1,4 +1,6 @@ -USING: arrays math.primes tools.test ; +USING: arrays math math.primes math.primes.miller-rabin +tools.test ; +IN: math.primes.tests { 1237 } [ 1234 next-prime ] unit-test { f t } [ 1234 prime? 1237 prime? ] unit-test @@ -7,3 +9,12 @@ USING: arrays math.primes tools.test ; { { 4999963 4999999 5000011 5000077 5000081 } } [ 4999962 5000082 primes-between >array ] unit-test + +[ 2 ] [ 1 next-prime ] unit-test +[ 3 ] [ 2 next-prime ] unit-test +[ 5 ] [ 3 next-prime ] unit-test +[ 101 ] [ 100 next-prime ] unit-test +[ t ] [ 2135623355842621559 miller-rabin ] unit-test +[ 100000000000031 ] [ 100000000000000 next-prime ] unit-test + +[ 49 ] [ 50 random-prime log2 ] unit-test diff --git a/basis/math/primes/primes.factor b/basis/math/primes/primes.factor index fa1cd5cb63..e3985fc600 100644 --- a/basis/math/primes/primes.factor +++ b/basis/math/primes/primes.factor @@ -1,8 +1,8 @@ ! Copyright (C) 2007-2009 Samuel Tardieu. ! See http://factorcode.org/license.txt for BSD license. -USING: combinators kernel math math.functions -math.primes.miller-rabin math.order math.primes.erato -math.ranges sequences ; +USING: combinators kernel math math.bitwise math.functions +math.order math.primes.erato math.primes.miller-rabin +math.ranges random sequences sets fry ; IN: math.primes } cond ; foldable : next-prime ( n -- p ) - next-odd [ dup really-prime? ] [ 2 + ] until ; foldable + dup 2 < [ + drop 2 + ] [ + next-odd [ dup really-prime? ] [ 2 + ] until + ] if ; foldable : primes-between ( low high -- seq ) [ dup 3 max dup even? [ 1 + ] when ] dip @@ -32,3 +36,34 @@ PRIVATE> : primes-upto ( n -- seq ) 2 swap primes-between ; : coprime? ( a b -- ? ) gcd nip 1 = ; foldable + +: random-prime ( numbits -- p ) + random-bits* next-prime ; + +: estimated-primes ( m -- n ) + dup log / ; foldable + +ERROR: no-relative-prime n ; + + [ 2 + (find-relative-prime) ] [ nip ] if ; + +PRIVATE> + +: find-relative-prime* ( n guess -- p ) + #! find a prime relative to n with initial guess + >odd (find-relative-prime) ; + +: find-relative-prime ( n -- p ) + dup random find-relative-prime* ; + +ERROR: too-few-primes n numbits ; + +: unique-primes ( n numbits -- seq ) + 2dup 2^ estimated-primes > [ too-few-primes ] when + 2dup '[ _ random-prime ] replicate + dup all-unique? [ 2nip ] [ drop unique-primes ] if ;