/******************************************************
素數檢測：
	Close["D:\\usbkey\\RSA\\tmp.txt"]
			Clear[i, result, v]
			For[i = 1; v = Read["D:\\usbkey\\RSA\\tmp.txt", {Number}], 
			  i < 300, i++, result =PrimeQ[v[[1]]]; Print[result]]
******************************************************/
#include "StdAfx.h"
#include ".\heck_prime.h"
#include <stdio.h>
#include <math.h>
#include ".\SmallPrime.h"

Check_Prime::Check_Prime(void)
{
}

Check_Prime::~Check_Prime(void)
{
}
void Check_Prime::generate_small_prime(unsigned int mymax)
{
	bool prime;
	unsigned int i,j,number=0;
	for(i=3;i<=mymax;i++)
	{
		prime=true;
		for(j=2;j<=sqrt((double)i);j++)
		{
			if(i%j==0)
			{
				prime=false;
				break;
			}
		}
		if(prime)
		{
			number++;
			printf(",%u",i);
		}
	}
	printf("\n%u\n",number);
	return;
}
//返回0表示是通過了小素數檢驗，其他數字是他的第一個因子
unsigned int Check_Prime::small_prime_check(unsigned int *small_prime)
{
	for(unsigned int i=1;i<small_prime[0];i++)
	{
		if(this->number.Mod((long)small_prime[i])==0)
		{
			return small_prime[i];
		}
	}
	return 0;
}
/*****************************************************************************
					Miller-Rabin素性檢驗算法
Input:一個奇整數n（n＞＝３）以及一個安全參數　t＞＝１
Output:n是合數或者可能是一個素數
１．讓Ｎ－１＝２＾Ｓ＊Ｒ　Ｒ必須為奇數
２．Ｉ＝１；Ｉ＜＝Ｔ；Ｉ＋＋
	２．１　隨機選取Ａ,（Ａ大于等于２，小于等于Ｎ－１）
	２．２　計算Ｙ＝Ａ＾Ｒmod N
	２．３　如果Ｙ不等于１且不等于Ｎ－１
			２．３．１　Ｊ＝１
			２．３．２　如果Ｊ小于等于Ｓ－１且Ｙ不等于Ｎ－１
						２．３．２．１　Ｙ＝Ｙ＾２　mod  Ｎ
						２．３．２．２　如果Ｙ不等于Ｎ－１　返回是合數
			２．３．３　如果Ｙ不等于Ｎ－１　返回是合數
３．返回可能是素數

通過測試，發現誤判的可能性非常小，已經達到了我們非常滿意的地步了

不過在馮登國和周玉潔兩位密碼學前輩的<公開密鑰密碼算法及其快速實現＞一書中
在這個地方是有問題的，我就是依據他們的算法來寫的，可是檢驗的時候發現有問題
都是我很欽佩的密碼學方面的專家

*****************************************************************************/
//返回0表示是合數，返回１表示可能是素數
unsigned int Check_Prime::prime_check(unsigned int iterator)
{
	int isprime=1;//默認為這個數位素數
	unsigned int small_check_result,weishu,i,j;
	CBigInt n_1(this->number.base),q(this->number.base),p(this->number.base),y;
	//如果數字本身為1，返回0
	if(this->number.Cmp(1)==0)
	{
		return 0;
	}
	small_check_result=this->small_prime_check(smallprime);
	//if(this->number.Cmp(small_check_result
	if(small_check_result!=0)
	{
		if(this->number.Cmp(small_check_result)==0)
			return 1;//
		else
			return 0;
	}
	//如果iterator==0說明由大數的位數載錯誤概率小于2E-80的條件下確定最優迭代次數
	if(iterator==0)
	{
		//確定被測驗的大數的位數，只能大概估計，100000000大約為27位，0x100000000為32
		if(this->number.base==100000000)
		{
			weishu=27*this->number.m_nLength;
		}
		else
		{
			weishu=32*this->number.m_nLength;
		}
		if	(weishu<	73)	iterator=37;
		else if (weishu< 105) iterator=32;
		else if (weishu< 137) iterator=25;
		else if (weishu< 197) iterator=19;
		else if (weishu< 220) iterator=15;
		else if (weishu< 235) iterator=13;
		else if (weishu< 253) iterator=12;
		else if (weishu< 275) iterator=11;
		else if (weishu< 300) iterator=10;
		else if (weishu< 332) iterator=9;
		else if (weishu< 375) iterator=8;
		else if (weishu< 433) iterator=7;
		else if (weishu< 514) iterator=6;
		else if (weishu< 638) iterator=5;
		else if (weishu< 847) iterator=4;
		else if (weishu<1275) iterator=3;
		else if (weishu< 2861) iterator=2;
		else iterator=1;
	}
	n_1.Mov(this->number);
	n_1.Mov(n_1.Sub(1));
	weishu=Fact_Two(n_1,q);//這個地方是值得重點檢查的地方

	//開始進行運算
	
	i=1;
	isprime=1;
	do
	{
		unsigned int tmp=smallprime[i++];
		p.Mov(tmp);
		y.Mov(p.Mon(q,this->number));//y=p^q mod n
		
		if((y.Cmp(1)!=0)&&(y.Cmp(n_1)!=0))//相等為０，大于為１，小于為－１；
		{
			j=1;
			//求平方，只要y不同于＋１和－１且weishu-1次迭代還沒有執行完
			while((y.Cmp(n_1)!=0)&&(j<weishu))
			{
				y.Mov(y.Mul(y));
				y.Mov(y.Mod(this->number));
				if(y.Cmp(1)==0)
				{
					isprime=0;
					break;
				}
				j++;
			}
			if(y.Cmp(n_1)!=0)
			{
				isprime=0;
			}
		}
	}while(((--iterator)>0)&&isprime);
	return isprime;
}

//分解n-1=2的k次方*q;
unsigned int Check_Prime::Fact_Two(CBigInt n,CBigInt& q)
{
	unsigned int k=0;
	while(n.Is_Even())
	{
		n.Mov(n.Div(2));
		k++;
	}
	q.Mov(n);
	return k;
};