php/Math   
Recreational Mathematics   
   home  |  library  |  contact
 Math Notes
 Math Programming [25]
 Regression [3]
 Data Mining [17]
 Notation [6]
 Linear Algebra [9]
 Stats & Prob [15]
 Math Cognition [5]
 Space & Physics [6]
 Formulas [5]
 Fun & Games [2]
 Haskell [1]
 Bayes Theory [1]
 Site News [0]
 Math Projects [5]
 Polynomials [1]
 Calculus [9]
 Number Theory [3]
 Optimization [2]
 Financial [1]

 Math Links
 PHP/ir
 Andrew Gelman
 Chance Wiki
 Daniel Lemire
 KD Knuggets
 Social Stats
 MySQL Performance
 Hunch.net
 Matthew Hurst
 JMLR
 JSS
 Hal Daume III
 Math Notes >> Number Theory

Primality testing [Number Theory
Posted on January 8, 2008 @ 03:27:14 AM by Paul Meagher

Primality testing refers to a test for whether a number is prime number or not. The is_prime.php script below offers two different approaches to primality testing:

<?php

/**
* Determines whether a number $n is prime by incrementing 
* a number $i and testing whether $n mod $i has a reminder 
* of 0.  If so, the number is not prime. 

* @see http://www.phpclasses.org/browse/package/3032.html
*/
function mod_prime($n) {
  
$sqrtn sqrt($n);
  for(
$i=2$i <= $sqrtn$i++ )
    if (
$n $i == 0) return false;
  return 
true ;
}

/** 
* Uses a clever regex method to determine if a number $n 
* is prime.
*
* @see http://www.noulakaz.net/weblog/2007/03/18/
*/
function regex_prime($n) {
  if (!
preg_match('/^1?$|^(11+?)\1+$/x'str_repeat('1'$n))) 
    return 
true;
  else 
    return 
false;
}

?>

The is_prime_timer.php script below tests the relative speeds of these two primality testing approaches:

<?php
/**
* Script to test relative speed of the two primality 
* testing methods.
*/

include "is_prime.php";

/**
* Timer convenience function.
*/
function microtime_float() {
  list(
$usec$sec) = explode(" "microtime());
  return ((float)
$usec + (float)$sec);
}

// 1. Time primality testing using mod_prime
$time_start microtime_float();
$n          2;
$end        5000;
$num_primes 0;
while (
$n $end) {
  if (
mod_prime($n))
    
$num_primes++; 
  
$n++;
}
$time_end microtime_float();
$time $time_end $time_start;
?>

<p>
<?php 
echo "The mod_prime function found $num_primes primes up to $end in $time seconds";
?>
</p>

<?php
// 2. Time primality testing using regex_prime
$time_start microtime_float();
$n          2;
$end        5000;
$num_primes 0;
while (
$n $end) {
  if (
regex_prime($n))
    
$num_primes++; 
  
$n++;
}
$time_end microtime_float();
$time $time_end $time_start;
?>

<p>
<?php 
echo "The regex_prime function found $num_primes primes up to $end in $time seconds";
?>
</p>

An execution of the is_prime_timer.php script produced the following output:

The mod_prime function found 669 primes up to 5000 in 0.033866882324219 seconds

The regex_prime function found 669 primes up to 5000 in 5.8958280086517 seconds

So it appears that while the regex approach may be clever it is not very fast.

FYI, the fastest primality testing program in existence right now is PRIMO. Further exploration of this topic would be a good way to explore number theory. Interestingly, there is a version of BASIC called UBASIC (Ultra high-precision BASIC for math) that specializes in number theory, primality testing, factoring, and large integers (up to 2600 digits).

Permalink 

Direct Search Factorization [Number Theory
Posted on January 4, 2008 @ 05:26:35 AM by Paul Meagher

Wikipedia defines prime factors as follows:

In number theory, the prime factors of a positive integer are the prime numbers that divide into that integer exactly, without leaving a remainder. The process of finding these numbers is called integer factorization, or prime factorization.

Direct Search Factorization is one method for finding prime factors. Consult Wikipedia's "Integer Factorization" entry for a discussion of other methods.

Below is a class implementing the Direct Search Factorization algorithm:

<?php

/**
* Factorization by trial division 

* Algorithm is described in:
*
* Donald E. Knuth, The Art of Computer Programming, 3rd edition, 
* Volume 2: Seminumerical Algorithms (Addison-Wesley, 1998), 
* section 4.5.4 (Factoring Into Primes), pp. 380-381. 

* @see http://www.merriampark.com/factor.htm
*/

include "SieveOfEratosthenes.php";

class 
DirectSearchFactorization {

  var 
$max 0;
  var 
$trialDivisors = array();

   
  
/** 
  * Constructor
  *
  * Generates trial divisors.  The trial divisors include all 
  * primes <= sqrt(n).  The trial divisors must also contain 
  * at least one value such that d(k) >= sqrt(n).
  */
  
public function DirectSearchFactorization($max=500000) {
    
    
$this->max $max;
    
    
$sqrtn = (int)ceil(sqrt($this->max));    
    
$sieve = new SieveOfEratosthenes;
    
$this->trialDivisors   $sieve->getPrimes($sqrtn);
    
$this->trialDivisors[] = $sqrtn;
    
  }

  
/**
  * Main method for computing the prime factors of some number n.
  */
  
public function getPrimeFactors($n) {
    
    
$factors = array();
  
    if (
$n 0) {
  
      if (
$n $this->max
        die(
"Illegal arguments: ".$n." > ".$this->max);
      
      
$k 0;
  
      while (
$n 1) {
          
        
$divisor   $this->trialDivisors[$k];
        
$quotient  $n $divisor;
        
$remainder $n $divisor;
          
        if (
$remainder == 0) {
          
$factors[] = $divisor;
          
$n         $quotient;
        } else {
          if (
$quotient $divisor
            
$k++;
          else {          
            
$factors[]=$n;
            break;
          }
        }
                
      }
      
      return 
$factors;
      
    } else 
      die(
"n must be greater than 0.");      
  }
  
}

?>

Notice that the class depends upon the SieveOfEratosthenes class discussed in my last blog.

Here is a simple test program that finds the prime factorization of 100:

<?php

include "DirectSearchFactorization.php";

$dsf     = new DirectSearchFactorization;
$factors $dsf->getPrimeFactors(100);

echo 
"<pre>";
print_r($factors);
echo 
"</pre>";

?>

The test program outputs the following prime factors:

Array
(
    [0] => 2
    [1] => 2
    [2] => 5
    [3] => 5
)

We can verify that 2 * 2 * 5 * 5 = 100 which means the class has found a prime factorization for 100.

Permalink 

Sieve Of Eratosthenes [Number Theory
Posted on January 3, 2008 @ 01:53:39 AM by Paul Meagher

The National Institute of Standards and Technology (NIST) defines the Sieve of Eratosthenes as follows:

Definition: An algorithm to find all prime numbers up to a certain N. Begin with an (unmarked) array of integers from 2 to N. The first unmarked integer, 2, is the first prime. Mark every multiple of this prime. Repeatedly take the next unmarked integer as the next prime and mark every multiple of the prime.

Here is a Wikpedia animation of the sieve process:

Here is a ported PHP class that implements the Sieve of Eratosthenes:

<?php
/**
* This class uses the sieve of Eratosthenes to
* find prime numbers
*
* @see http://www.merriampark.com/factor.htm
*/

class SieveOfEratosthenes {

  function 
getPrimes($upTo) {

    if (
$upTo 0) {
  
      
$size   $upTo 1;
      
$flags  = array();
      
$primes = array();
      
$limit  sqrt($size);
  
       
// Set flags
       
for ($i=2$i $size$i++) 
         
$flags[$i] = true;
  
       
// Cross out multiples of 2
       
$j 2;
       for (
$i=$j+$j$i $size$i=$i+$j
         
$flags[$i] = false;
  
       
// Cross out multiples of odd numbers
       
for ($j=3$j <= $limit$j=$j+2) {
         if (
$flags[$j]) {
           for (
$i=$j+$j$i $size$i=$i+$j
             
$flags[$i] = false;
         }
       }
  
       
// Build list of primes from what is left
       
for ($i=2$i $size$i++) {
         if (
$flags[$i]) 
           
$primes[] = $i;
       }
  
       return 
$primes;

    } else     
      die(
"n must be greater than 0.");
      
  }

}

?>

Here is a class testing program:

<?php

include "SieveOfEratosthenes.php";

$sieve  = new SieveOfEratosthenes;
$primes $sieve->getPrimes(120);

echo 
"<pre>";
print_r($primes);
echo 
"</pre>";

?>

Here is the output of the test program:

Array
(
    [0] => 2
    [1] => 3
    [2] => 5
    [3] => 7
    [4] => 11
    [5] => 13
    [6] => 17
    [7] => 19
    [8] => 23
    [9] => 29
    [10] => 31
    [11] => 37
    [12] => 41
    [13] => 43
    [14] => 47
    [15] => 53
    [16] => 59
    [17] => 61
    [18] => 67
    [19] => 71
    [20] => 73
    [21] => 79
    [22] => 83
    [23] => 89
    [24] => 97
    [25] => 101
    [26] => 103
    [27] => 107
    [28] => 109
    [29] => 113
)

I have some more prime number theory code lying around that I will post in the coming days.

Permalink 

 Archive 
 


php/Math Project
© 2011. All rights reserved.