我发现的大数表示法，请各位评价下

zeroieme · 发表于 2011-3-30 13:41:14

马上注册，结交更多好友，享用更多功能，让你轻松玩转社区。

您需要登录才可以下载或查看，没有账号？欢迎注册

×

当单个数字不足够表示我们所需精度时，基本都是利用进位制多个单精度数连接。$\sum _{i=0}^n a_iM^i$

假如是M进位，通常是设置每位$[0,M-1]$，再加设一个符号位。而我发现设置每位$[-\frac{M}{2},\frac{M}{2}-1]$计算更便利。（设M是偶数）

1、现在已经是64位机时代，1个bit的符号位或者浪费，或者令首位与后面其他位不同，而增加计算复杂性。
2、加减不用对比绝对值，可直接计算。按多项式原理，乘除法也能直接计算。
3、最重要的是降低了进位的消耗。
$[0,M-1]$的加法每位得数平均值$\frac{M-1}{2}$，$[-\frac{M}{2},\frac{M}{2}-1]$的加法平均值$-\frac{1}{M}$，进位几率是常规的1/2。
乘法差异更大，$[-\frac{M}{2},\frac{M}{2}-1]$的乘法每位得数平均值$\frac{1}{4}$。而且大数乘法就是卷积，常规表示$[0,M-1]$中段会产生$O(M^2n)$的大结果。用每位$[-\frac{M}{2},\frac{M}{2}-1]$，因为各位正负参半，积也如此，因此它们的和能相互抵消，不会过度增长。
4、除法竖式法直除不用试除检查余数。不过我知道牛顿法以后没用竖式法了。

这是我发现的方法，我只在C上用带符号运算实现了，数论变换做乘法实现了2¹¹位。没有汇编优化，不知道带符号运算与无符号运算速度有多大差异。
请各位大大评价一下。或者有前人做过类似工作让我参考一下。

gxqcn · 发表于 2011-3-30 14:11:49

乘法卷积要考虑极端情形（而不是平均情况），比如每位的值的都是$-M/2$，此时与通常的表达形式没有任何优势。

该计数法在输入输出时的转换将成问题。
而且算法中对进位/借位的处理将麻烦一些。

zeroieme · 发表于 2011-3-30 14:22:03

都是 -M/2，平方是M^2/4，依然比(M-1)^2有优势。

由于依然是M进制，转换就是超M/2进1。是O(n)的复杂度。
由于每位都是带符号数，直接统一为进位，借位归结为进负数。

gxqcn · 发表于 2011-3-30 14:40:11

M^2/4 与 (M-1)^2 仅1:4的关系，几乎没什么优势。
而表达同样大的整数，你的计数法则可能得多用一“位”。

感觉这种计数法对固定长度的大整数计算还可以，
对于动态的大整数计算，比如一个小规模整数减去一个超大规模整数，新算法并没带来什么好处。

如果当$M$是奇数时，某些公式可能看起来更美妙一些。

liangbch · 发表于 2011-3-30 15:59:52

ooura fft代码中求圆周率的运算采用类似的表示法，使得作FFT运算时，误差更小。不同意大数运算就是做卷积，只有大数的规模很大时，卷积运算才有优势，在数的规模相对较小时，3,4 Toom-Cook算法仍然是最佳选择。将每个单元从 0 --> M-1 转为-M/2 --> M/2-1花不了多少时间，当需要FFT运算时，将每个单元转为-M/2 --> M/2-1 即可。

zeroieme · 发表于 2011-3-30 16:49:24

本帖最后由 zeroieme 于 2011-3-30 16:52 编辑

同样“位数”，[-M/2,M/2-1]能比[0,M-1]多表示2^(L-1)倍，假设每位字长L。

卷积是一种向量运算定义，FFT是计算卷积的方法之一。恰好大数乘法的定义和卷积一致。
倒过来看， Toom-Cook3,4也能作为小向量的卷积速算法。

极限状态几乎没什么优势，起码没有更慢。但这个极限的概率是多少？极限没什么优势，平均有上明显优势。

算了下每位乘法结果的分布状态
[0,M-1]
$E=\frac{1}{m^2}(\sum _{i=0}^{m-1} \sum _{j=0}^{m-1} i j)=\frac{1}{4} (-1+m)^2$
$D=\frac{1}{m^2}(\sum _{i=0}^{m-1} \sum _{j=0}^{m-1} (i j-E)^2)=\frac{1}{144} (-5+12 m-2 m^2-12 m^3+7 m^4)$

[-M/2,M/2-1]
$E=\frac{1}{m^2}(\sum _{i=-\frac{m}{2}}^{\frac{m-1}{2}} \sum _{j=-\frac{m}{2}}^{m\frac{m-1}{2}} i j)=\frac{1}{4}$
$D=\frac{1}{m^2}(\sum _{i=-\frac{m}{2}}^{\frac{m-1}{2}} \sum _{j=-\frac{m}{2}}^{m\frac{m-1}{2}} (i j-E)^2)=\frac{1}{144} (-5+4 m^2+m^4)$
奇怪的Mathematica输出

showjim · 发表于 2011-4-1 09:03:16

也许是我看不懂，怎么看都是把简单的事搞复杂了，一点优势也没有

gxqcn · 发表于 2011-4-1 10:10:27

无符号数相加时，需要考虑上溢出问题；
有符号数相加时，还得考虑下溢出问题。
而比较、进借位的判断，是比较耗clock的。

G-Spider · 发表于 2011-4-1 12:13:41

32位机取31位或32位还是多一些。32位进位有些麻烦，要经常看CF等标志，这又是汇编擅长的，好处之一是表示的位数多了，并且隐含了取模运算。符号位还是单独抽出来比较妥当，有符号数相加时，既要考虑上溢又要考虑下溢。

liangbch · 发表于 2011-4-1 12:57:20

8# gxqcn
　　我来谈一下关于分支预测的话题。

　　当代CPU均采用pipeline来实现多级流水线，在当前指令的执行阶段，其后的指令，有的已经完成译码，有的已经完成取指等操作，所以顺序执行的指令速度很快，但一旦遇到条件转移，事情就比较复杂了。

　　早期的CPU(如8086)没有分支预测部件，所以在执行条件转移指令时，如不跳转则速度很快，如跳转则速度很慢，这是因为，当执行跳转指令时，pipeline中的所有内容必须被清除，已经取到的指令，已经完成的译码等操作都白干了。如果你查看8086参考手册，你会发现，对于条件转移指令，如不跳转则仅需4个周期，如跳转则需要16个周期，见本站链接：资料：Intel 和 AMD CPU 资料。

　　当前的CPU都有分支预测部件，CPU在执行每个分支指令(比较跳转指令对)时，会保存这个指令的跳转情况，当下次执行这条指令时，根据这个指令的跳转情况历史记录，执行跳转概率比较高的那个分支，这个分支的后续指令会提前进入pipeline，这就是分支预测。如果预测正确，则跳转指令的执行速度很快，和普通的指令速度相当。如果预测失败，则必须付出相当的代价。对于当前的CPU，依赖于具体的程序逻辑，分支预测的正确率会有所不同，不过平均起来，分支预测的正确率很高，按照intel的文档，当前的CPU分支预测的正确率一般可达90%以上。

　　分支有数据依赖分支和非数据依赖分支，对于后者，分支预测的正确率非常高，速度很快，但对于前者速度就相当的慢了。下面的举2个例子来说明。

例子1, 非数据依赖分支

C代码：

int sum(int *p, int len)
{
int s=0;
int i;
for (i=0;i<len;i++)
s+=len;
return s;
}

复制代码

下面是我写的对应的汇编子程序

_sum PROC NEAR
push ebx
mov ebx, DWORD PTR [esp+4+4] //ebx为p
mov edx, DWORD PTR [esp+4+8] //edx为len
xor eax, eax //eax 为返回值是s
test edx, edx
jle SHORT this_exit //len<0, 直接返回0, eax为返回值，这里eax=0
xor ecx, ecx //ecx 为循环变量i,初值为0
//循环部分
Loop_start:
mov edx, DWORD PTR [ebx+ecx*4] //edx 为p[ i ]
add eax, edx //s+= p[ i ]
inc ecx //i++;
cmp ecx, DWORD PTR [esp+4+8] // i<len?
jl Loop_start // i<len 则继续循环
this_exit:
pop ebx
ret 0
_sum ENDP

复制代码

在这个子程序中，共有2处比较跳转指令对
第1处是：
　　test  edx, edx
　　jle SHORT this_exit       //len<0, 直接返回0, eax为返回值，这
第2处是：
　　cmp ecx, DWORD PTR [esp+4+8] // i<len?
　　jl  Loop_start // i<len 则继续循环

　　我们重点说一下第2处，对应的c代码的语义为"如果i<len，则循环继续"。这个分支，共须执行len次，只有最有一次不跳转，其余的len-1均需要跳转。对于这样的分支，分支预测部件的历史记录显示绝大多数是跳转的，所有分支预测部件总是能成功预测。
　　这样的循环控制分支，和数据无关，分支预测的正确率很高，速度很快。

例子2. 数据依赖分支
　　下面的函数计算 array1[ i ]+array2[ i ] (i=0..len)并将结果存入array1[ i ]
　　数组 array1和array2的每个元素均为0..9999

int arrayAdd(int array1[], int array2, int len)
{
int i,sum,carry;
carry=0;
for (i=0;i<len;i++)
{
sum= array1[ i ]+ array2[ i ] + carry;
carry=0;
if (sum>=10000)
{
carry=1;
sum-=10000;
}
array1[ i ]=sum;
}
return carry;
}

复制代码

我写的等价的汇编代码

_arrayAdd PROC NEAR
push esi
push edi
xor eax, eax ; eax：返回值carry
mov edx, DWORD PTR [esp+8+12] //[esp+8+12]： len
test edx, edx
jle SHORT this_exit //len<=0，则直接返回0
mov esi, DWORD PTR [esp+8+4] //esi: array1
mov edi, DWORD PTR [esp+8+8] //edi: array2
xor ecx,ecx //ecx=0: ecx：循环变量i，
xor edx,edx //edx=0：edx：carry
loop_start:
mov eax, dword ptr [esi+ecx*4] //eax =array1[ i ]
add eax, dword ptr [edi+ecx*4] //eax对应于sum，sum= array1[ i ]+array2[ i ]
add eax, edx //eax对应与sum, sum= array1[ i ]+array2[ i ]+carry
xor edx, edx //carry=0
cmp eax, 10000
jl SHORT next_comp
mov edx, 1 //carry=1
sub eax, 10000 //sum-=10000
next_comp:
mov DWORD PTR [esi+ecx*4], eax //array1[ i ]=sum
add ecx, 1 //i++;
cmp ecx, DWORD PTR [esp+8+12] // i<len?
jl SHORT loop_start //i<len 则继续循环
mov eax,edx //eax:为返回值
pop edi
pop esi
this_exit:
ret 0
_arrayAdd ENDP

复制代码

在这个函数/过程中，共有3处比较跳转指令对
第2处是：
　　cmp eax, 10000
　　jl SHORT next_comp

第3处是：
　　cmp ecx, DWORD PTR [esp+8+12] // i<len?
　　jl SHORT loop_start //i<len 则继续循环

第3处用于循环控制，和第1个例子一样，是一个非数据依赖分支，cpu分支预测的正确率很高，速度很快。
但是，第2处比较跳转指令对是一个数据依赖分支，两个0-9999之间的数相加，需要进位和不进位的概率各为50%，CPU分支预测部件无法知道下一次这个分支指令的哪一个分支将被执行，因此无论执行哪一个分支，都有50%的情况预测失败，cpu在执行这个跳转指令时，速度很慢。

账号		自动登录	找回密码
密码			欢迎注册

[原创] 我发现的大数表示法，请各位评价下

马上注册，结交更多好友，享用更多功能，让你轻松玩转社区。