主题:【求助】由一个编程考试想到的 -- 东方射日
前两天有一个remote coding test,题目是:
Find the first 10 consecutive digits of the mathematical constant "e" for which every digit from 0-9 is present only once.
看这个题目,我马上上网搜到e的定义是自然数阶乘的倒数和的极限。有了这个定义,我大约花了两个小时编程,调试和优化,后来有了代码和结果。(有兴趣的可以给我发短消息,我可以给你发我的代码,看看还有没有优化空间)
// My result:
// The first place find 0 - 9 is begin at 1729 :
// 0719425863
我试了下效率还算满意,e计算到小数点后10000位在我那台烂机器(p4 2G, 512M : 不许笑!!)上也只用1秒左右。
发出去后,突然想到是不是可以也计算一下其他的无理数,例如根号2。那么应该有以下函数来精确计算一个正整数的平方根:
string sqrt(int n, int num)
// 返回n的平方根精确到小数点后num位
不过,平方根的计算我就捉瞎了,搜一下,找到例程:
double Sqrt(double a,double p)//a是被开平方根数,p是所求精度
{
double x=1.0;
double check;
do{
x=(a/x+x)/2.0;
check=x*x-a;
}while(check<-p || check>p);
return x;
}
有了这个迭代就好办了,反正用字符串来精确计算的代码上面已经有了。不过我的问题是:
以上红字部分:x=(a/x+x)/2.0;
是怎么来的?为什么?
立方根又如何计算了?
x=(a/x+x)/2.0;
按这个迭代,在当x=sqrt(a) 时应当有:
sqrt(a) = (a/sqrt(a) + sqrt(a))/2
显然成立。
同样对立方根也有类似关系
所以将2换为3显然不对
x = (a/(x*x) + x)/2 还有可能
具体的您可以看看这个链接
http://mitpress.mit.edu/sicp/full-text/sicp/book/node12.html
这个算法是怎么来的?
这里有一些相关的讨论,可以参考
http://topic.csdn.net/u/20080908/18/20339824-3be7-4ba6-baed-6641240566a9.html
不知道5000位数字中,会不会出现所需要的结果?
好像这样的数字在加密安全、随机数方面有很多应用
对于 任意一方程 f(x)=0;, 求 x 值,
x[n+1]=x[n]-f(x[n])/f'(x[n])
对于 f(x)=x^2-a=0, (即 x^2=a);
x[n+1]=x[n]- (x[n]^2-a)/(2x[n])
即
x[n+1]=(a+x[n]^2)/(2x[n])
如果是 f(x)=x^3-a, (即 x^3=a)牛顿法迭代为:
x[n+1]=(2x[n]^3+a)/(3x[n]^2)
牛顿厉害!
讨论见
http://en.wikipedia.org/wiki/Newton%27s_method#History
比如三分之一得到的结果是3.3333333333333331482961625624739099293947219848632812500000000000000000000000000000000000000000000000e-01
如果你算到上千位那结果可能根本就不正确。希望得到你的回复。
要想高精度计算,绝对不能用浮点数。由于计算机是二进制的,在电脑中的浮点数总有一个误差。例如你所说的1/3就是一个例子。
这个误差会有很多问题,例如:
float a = 1.0f / 3.0f;
a /= 5.0f;
float b = 1.0f / 5.0f;
b /= 3.0f;
这里你可以试试比较a和b,应当得到是不等的。
所以通常我们比较两个浮点数是否相等不是直接用“==”而是比较两者差值小于一个我们可以接受的误差即认为相等。
如果如我说说的进行高精度计算只能自己用一个数组进行管理。
例如我们要计算到小数点后10000位,就建立一个至少10000位的数组,当然你最好建立一个类:
class Digit
{
public:
Digit();
~Digit();
friend Digit & operator + (const Digit & a, const Digit & b);
// 你自己可以添加其他诸如 - * / > < ==等运算符的实现
private:
int mInterger;
char mDecimal[10000]; // 当然最好在构造函数中指定位数,或提供一个成员函数来设定或变更精确的位数
}
下面,我把e计算的源代码贴出来,请各位指教,看看是否还有性能提高的余地。
不过该程序是测试时写的,如果想作为库还是不合格,因为 result在内部分配的内存,不符合自动释放或谁分配谁释放的原则,存在内存泄露的危险
// num:小数点后位数,为保证精度,调用时可以多一个经验值,例如10位。比如需要精确到小数点后10000位,就设num位10010
//result:返回结果,第一位为整数部分,其余为小数,每字符为一位,为0-9;如果需要直接可打印的结果,可以将各字符加上‘0’
void calculateE(int num, char * &result)
{
// e = 1 + 1 + 1/2! + 1/ 3! + 1/4! + ..... + 1/n!
if (result) SAFE_DELETE(result);
result = new char[num];
memset(result ,0,num);
assert(result);
char * tmp;
char * lastTmp;
char * tempResult[2];
tempResult[0] = new char[num];
assert(tempResult[0]);
tempResult[1] = new char[num];
assert(tempResult[1]);
tmp = tempResult[0];
lastTmp = tempResult[1];
int tmpIndex = 0;
memset(lastTmp,0,num);
lastTmp[0] = 1;
for (int cal = 2;cal <= num; ++cal)
{
memset(tmp,0,num);
int firstNoZero = 0;
static int div = 1;
for (;div <= cal; ++div)
{
int quo = 0;
int rem = 0;
for (int j = 0; j< num; ++j)
{
quo = rem * 10 + lastTmp[j];
rem = quo % div;
quo = quo / div;
tmp[j] = quo;
if (firstNoZero == 0 && div == cal && tmp[j] != 0)
{
firstNoZero = j;
}
}
}
div = cal;
int flag = 0;
if (firstNoZero == 0)
{
break; // no need to calculate anymore 此时在需要的精度内全为零,不需要继续计算了
}
for (int i = num - 1; i >= firstNoZero - 1; --i)
{
int tmp_sum = tmp[i] + result[i] + flag;
flag = (tmp_sum >= 10 ? 1 : 0);
result[i] = (flag ? tmp_sum - 10 : tmp_sum);
}
lastTmp = tmp;
tmpIndex = (tmpIndex ? 0 : 1);
tmp = tempResult[tmpIndex];
}
delete[] tempResult[0];
delete[] tempResult[1];
result[0] += 2;
}