Thread: Math Programming in C++

  1. #1
    Registered User
    Join Date
    Oct 2005
    Posts
    1

    Math Programming in C++

    have been doing a little math programming over the past few weeks for science fair and my object is to make an efficient program that can find very large mersenne primes, I am currently using Lucas Lehmer's formula and a small statement that finds the first 500,000 prime number and divides the Mersenne prime by those numbers every time to help speed it up. My goal of this program is to make a type of filter that take out uneeded composite number before lucas takes over since it is kind of slow compared to some probablistic tests. One test I am really wanting to integrate is the Rabin-Miller strong pseudoprime test but I am unable to comperehend it, can someone please put it into C++ or a way in which I can understand it? Any suggestions?

    Code:
    /*
      Name: MPrime Finder
      Copyright: 2005 
      Author: Devon Taylor
      Date: 13/10/05 22:07
      Description: Uses high precision gmp library and lucas lehmers 
      determinant formula to find extremely large Mersenne Primes.
      */
    #include <windows.h>
    #include <fstream> 
    #include <iostream> 
    #include <stdlib.h>
    #include <math.h>
    #include <stdio.h>
    #include <gmp.h>                              //Import needed includes
    using namespace std;
    int main(int argc, char* argv[]){                        //Initiate program
    cout<<"Program is making preparations, one moment...\n";
    long int n_check;
    long int co=0;
    long int prime_t;
    long int test_primes[500001];
    mpz_t d1;
    mpz_init_set_si(d1,1);
    long int count1=0;
    long int d11; 
    while(count1<500000){
                         mpz_nextprime(d1,d1);
                         d11 =mpz_get_si(d1); 
                         test_primes[count1]=d11;
                         count1=count1+1;
    }
    long int input_low;
    long int input_high;                                     //Initiate Long Integer variables
    long int color;
            color = atoi(argv[1]);                           //Set color                          
            HANDLE console;                                  //Initiate console handle     
            console = GetStdHandle(STD_OUTPUT_HANDLE);
    SetConsoleTextAttribute(console, 3);
    cout<<"Program Created By Devon Taylor in 2005 \n";
    SetConsoleTextAttribute(console, 6);
    cout<<"Enter an odd # > 3 and < ten mill to start calc. at: ";  //Ask for Lower Boundry
    SetConsoleTextAttribute(console, 2);
    cin>>input_low;
    SetConsoleTextAttribute(console, 6);
    cout<<"Enter an odd # > 3 and < ten mill to stop calc. at: ";   //Ask for High Boundry
    SetConsoleTextAttribute(console, 2);
    cin>>input_high;
    SetConsoleTextAttribute(console, 8);
    cout<<"Calculating.........";
    FILE *fp; 
    ofstream a_file ( "Mersenne_primes.txt", ios::app );         //Create stream for a file
                                      fp = fopen("Mersenne_primes.txt", "a");                      //Open the file
                                      //Create file pointer
    mpz_t test_prime2;
    mpz_t p;                                    //Initiate High Precision Variables  
    mpz_t b;                              
    mpz_t c;                             
    mpz_t dd;                                 
    mpz_t x; 
    mpz_t two;
    mpz_t r;  
    mpz_t k;  
    mpz_t w;
    mpz_t y;
    mpz_t m;
    mpz_t s;
    mpz_t kr;
    mpz_t test;
    mpz_t one;
    mpz_t zero;
    mpz_t limit;
    mpz_t count;
    mpz_init(limit);
    mpz_init(kr);
    mpz_init2(m,10000000);                                    
    mpz_init2(y,10000000);
    mpz_init(test_prime2); 
    mpz_init2(s,1000000000);
    mpz_init(c);
    mpz_init(r);
    mpz_init_set_si(zero,0);
    mpz_init_set_si(one,1);
    mpz_set_si(s,4);
    mpz_init_set_ui(dd,input_high);
    mpz_init_set_si(p,input_low); 
    mpz_init_set_si(two,2);
    mpz_init_set_ui(c,3);
    mpz_init_set_si(w,2);
    mpz_init_set_si(k,3);
    mpz_init_set_si(b,3); 
    mpz_init_set_si(x,3);
    mpz_init_set_si(test,3);
    mpz_init_set_ui(count,1);
    int a=2;                                    //Initiate Variables                                  //|              //|
    int mt;                                   
    int t;  
    int z;
    int s_limit=mpz_cmp(p,dd);
    int div_three;
    int div_five;
    int div_seven;
    int h;
    int l;
    int cvd_prm; 
    int cvd_p;
    int is_prime;
    int s_check;
    while (s_limit<0){                          //While the program stays within its boundries
                     mpz_add_ui(p,p,2);                          //Add 2 to the Mersenne Prime Exponent
                     s_limit=mpz_cmp(p,dd);                      //Check to see if number is within boundries
                     t=0;                                        //Set exit variable to 0
                     mpz_set_si(b,3);                            //Set b to 3
                     mpz_sqrt(c,p);                          //Divide the Prime Exponent by 3
    while(t==0){                                //While exit statement is not true
               z=mpz_divisible_p(p,b);                     //Check if Prime Exponent is divisible by b
               mpz_add_ui(b,b,2);                          //If not divisible add 2 to b and try again
               mt=mpz_cmp(c,b);                            //See if b becomes greater than c
          if (z>0){                                   //If p is ever divisible, break while statement
                  t=1;                                        //Set exit statement to true
                  } 
          if(mt<0){                                   //If the Prime Exponent is proven to be prime 
            prime_t =mpz_get_si(p); 
            mpz_set_si(x,3);                    
            l =mpz_get_si(p);                  
            mpz_mul_2exp(y,one,l);       
            mpz_sub_ui(y,y,1);           
            mpz_set_si(c,500000);  
            co=0;
              while(t==0  && co<500000){           //Test for divisibility of prime #'s under 500,000        
                    mpz_set_ui(x,test_primes[co]);
                    co=co+1;
                    z=mpz_divisible_p(y,x);      
                    h=mpz_cmp(c,x);        
                    if(z>0){               
                           t=1;                   
                           SetConsoleTextAttribute(console, 5);
                           cout<<"Testing- "<<l<<", ";  
                    } 
                    if(h<0){                  
                           SetConsoleTextAttribute(console, 5);
                           cout<<" Testing ";
                           gmp_printf("%Zd",p);                              //State what's being tested
                           cout<<", ";   
                           cvd_prm =mpz_get_si(p);                           //Convert Prime Exponent to an integer
                           mpz_mul_2exp(m,one,cvd_prm);                        //Put 2 to the power of the Prime Exponent
                           mpz_sub_ui(m,m,1);                                //Subtract 1 from the Mersenne Number from above
                           mpz_set_si(s,4);                                  //Set s to 4
                           mpz_sub_ui(limit,p,1);                            //Set limit for Lucas test
                           mpz_set_ui(count,1);                              //Start process count at 1
                           is_prime=mpz_cmp(limit,count);                // Check if count is over limit
            while(t==0){                                  //While exit statement is false                                
                       mpz_mul(s,s,s);                   //Multiply s to itself
                       mpz_sub_ui(s,s,2);                //Subtract 2 from s
                       mpz_mod(s,s,m);                   //Get modulous from s and m
                       mpz_add_ui(count,count,1);        //Add 1 to the count
                       is_prime=mpz_cmp(limit,count);    //Check if count exceeds the limit
                       s_check=mpz_cmp(s,zero);          //Compare s to zero        
                       if(is_prime<0){                   //If the limit is reached its not prime, exit
                                     t=1;     
                                     }
                
                       if (s_check==0){                  //If s is equal to zero then the number is a Mersenne Prime
                                      cvd_p =mpz_get_si(p);           //Convert p to cvd_p
                                      mpz_mul_2exp(r,one,cvd_p);          //Store the mersenne prime 
                                      mpz_sub_ui(r,r,cvd_p);          //Subtract 1         
    
                                      mpz_out_str(fp,0,r);                                 //Write the Prime to the file
                                      a_file<<" & ";
                                      a_file<<cvd_p;
                                      a_file<<", ";
    
                                      SetConsoleTextAttribute(console, 7);
                                      gmp_printf("%Zd",p);                                //Print the prime_exp to the screen
                                      cout<<" is a Mersenne Prime, ";
                                      t=1;                                               //Exit to test a new #
                                      }
                  
               }
    }
    }
    }
    }
    }
    fclose(fp); 
    cout<<".....Finished Calculating......";
    system("pause");
    }
    Last edited by spiderman45144; 10-20-2005 at 08:44 AM.

Popular pages Recent additions subscribe to a feed

Similar Threads

  1. Math
    By knightjp in forum A Brief History of Cprogramming.com
    Replies: 16
    Last Post: 04-01-2009, 05:36 PM
  2. Help with C++ Math
    By aonic in forum C++ Programming
    Replies: 4
    Last Post: 01-29-2005, 04:40 AM
  3. Basic Math Problem. Undefined Math Functions
    By gsoft in forum C Programming
    Replies: 1
    Last Post: 12-28-2004, 03:14 AM
  4. Math Header?
    By Rune Hunter in forum C++ Programming
    Replies: 26
    Last Post: 09-17-2004, 06:39 AM
  5. toughest math course
    By axon in forum A Brief History of Cprogramming.com
    Replies: 12
    Last Post: 10-28-2003, 10:06 PM