|
-
May 13th, 2009, 11:04 PM
#1
Miller Rabin Primality Test Wrong Result
Hello to all, i try to code miller rabin algorithm with C++ but the program doesn't did what it suppose to be.
What wrong with my program ?
Code:
#include <bitset>
#include <vector>
#include <limits>
#include <algorithm>
#include <functional>
#include <utility>
#include <ctime>
using namespace std;
const unsigned short NITE = 6;
typedef unsigned long long ulong;
void Userinput(ulong&);
bool MillerRabinPrimeTest(const ulong&);
ulong ComputeModularExponentiation(const ulong&, const ulong&, const ulong&);
ulong ComputeExponentiation(const ulong&, const ulong&);
// =========================================
int main()
{
ulong number(0);
while (1)
{
Userinput(number);
cout << MillerRabinPrimeTest(number);
}
return 0;
}
// =========================================
void Userinput(ulong& number)
{
cout << "Enter Number : ";
cin >> number;
while (cin.fail())
{
cin.clear();
cin.ignore(numeric_limits<int>::max(), '\n');
cout << "\nEnter Number : ";
cin >> number;
}
}
// =========================================
bool MillerRabinPrimeTest(const ulong& number)
{
ulong a(0), x(0), y(0), tempNumber(0), s(1);
vector<int> primeFactor(10);
bool IsPrime(false);
tempNumber = number - 1;
primeFactor[0] = 2;
primeFactor[1] = 3;
primeFactor[2] = 5;
primeFactor[3] = 7;
primeFactor[4] = 11;
primeFactor[5] = 13;
primeFactor[6] = 19;
primeFactor[7] = 23;
primeFactor[8] = 29;
primeFactor[9] = 31;
if (number > 2)
{
// Write them as 2 ^ s * d
while (tempNumber % 2 == 0)
{
// tempNumber is d
tempNumber /= 2;
}
// compute s
while (ComputeExponentiation(2, s) * tempNumber != number - 1)
{
++s;
}
for (int loop = 0;loop < NITE;++loop)
{
srand((unsigned int)time(0));
// rand() % upperBound + startnumber
a = rand() % (number - 2) + 2;
// a = primeFactor[loop];
x = ComputeModularExponentiation(a, tempNumber, number);
for (int r = 0; r < s-1;++r)
{
y = (x * x) % number;
if (x == 1 || y == -1)
{
return true;
}
}
}
}
return IsPrime;
}
// =========================================
ulong ComputeModularExponentiation(
const ulong& a, const ulong& d,
const ulong& n)
{
enum {NBITS = numeric_limits<ulong>::digits };
string bitExponent = bitset<NBITS>((unsigned long)d).to_string();
typedef string::size_type strType;
strType MSSB = bitExponent.find_first_of('1');
ulong result(a % n);
MSSB += 1;
while (MSSB < NBITS)
{
result *= result;
if (bitExponent[MSSB] == '1')
{
result = result * a % n;
}
++MSSB;
}
return result;
}
// =========================================
ulong ComputeExponentiation(const ulong& base, const ulong& exponent)
{
enum {NBITS = numeric_limits<ulong>::digits };
string bitExponent = bitset<NBITS>((unsigned long)exponent).to_string();
typedef string::size_type strType;
strType MSSB = bitExponent.find_first_of('1');
ulong result(base);
MSSB += 1;
while (MSSB < NBITS)
{
result *= result;
if (bitExponent[MSSB] == '1')
{
result *= base;
}
++MSSB;
}
return result;
}
Please explain where my program goes wrong.
Thanks.
Thanks for your help.
Posting Permissions
- You may not post new threads
- You may not post replies
- You may not post attachments
- You may not edit your posts
-
Forum Rules
|
Click Here to Expand Forum to Full Width
|