num_theory/miller_rabin/index.php

<?php
require_once("style/std/stdwpstyle.inc");
//----------------------------------------------------------------------------------------------------
//Copyright (c) 2003, 2018 David T. Ashley.
//
//Permission is hereby granted, free of charge, to any person obtaining a copy
//of this software and associated documentation files (the "Software"), to deal
//in the Software without restriction, including without limitation the rights
//to use, copy, modify, merge, publish, distribute, sublicense, and/or sell
//copies of the Software, and to permit persons to whom the Software is
//furnished to do so, subject to the following conditions:
//
//The above copyright notice and this permission notice shall be included in all
//copies or substantial portions of the Software.
//
//THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
//IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
//FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
//AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
//LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,
//OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE
//SOFTWARE.
//----------------------------------------------------------------------------------------------------
//
//Returns TRUE if the passed argument is a digit, or FALSE otherwise.  If the
//passed string has length > 1, only the least significant digit is
//considered.
function miller_rabin_isdigit($Digit)
   {
   $Digit = (string)$Digit;

   //If the string has zero length, this is an error condition and return 0.
   $Length = strlen($Digit);  

   if ($Length==0)
      {
      return(0);
      }

   //Obtain the trailing character of the string.
   $Trailing_char = Substr($Digit, $Length-1, 1);

   //Do a switch and return an integer based on the character.
   switch($Trailing_char)
      {
      case "9":
         return 1;
         break;
      case "8":
         return 1;
         break;
      case "7":
         return 1;
         break;
      case "6":
         return 1;
         break;
      case "5":
         return 1;
         break;
      case "4":
         return 1;
         break;
      case "3":
         return 1;
         break;
      case "2":
         return 1;
         break;
      case "1":
         return 1;
         break;
      case "0":
         return 1;
         break;
      default:
         return 0;
         break;
      }
   }
//
//
//Returns TRUE if the passed argument is a simple integer (no error   
//codes, no commas), or FALSE otherwise.
function miller_rabin_issimpleint($Intnum)
   {
   $Intnum = (string)$Intnum;

   //If the string has zero length, this is an error condition and return 0.
   $Length = strlen($Intnum);  

   if ($Length==0)
      {
      return(0);
      }

   //If the first digit is "0", this must be the only character in the string,
   //otherwise there is a leading-zero violation.
   if (Substr($Intnum, 0, 1)=="0")
      {
      if ($Length==1)
         {
         return 1;
         }
      else
         {
         return 0;
         }
      }

   //If the first character is a "-" sign, it may not be the only character,
   //and it may not precede a "0" (which would again be a leading 0 violation).
   //If it is a "-", we can exclude "0" and parse as a positive.
   if (Substr($Intnum, 0, 1)=="-")
      {
      if ($Length==1)
         {
         return 0;
         }
      else
         {
         if (Substr($Intnum, 1, 1) == "0")
            {
            //This is illegal (a "-" followed by a "0").
            return 0;
            }
         else
            {
            //We can just remove the "-" and parse as a positive.
            $Intnum = Substr($Intnum, 1, $Length - 1);
            $Length--;
            }
         }
      }

   //It is known at this point that the length is at least 1 and the leading 
   //character is not zero.  All we need to do at this point is check that each 
   //and every character is a digit.
   for ($index=0; $index < $Length; $index++)
      {
      if (!miller_rabin_isdigit(Substr($Intnum, $index, 1)))
         return 0;
      }

   return 1;
   }
//
//
//Returns the number of digits in the argument, which must be a valid simple
//integer.
function miller_rabin_ndigits($Intnum)
   {
   $Intnum = (string)$Intnum;

   if (Substr($Intnum, 0, 1) == "-")
      {
      $rv = strlen($Intnum) - 1;
      }
   else
      {
      $rv = strlen($Intnum);
      }

   return($rv);
   }
//
//
//Given a string, converts the final character to an integer.  If the character
//is not a digit, returns 0.
function miller_rabin_digittoint($Digit)
   {
   //echo "Digit is: ";
   //echo $Digit;
   //echo "<br>";

   $Digit = (string)$Digit;

   //echo "Digit is: ";
   //echo $Digit;
   //echo "<br>";
   
   //If the string has zero length, this is an error condition and return 0.
   $Length = strlen($Digit);  

   //echo "Length is: ";
   //echo $Length;
   //echo "<br>";

   if ($Length==0)
      {
      return(0);
      }

   //Obtain the trailing character of the string.
   $Trailing_char = Substr($Digit, $Length-1, 1);

   //echo "Trailing char is: ";
   //echo $Trailing_char;
   //echo "<br>";

   //Do a switch and return an integer based on the character.
   //The default case gives 0.

   switch($Trailing_char)
      {
      case "9":
         return 9;
         break;
      case "8":
         return 8;
         break;
      case "7":
         return 7;
         break;
      case "6":
         return 6;
         break;
      case "5":
         return 5;
         break;
      case "4":
         return 4;
         break;
      case "3":
         return 3;
         break;
      case "2":
         return 2;
         break;
      case "1":
         return 1;
         break;
      default:
         return 0;
         break;
      }
   }

//Decommanates an integer and removes all other extraneous characters without respect for 
//whether the commas are in the right place or the integer is syntactically valid.
function miller_rabin_decommanate(&$arg)
   {
   $temp = "";

   $len = strlen($arg);

   for ($i=0; $i<$len; $i++)
      {
      $c = (string)Substr($arg, $i, 1);
   
      if (miller_rabin_isdigit($c) || ($c == "-"))
         $temp = $temp . $c;
      }

   $arg = (string)$temp;
   } 

//Commanates an integer which is assumed to by syntactically correct.
function miller_rabin_commanate(&$arg)
   {
   miller_rabin_decommanate($arg);
 
   $temp = "";
  
   $len = strlen($arg);

   for ($i=0; $i<$len; $i++)
      {
      $c = (string)Substr($arg, $len - $i - 1, 1);

      $temp = $c . $temp;
      
      if ($i != ($len-1))
         {
         if (Substr($arg, $len-$i-2, 1) != "-")
            {
            if (($i % 3) == 2)
               $temp = "," . $temp;
            }
         }
      }

   $arg = (string)$temp;
   } 

function do_header(&$style, $miller_rabin_nrounds)
   {
   $style->static_page_header_title_std("Probable Primality Test (Miller-Rabin)",
                        "Miller-Rabin Probable Primality Test",
                        "This utility uses " . $miller_rabin_nrounds . " applications of the " .
                        "<a href=\"http://mathworld.wolfram.com/Rabin-MillerStrongPseudoprimeTest.html\">" .
                        "Miller-Rabin test</a> to determine " .
                        "with high probabilistic certainty whether an integer " .
                        "is " .
                        "<a href=\"http://mathworld.wolfram.com/PrimeNumber.html\">prime</a> or " .
                        "<a href=\"http://mathworld.wolfram.com/CompositeNumber.html\">composite</a>." .
                        "&nbsp;  This utility does <i>not</i> " .
                        "factor composite numbers.");
   }  

function do_footer(&$style)
   {
//   $style->static_page_hrule_std();
//   echo "<p align=\"center\">This utility is powered by the ";
//   echo "<a href=\"http://www.swox.com/gmp\"><img hspace=\"5\" vspace=\"5\" src=\"../genimages/gmplogo1.png\" width=\"232\" height=\"104\" align=\"absmiddle\"></a>.&nbsp; All source code is available ";
//   echo "<a href=\"../howtos/nth_web_gmp_src_code_dist/obtain_all_source.htm\">here</a>.</p>\n";
   $style->static_page_footer_std();
   }

function do_form($N_refresh, $Buttontext, $max_digits)
   {
   echo "<form method=get action=\"index.php\" width=\"100%\">\n";
   echo "<table align=\"center\">\n";
   echo "<tr>\n";
   echo "  <td width=\"10%\">\n";
   echo "   <p align=\"right\"><b>N</b> (" . number_format($max_digits) . "-digit maximum):</td>\n";
   echo "    <td width=\"18%\"><p align=\"left\"><input type=\"text\" ";
   if (strlen($N_refresh))
     {
     echo "value=\"" . $N_refresh . "\" ";
     }
   echo "name=\"N\" size=\"50\"></p></td>\n";
   echo "  </tr>\n";
   echo "  <tr>\n";
   echo "    <td width=\"20%\" colspan=\"2\">\n";
   echo "      <p align=\"center\"><input type=\"submit\" value=\"";
   echo $Buttontext;
   echo "\" name=\"B1\" style=\"margin-top: 12\"></td>\n";
   echo "  </tr>\n";
   echo "</table>\n";
   echo "</form>\n";
}

function do_err_msg($Err_msg)
   {
   echo "<p align=\"center\"><b><font color=\"#FF0000\">";
   echo $Err_msg;
   echo "</font></b></p>\n";
   }


//Main script begins here.
//
//
//Let's agree on a maximum for number of digits in the number to test.  It is known that
//Apache won't process URL's larger than 1024, so let's be conservative and say 500
//digits.
$max_digits = 500;
//
//Let's agree on a number of rounds for the Miller-Rabin test.  This influences the
//probabilities.
$miller_rabin_nrounds = 25;
//
//We need to establish the probability that a number which passes the rounds specified
//above is still in fact composite.  This values need to be calculated by hand whenever
//the value above is changed.
$miller_rabin_p_real                 =  "0.99999999999999911182";
$miller_rabin_p_percent              = "99.999999999999911182";
$miller_rabin_residual_real          =  "0.00000000000000088818";
$miller_rabin_residual_percent       =  "0.000000000000088818";
$miller_rabin_residual_lt_reciprocal =  "1,125,899,906,842,624";
//
$style = new StdWpStyle;
//Assign the current style in force.  Also, starts the CPU
//usage clock.

//Do the header unconditionally.  The header is always used on this page.
do_header($style, $miller_rabin_nrounds);
//
if (isset($_GET['N']))
   $N = $_GET['N'];
//
//If N was supplied, decommanate it and be sure it is a string.
if (isset($N))
   miller_rabin_decommanate($N);
//
//There are a few cases to break into here, and this affects the overall layout of the page.
//The header and footer are fixed, but the "guts" will change.
if (!isset($N))
   {
   //In this case, we are probably visiting the form for the first time, and have not done a submit.
   //Just display the form itself.
   do_form("", "Is This Prime?", $max_digits);
   }
elseif (!miller_rabin_issimpleint($N))
   {
   do_err_msg("N is a malformed integer. &nbsp;Choose a new value for N and try again.");
   $style->static_page_hrule_std();  
   miller_rabin_commanate($N);
   do_form($N, "Is This Prime?", $max_digits);
   }  
elseif ($N == "1")
   {
   do_err_msg("It is not a valid question to ask if 1 is prime.&nbsp; It is agreed " .
              "in number theory (largely to minimize exception cases in proofs and results) " .
              "that 1 is neither prime nor composite.&nbsp; " .
              "Choose a new value for N and try again.");
   $style->static_page_hrule_std();  
   miller_rabin_commanate($N);
   do_form($N, "Is This Prime?", $max_digits);
   }  
elseif (strlen($N) > $max_digits)
   {
   do_err_msg("N is too long (a " . number_format($max_digits) . "-digit maximum is " .
               "in effect).&nbsp; Choose a shorter value and try again.");
   $style->static_page_hrule_std();  
   miller_rabin_commanate($N);
   do_form($N, "Is This Prime?", $max_digits);
   }
else
   {
   //N seems to be a valid integer.
   //
   //If N is negative, make it positive.
   if (Substr($N, 0, 1)=="-")
     $N = Substr($N, 1, strlen($N)-1);

   //If N is "0", make it "2".
   if ($N=="0")
     $N = (string)"2";

   //We can now run the external program to do the calculation.
   $cgi_command_string = "/hl/cgibin/aux_progs/dtats_cgi_aux_arith_large gmp_prob_prime " . $N . " " . $miller_rabin_nrounds;

   //echo "CGI command string is : " . $cgi_command_string . "<br>\n";

   exec($cgi_command_string, $cgi_result);
   $cgi_result_nelem = count($cgi_result);

   //echo "Result is : " . $cgi_result_nelem . "<br>\n";

   //for ($i=0; $i<$cgi_result_nelem; $i++)
   //{
   //  echo $cgi_result[$i] . "<br>\n";
   //}

   //We need to perform some sanity checks on the CGI output to be sure it
   //is what we want.  If it seems wrong, we need to generate an error message
   //and not try to display the output.
   $cgi_output_is_sane = 1;
   if ($cgi_result_nelem == 0)
      $cgi_output_is_sane = 0;
   if ($cgi_result_nelem > 0)
      {
      if ($cgi_result_nelem != 5)
         $cgi_output_is_sane = 0;
      if ($cgi_result[0] != "S")
         $cgi_output_is_sane = 0;
      if (($cgi_result[3] != "0") && ($cgi_result[3] != "1") && ($cgi_result[3] != "2"))
         $cgi_output_is_sane = 0;
      if ($cgi_result[4] != "S")
         $cgi_output_is_sane = 0;
      }

   if (!$cgi_output_is_sane)
      {
      do_err_msg("An unspecified error occurred when interacting with the CGI-BIN program.&nbsp; Please advise the webmaster.");
      $style->static_page_hrule_std();  
      miller_rabin_commanate($N);
      do_form($N, "Is This Prime", $max_digits);
      }
   else
      {
      miller_rabin_commanate($cgi_result[1]);

      if ($cgi_result[3] == 0)
         {
         echo "<p align=\"center\"><font size=\"4\">The integer ";
         echo $cgi_result[1];
         echo " is definitely (with 100% certainty) composite.&nbsp; Note, however, that although ";
         echo "it was relatively computationally inexpensive to determine that this integer is composite, ";
         echo "factoring it may be computationally <i>much</i> more expensive.</font></p>\n";
         }
      elseif ($cgi_result[3] == 1)
         {
         echo "<p align=\"center\"><font size=\"4\">Based on " . $miller_rabin_nrounds . " applications of ";
         echo "the Miller-Rabin test, the integer ";
         echo $cgi_result[1];
         echo " is probably (with " . $miller_rabin_p_percent . "% certainty) ";
         echo "prime.&nbsp; (This means that there is no more than 1 chance in ";
         echo $miller_rabin_residual_lt_reciprocal . " that this integer is ";
         echo "composite.)</font></p>\n";
         }
      else
         {
         echo "<p align=\"center\"><font size=\"4\">The integer ";
         echo $cgi_result[1];
         echo " is prime (with 100% certainty).&nbsp; This means that the ";
         echo "GMP";
         echo " library routine was able to attempt trial division with divisors ";
         echo "up to <i>floor</i>(<i>sqrt</i>(N)).</font></p>\n";
         }

      $style->static_page_hrule_std();  
?>
<p align="left"><b>Notes on the algorithm:</b>
<ul>
    <li>It is agreed in number theory that 
        <a href="http://www.mth.uct.ac.za/~digest/faq/prime.html">'1' is neither prime nor 
        composite</a>.&nbsp; This definition
        is largely <a href="http://dictionary.reference.com/search?q=arbitrary">arbitrary</a>, 
        and there is some discussion in Usenet newsgroups
        such as <i>sci.math</i> from time to time about this topic.&nbsp; It substantially shortens
        and neatens many results and proofs to define '1' in this way.&nbsp; For this reason,
        this utility refuses to <a href="http://dictionary.reference.com/search?q=opine">opine</a> on '1'. 
        </li>
    <li>The <a href="http://mathworld.wolfram.com/Rabin-MillerStrongPseudoprimeTest.html">Miller-Rabin 
        test</a> and other 
        <a href="http://mathworld.wolfram.com/StrongPseudoprime.html">strong pseudoprime</a> tests 
        are part of the backbone of popular 
        <a href="http://mathworld.wolfram.com/Public-KeyCryptography.html">public-key cryptography</a>,
        specifically 
        <a href="http://mathworld.wolfram.com/RSAEncryption.html">RSA public-key cryptography</a>,
        which relies on the fact that it is much easier to create a large composite number with
        large prime components than to factor it.        
        </li>
    <li>This utility uses the <i>mpz_prob_prime_p()</i> function from
        the <a href="http://www.swox.com/gmp">GMP</a>
        to make its primality assessment.&nbsp; This function has not been proofread
        to verify how it selects bases and applies the Miller-Rabin test.&nbsp; It is assumed 
        that the standard upper bound of 0.25^N can be applied where N is the <i>reps</i>
        parameter to the function (in the case of this utility, <i>reps</i> is  
        <?php echo " " . $miller_rabin_nrounds ; ?>). 
        </li>
    <li>This web site also contains a <a href="pfact18digit.php">factorization
        utility</a> which will factor integers subject to size and CPU time constraints.                
</ul>
<?php
      $style->static_page_hrule_std();  

      miller_rabin_commanate($N);
      do_form($N, "Is This Prime?", $max_digits);
      }
   }

//Now do the footer unconditionally.  The footer is always used on this web page.
do_footer($style);  
?>
1	dashley	23	<?php
2	dashley	190	require_once("style/std/stdwpstyle.inc");
3			//----------------------------------------------------------------------------------------------------
4			//Copyright (c) 2003, 2018 David T. Ashley.
5	dashley	23	//
6	dashley	190	//Permission is hereby granted, free of charge, to any person obtaining a copy
7			//of this software and associated documentation files (the "Software"), to deal
8			//in the Software without restriction, including without limitation the rights
9			//to use, copy, modify, merge, publish, distribute, sublicense, and/or sell
10			//copies of the Software, and to permit persons to whom the Software is
11			//furnished to do so, subject to the following conditions:
12	dashley	23	//
13	dashley	190	//The above copyright notice and this permission notice shall be included in all
14			//copies or substantial portions of the Software.
15	dashley	23	//
16	dashley	190	//THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
17			//IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
18			//FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
19			//AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
20			//LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,
21			//OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE
22			//SOFTWARE.
23			//----------------------------------------------------------------------------------------------------
24			//
25	dashley	23	//Returns TRUE if the passed argument is a digit, or FALSE otherwise. If the
26			//passed string has length > 1, only the least significant digit is
27			//considered.
28			function miller_rabin_isdigit($Digit)
29			{
30			$Digit = (string)$Digit;
31
32			//If the string has zero length, this is an error condition and return 0.
33			$Length = strlen($Digit);
34
35			if ($Length==0)
36			{
37			return(0);
38			}
39
40			//Obtain the trailing character of the string.
41			$Trailing_char = Substr($Digit, $Length-1, 1);
42
43			//Do a switch and return an integer based on the character.
44			switch($Trailing_char)
45			{
46			case "9":
47			return 1;
48			break;
49			case "8":
50			return 1;
51			break;
52			case "7":
53			return 1;
54			break;
55			case "6":
56			return 1;
57			break;
58			case "5":
59			return 1;
60			break;
61			case "4":
62			return 1;
63			break;
64			case "3":
65			return 1;
66			break;
67			case "2":
68			return 1;
69			break;
70			case "1":
71			return 1;
72			break;
73			case "0":
74			return 1;
75			break;
76			default:
77			return 0;
78			break;
79			}
80			}
81			//
82			//
83			//Returns TRUE if the passed argument is a simple integer (no error
84			//codes, no commas), or FALSE otherwise.
85			function miller_rabin_issimpleint($Intnum)
86			{
87			$Intnum = (string)$Intnum;
88
89			//If the string has zero length, this is an error condition and return 0.
90			$Length = strlen($Intnum);
91
92			if ($Length==0)
93			{
94			return(0);
95			}
96
97			//If the first digit is "0", this must be the only character in the string,
98			//otherwise there is a leading-zero violation.
99			if (Substr($Intnum, 0, 1)=="0")
100			{
101			if ($Length==1)
102			{
103			return 1;
104			}
105			else
106			{
107			return 0;
108			}
109			}
110
111			//If the first character is a "-" sign, it may not be the only character,
112			//and it may not precede a "0" (which would again be a leading 0 violation).
113			//If it is a "-", we can exclude "0" and parse as a positive.
114			if (Substr($Intnum, 0, 1)=="-")
115			{
116			if ($Length==1)
117			{
118			return 0;
119			}
120			else
121			{
122			if (Substr($Intnum, 1, 1) == "0")
123			{
124			//This is illegal (a "-" followed by a "0").
125			return 0;
126			}
127			else
128			{
129			//We can just remove the "-" and parse as a positive.
130			$Intnum = Substr($Intnum, 1, $Length - 1);
131			$Length--;
132			}
133			}
134			}
135
136			//It is known at this point that the length is at least 1 and the leading
137			//character is not zero. All we need to do at this point is check that each
138			//and every character is a digit.
139			for ($index=0; $index < $Length; $index++)
140			{
141			if (!miller_rabin_isdigit(Substr($Intnum, $index, 1)))
142			return 0;
143			}
144
145			return 1;
146			}
147			//
148			//
149			//Returns the number of digits in the argument, which must be a valid simple
150			//integer.
151			function miller_rabin_ndigits($Intnum)
152			{
153			$Intnum = (string)$Intnum;
154
155			if (Substr($Intnum, 0, 1) == "-")
156			{
157			$rv = strlen($Intnum) - 1;
158			}
159			else
160			{
161			$rv = strlen($Intnum);
162			}
163
164			return($rv);
165			}
166			//
167			//
168			//Given a string, converts the final character to an integer. If the character
169			//is not a digit, returns 0.
170			function miller_rabin_digittoint($Digit)
171			{
172			//echo "Digit is: ";
173			//echo $Digit;
174			//echo "<br>";
175
176			$Digit = (string)$Digit;
177
178			//echo "Digit is: ";
179			//echo $Digit;
180			//echo "<br>";
181
182			//If the string has zero length, this is an error condition and return 0.
183			$Length = strlen($Digit);
184
185			//echo "Length is: ";
186			//echo $Length;
187			//echo "<br>";
188
189			if ($Length==0)
190			{
191			return(0);
192			}
193
194			//Obtain the trailing character of the string.
195			$Trailing_char = Substr($Digit, $Length-1, 1);
196
197			//echo "Trailing char is: ";
198			//echo $Trailing_char;
199			//echo "<br>";
200
201			//Do a switch and return an integer based on the character.
202			//The default case gives 0.
203
204			switch($Trailing_char)
205			{
206			case "9":
207			return 9;
208			break;
209			case "8":
210			return 8;
211			break;
212			case "7":
213			return 7;
214			break;
215			case "6":
216			return 6;
217			break;
218			case "5":
219			return 5;
220			break;
221			case "4":
222			return 4;
223			break;
224			case "3":
225			return 3;
226			break;
227			case "2":
228			return 2;
229			break;
230			case "1":
231			return 1;
232			break;
233			default:
234			return 0;
235			break;
236			}
237			}
238
239	dashley	190	//Decommanates an integer and removes all other extraneous characters without respect for
240	dashley	23	//whether the commas are in the right place or the integer is syntactically valid.
241			function miller_rabin_decommanate(&$arg)
242			{
243			$temp = "";
244
245			$len = strlen($arg);
246
247			for ($i=0; $i<$len; $i++)
248			{
249			$c = (string)Substr($arg, $i, 1);
250
251			if (miller_rabin_isdigit($c) \|\| ($c == "-"))
252			$temp = $temp . $c;
253			}
254
255			$arg = (string)$temp;
256			}
257
258			//Commanates an integer which is assumed to by syntactically correct.
259			function miller_rabin_commanate(&$arg)
260			{
261			miller_rabin_decommanate($arg);
262
263			$temp = "";
264
265			$len = strlen($arg);
266
267			for ($i=0; $i<$len; $i++)
268			{
269			$c = (string)Substr($arg, $len - $i - 1, 1);
270
271			$temp = $c . $temp;
272
273			if ($i != ($len-1))
274			{
275			if (Substr($arg, $len-$i-2, 1) != "-")
276			{
277			if (($i % 3) == 2)
278			$temp = "," . $temp;
279			}
280			}
281			}
282
283			$arg = (string)$temp;
284			}
285
286			function do_header(&$style, $miller_rabin_nrounds)
287			{
288	dashley	190	$style->static_page_header_title_std("Probable Primality Test (Miller-Rabin)",
289	dashley	23	"Miller-Rabin Probable Primality Test",
290			"This utility uses " . $miller_rabin_nrounds . " applications of the " .
291			"<a href=\"http://mathworld.wolfram.com/Rabin-MillerStrongPseudoprimeTest.html\">" .
292			"Miller-Rabin test</a> to determine " .
293			"with high probabilistic certainty whether an integer " .
294			"is " .
295			"<a href=\"http://mathworld.wolfram.com/PrimeNumber.html\">prime</a> or " .
296			"<a href=\"http://mathworld.wolfram.com/CompositeNumber.html\">composite</a>." .
297			"  This utility does <i>not</i> " .
298			"factor composite numbers.");
299			}
300
301			function do_footer(&$style)
302			{
303	dashley	190	// $style->static_page_hrule_std();
304			// echo "<p align=\"center\">This utility is powered by the ";
305			// echo "<a href=\"http://www.swox.com/gmp\"><img hspace=\"5\" vspace=\"5\" src=\"../genimages/gmplogo1.png\" width=\"232\" height=\"104\" align=\"absmiddle\"></a>.  All source code is available ";
306			// echo "<a href=\"../howtos/nth_web_gmp_src_code_dist/obtain_all_source.htm\">here</a>.</p>\n";
307			$style->static_page_footer_std();
308	dashley	23	}
309
310			function do_form($N_refresh, $Buttontext, $max_digits)
311			{
312	dashley	190	echo "<form method=get action=\"index.php\" width=\"100%\">\n";
313	dashley	23	echo "<table align=\"center\">\n";
314			echo "<tr>\n";
315			echo " <td width=\"10%\">\n";
316			echo " <p align=\"right\"><b>N</b> (" . number_format($max_digits) . "-digit maximum):</td>\n";
317			echo " <td width=\"18%\"><p align=\"left\"><input type=\"text\" ";
318			if (strlen($N_refresh))
319			{
320			echo "value=\"" . $N_refresh . "\" ";
321			}
322			echo "name=\"N\" size=\"50\"></p></td>\n";
323			echo " </tr>\n";
324			echo " <tr>\n";
325			echo " <td width=\"20%\" colspan=\"2\">\n";
326			echo " <p align=\"center\"><input type=\"submit\" value=\"";
327			echo $Buttontext;
328			echo "\" name=\"B1\" style=\"margin-top: 12\"></td>\n";
329			echo " </tr>\n";
330			echo "</table>\n";
331			echo "</form>\n";
332			}
333
334			function do_err_msg($Err_msg)
335			{
336			echo "<p align=\"center\"><b><font color=\"#FF0000\">";
337			echo $Err_msg;
338			echo "</font></b></p>\n";
339			}
340
341
342			//Main script begins here.
343			//
344			//
345			//Let's agree on a maximum for number of digits in the number to test. It is known that
346			//Apache won't process URL's larger than 1024, so let's be conservative and say 500
347			//digits.
348			$max_digits = 500;
349			//
350			//Let's agree on a number of rounds for the Miller-Rabin test. This influences the
351			//probabilities.
352			$miller_rabin_nrounds = 25;
353			//
354			//We need to establish the probability that a number which passes the rounds specified
355			//above is still in fact composite. This values need to be calculated by hand whenever
356			//the value above is changed.
357			$miller_rabin_p_real = "0.99999999999999911182";
358			$miller_rabin_p_percent = "99.999999999999911182";
359			$miller_rabin_residual_real = "0.00000000000000088818";
360			$miller_rabin_residual_percent = "0.000000000000088818";
361			$miller_rabin_residual_lt_reciprocal = "1,125,899,906,842,624";
362			//
363	dashley	190	$style = new StdWpStyle;
364	dashley	23	//Assign the current style in force. Also, starts the CPU
365			//usage clock.
366
367			//Do the header unconditionally. The header is always used on this page.
368			do_header($style, $miller_rabin_nrounds);
369			//
370	dashley	190	if (isset($_GET['N']))
371			$N = $_GET['N'];
372	dashley	23	//
373			//If N was supplied, decommanate it and be sure it is a string.
374			if (isset($N))
375			miller_rabin_decommanate($N);
376			//
377			//There are a few cases to break into here, and this affects the overall layout of the page.
378			//The header and footer are fixed, but the "guts" will change.
379			if (!isset($N))
380			{
381			//In this case, we are probably visiting the form for the first time, and have not done a submit.
382			//Just display the form itself.
383			do_form("", "Is This Prime?", $max_digits);
384			}
385			elseif (!miller_rabin_issimpleint($N))
386			{
387			do_err_msg("N is a malformed integer.  Choose a new value for N and try again.");
388	dashley	190	$style->static_page_hrule_std();
389	dashley	23	miller_rabin_commanate($N);
390			do_form($N, "Is This Prime?", $max_digits);
391			}
392			elseif ($N == "1")
393			{
394			do_err_msg("It is not a valid question to ask if 1 is prime.  It is agreed " .
395			"in number theory (largely to minimize exception cases in proofs and results) " .
396			"that 1 is neither prime nor composite.  " .
397			"Choose a new value for N and try again.");
398	dashley	190	$style->static_page_hrule_std();
399	dashley	23	miller_rabin_commanate($N);
400			do_form($N, "Is This Prime?", $max_digits);
401			}
402			elseif (strlen($N) > $max_digits)
403			{
404			do_err_msg("N is too long (a " . number_format($max_digits) . "-digit maximum is " .
405			"in effect).  Choose a shorter value and try again.");
406	dashley	190	$style->static_page_hrule_std();
407	dashley	23	miller_rabin_commanate($N);
408			do_form($N, "Is This Prime?", $max_digits);
409			}
410			else
411			{
412			//N seems to be a valid integer.
413			//
414			//If N is negative, make it positive.
415			if (Substr($N, 0, 1)=="-")
416			$N = Substr($N, 1, strlen($N)-1);
417
418			//If N is "0", make it "2".
419			if ($N=="0")
420			$N = (string)"2";
421
422			//We can now run the external program to do the calculation.
423	dashley	190	$cgi_command_string = "/hl/cgibin/aux_progs/dtats_cgi_aux_arith_large gmp_prob_prime " . $N . " " . $miller_rabin_nrounds;
424	dashley	23
425			//echo "CGI command string is : " . $cgi_command_string . "<br>\n";
426
427			exec($cgi_command_string, $cgi_result);
428			$cgi_result_nelem = count($cgi_result);
429
430			//echo "Result is : " . $cgi_result_nelem . "<br>\n";
431
432			//for ($i=0; $i<$cgi_result_nelem; $i++)
433			//{
434			// echo $cgi_result[$i] . "<br>\n";
435			//}
436
437			//We need to perform some sanity checks on the CGI output to be sure it
438			//is what we want. If it seems wrong, we need to generate an error message
439			//and not try to display the output.
440			$cgi_output_is_sane = 1;
441			if ($cgi_result_nelem == 0)
442			$cgi_output_is_sane = 0;
443			if ($cgi_result_nelem > 0)
444			{
445			if ($cgi_result_nelem != 5)
446			$cgi_output_is_sane = 0;
447			if ($cgi_result[0] != "S")
448			$cgi_output_is_sane = 0;
449			if (($cgi_result[3] != "0") && ($cgi_result[3] != "1") && ($cgi_result[3] != "2"))
450			$cgi_output_is_sane = 0;
451			if ($cgi_result[4] != "S")
452			$cgi_output_is_sane = 0;
453			}
454
455			if (!$cgi_output_is_sane)
456			{
457			do_err_msg("An unspecified error occurred when interacting with the CGI-BIN program.  Please advise the webmaster.");
458	dashley	190	$style->static_page_hrule_std();
459	dashley	23	miller_rabin_commanate($N);
460			do_form($N, "Is This Prime", $max_digits);
461			}
462			else
463			{
464			miller_rabin_commanate($cgi_result[1]);
465
466			if ($cgi_result[3] == 0)
467			{
468			echo "<p align=\"center\"><font size=\"4\">The integer ";
469			echo $cgi_result[1];
470			echo " is definitely (with 100% certainty) composite.  Note, however, that although ";
471			echo "it was relatively computationally inexpensive to determine that this integer is composite, ";
472			echo "factoring it may be computationally <i>much</i> more expensive.</font></p>\n";
473			}
474			elseif ($cgi_result[3] == 1)
475			{
476			echo "<p align=\"center\"><font size=\"4\">Based on " . $miller_rabin_nrounds . " applications of ";
477			echo "the Miller-Rabin test, the integer ";
478			echo $cgi_result[1];
479			echo " is probably (with " . $miller_rabin_p_percent . "% certainty) ";
480			echo "prime.  (This means that there is no more than 1 chance in ";
481			echo $miller_rabin_residual_lt_reciprocal . " that this integer is ";
482			echo "composite.)</font></p>\n";
483			}
484			else
485			{
486			echo "<p align=\"center\"><font size=\"4\">The integer ";
487			echo $cgi_result[1];
488			echo " is prime (with 100% certainty).  This means that the ";
489			echo "GMP";
490			echo " library routine was able to attempt trial division with divisors ";
491			echo "up to <i>floor</i>(<i>sqrt</i>(N)).</font></p>\n";
492			}
493
494	dashley	190	$style->static_page_hrule_std();
495	dashley	23	?>
496			<p align="left"><b>Notes on the algorithm:</b>
497			<ul>
498			<li>It is agreed in number theory that
499			<a href="http://www.mth.uct.ac.za/~digest/faq/prime.html">'1' is neither prime nor
500			composite</a>.  This definition
501			is largely <a href="http://dictionary.reference.com/search?q=arbitrary">arbitrary</a>,
502			and there is some discussion in Usenet newsgroups
503			such as <i>sci.math</i> from time to time about this topic.  It substantially shortens
504			and neatens many results and proofs to define '1' in this way.  For this reason,
505			this utility refuses to <a href="http://dictionary.reference.com/search?q=opine">opine</a> on '1'.
506			</li>
507			<li>The <a href="http://mathworld.wolfram.com/Rabin-MillerStrongPseudoprimeTest.html">Miller-Rabin
508			test</a> and other
509			<a href="http://mathworld.wolfram.com/StrongPseudoprime.html">strong pseudoprime</a> tests
510			are part of the backbone of popular
511			<a href="http://mathworld.wolfram.com/Public-KeyCryptography.html">public-key cryptography</a>,
512			specifically
513			<a href="http://mathworld.wolfram.com/RSAEncryption.html">RSA public-key cryptography</a>,
514			which relies on the fact that it is much easier to create a large composite number with
515			large prime components than to factor it.
516			</li>
517			<li>This utility uses the <i>mpz_prob_prime_p()</i> function from
518			the <a href="http://www.swox.com/gmp">GMP</a>
519			to make its primality assessment.  This function has not been proofread
520			to verify how it selects bases and applies the Miller-Rabin test.  It is assumed
521			that the standard upper bound of 0.25^N can be applied where N is the <i>reps</i>
522			parameter to the function (in the case of this utility, <i>reps</i> is
523			<?php echo " " . $miller_rabin_nrounds ; ?>).
524			</li>
525			<li>This web site also contains a <a href="pfact18digit.php">factorization
526			utility</a> which will factor integers subject to size and CPU time constraints.
527			</ul>
528			<?php
529	dashley	190	$style->static_page_hrule_std();
530	dashley	23
531			miller_rabin_commanate($N);
532			do_form($N, "Is This Prime?", $max_digits);
533			}
534			}
535
536			//Now do the footer unconditionally. The footer is always used on this web page.
537			do_footer($style);
538			?>