prime.cpp

用概率算法实现求素数10000以内的素数的问题

// prime.cpp : 定义控制台应用程序的入口点。
//

#include "stdafx.h"

#include <stdio.h>
#include <stdlib.h>
#include <math.h>
#include <time.h>

int primeMb[5000];

//产生2~n-2的一个随机数
int uniform(int MIN,int MAX)
{
	int random = (rand()%(MAX-MIN+1))+MIN;  // 
	
	return random;
}

// a^t mod n
int ModuleExp(int a, int i, int n)
{
	int s = 1;

	while( i > 0)
	{
		if(i%2 == 1)  s = s*a % n;
		a = a*a %n;
		i = i / 2;
	}

	return s;
}


//返回真说明n为强伪素数
bool Btest(int a, int n)
{
	int i,x;
	int s = 0;
	int t = n-1;

	do{  // n-1= 2^s * t
		s++;
		t = t/2;
	}while(t%2 == 1);

	x = ModuleExp(a,t,n);

	if( x==1 || x==n-1) return true;
	
	for(i=1; i < s; i++)
	{
		x = x*x % n;
		if( x==n-1) return true;
	}

	return false;
}


// 返回真时表示强伪素数，返回假表示合数
bool MillRab(int n)
{
	int a;

	a = uniform(2,n-2);
	return Btest(a,n);

}//该算法是0.75正确的，偏假的


//重复调用MillRab k次，减少错误的概率
bool ReapeatMillRab(int n,int k)
{
	
	for(int i = 0; i <= k; i++)
	{
		if( !MillRab(n) )
			return false;
	}
	return true;

}

//main
int _tmain(int argc, _TCHAR* argv[])
{
	int n;

	primeMb[1] = 2;
	primeMb[2] = 3;
	primeMb[0] = 2; //count
    n = 5;

	srand((unsigned)time(NULL));

	do{
		if(ReapeatMillRab(n,(int)log(n)))
		{
			primeMb[0]++;
			primeMb[primeMb[0]] = n;
		}
		n = n + 2;
		
	}while(n < 10000);

	for(int k = 1; k <= primeMb[0]; k++)
	{
		printf("%d ",primeMb[k]);
	}

	return 0;
}