//$Header$ // //******************************************************************************** //Copyright (C) 2003 David T. Ashley //******************************************************************************** //This program or source file is free software; you can redistribute it and/or //modify it under the terms of the GNU General Public License as published by //the Free Software Foundation; either version 2 of the License, or (at your //option) any later version. // //This program or source file is distributed in the hope that it will //be useful, but WITHOUT ANY WARRANTY; without even the implied warranty of //MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the //GNU General Public License for more details. // //You may have received a copy of the GNU General Public License //along with this program; if not, write to the Free Software //Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA //******************************************************************************** // //This module attempts to factor an integer of up to 18 decimal //digits using a sieve method and subject to elapsed time constraints. //The constraint of "up to 18 digits" //is because that way trial divisors can stay in one unsigned //integer, which speeds the division. Additionally, there isn't //much use in trying to factor larger integers on a web page. //This subfunction is "parlor trivia" grade only. // //INPUT PARAMETERS //---------------- //This subfunction accepts the following parameters, in order. // // (a) The number to factor. // // (b) The number of Miller-Rabin rounds to use when looking // at the probable primality of the original argument or of // residues. // // (c) The maximum time that should be allowed to elapse before the // program terminates (perhaps without finding factors). // This constraint is to protect the server capacity. (In the long // term, however, the web page may need to be removed--if too many // people use it, it will bring any server to a crawl.) // //OUTPUT RESULTS //-------------- //The notation below gives the outputs of the program. In some cases, [i] notation //is used to indicate line numbers. // //[01] An overall success or failure code for the operation. Valid responses // are: // S : Success (for this subfunction, the only possible outcome). // //[02] The fully normalized first integer entered. This means it has been // stripped of all weird characters, etc., and also perhaps assigned // a default value if the original parameter wasn't acceptable. // This can be used by the PHP script to repopulate form boxes. //[03] The fully normalized second integer entered. //[04] The fully normalized second integer entered. //[05] The type of factor specified on the next line. // Possibilities are: // "P" : The factor is definitely prime. // "C" : The factor is definitely composite. // "p" : The factor is probably prime, established by // Miller-Rabin. // // Note that the "C" code can only appear if the utility has to give // up because it runs out of allowed time. This means that the input // number has not been fully factored or perhaps not factored at all. // The "C" code can only appear for the last factor or the only factor. // // The "p" code can also only appear for the last factor. This occurs when // either the original input argument is prime as established by // Miller-Rabin or else a division leaves a result that is similarly // established by Miller-Rabin. //[06] The factor itself. //[07] Its multiplicity. //[..] Lines 5-7 are repeated for as many factors as are located. //[..] The second-to-last line will contain "X". //[..] The last line will contain an "S". // //Note the following: // (a) Any valid output will contain at least 9 lines. // (b) Any valid output will have a number of lines divisible by 3. // (c) The number of factors found is (nlines - 6)/3 or alternately // nlines/3 - 2. // (d) What happened can be determined using the number of lines plus // the code of the last factor. // //The return value (exit code) from this subfunction is always 0. // #define MODULE_SUBFUNC_PFACT_18 #include #include #include #include #include #include #include #include #include "auxfuncs.h" #include "subfunc_pfact_18.h" #include "sieve_eratosthenes.h" int SUBFUNC_PFACT_18_main(int argc, char *argv[]) { //The time snapshot against which we compare to see if we're over //time budget. time_t time_snapshot; //Normalized first and second parameters (the integers). char *arg1; char *arg2; char *arg3; //Number of Miller-Rabin iterations to establish probable primality. int miller_rabin_iterations; //Maximum elapsed time allowed. int max_time; //Temporary iteration integers. int i; //Mask counter. Only check for timeout periodically, as calling the OS to get time //is presumed expensive. int mask_counter; //The current trial divisor. unsigned long trial_divisor, new_trial_divisor; //The current sieve table index. int sieve_table_index; //An exit flag kept to remember if we should bail the loop. int exit_flag; //The number to factor. mpz_t number_to_factor; //The square root limit. We only need to go that far. mpz_t square_root_limit; //The multiplicity of any factors we find. int multiplicity; //Initialize all of the GMP variables. mpz_init(number_to_factor); mpz_init(square_root_limit); //Grab a time snapshot. time_snapshot = time(NULL); //There must be an acceptable number of command-line arguments. If not, //abort the progam with phony data. if (argc != 5) { printf("S\n2\n25\n2\nC\n2\n2\nX\nS\n"); return(0); } //Copy the command-line arguments to a safe place where we can manipulate them. //Leave 2 characters of space in case we assign a "0". arg1 = (char *)malloc((AUXFUNCS_size_t_max(1, strlen(argv[2])) + 1) * sizeof(char)); arg2 = (char *)malloc((AUXFUNCS_size_t_max(1, strlen(argv[3])) + 1) * sizeof(char)); arg3 = (char *)malloc((AUXFUNCS_size_t_max(1, strlen(argv[4])) + 1) * sizeof(char)); if ((arg1 == NULL) || (arg2 == NULL) || (arg3 == NULL)) { printf("S\n2\n25\n2\nC\n2\n2\nX\nS\n"); return(0); } strcpy(arg1, argv[2]); strcpy(arg2, argv[3]); strcpy(arg3, argv[4]); //Strip all of the non-digit trash out of the arguments. AUXFUNCS_remove_non_digits(arg1); AUXFUNCS_remove_non_digits(arg2); AUXFUNCS_remove_non_digits(arg3); //Remove all leading zeros from arguments. AUXFUNCS_remove_leading_zeros(arg1); AUXFUNCS_remove_leading_zeros(arg2); AUXFUNCS_remove_leading_zeros(arg3); //If an argument is zero length, fill it in with 0. if (!strlen(arg1)) strcpy(arg1, "0"); if (!strlen(arg2)) strcpy(arg2, "0"); if (!strlen(arg3)) strcpy(arg3, "0"); //We are not allowed to have 0's in this function, so //abort on zeros. Also, we can't have 1 for a number to check, //as 1 can't be factored. if ((!strcmp(arg1, "0")) || (!strcmp(arg2, "0")) || (!strcmp(arg3, "0")) || (!strcmp(arg1, "1"))) { printf("S\n2\n25\n2\nC\n2\n2\nX\nS\n"); return(0); } //If the number to factor exceeds 18 digits, abort. Anything that has //18 or fewer digits is allowed. if (strlen(arg1) >18) { printf("S\n2\n25\n2\nC\n2\n2\nX\nS\n"); return(0); } //Assign the number to factor. This is definitely a //GMP type. mpz_set_str(number_to_factor, arg1, 10); //Assign the number of Miller-Rabin repetitons to use, subject to an //absolute maximum of 1000 and minimum of 1. miller_rabin_iterations = 25; if (strlen(arg2) > 3) { miller_rabin_iterations = 1000; } else { sscanf(arg2, "%d", &miller_rabin_iterations); if (miller_rabin_iterations < 1) miller_rabin_iterations = 1; } //Assign the maximum time allowed, subject to a maximum of 1000 and minimum of 3. if (strlen(arg3) > 3) { max_time = 1000; } else { sscanf(arg3, "%d", &max_time); if (max_time < 3) max_time = 3; } //Output the header information before beginning the search. printf("S\n"); mpz_out_str(stdout, 10, number_to_factor); printf("\n"); printf("%d\n", miller_rabin_iterations); printf("%d\n", max_time); //Loop through the list of divisors we should try before starting the sieve. //This won't possibly exceed our time, so do it. for (i=0; i max_time) exit_flag = 1; //Does the Miller-Rabin test say that the remaining number //is with near perfect certainty prime? if (mpz_probab_prime_p(number_to_factor, miller_rabin_iterations) >= 1) exit_flag = 1; //Is the remaining number "1", indicating we've factored it fully? if (!mpz_cmp_ui(number_to_factor, 1)) exit_flag = 1; //Set multiplicity to zero again so don't re-enter until successful //division again. multiplicity = 0; } //Advance to our next trial divisor. new_trial_divisor = trial_divisor + SIEVE_ERATOSTHENES_sieve[sieve_table_index]; //If the new is < the old, we've rolled over. This means an exit is necessary. if (new_trial_divisor < trial_divisor) exit_flag = 1; //Advance the sieve index. sieve_table_index = (sieve_table_index + 1) % SIEVE_ERATOSTHENES_N_SIEVE; //This is our only chance to check for termination conditions that //don't come about from a successful division. But we don't do //this often. mask_counter++; if (!(mask_counter & 0xFFFFF)) { //First, have we exceeded the square root bound? if (mpz_cmp_ui(square_root_limit, trial_divisor) <= 0) exit_flag = 1; //Second, are we over the time budget? if ((time(NULL) - time_snapshot) > max_time) exit_flag = 1; } } //If we've made it out of the loop, there are a variety of reasons for that. //Find the right one and close up. if (!mpz_cmp_ui(number_to_factor, 1)) { //We divided the number down to 1. There is nothing further to do. } else if (mpz_cmp_ui(square_root_limit, trial_divisor) <= 0) { //We are at or over the square root limit. The remaining number is definitely prime. printf("P\n"); mpz_out_str(stdout, 10, number_to_factor); printf("\n"); printf("1\n"); } else if (mpz_probab_prime_p(number_to_factor, miller_rabin_iterations) == 0) { //Miller-Rabin says the remaining number is definitely composite. printf("C\n"); mpz_out_str(stdout, 10, number_to_factor); printf("\n"); printf("1\n"); } else if (mpz_probab_prime_p(number_to_factor, miller_rabin_iterations) == 1) { //Miller-Rabin says the remaining number is probably prime. printf("p\n"); mpz_out_str(stdout, 10, number_to_factor); printf("\n"); printf("1\n"); } else if (mpz_probab_prime_p(number_to_factor, miller_rabin_iterations) == 2) { //Miller-Rabin says the remaining number is definitely prime. printf("P\n"); mpz_out_str(stdout, 10, number_to_factor); printf("\n"); printf("1\n"); } //Output the invariant footer information. done: printf("X\n"); printf("S\n"); //Always return 0. return(0); } //******************************************************************************** // $Log: subfunc_pfact_18.c,v $ // Revision 1.7 2003/07/01 03:46:58 dtashley // Edits towards working continued fraction best rational approximation // functionality. // // Revision 1.6 2003/04/17 20:02:05 dtashley // License text for the GPL added. All source files are now under the GPL, // after some discussion on the GMP list. // // Revision 1.5 2003/04/16 07:22:37 dtashley // All checks completed. // // Revision 1.4 2003/04/16 07:02:06 dtashley // Spelling of Greek name corrected to Eratosthenes from incorrect Erastothenes. // // Revision 1.3 2003/04/16 06:49:21 dtashley // Edits. // // Revision 1.2 2003/04/16 03:25:15 dtashley // Seems to be working correctly. Only a careful proofreading and some // testing remain. // // Revision 1.1 2003/04/15 23:54:58 dtashley // Initial checkin. //******************************************************************************** // End of SUBFUNC_PFACT_18.C. //********************************************************************************