由RSA引出的问题。（欢迎讨论）

[z07165864809]

起因：密码学老师让写一个RSA实现的算法．

经过：以前我就写过了，这次，我想写的高深些，使用最好的办法．

技术上的问题：
1．大数存储：
　　这样应该没有争论的，大家使用的都是３２位整数数组，以２＾３２为进位制．存储方式是由低到高的存储．就是0x aaaaaaaa ffffffff bbbbbbbb 在数组int a[]中表示为：a[0]=0xbbbbbbbb;a[1]=0xffffffff;a[2]=0xbbbbbbbb.

2．加法和减法使用汇编最快了，这个也没有问题．
３．rsa中的各种算法短时间里面也没有办法改进了。这些都不管他们了．

４．有改进的地方只有乘法和除法的实现了。！！！！！！

下面要和大家讨论的就是乘法和除法的快速计算的问题．

乘法：

通用的算法有３个．
一：传统乘法。对少于３２*４０位比较有效．
二：Karatsuba乘法。对于少于32*600位比较有效。
三：FFT乘法。100万位的时候，没有算法比它快。

[z07165864809]

RSA使用的是512-2048位的密码，所以，使用传统的乘法就可以了。
这里有一个问题。就是，怎么计算a=a^2??
如果使用一个个的乘（就是计算a*a)，要o(n^2) （n为a的数据数组长度）。
有没有简单的办法呢？因为这个在计算a^b mod c时很常见。

第二个问题出现在除法上面。
计算a=b/c，使用的办法就是，先计算2^n/c（这里的n是2*BitLength(c-1),Newton方式)；
再计算d=b*2^n/c；
最后再a=d>>n；
这个算法是常见的一种。也是giantint（http://www.perfsci.com/index.html#news）大数计算库使用的办法之一。
这里有2个问题。
一：在giantint的介绍中，它说它使用了基于Barrett's method的大数除法，可是，这个到底是什么方式的？？？
二：我以前使用的是试除法。就是用最高位除以最高位（比如，a[20]/b[20]),结果是i的话，让a=a-b*(i-1),再试除。直到，a<b为止。
问题是，到底这个办法和上面的求倒的办法，谁更好呢？有没有谁能给出点数学证据来？

[z07165864809]

使用Barrett's method的计算方法书上其实有描述。
只是，其他的问题还是想和大家讨论讨论。

[NowCan]

研究的比较深入了。

[z07165864809]

上面有一个笔误。
“就是0x aaaaaaaa ffffffff bbbbbbbb 在数组int a[]中表示为：a[0]=0xbbbbbbbb;a[1]=0xffffffff;a[2]=0xbbbbbbbb. ”
这里的a[2]应该是=0xaaaaaaaa.

下面我把具体怎么实现除法的方式写出来，也让不知道的兄弟来参加讨论。
比如，10进制下，我们计算755/31（注意，这里我们只要quotient(商），remainder（余数）可以用减法得到）。
方法一（Newton迭代）：

先求10000/31。 （这里是不是可以改进成1000/31？）
等于求函数f(x)=10000/x-31的零点。
知道Newton迭代公式：x(n+1)=2*x(n)-31*x(n)^2/10000.
令x(0)=100；
x(1)=2*x(0)-31*x(0)^2/10000=169;
x(2)=2*x(1)-31*x(1)^2/10000=250;
x(3)=2*x(2)-31*x(2)^2/10000=307;
x(4)=2*x(3)-31*x(3)^2/10000=322;
到这里，就不用再算下去了。看的出，322*31=9982最接近10000了。
计算755*322/10000;
结果是a=755*322/10000=24;

我们分析一下上面的复杂度，这里，假设加法、减法常数，乘法不是。
计算10000/31的时候，进行了3*4次=12次乘法。而除法中，除以10000可以简单的移位来实现，所以时间可以是常数。最后，计算755*322进行了一次乘法。共13次乘法。推而广之，可以估计，除数有n位的除法，可以在2*n+1次乘法中完成。

[z07165864809]

方法二（试除法）：

d=755;
a=d;
b=31;
c=0;
使用最高位相除的办法。（因为我们假设处理器可以快速的处理2个小数相除）那么有：
c=7/(3+1)*10=10; （时间为0）
a=d-b*c=445; （一个乘法）
c=c+4/(3+1)*10=20; （时间为0）
a=d-b*c=135; （一个乘法）
c=c+13/(3+1)=23; （时间为0）
a=d-b*c=42; （一个乘法）
c=c+4/(3+1)=24; （时间为0）
a=d-b*c=11; （一个乘法）
因为这里有a<b;所以，计算结束。

好象只要4个乘法呢！是不是比上面的快？？
我很有疑惑的是，为什么象gaintint这样的高级库不用这样的方式呢？？

[z07165864809]

呵呵，这些不是最好的办法。
已经公布的最好的办法比这些要好。
可是，他们对小数（100万位以下的数相乘没有用处啊！）
我想看看有没有自己设计出了高水平的算法的兄弟。

难道这里没有高手了么？？？？？？？难道你们写的RSA算法的程序都是照抄的别人的么？？？？

[BlueSky2008]

幂模不是直接算得吧，而应该使用"平方-和-乘法"，
至于平方，用平方和公式(a+b)^2 = a^2+b^2+2*ab不是可以更快一点吗？

[z07165864809]

研究了几天《半数值算法》后，终于有一些了解了。
最后的结论还是：
在2048位以下，使用试除法是最快的方式。（经过改进的算法，可以在2*n*n个乘法完成）

可是平方还是没有好的办法。
期待大家的改进。

[liangbch]

对于第一个问题，我的看法是：
Karatsuba乘法大概就是分治法吧，其基本原理是：2n位 *2n位 的乘法运算可以转化为3次 n位 * n位的乘法运算。在我的阶乘计算器3.0版本中，就使用了这个算法，为了提高大数的显示速度，在我的程序中，大数的表示是1000000000进制的，经过测定，我的程序中，在PIII/赛扬533-733 CPU 运行，当乘数的长度>64个DWORD(32位整数),分治法开始快于普通乘法。

对于第二个问题，我的看法如下：
1。试商法的复杂度分析，假设被除数为2n位，除数为n位，则商为n位（一般的n位数 * n位数的结果为 2n-1 或者 2n位)，在整个除法过程中，一般 需要重复3个过程：
1。试商，复杂度可忽略：
2。计算部分积，需要n次 1位数乘法(一位数一个基本数据类型，目前的pc位一个整数，32bit,下同)
3。余数 减 部分积，需要n+1次 一位数减法。
4。这只是一般的过程，有时还需要调整 新试的商，因为单从 被除数 和除数 的前几位来试商，可能会使试商偏大或者 偏小，当部分积计算完成时，如果部分积 大于被除数，则试商偏大，需要部分积减去被除数；当部分积计算完成时，如果被除数减去部分积的结果仍然大于除数，在试商偏小，需要部分积加上除数。
如果这4个过程中，1和3，4可以忽略，因此需要n*n次乘法。

2。Newton迭代 的复杂度分析：
在通过newTon 迭代 求 除法时：我见过的算法是这样的， 首先求除数的倒数，使用迭代法：每次迭代可以将计算精度提高1倍，例如，初值的精度为1位，经过5次迭代，精度依次达到 2，4，8，16，32位，假设结果为n位，则这一过程 需要 4/3 * n^2 次乘法，其次求倒数和被除数的积，这一步需要n^2 次乘法，总共需要 7/4* n^2 乘法，因此单从乘法次数看，的确比试商法慢，但是，迭代法使用大数乘法，可以使用分治法来提高速度，而试商法 只能使用 1位 乘 多位数的普通乘法，不能使用分治法，所以从这个角度来看，Newton迭代法优于试商法。
3。第三种方案，2者的结合（我自己想出来的)，例如 计算 a / b, A.先用迭代法求b的倒数c，只求8/16/32...位，B.用 a * c 来试商，这样一次可以8/16/32位部分商 C.再用部分商和除数相乘得到部分积，D)用余数 和 C 的计算结果想减求新的余数，重复B,C,D 这3部，直到余数小于除数为之。 这个算法结合了前2种优点，在计算乘法时，可以使用分治法，当计算基于十进制表示的小数（浮点数，非整数)时很有用。

对于你，我的看法时，1. 太小气，这么大个问题只给20分。 2. 太轻狂“难道这里没有高手了么？？？？？？？“，连用7个问号。

[z07165864809]

终于引起（宝宝）老兄的不满了？哈哈。

首先，对于给多少分，我没有概念，因为我自己虽然注册了这么长的时间，都只是看看，没有问过问题。

说“难道这里没有高手了么？？？？？？”的目的是为了引起别人的注意。

不知道宝宝兄写过RSA没有？
因为它的位数小，所以，２分法或FFT都没有用。
对于除法，用试的办法，我已经有比较好的解决方法了。使用的不是直接除法（Knuth的书上的），是改进的方式，很复杂，难以表达出来。
今天我将重点研究宝宝兄的算法。或许有用的。

还有啊！我想高手不应该只会用别人的算法，要自己有更好的解决方案才行．
很感谢宝宝兄的第３种方式．

但愿宝宝兄可以帮我解决快速的平方的问题。
可能这个才是关键的，因为平方在RSA中出现的次数太多了。

注：已经研究过若干算法了。没有一个有用的。

[z07165864809]

宝宝兄：

思考了半天，发现你的第三中方式在这里没有用处，它在这里不是集中了优点，而是集中了缺点。
因为位数小，所以只能用传统的乘法，那么newton迭代就失去了意义。

而，试二位和二次试一位是一样的，他们的乘法都是2*n次．

[Xcoder]

faint,偶做大数除法都是用减法做的

[z07165864809]

倒，要是减法的话，
999999999999999999999999999999999999999999999999999999/1
你要做到什么时候？？

[NowCan]

诶，数学功底不够，完全看不懂。
我还是先看看别人写的RSA吧。
MIRACL和GMP里都有吧。

[liangbch]

同意 BlueSky2008 的方法，，用平方和公式(a+b)^2 = a^2+b^2+2*ab，
另f(n)表示 n位数平方的时间复杂度，
将n位数一分为2，前n/2为数表示为a, 后 n/2为数表示为b,

则f(n)=f(n/2) + f(n/2) + 1/4 *n^2
计算a的复杂度 计算b的复杂度+ 计算2*a*b的复杂度

注：乘以2只需要做加法。

如此，在计算n位数*n位数的平方时，可以将n位数一分为2，分别计算a^2(递归调用)，b^2(递归调用)，a*b 这3部分，然后错位相加即可。下面用例子所一个算法分析，公式你去推到吧。
1。将 f(2)= 2位 * 2位=4次 乘法 作为1个基本单位，计为x

递归的复杂度 直接计算的复杂度
2。则 f(4)= 2*f(2) + 2*2 = 3x 4x
f(8)=2*f(4) + 4*4 = 6x + 4x=10x 16x
f(16)=2*f(8) + 8*8 = 20x + 16x=36x 64x
f(32)=2*f(16) + 16*16 = 72x + 64x=136x 256x
f(64)=2*f(32) + 32*32 = 272x + 256x=528x 1024x
f(128)=2*f(64)+ 64*64 = 1056x + 1024x=2080x 4096x
f(256)=2*f(128)+128*128=4160x + 4096x=8256x 16384x
直接计算f(256)=256*256=65536,递归计算需要8256*4=33024,速度提高了近1倍,这只是理论值，在实际计算中，由于要做 加法运算，函数的递归调用等，实际的复杂度要比这个值高。

[z07165864809]

是啊！上面的办法可以加快1倍的。

已经有简单的方法了。

比如，abcd * abcd

可以的方式是：

a*a --->结果的最低1、2位。
2*a*b--->到结果的2、3位。
2*a*c+b^2--->到结果的3、4位。
2*(a*d+b*c)--->到结果的4、5位。


有鉴于这个的难度太大了。
所以，在思考了7天后，研究了若干本牛鬼的书和大数库以后，我已经失去了信心了，我打算再等二天就结了。


--------因为，大家好象都没有好的办法。可能只能到数学上面找办法了，而不是算法上面。

谢谢大家的支持。