vlong.cpp

椭圆曲线密码C实现的

#include "stdafx.h"
#include "vlong.h"

#define PTLSIZE 550
const static int pt[PTLSIZE] =
{   3,    5,    7,    11,   13,   17,   19,   23,   29,   31,
    37,   41,   43,   47,   53,   59,   61,   67,   71,   73,
	79,   83,   89,   97,   101,  103,  107,  109,  113,  127, 
	131,  137,  139,  149,  151,  157,  163,  167,  173,  179, 
	181,  191,  193,  197,  199,  211,  223,  227,  229,  233, 
	239,  241,  251,  257,  263,  269,  271,  277,  281,  283, 
	293,  307,  311,  313,  317,  331,  337,  347,  349,  353, 
	359,  367,  373,  379,  383,  389,  397,  401,  409,  419, 
	421,  431,  433,  439,  443,  449,  457,  461,  463,  467, 
	479,  487,  491,  499,  503,  509,  521,  523,  541,  547, 
	557,  563,  569,  571,  577,  587,  593,  599,  601,  607, 
	613,  617,  619,  631,  641,  643,  647,  653,  659,  661, 
	673,  677,  683,  691,  701,  709,  719,  727,  733,  739, 
	743,  751,  757,  761,  769,  773,  787,  797,  809,  811, 
	821,  823,  827,  829,  839,  853,  857,  859,  863,  877,
	881,  883,  887,  907,  911,  919,  929,  937,  941,  947, 
	953,  967,  971,  977,  983,  991,  997,  1009, 1013, 1019, 
	1021, 1031, 1033, 1039, 1049, 1051, 1061, 1063, 1069, 1087,
	1091, 1093, 1097, 1103, 1109, 1117, 1123, 1129, 1151, 1153, 
	1163, 1171, 1181, 1187, 1193, 1201, 1213, 1217, 1223, 1229, 
	1231, 1237, 1249, 1259, 1277, 1279, 1283, 1289, 1291, 1297, 
	1301, 1303, 1307, 1319, 1321, 1327, 1361, 1367, 1373, 1381,
	1399, 1409, 1423, 1427, 1429, 1433, 1439, 1447, 1451, 1453, 
	1459, 1471, 1481, 1483, 1487, 1489, 1493, 1499, 1511, 1523,
	1531, 1543, 1549, 1553, 1559, 1567, 1571, 1579, 1583, 1597, 
	1601, 1607, 1609, 1613, 1619, 1621, 1627, 1637, 1657, 1663, 
	1667, 1669, 1693, 1697, 1699, 1709, 1721, 1723, 1733, 1741, 
	1747, 1753, 1759, 1777, 1783, 1787, 1789, 1801, 1811, 1823, 
	1831, 1847, 1861, 1867, 1871, 1873, 1877, 1879, 1889, 1901, 
	1907, 1913, 1931, 1933, 1949, 1951, 1973, 1979, 1987, 1993, 
	1997, 1999, 2003, 2011, 2017, 2027, 2029, 2039, 2053, 2063,
	2069, 2081, 2083, 2087, 2089, 2099, 2111, 2113, 2129, 2131, 
	2137, 2141, 2143, 2153, 2161, 2179, 2203, 2207, 2213, 2221, 
	2237, 2239, 2243, 2251, 2267, 2269, 2273, 2281, 2287, 2293,
	2297, 2309, 2311, 2333, 2339, 2341, 2347, 2351, 2357, 2371,
	2377, 2381, 2383, 2389, 2393, 2399, 2411, 2417, 2423, 2437, 
	2441, 2447, 2459, 2467, 2473, 2477, 2503, 2521, 2531, 2539, 
	2543, 2549, 2551, 2557, 2579, 2591, 2593, 2609, 2617, 2621, 
	2633, 2647, 2657, 2659, 2663, 2671, 2677, 2683, 2687, 2689, 
	2693, 2699, 2707, 2711, 2713, 2719, 2729, 2731, 2741, 2749, 
	2753, 2767, 2777, 2789, 2791, 2797, 2801, 2803, 2819, 2833, 
	2837, 2843, 2851, 2857, 2861, 2879, 2887, 2897, 2903, 2909,
	2917, 2927, 2939, 2953, 2957, 2963, 2969, 2971, 2999, 3001,
	3011, 3019, 3023, 3037, 3041, 3049, 3061, 3067, 3079, 3083,
	3089, 3109, 3119, 3121, 3137, 3163, 3167, 3169, 3181, 3187, 
	3191, 3203, 3209, 3217, 3221, 3229, 3251, 3253, 3257, 3259, 
	3271, 3299, 3301, 3307, 3313, 3319, 3323, 3329, 3331, 3343,
	3347, 3359, 3361, 3371, 3373, 3389, 3391, 3407, 3413, 3433, 
	3449, 3457, 3461, 3463, 3467, 3469, 3491, 3499, 3511, 3517, 
	3527, 3529, 3533, 3539, 3541, 3547, 3557, 3559, 3571, 3581,
    3583, 3593, 3607, 3613, 3617, 3623, 3631, 3637, 3643, 3659, 
	3671, 3673, 3677, 3691, 3697, 3701, 3709, 3719, 3727, 3733, 
	3739, 3761, 3767, 3769, 3779, 3793, 3797, 3803, 3821, 3823, 
	3833, 3847, 3851, 3853, 3863, 3877, 3881, 3889, 3907, 3911, 
	3917, 3919, 3923, 3929, 3931, 3943, 3947, 3967, 3989, 4001
};
//利用费儿马定理证明一个素数
static int is_probable_prime( const vlong &p )
{
	// Test based on Fermats theorem a**(p-1) = 1 mod p for prime p
	// For 1000 bit numbers this can take quite a while
	const rep = 4;	
	vlong mod;
	const unsigned any[rep] = { 2,3,5,7 };
	for ( unsigned i=0; i<rep; i+=1 )
	{
		::modexp(any[i],p-1,p,mod);
		if( mod != vlong(1) )
			return 0;
	}
	return 1;
}

vlong::vlong ( int	x )
{
	if(x<0){	
		m_nNegative = 1;
	}else{
		m_nNegative = 0;
	}
	m_pValue= new uvlong(x);
	m_clint = m_pValue->m_clint;

}

vlong::vlong ( const vlong& x ) // copy constructor
{
  m_nNegative	=	x.m_nNegative;
  m_pValue		=	x.m_pValue;
  m_pValue->m_share += 1;
  m_clint = m_pValue->m_clint;
}

vlong::~vlong()
{
	if(m_pValue->m_share) { m_pValue->m_share -= 1;}else{	delete m_pValue;}

}

vlong& vlong::operator =(const vlong& x)
{
	copy_cl(m_clint,x.m_clint);
	m_nNegative = x.m_nNegative;
	return *this;
}

vlong& vlong::operator +=(const vlong& x)
{
	if ( m_nNegative == x.m_nNegative )
	{
		add_cl(m_clint, x.m_clint, m_clint );
	}else{ 	
		if(comp_cl(m_clint,x.m_clint) >= 0)
		{
			sub_cl( m_clint, x.m_clint, m_clint );
		}else{
			ULINT	tmp;
			copy_cl( tmp, x.m_clint);
			sub_cl( tmp, m_clint, m_clint );
			m_nNegative = x.m_nNegative; 
		}
	}
	
	return *this;
}

vlong& vlong::operator -=(const vlong& x)
{
	if ( m_nNegative != x.m_nNegative )
	{
		add_cl(m_clint, x.m_clint, m_clint );
	}else{
		if ( comp_cl(m_clint,x.m_clint) >= 0 )
		{
			sub_cl( m_clint, x.m_clint, m_clint );
		}else{
			ULINT	tmp;
			copy_cl( tmp, x.m_clint);
			sub_cl( tmp, m_clint, m_clint );
			m_nNegative = 1 - m_nNegative; 			
		}
	}
	return *this;
}

vlong	operator +( const vlong& x, const vlong& y )
{
	vlong	result = x;
	result	+= y;
	return result;
}


vlong	operator -( const vlong& x, const vlong& y )
{
	vlong	result = x;
	result	-= y;
	return result;
}

vlong	operator *( const vlong& x, const vlong& y )
{
	vlong result;
	int a = x.m_pValue->GetBits();
	int b = y.m_pValue->GetBits();
	
	mul_cl(x.m_clint,y.m_clint,result.m_clint,a+b);
	result.m_nNegative = x.m_nNegative ^ y.m_nNegative;
	return result;
}

vlong	operator /( const vlong& x, const vlong& y )
{
	vlong result;
	ULINT rem; 

	div_cl( x.m_clint, y.m_clint, result.m_clint, rem);
	result.m_nNegative = x.m_nNegative ^ y.m_nNegative;
	return result;
}


vlong	operator %( const vlong& x, const vlong& y )
{
	vlong	result;
	ULINT	divid;
	if(y.m_nNegative == 0)
	{
		div_cl( x.m_clint, y.m_clint, divid,result.m_clint );
		result.m_nNegative = x.m_nNegative; // not sure about this?
	}
	
	return result;
}

vlong::operator unsigned()
{
	int x= ( m_clint[0] > 2) ? 2 : m_clint[0];
	
	int n=0;
	for(int i = 1 ; i <= x ; i++ )
	{
		n += (m_clint[i]<<(BPU*(i-1)));
	}
	
	if(m_nNegative>0)
		n = -n;
	
	return n;
}

int comp_vl( const vlong  a,const vlong   b)
{
	if(a.m_nNegative < b.m_nNegative) { return +1; }
	if(a.m_nNegative > b.m_nNegative) { return -1; }
	int x = comp_cl(a.m_clint, b.m_clint );
	if(x == 0){ return 0; }
	
	if(( a.m_nNegative > 0  && x > 0 ) || ( a.m_nNegative == 0 && x < 0 )){
		return -1;
	}else{
		return +1;
	}
}

bool  modinv1(const vlong &e, const vlong &m, vlong &u, vlong &v )
{
	//u*a+v*b=(a,b)
	vlong	a,b;
	vlong	s,t,zero,g,w;//u,v,
	vlong	q,r;
	
	zero = ushort(0);
	
	if(m>e)
	{
		a = m;
		b = e;
	}else{
		if(m<e){
			a = e;
			b = m;
		}else{
			return false;
		}
	}
	
	s = 1; u = 0;
	t = 0; v = 1;
	do{
		//a = q*b+r
		if(b == vlong(0) )
			return false;
		q = a/b;
		r = a%b;
		
        if(r !=vlong(0) )
		{
			a = b;		b = r;
			w = u;		u = s-q*u; 
			s = w;		w = v;
			v = t-q*v;	t = w;
		}else{
			g = b;
			break;
		}
	}while(1);
	
	if(u<vlong(0)) { u = (m>e) ? u+e : u+m; }
	if(v<vlong(0)) { v = (m>e) ? v+m : v+e; }
	
	
	return true;
}

bool modinv2(const vlong &e, const vlong &m, vlong &d)
{
	vlong j,b,c,x,y;
	c.m_pValue->gcd(*e.m_pValue,*m.m_pValue) ;
	if( c != vlong(1))
		return false;

	b=m;j=1;c=e;d = 0;
	
	vlong t(0);
	while ( c != t )
	{
		if(c == vlong(0))
			return false;
		x = b / c;
		y = b - x*c;
		b = c;
		c = y;
		y = j;
		j = d - j*x;
		d = y;
	}
	if ( d < t )
		d += m;

	if( ((e*d)%m)!= vlong(1))//(e*d)%m
		return false;

	return true;
}

//产生随机大数
bool    randval(vlong &out,unsigned int bitlen)//len 输入随即数得二进制位数
{
	if(bitlen <= 0 )
		return false;

	out = 0;
	srand(GetTickCount());
	
	unsigned  len	= bitlen/16;
	unsigned  j		= bitlen%16;  
	
	for(unsigned i = 1; i <=  len; i ++ )
	{
		out.m_clint[i] = rand()%0x010000;
		out.m_clint[0] = i;
	}
	
	if(j == 0)	
	{
		out.m_clint[i-1] |= 1 << 0xf;
	}else{	
		out.m_clint[i] = rand()%0x010000;
		out.m_clint[i] &= (0xffff>>(16-j));
		out.m_clint[i] |= 1 << (j-1);
		out.m_clint[0] += 1;  
	}
	out.m_clint[1] |= 0x01;

	return true;
}
/*
bool  findprime(vlong &p,unsigned int bitlen)
{
//	if( bitlen<16 || bitlen>(CLINTMAXDIGIT*BPU)/4)
//		MessageBox(NULL,"素数长度不在合法范围之内!","错误" ,MB_ICONERROR | MB_OK);
	
	vlong	p_1,p_1_div_2,lp,a,tmp;
	int		i,k=0;
	int		end = PTLSIZE;
	// 产生一个随机大数(奇数) p

Next:	++k;
		if(!randval(p,bitlen))
			return false;

		// 再用 Lehmann 方法测试
		p_1	 =	p; 
		p_1	-=	1;
		p_1_div_2 = p_1/vlong(2) ;

		for( i = 0; i < 15; ++i )// 产生一个小随机大数 A
		{   
			while(1)
			{
				a = (rand() & 0x0fffffff);
				if(a < p_1)
					break;
			}
			tmp.modexp(a,p_1_div_2,p);			
			//思考：当p为素数时，tmp的值为1或p_1??
			if(( tmp != vlong(1) ) && (tmp != p_1))
			{   
				goto Next;
			}	
		}

		return true;

}
*/
bool  findprime(vlong &p,unsigned int bitlen)
{
	int i,SS = PTLSIZE*5,tested = 0; 
	unsigned char * b = new unsigned char[SS];
	vlong pl,r;
	
	if(!randval(p,bitlen))
		return false;
    while (1)
	{
		for (i=0;i<SS;i+=1)
			b[i] = 1;
		for (i=0;i<PTLSIZE;i+=1)
		{
			pl = pt[i];
			r = p % vlong(pl); 
			
			if (r) r = pl - r;
			
			while ( r < vlong(SS) )
			{
				b[r] = 0;
				r += pl;
			}
		}
		
		for (i=0;i<SS;i+=1)
		{
			if ( b[i] )
			{
				tested += 1;
				if ( is_probable_prime(p) )
				{
					goto end;
				}
			}
			p += 1;
		}
	}
end :
	delete []b;
	return true;

}

bool modexp( const vlong & m, const vlong & e, const vlong & mod , vlong & result)
{
	if(result.m_pValue != NULL)
		result.m_pValue->modexp(*m.m_pValue, *e.m_pValue, *mod.m_pValue);

	return true;
}