当我看到这个问题时,我决定解决它,但那时我非常忙碌。
上个周末我终于有了一些空闲时间,所以考虑了一下我的待处理挑战。
首先,我建议您考虑上面的回答。我从未使用过GMP库,但我确信它比手工编写的代码更好。
此外,您可能有兴趣分析bc计算器的代码;它可以处理大数,我用它来测试自己的代码。
好的,如果您仍然有兴趣自己编写代码(只支持C语言和标准C库),我可以给您一些东西。
在所有之前,先讲一点理论。在基本的数字理论(模运算级别)中,有一个算法启发了我得出一个解决方案;乘法和幂算法用于求解a^N对m取模的值:
Result := 1;
for i := k until i = 0
if n_i = 1 then Result := (Result * a) mod m;
if i != 0 then Result := (Result * Result) mod m;
end for;
其中k是N的二进制表示中数字位数减一,n_i是第i个二进制数字。例如(N为指数):
N = 44 -> 1 0 1 1 0 0
k = 5
n_5 = 1
n_4 = 0
n_3 = 1
n_2 = 1
n_1 = 0
n_0 = 0
当我们进行模块操作时,例如整数除法,我们可能会丢失部分数字,因此我们只需要修改算法以避免丢失相关数据。
以下是我的代码(请注意,这是一个特定情况下的代码,强烈依赖于我的计算机体系结构。基本上我在使用C语言的数据长度,所以请小心,因为我的数据长度可能不同):
#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include <math.h>
enum { SHF = 31, BMASK = 0x1 << SHF, MODULE = 1000000000UL, LIMIT = 1024 };
unsigned int scaleBigNum(const unsigned short scale, const unsigned int lim, unsigned int *num);
unsigned int pow2BigNum(const unsigned int lim, unsigned int *nsrc, unsigned int *ndst);
unsigned int addBigNum(const unsigned int lim1, unsigned int *num1, const unsigned int lim2, unsigned int *num2);
unsigned int bigNum(const unsigned short int base, const unsigned int exp, unsigned int **num);
int main(void)
{
unsigned int *num, lim;
unsigned int *np, nplim;
int i, j;
for(i = 1; i < LIMIT; ++i)
{
lim = bigNum(i, i, &num);
printf("%i^%i == ", i, i);
for(j = lim - 1; j > -1; --j)
printf("%09u", num[j]);
printf("\n");
free(num);
}
return 0;
}
unsigned int bigNum(const unsigned short int base, const unsigned int exp, unsigned int **num)
{
unsigned int m, lim, mem;
unsigned int *v, *w, *k;
mem = ( (unsigned int) (exp * log10( (float) base ) / 9 ) ) + 3;
v = (unsigned int *) malloc( mem * sizeof(unsigned int) );
w = (unsigned int *) malloc( mem * sizeof(unsigned int) );
for(m = BMASK; ( (m & exp) == 0 ) && m; m >>= 1 ) ;
v[0] = (m) ? 1 : 0;
for(lim = 1; m > 1; m >>= 1)
{
if( exp & m )
lim = scaleBigNum(base, lim, v);
lim = pow2BigNum(lim, v, w);
k = v;
v = w;
w = k;
}
if(exp & 0x1)
lim = scaleBigNum(base, lim, v);
free(w);
*num = v;
return lim;
}
unsigned int scaleBigNum(const unsigned short scale, const unsigned int lim, unsigned int *num)
{
unsigned int i;
unsigned long long int n, t;
for(n = 0, t = 0, i = 0; i < lim; ++i)
{
t = (n / MODULE);
n = ( (unsigned long long int) scale * num[i] );
num[i] = (n % MODULE) + t;
}
num[i] = (n / MODULE);
return ( (num[i]) ? lim + 1 : lim );
}
unsigned int pow2BigNum(const unsigned int lim, unsigned int *src, unsigned int *dst)
{
unsigned int i, j;
unsigned long long int n, t;
unsigned int k, c;
for(c = 0, dst[0] = 0, i = 0; i < lim; ++i)
{
for(j = i, n = 0; j < lim; ++j)
{
n = ( (unsigned long long int) src[i] * src[j] );
k = i + j;
if(i != j)
{
t = 2 * (n % MODULE);
n = 2 * (n / MODULE);
dst[k] = ( (k > c) ? ((c = k), 0) : dst[k] ) + (t % MODULE);
++k;
dst[k] = ( (k > c) ? ((c = k), 0) : dst[k] ) + ( (t / MODULE) + (n % MODULE) );
++k;
dst[k] = ( (k > c) ? ((c = k), 0) : dst[k] ) + (n / MODULE);
}
else
{
dst[k] = ( (k > c) ? ((c = k), 0) : dst[k] ) + (n % MODULE);
++k;
dst[k] = ( (k > c) ? ((c = k), 0) : dst[k] ) + (n / MODULE);
}
for(k = i + j; k < (lim + j); ++k)
{
dst[k + 1] += (dst[k] / MODULE);
dst[k] %= MODULE;
}
}
}
i = lim << 1;
return ((dst[i - 1]) ? i : i - 1);
}
unsigned int addBigNum(const unsigned int lim1, unsigned int *num1, const unsigned int lim2, unsigned int *num2)
{
unsigned long long int n;
unsigned int i;
if(lim1 > lim2)
return 0;
for(num2[lim2] = 0, n = 0, i = 0; i < lim1; ++i)
{
n = num2[i] + num1[i] + (n / MODULE);
num2[i] = n % MODULE;
}
for(n /= MODULE; n; ++i)
{
num2[i] += n;
n = (num2[i] / MODULE);
}
return (lim2 > i) ? lim2 : i;
}
编译方法:
gcc -o bgn <name>.c -Wall -O3 -lm //Math library if you wants to use log func
要检查结果,请将直接输出作为bc的输入。以下是一个简单的shell脚本:
#!/bin/bash
select S in ` awk -F '==' '{print $1 " == " $2 }' | bc`;
do
0;
done;
echo "Test Finished!";
我们拥有一个无符号整型的数组(4个字节),在数组的每个整数中存储一个九位数字(%1000000000UL); 因此num [0]将具有前9位数字,num [1]将具有第10到18位数字,num [2] ... 我们使用传统内存进行工作,但可以通过动态内存进行改进。好的,但是这个数组应该有多长?(或者需要分配多少内存?)使用bc计算器(bc -l带mathlib)我们可以确定一个数字的位数:
l(a^N) / l(10) // Natural logarith to Logarithm base 10
如果我们知道数字,就知道我们需要多少个整数:
( l(a^N) / (9 * l(10)) ) + 1 // Truncate result
如果您处理的数值为(2^k)^N这样的式子,您可以使用这个表达式进行对数运算:
如果用LaTeX表示为$\log_2{(2^k)^N}=kN$
( k*N*l(2)/(9*l(10)) ) + 1 // Truncate result
确定整数数组的确切长度。例如:
256^800 = 2^(8*800) ---> l(2^(8*800))/(9*l(10)) + 1 = 8*800*l(2)/(9*l(10)) + 1
值为 1000000000UL(10^9)的常量非常重要。像10000000000UL(10^10)这样的常量不起作用,因为可能会产生未检测到的溢出(尝试使用数字16^16和10^10常量会发生什么),而更小的常量,如1000000000UL(10^8),是正确的,但我们需要保留更多的内存并执行更多步骤。 10^9是32位无符号整数和64位无符号长整数的关键常量。
代码有两部分,乘法(易)和2的幂次方(更难)。乘法只是乘法并扩展整数溢出。它采用数学中的结合律原理来做完全相反的原理,因此如果k(A + B + C),我们想要kA + kB + kC,其中数字将为k * A * 10 ^ 18 + k * B * 10 ^ 9 + kC。显然,kC操作可能会生成一个大于999999999的数字,但永远不会超过0xFFFFFFFFFFFFFFFF。在乘法中,大于64位的数字永远不会出现,因为C是32位无符号整数,而k是16位无符号短整数。在最坏的情况下,我们将拥有这个数字:
k = 0x FF FF;
C = 0x 3B 9A C9 FF; // 999999999
n = k*C = 0x 3B 9A | 8E 64 36 01;
n
n / 1000000000 = 0x FF FE;
在Mul k
B之后,我们需要从C的上一次乘法中添加0xFFFE(B = kB +(C / module)),以此类推(我们有18位算术偏移量,足以保证正确的值)。
Power更复杂,但本质上是相同的问题(乘法和加法),因此我提供了一些关于代码功率的技巧:
- 数据类型非常重要
- 如果尝试将无符号整数与无符号整数相乘,则会得到另一个无符号整数。使用显式转换来获取unsigned long long int,并且不会丢失数据。
- 总是使用unsigned修饰符,不要忘记它!
- 2的幂可以直接修改当前索引的2个索引
- gdb是你的朋友
我开发了另一种将大数相加的方法。这些我没有进行太多测试,但我认为它运行良好。如果它有错误,请不要对我太过苛刻。
...就这样!
PD1:在
中开发。
Intel(R) Pentium(R) 4 CPU 1.70GHz
Data length:
unsigned short: 2
unsigned int: 4
unsigned long int: 4
unsigned long long int: 8
数字如 256^1024 要花费:
real 0m0.059s
user 0m0.033s
sys 0m0.000s
一个循环计算i^i,其中i的范围为1到1024:
real 0m40.716s
user 0m14.952s
sys 0m0.067s
对于像65355^65355这样的数字,计算所需时间非常长。
PD2: 我的回复很晚但我希望我的代码能够有用。
PD3: 抱歉,用英语解释是我最大的缺陷之一!
最后更新:我刚刚想到了一个新的实现方式,使用相同的算法来改善响应并减少内存使用量(我们可以使用unsigned int的完整位数)。秘密:n^2 = n * n = n * (n - 1 + 1) = n * (n - 1) + n。
(我不会写这个新代码,但如果有人感兴趣,可能会在考试后...)