Primes are positive integers >1 that can only be divided (without remainder) by 1 and themselves. Examples from [RS]'s collection: 1234567891 987654321987654329 9223372036854775783 Two more primes, which may be algorithmically useful since they are close to pow(2,15) and pow(2,16): 32003 65521 Here's a very simple prime checker, which returns 1 for primes, else the first found divisor (and the quotient): proc primetest n { set max [expr wide(sqrt($n))] if {$n%2==0} {return [list 2 [expr $n/2]]} for {set i 3} {$i<=$max} {incr i 2} { if {$n%$i==0} {return [list $i [expr $n/$i]]} } return 1 } [LV] 2009-May-4 One mathematical nit. I seem to recall that the number '''1''' is not considered a prime because it has only one factor - itself, whereas the definition for a prime number is that the number has only 2 distinct factors - 1 and itself. the primetest function returns a '''1''' when called with a 1. I have no idea what primetest should return in this case - perhaps a [[list 1 1]] ? [AM] commented in the [Tcl chatroom]: "The trick: divide by the first ten or twenty primes, after that: numbers of the form 6n+1 and 6n+5 may act as divisors. If you use only numbers of the form I mentioned, you gain about 50% with no loss of information." This uses the first part of the proposal: proc primetest2 n { set max [expr int(sqrt($n))] foreach i { 2 3 5 7 11 13 17 19 23 29 31 37 41 43 47 53 59 61 67 71 73 79 83 89 97 101 103 107 109 113 127 131 137 139 149 151 157 163 167 173 179 181 191 193 197 199 } { if {$i>$max} break if {$n%$i == 0} {return [list $i [expr $n/$i]]} } for {incr i 2} {$i<$max} {incr i 2} { if {$n%$i == 0} {return [list $i [expr $n/$i]]} } return 1 } ;# RS [CL]: A comment on prime calculations: for modest sizes--up to a few tens of decimal digits--I've always found it much the most satisfying to cache a list of known primes. I'm working in a stateful environment. My basic algorithm: test a trial number for quotients against my list of known primes, up to the square root of the trial number. Extend the list of primes, if it doesn't reach that high. It's a compact algorithm, it takes little space to keep the list, and I get excellent speed, especially for the frequent case where I want to know the primality of several different integers. CL offers an example implementation [[ugh: need to fix corrupted indentation]]: proc initialize_primes {} { namespace eval primes { variable list {} proc is_prime candidate { variable list foreach prime $list { # If a candidate has a factor, then # it's composite. if {![expr $candidate % $prime]} { return false } # If a candidate has no prime factor # smaller than or equal to its # square root, then it has no prime # factor, and therefore no factor, # and therefore is itself prime. if {[expr $prime * $prime > $candidate]} { return true } } # Use the list of primes we know to test the # primality of the candidate. If the tests # are inconclusive, extend the list, until # the question is settled. In any case, be # careful to maintain the list of known primes; # this is valuable in consideration of new # candidates. while true { set largest [get_next_prime] if {[expr $largest * $largest >= $candidate]} { return [is_prime $candidate] } } } proc get_next_prime {} { variable list if {{} == $list} { set candidate 2 } else { set candidate [lindex $list end] while true { incr candidate if [is_prime $candidate] break } } lappend list $candidate return $candidate } } } initialize_primes # Print the first thousand primes. for {set i 1} {$i <= 1000} {incr i} { puts [primes::get_next_prime] } [DPE] '8 Aug 2006' I've optimised the above code by simple changes such as removing [expr] where [if] is used. On my machine it runs approximately 20x faster (e.g. primes::is_prime 1234567891 takes 0.17s vs 4.08s). [glennj] 2009-04-20 small optimization in get_next_prime: incr by 2 as even numbers are obviously not prime. Also, add "restart" proc to complement "reset" without throwing away the known primes. namespace eval primes {} proc primes::reset {} { variable list [list] variable current_index end } namespace eval primes {reset} proc primes::restart {} { variable list variable current_index if {[llength $list] > 0} { set current_index 0 } } proc primes::is_prime {candidate} { variable list foreach prime $list { # If a candidate has a factor, then # it's composite. if {$candidate % $prime == 0} { return false } # If a candidate has no prime factor # smaller than or equal to its # square root, then it has no prime # factor, and therefore no factor, # and therefore is itself prime. if {$prime * $prime > $candidate} { return true } } # Use the list of primes we know to test the # primality of the candidate. If the tests # are inconclusive, extend the list, until # the question is settled. In any case, be # careful to maintain the list of known primes; # this is valuable in consideration of new # candidates. while true { set largest [get_next_prime] if {$largest * $largest >= $candidate} { return [is_prime $candidate] } } } proc primes::get_next_prime {} { variable list variable current_index if {$current_index ne "end"} { set p [lindex $list $current_index] if {[incr current_index] == [llength $list]} { set current_index end } return $p } switch -exact -- [llength $list] { 0 {set candidate 2} 1 {set candidate 3} default { set candidate [lindex $list end] while true { incr candidate 2 if {[is_prime $candidate]} break } } } lappend list $candidate return $candidate } # For testing # Print the first thousand primes. proc printAllPrimes {num} { for {set i 1} {$i <= $num} {incr i} { puts [primes::get_next_prime] } } printAllPrimes 1000 ---- The first 6 million primes: http://www.geocities.com/ResearchTriangle/Thinktank/2434/prime/primenumbers.html The largest known prime as of 2006 is the Mersenne prime 2**30402457 - 1 [http://www.mersenne.org] A prime tester on the Web at http://alpertron.com.ar/ecm.htm ---- Another useful prime-related concept is that of ''co-primality''. Two numbers are co-prime if the greatest common divisor is 1. Note that if ''P'' is a prime number, then for any ''Q'', ''P'' and ''Q'' are co-prime. ---- [RS] wrote this cute recursive prime tester for breakfast fun: proc prime i {prime0 $i [expr {int(sqrt($i))}]} proc prime0 {i div} { expr {!($i%$div)? 0: $div<=2? 1: [prime0 $i [incr div -1]]} } ---- [DL] wrote a Tcl script to let kids experiment with [Eratosthenes Sieve]s. ---- '''Complex factorization of primes''': Carl Friedrich Gauss knew it some 200 years ago, [RS] just read Easter 2005 that some (very roughly, half of the) primes can be factorized into a pair of complex numbers, where the second is the conjugate of the first, e.g. for 2: : (1+''i'')*(1-''i'') = 1 + ''i'' - ''i''*(1+''i'') = 1 + ''i'' - ''i'' - ''i''² = 2 I first tested this using [Complex math with TOOT], but as the product of a complex number a+''i''b and its conjugate a-''i''b is just a² + b², the test can be done easier as the sum of the two squares. The following proc returns the first list of {a b} found (the second is always {b a}), so that a²+b² == x, or an empty list when nothing was found: proc gaussfactor x { for {set real 1} {$real <=sqrt($x)} {incr real} { for {set im 1} {$im<=sqrt($x)} {incr im} { if {$real*$real+$im*$im == $x} { return [list $real $im] } } } } Here's the first few primes tested - five of ten have a complex factorization: % gaussfactor 2 1 1 % gaussfactor 3 % gaussfactor 5 1 2 % gaussfactor 7 % gaussfactor 11 % gaussfactor 13 2 3 % gaussfactor 17 1 4 % gaussfactor 19 % gaussfactor 23 % gaussfactor 29 2 5 [AM] (29 march 2005) I am not sure about the details but it seems that certain problems regarding prime numbers are easier to solve using this complex factorisation than if you stick to ordinary integers. Complex numbers with integer real and imaginary fractions are also known as Gaussian integers. [SMH] (05 Jan 2006) The inner loop (for im) is not necessary. For any value of `$real` only one `$im` is possible. proc gaussfactor2 x { set sx [expr {int(sqrt($x))}] for {set real 1} {$real < $sx} {incr real} { set d [expr {$x - $real * $real}] set im [expr {int(sqrt($d))}] if { $im * $im == $d } { return [list $real $im] } } return "" } ---- A SourceForge bug [http://sourceforge.net/tracker/index.php?func=detail&aid=1513763&group_id=10894&atid=110894%|% 1513763%|%] brought this very cute, [regexp]-based prime finder (based on an old [perl] hack): proc isprime x { expr {$x>1 && ![regexp {^(oo+?)\1+$} [string repeat o $x]]} } ---- [Twylite] If you need primes for cryptographic purposes (e.g. RSA keys), consider using the prime-finding algorithms in http://csrc.nist.gov/publications/PubsFIPS.html%|%FIPS 186-3%|%. It is often preferable (and quicker) to find a prime p such that (p-1) and (p+1) have large prime factors, as described in that standard. To speed up the process use trial division by small primes to sieve the candidates before performing the Miller-Rabin test. ---- !!!!!! %| [Category Concept] | [Category Mathematics] |% [Arts and crafts of Tcl-Tk programming] !!!!!!