我想的圆周率计算法，１万位十进制精度计算时间在数秒至数十秒左右

abiao

我用自制的长浮点数计算ｐｉ，１万位十进制精度计算时间在数秒至数十秒左右，视机器性能而定， 基本算法是 x=3或者3.1或者3.14等等; 循环(x+=sin(x)); 每计算一次sin可以获得一倍的精度，sin是用级数叠加的，每叠加一次需要计算4次乘法，一次除法，一次加法，获得２位左右的10进制精度 加乘除法是用最常规的方法计算的，没有使用傅立叶变换，只是在汇编级进行了优化，利用了x86的指令集的特性 基本的数据结构是使用一个无符号整数数组，每个都是一个计算单位,(2^32=4G)进制，表示时再转回10进制，一个符号标志，一个指数标志 性能在使用常规方法计算中应该算快的了吧 这个计算方法是在我算sin(3.14)=0.00159xxxx时想到的

gxqcn

我估计这样的算法流程是比较慢的。 因为： 1、计算pi的数学公式不够高效； 2、未使用高级大数乘法算法（FFT之类）； 3、先2进制，再转10进制，这个转换是非常耗时的，如果处理不好的话。

mathe

最主要是sin函数本身效率不会很高。 至于说高级大数乘法运算，是具体如何做乘法的问题，倒是同楼主这里的算法没有关系。 二进制到十进制转化倒还好，这个转化过程如果实现得好，相对于计算过程花费的时间是很少的。我记得最近计算派的记录的那个链接里面就是采用2进制计算的，然后最后转化为10进制的

liangbch

计算sin(x)本身就很麻烦，这样会很慢。常规算法中，Martin公式应该比你的算法更快，而且更简单。

abiao

算pi的主要目的是测试三角函数的正确性

abiao

另外Martin公式太复杂，我的这个比较简洁，虽然说sin要计算多次，不过前几次的时间可以忽略不计的

〇〇

源代码呢

abiao

源代码比较多，贴出乘法的吧 void mul(int Precision, unsigned int* Num_Array, int* pSign, int* pExponent, unsigned int* num_Num_Array, int num_Sign, int num_Exponent, unsigned int* buf_Num_Array) { int i; int cutOff = 0; for (i = Precision - 1; i >=0; i--) { buf_Num_Array[i] = Num_Array[i]; if(cutOff == 0 && Num_Array[i] != 0) { cutOff = Precision - i; } Num_Array[i] = 0; } cutOff--; int temp_Exponent = (*pExponent) + num_Exponent; *pExponent = 0; *pSign = (*pSign) * num_Sign; for (int idxMain = 0; idxMain < Precision; idxMain++) { int n = num_Num_Array[idxMain]; if (n!=0) { int shift = idxMain + (*pExponent);//(*pExponent)修正量 int cnt = Precision - (cutOff > shift ? cutOff : shift); int carry = 0; if (cnt>0) { unsigned int* p0 = &Num_Array[shift - 1]; unsigned int* p1 = &buf_Num_Array[-1]; _asm { mov ecx,cnt mov edi,p0 mov esi,p1 mov ebx,0 lp: mov eax,n cdq mul DWORD PTR [esi+ecx*4] // low:eax,high:edx add eax,[edi+ecx*4] adc edx,0 add eax,ebx adc edx,0 mov [edi+ecx*4],eax mov ebx,edx loop lp mov carry,ebx } } //继续进位 if (shift>0) { unsigned int* p0 = &Num_Array[shift - 1]; _asm { mov ecx,shift mov edi,p0 mov edx,carry lp2: mov eax,[edi] add eax,edx mov [edi],eax mov edx,0 jnc ext mov edx,1 sub edi,4 loop lp2 ext: mov carry,edx } } if (carry!=0) { (*pExponent)++; for (i = Precision - 1; i > 0; i--) { Num_Array[i] = Num_Array[i - 1]; } Num_Array[0] = carry; } } } (*pExponent) += temp_Exponent - 1; }

abiao

我用了Martin公式(pi=16*arctan(1/5)-4*arctan(1/239)) 发现比我的算法还慢 因为sin的收敛的速度比atan快，叠加的次数少 sin x = x-x^3/3!+x^5/5!-...(-1)k-1*x2k-1/(2k-1)!+... (-∞

liangbch

我试了下，楼主的算法确实可行。 就楼主的算法而言，每次迭代精度加倍，应该是比较高效的算法了。麻烦的在于计算sin(x) 计算sin(x)是不需要计算出全部精度的，比如当x=3时，计算sin(x)只需要精确2位小数，当x=3.14时，计算sin(x)只需要精确4位小数，当x越接近于pi,需要的sin(x)的精度越高。 但是，这里需要解决一个核心问题，如何计算sin(x)，当x接近于pi时，用级数展开式计算sin(x)时收敛的很慢。