现在的位置: 首页 > 综合 > 正文

圆周率的精确计算方法

2013年05月05日 ⁄ 综合 ⁄ 共 1072字 ⁄ 字号 评论关闭
精确计算圆周率的原理实际上利用Machin公式最有效。

Machin公式:




这个公式由英国天文学教授约翰·马青于1706年发现。

马青公式每计算一项可以得到1.4位的十进制精度。因为它的计算过程中被乘数和被除数都不大于长整数,所以可以很容易地

在计算机上编程实现。还有很多类似于马青公式的反正切公式。在所有这些公式中,马青公式似乎是最快的了。


虽然如此,如果要计算更多的位数,就把两大数的乘除运算用FFT实现即可。

#include <iostream>
#include <string.h>
#include <stdio.h>

using namespace std;

int main()
{
    int const N=7200;
    int const M=10000;
    int const B=10000;
    int const L=4;
    // Machin公式 计算PI到一万位
    int s[M/L];
    int r1[N]={0},r2[N]={0},d1[N]={0},d2;
    int r3[N]={0},r4[N]={0},d3[N]={0},d4;
    int i,k,t,p=0,mp=M/L/20;

    r1[0]=1;
    r1[1]=3;
    r3[0]=4;

    printf("正在计算,请等待!\n");

    for(k=0;k<M/L;++k)
    {
        t=r1[0]*B;
        d1[0]=t/0x5;
        r1[0]=t%0x5;

        t=r3[0]*B;
        d3[0]=t/0xEF;
        r3[0]=t%0xEF;
        s[k]=d1[0]-d3[0];
        int tag=0;
        for(i=1;i<N;++i)
        {
            t=r1[i]*B+d1[i-1];
            d1[i]=t/0x19;
            r1[i]=t%0x19;
            t=r2[i]*B+d1[i];
            d2=t/(2*i+1);
            r2[i]=t%(2*i+1);

            t=r3[i]*B+d3[i-1];
            d3[i]=t/0xDF21;
            r3[i]=t%0xDF21;
            t=r4[i]*B+d3[i];
            d4=t/(2*i+1);
            r4[i]=t%(2*i+1);
            if(tag)
            {
                s[k]+=(d2-d4);
                tag=0;
            }
            else
            {
                s[k]+=(d4-d2);
                tag=1;
            }
        }
        if(p==mp)
        {
            printf(">");
            p=0;
        }
        else
            p++;
    }
    for(i=M/L-1;i>=0;i--)
    {
        while(s[i]>=B)
        {
            s[i-1]++;
            s[i]-=B;
        }
        while(s[i]<0)
        {
            s[i-1]--;
            s[i]+=B;
        }
    }
    printf("\nPI=3.\n");
    for(i=0;i<M/L;++i)
        printf("%04d",s[i]);
    return 0;
}

抱歉!评论已关闭.