RSA算法实现与蒙哥马利算法(转)

,需动用价值数千万美金的大型计算机系统并耗费一年的时间。 RSA 密钥的选取和加解密过程都非常简洁,在算法上主要要实现四个问题: 1、如何处理大数运算 2、如何求解同余方程XY % M = 1 3、如何快速进行模幂运算 4、如何获取大素数 实际上,在实现RSA 算法的过程中大家会发现后三个问题 不是各自独立的,它们互有关联,环环相套,相信届时你会意识到:RSA算法是一种“优美”的算法! 大数存储: RSA 依赖大数运算,目前主流RSA 算法都建立在1024位的大数运算之上。而大多数的编译器只能支持到64位的整数运算,即我们在运算中所使用的整数必须小于等于64位,即:0xffffffffffffffff,也就是15,这远远达不到RSA 的需要,于是需要专门建立大数运算库来解决这一问题。 最简单的办法是将大数当作数组进行处理,数组的各元素也就是大数每一位上的数字,通常采用最容易理解的十进制数字0—9。然后对“数字数组”编写加减乘除函数。但是这样做效率很低,因为二进制为1024位的大数在十进制下也有三百多位,对于任何一种运算,都需要在两个有数百个元素的数组空间上多次重循环,还要许多额外的空间存放计算的进退位标志及中间结果。另外,对于某些特殊的运算而言,采用二进制会使计算过程大大简化,而这种大数表示方法转化成二进制显然非常麻烦,所以在某些实例中则干脆采用了二进制数组的方法来记录大数,当然这样效率就更低了。 一个有效的改进方法是将大数表示为一个n 进制数组,对于目前的32位系统而言n 可以取值为2 的32次方,即 0x100000000,假如将一个二进制为1024位的大数 转化成0x10000000进制,就变成了32位,而每一位的取值范围不再是二进制的0—1或十进制的0—9,而是0-0xffffffff,我们正好可以用一个32位的DWORD (如:无符号长整数, unsigned long)类型来表示该值。所以1024位的大数就变成一个含有32个元素的DWORD数组,而针对DWORD 数组进行各种运算所需的循环规模至多32次而已。而且 0x100000000 进制与二进制,对于计算机来说,几乎是一回事,转换非常容易。 例如大数15,等于0xffffffff ffffffff,就相当于十进制的99:有两位,每位都是0xffffffff。而16等于0x00000001 00000000 00000000,就相当于十进制的100:有三位,第一位是1 ,其它两位都是0 ,如此等等。在实际应用中,“数字数组”的排列顺序采用低位在前高位在后的方式,这样,大数A 就可以方便地用数学表达式来表示其值: A=Sum[i=0 to n](A[i]*r**i),r=0x100000000,0<=A=B: A=Sum[i=0 to p](A[i]*r**i) B=Sum[i=0 to q](B[i]*r**i) C=Sum[i=0 to n](C[i]*r**i) r=0x100000000(32位),0<=A[i],B[i],C[i]=q(A[i],B[i],C[i]都是32位的数) 则当C=A+B、C=A-B、C=A*B时,我们都可以通过计算出C来获得C: 1.加法 C=A+B,显然C[i]不总是等于A[i]+B[i],因为A[i]+B[i]可能 >0xffffffff,而C[i]必须<=0xffffffff,这时就需要进位,当然在计算C[i-1]时也可能产生了进位,所以计算C[i]时还要加上上次的进位值。如果用一个64位变量result来记录和(64位是为乘法准备的,实际加减法只要33位即可),另一个32位变量carry来记录进位(为什么要32位?为乘法准备的,实际加减法进位只有1),则有: carry=0; for(i=0;i<=p;i++) {//i从0到p 因为A>B result=A[i]+B[i]+carry; C[i]=result%0x100000000 ; //从这里看result应该大于64位,至少65位 carry=result/0x100000000; } if(carry=0) n=p; else n=p+1; 2.减法 C=A-B,同理C[i]不总是等于A[i]-B[i],因为A[i]-B[i]可能<0,而C[i]必须>=0,这时就需要借位,同样在计算C[i-1]时也可能产生了借位,所以计算C时还要减去上次的借位值: carry=0 for(i=0;i<=p;i++) { //i从0到p 因为A>B if((A[i]-B[i]-carry)>=0){ C[i]=A[i]-B[i]-carry; carry=0; } else{ C[i]=0x100000000+A[i]-B[i]-carry; carry=1; } } n=p; while (C[n]==0) n=n-1; //将前边的0去掉 3.乘法 C=A*B,首先我们需要观察日常做乘法所用的“竖式计算”过程: A3 A2 A1 A0 * B2 B1 B0 ------------------------------------------ = A3B0 A2B0 A1B0 A0B0 + A3B1 A2B1 A1B1 A0B1 + A3B2 A2B2 A1B2 A0B2 ------------------------------------------ = C5 C4 C3 C2 C1 C0 可以归纳出:C[i]=Sum[j=0 to q](A[i-j]*B[j]) (注意是C[i]), 其中i-j必须>=0且<=p。 当然这一结论没有考虑进位,虽然计算A[i-j]*B[j]和Sum的时候都可能发生进位,显然这两种原因产生的进位可以累加成一个进位值。最终可用如下算法完成乘法: C = Sum[i= 0 to n](C[i]*r**i) = Sum[i= 0 to n] ( Sum[j=0 to q](A[i-j]*B[j]) *r**i).(这里的n=p+q-1,但当第n位的运算有进位时n应加1) C也可以表示成C= Sum[i= 0 to q](A*B[i] *r**i) n=p+q-1 carry=0 for(i=0;i<=n;i++){ result=carry; for(j=0;j<=q;j++){ if (0<=i-j<=p ){ result=result+A[i-j]*B[j]; C[i]=result%0x100000000; carry=result/0x100000000; } } } if(carry!=0) { n=n+1; C[n]=carry } 4.除法 对于C=A/B,由于无法将B 对A“试商”,我们只能转换成B[q]对A[p]的试商来得到一个近似值,所以无法直接通过计算C来获得C,只能一步步逼近C。由于: B*(A[p]/B[q]-1)*0x100000000**(p-q) (b)11 x - 5 y = 1 11%5 =1 -> (c)x - 5 y = 1 令y=0 代入(c)得x=1 令x=1 代入(b)得y=2 令y=2 代入(a)得x=9 同理可使用递归算法求得任意ax-by=1(a、b互质)的解。实际上通过分析归纳将递归算法转换成非递归算法就变成了大衍求一术: x=0,y=1 WHILE a!=0 i=y y=x-y*a/b x=i i=a a=b%a b=i IF x<0 x=x+b RETURN x 模幂运算 模幂运算是RSA 的核心算法,最直接地

决定

郑伟家庭教育讲座全集个人独资股东决定成立安全领导小组关于成立临时党支部关于注销分公司决定

了RSA 算法的性能。针对快速模幂运算这一课题,西方现代数学家提出了大量的解决

方案

气瓶 现场处置方案 .pdf气瓶 现场处置方案 .doc见习基地管理方案.doc关于群访事件的化解方案建筑工地扬尘治理专项方案下载

,通常都是先将幂模运算转化为乘模运算。例如求D=C**15 % N,由于:a*b % n = (a % n)*(b % n) % n,所以: C1 =C*C % N =C**2 % N C2 =C1*C % N =C**3 % N C3 =C2*C2 % N =C**6 % N C4 =C3*C % N =C**7 % N C5 =C4*C4 % N =C**14 % N C6 =C5*C % N =C**15 % N 即:对于E=15的幂模运算可分解为6 个乘模运算,归纳分析以上方法可以发现对于任意E,都可采用以下算法计算 D=C**E % N: D=1 WHILE E>=0 IF E%2=0 C=C*C % N E=E/2 ELSE D=D*C % N E=E-1 RETURN D 继续分析会发现,要知道E 何时能整除2,并不需要反复进行减一或除二的操作,只需验证E 的二进制各位是0 还是1 就可以了,从左至右或从右至左验证都可以,从左至右会更简洁,设E=Sum[i=0 to n](E*2**i),0<=E<=1,则: D=1 FOR i=n TO 0 D=D*D % N IF E[i]=1 D=D*C % N RETURN D 这样,模幂运算就转化成了一系列的模乘运算。 模乘运算 对于乘模运算A*B%N,如果A、B都是1024位的大数,先计算A*B,再% N,就会产生2048位的中间结果,如果不采用动态内存分配技术就必须将大数定义中的数组空间增加 一倍,这样会造成大量的浪费,因为在绝大多数情况下不会 用到那额外的一倍空间,而采用动态内存分配技术会使大数存储失去连续性而使运算过程中的循 环操作变得非常繁琐。所以模乘运算的首要原则就是要避免直接计算A*B。 设A=Sum[i=0 to k](A[i]*r**i),r=0x10000000,0<=A=N C'=C'-N RETURN C' 由于在RSA中A、B总是小于N,又0<=A,C'[0]<=1,所以: C' = (C'+A*B+C'[0]*N)/2 C' < (C'+2N)/2 2C' < C'+2N C' < 2N 既然C'总是小于2N,所以求C' %N 就可以很简单地在结束循环后用一次减法来完成,即在求A*B*2**(-k) %N的过程中不用反复求模,达到了我们避免做除法的目 的。当然,这一算法求得的是A*B*2**(-k) %N,而不是我们最初需要的A*B %N。但是利用A*B*2**(-k)我们同样可以求得A**E %N。 设R=2**k %N,R'=2**(-k) %N,E=Sum[i=0 to n](E[i]*2**i): A'=A*R %N //这一步是怎么求的? X=A' FOR i FROM n TO 0 X=X*X*R' %N IF E[i]=1 X=X*A'*R' %N X=X*1*R' %N RETURN X 最初: X = A*R %N, 开始循环时: X = X*X*R' %N = A*R*A*R*R' %N = A**2*R %N 反复循环之后: X = A**E*R %N 最后: X = X*1*R' %N = A**E*R*R' %N = A**E %N 如此,我们最终实现了不含除法的模幂算法,这就是著名的蒙哥马利算法,而X*Y*R' %N 则被称为“蒙哥马利模乘”。以上讨论的是蒙哥马利模乘最简单,最容 易理解的二进制形式。蒙哥马利算法的核心思想在于将求 A*B %N转化为不需要反复取模的A*B*R' %N,(移位即可,因为R是2^K,总之R是与进制相关的数),但是利用二进制算法求1024位的A*B*R' %N,需要循环1024次之多,我么必然希望找到更有效的计算A*B*R' %N的算法。 考虑将A表示为任意的r进制: A = Sum[i=0 to k](A*r**i) 0<=A<=r 我们需要得到的蒙哥马利乘积为: C'= A*B*R' %N R'=r**(-k) 则以下算法只能得到C'的近似值 C'=0 FOR i FROM 0 TO k C'=C'+A*B C'=C'/r IF C'>=N C'=C'-N RETURN C' 因为在循环中每次C'=C'/r 时,都可能有余数被舍弃。假如我们能够找到一个系数q,使得(C' + A*B + q*N) %r =0,并将算法修改为: C'=0 FOR i FROM 0 TO k C'=C'+A*B+q*N C'=C'/r IF C'>=N C'=C'-N RETURN C' 则C'的最终返回值就是A*B*R' %N的精确值,所以关键在于求q。由于: (C' + A*B + q*N) %r =0 ==> (C' %r + A*B %r + q*N %r) %r =0 ==> (C'[0] + A*B[0] + q*N[0]) %r =0 若令N[0]*N[0]' %r =1,q=(C'[0]+A*B[0])*(r-N[0]') %r,则: (C'[0] + A*B[0] + q*N[0]) %r = (C'[0]+A*B[0] - (C'[0]+A*B[0])*N[0]'*N[0]) %r) %r = 0 于是我们可以得出r为任何值的蒙哥马利算法: m=r-N[0]' C'=0 FOR i FROM 0 TO k q=(C'[0]+A*B[0])*m %r C'=(C'+A*B+q*N)/r IF C'>=N C'=C'-N RETURN C' 如果令r=0x100000000,则%r 和/r 运算都会变得非常容易,在1024位的运算中,循环次数k 不大于32,整个运算过程中最大的中间变量C'=(C'+A*B+q*N) < 2*r*N < 1057位,算法效率就相当高了。唯一的额外负担是需要计算N[0]',使N[0]*N[0]' %r =1,而这一问题前面已经用欧几里德算法解决过了,而且在模幂运算转化成反复模乘运算时,N是固定值,所以N[0]'只需要计算一次,负担并不大。 素数测试方法 数论学家利用费马小定理研究出了多种素数测试方法,目前最快的算法是拉宾米勒测试算法,其过程如下: (1)计算奇数M,使得N=(2**r)*M+1 (2)选择随机数A