一日在一BC++编程论坛上看到一篇用BC++计算π值的文章，他的代码基本上就是下面的第一段。我看了一下，觉得这段代码还有改进的必要，因为计算π值最关键的是速度问题。在这里我想用VC++6来说明一下我的改进方法，虽然速度快了不少，但还有没有可改进的地文方？怎么改？希望高手给予技术指导。

/***************************************************************************************************************
我加入下面几个　define　是为了改动数据类型方便。
用
　　　#define PART 12
　　　#define DATA 1000000000000
　　　#define TYPE UINT64
定义的方式计算竟然没有用
　　　#define PART 4
　　　#define DATA 10000
　　　#define TYPE UINT
定义的方式计算快，是因为我的系统不是64位的吗？？？
****************************************************************************************************************/
#define JISIWEI 3

#define PART 4
#define DATA 10000
#define TYPE UINT

//#define PART 12
//#define DATA 1000000000000
//#define TYPE UINT64

/*******************************************************************************************************
基本原理：
　　x=z=2;a=1;b=3;
　　while(z*=a/b)
　　{  x+=z;a++;b+=2;
　　}
　　计算结果在 x 中（相关知识去看高数^O^）

　　一个 double 型数值精度很有限，只有12位，要计算更多的位可以这样考：在每个字节中存一位结果，要计算多少位，就用多少个字节来存放结果，然后模拟人手工计算乘、除、加的过程。这样算的好处是容易理解，精度要多高都可，但至命弱点是速度太慢。

未优化算法：
*******************************************************************************************************/
void CPIDlg::OnButton1()
{const int DISPCNT=GetDlgItemInt(IDC_WEISHU,NULL,TRUE)+2;
const int ARRSIZE=DISPCNT+DISPCNT/50+JISIWEI; //计算到ARRSIZE位数，最后若干位不准，所以只显示前DISPCNT位数
  char *x=new char[ARRSIZE], *z=new char[ARRSIZE]; //x[0] x[1] . x[2] x[3] x[4] .... x[ARRSIZE-1]
  int a=1, b=3, c, d,i, Run=1;
  memset(x,0,ARRSIZE);
  memset(z,0,ARRSIZE);

  x[1] = z[1] = 2; // x[1]和z[1]表示右起第一位数，也就是最高位
  DWORD time1=GetTickCount();
   while(Run && (a<200000000))
  { // z *= a/b;把此计算分成两步：先算z*=a，再算z/=b
    d = 0; //z*=a的计算，z[i]中每字节存一位数字，从低位向高位乘。
    for(i=ARRSIZE-1; i>0; i--)
    { c = z[i]*a + d; //当前位乘a再加上低位的进位d
      z[i] = c % 10; //个位放在z[i]中
      d = c / 10; //进位放在d中
    }

    d = 0; //z/=b 的计算，z[i]中每字节存一位数字，从高位向低除。
    for(i=0; i<ARRSIZE; i++)
    { c = z[i]+d*10; //上一位除后的余数乘１０加当前位得到被除数
      z[i] = c / b; // z[i]中放商
      d = c % b; // d中放余数
    }
  
    Run = 0; //x+=z的计算，从低位向高位计算
    for(i=ARRSIZE-1; i>0; i--)
    { c = x[i] + z[i]; //加后和放在c中
      x[i] = c%10; //和的个位放在当前位置上
      x[i-1] += c/10; //和的十位与高一位置上的数相加
      Run |= z[i]; //只要z不为零则继续计算
    }
    a++;
    b+=2;
  }
   CString ss;
   ss.Format("未优化算法\n用时%g秒",(float)(::GetTickCount()-time1)/1000);
   SetDlgItemText(IDC_BUTTON1,ss);
   ss.Format("(未优化)%d",a-1);
   SetDlgItemText(IDC_TIMES,ss);
   for(x[0]=x[1]+0x30,x[1]='.',i=2; i<DISPCNT; i++)
                x[i]+=0x30;
   x[i]=x[i+1]=0;
        SetDlgItemText(IDC_OUTPUT,x);
        delete[] x;delete[] z;
}

/*******************************************************************************************************
　　上面的方法中有这样一个缺点：当 z 值很小时，z 值的前面许多位（高位）已经为0了，但是我们还在按步就班地做乘、除、加这样的无用功。对已经为0的高位，我们完全可以不用再管了，因为 z 是越来越小的。这样做可以省一半的时间。

优化算法：
*******************************************************************************************************/
void CPIDlg::OnButton2()
{ const int DISPCNT=GetDlgItemInt(IDC_WEISHU,NULL,TRUE)+2;
const int ARRSIZE=DISPCNT+DISPCNT/50+JISIWEI; //计算到ARRSIZE位数，最后若干位不准，所以只显示前DISPCNT位数
  char c,*x=new char[ARRSIZE], *z=new char[ARRSIZE+1],*ix,*iz,*ie,*i0; //x[0] x[1] . x[2] x[3] x[4] .... x[ARRSIZE-1]
  int a=1, b=3, d;
  memset(x,0,ARRSIZE);
  memset(z,0,ARRSIZE);
  *(x+1) = *(i0=z+1) = 2; // x[1]和z[1]表示右起第一位数，也就是最高位
  z[ARRSIZE]=10;
  DWORD time1=GetTickCount();
  do
  {        // z *= a/b;把此计算分成两步：先算z*=a，再算z/=b
        //z*=a的计算，z[i]中每字节存一位数字，从低位向高位乘。
        for(d=0,iz=z+ARRSIZE-1;i0<=iz;iz--)
    { *iz = (d += (*iz)*a) % 10; //当前位乘a再加上低位的进位d//个位放在z[i]中
      d /= 10; //进位放在d中
    }
        while(d)
        { *iz-- = d % 10; //个位放在z[i]中
     d /= 10; //进位放在d中
        }

    //d = 0; //z/=b 的计算，z[i]中每字节存一位数字，从高位向低除。
    for(i0=++iz,ie=z+ARRSIZE; iz<ie; iz++)
    { //上一位除后的余数乘１０加当前位得到被除数
      *iz = (d += *iz) / b; // z[i]中放商
      d = (d % b)*10; // d中放余数
    }

        while(*i0==0)i0++;//在z[i]中找到首位非０数//若z[i]中每位都为０则结束
        //if(i0>iz)break;
       
        c=0; //x+=z的计算，从低位向高位计算
    for(ix=x+ARRSIZE-1,ie=iz=z+ARRSIZE-1; i0<=iz; iz--,ix--)
    { *ix=(c+=*iz+*ix)%10; //加后和放在c中//和的个位放在当前位置上
      c/=10; //c中为进位数１或０
    }
        while((*ix+=c)>9)*ix--=0;

    a++;b+=2;
  }while(i0<=ie);
   CString ss;
   ss.Format("优化算法\n用时%g秒",(float)(::GetTickCount()-time1)/1000);
   SetDlgItemText(IDC_BUTTON2,ss);
   ss.Format("(优化)%d",a-1);
   SetDlgItemText(IDC_TIMES,ss);
   for(*x=x[1]+0x30,x[1]='.',ix=x+2,ie=x+DISPCNT; ix<ie; ix++)
        *ix+=0x30;
   *(int*)ix=0;
   SetDlgItemText(IDC_OUTPUT,x);
   delete[] x;delete[] z;
}

/*******************************************************************************************************
　　上面的方法还可以再改进：用一个字节存一位结果，不但浪费内存，而且速度也慢。因为“字节运算”一次，只算了一位，如果计算一次能多算几位，就可以省更多的时间了。
　　我的系统是32位的，每次可以算32的二进制数，也就是四字节，最大为0xffffffff,也就是4294967295。这样我们可以用象DWORD这样的数据类型保存计算结果，一次可以多算几位。这里我让每四字节保存四位计算结果，你也许要问
我：为什么不让它保存六位、七位甚至是九位的结果呢？哈哈，因为我知道当计算比较高时，算到一定时候 a 会很大，那样的话 z[i]*a 的结果就会大过4294967295而溢出出错。每次计算四位，虽然未节省内存，时间约为原来的四分之一。

改良算法：（精度不能太高~_~）
*******************************************************************************************************/
void CPIDlg::OnButton3()
{ const TYPE DISPCNT=GetDlgItemInt(IDC_WEISHU,NULL,TRUE)+2;
  const TYPE ARRSIZE=(DISPCNT+DISPCNT/100+JISIWEI+PART-1)/PART+1; //计算到ARRSIZE位数，最后若干位不准，所以只显示前DISPCNT位

数
  TYPE c,*x=new TYPE[ARRSIZE], *z=new TYPE[ARRSIZE+1],*ix,*iz,*ie,*i0; //x[0] x[1] . x[2] x[3] x[4] .... x[ARRSIZE-1]
  TYPE a=1, b=3, d;
　CString ss;
  memset(x,0,ARRSIZE*sizeof(TYPE));
  memset(z,0,(ARRSIZE+1)*sizeof(TYPE));
  *(x+1) = *(i0=z+1) = 2*DATA; // x[1]和z[1]表示右起第一位数，也就是最高位
  z[ARRSIZE]=10;
  DWORD time1=GetTickCount();
  do
  {        // z *= a/b;把此计算分成两步：先算z*=a，再算z/=b
        //z*=a的计算，z[i]中每字节存一位数字，从低位向高位乘。
        for(d=0,iz=z+ARRSIZE-1;i0<=iz;iz--)
    { *iz = (d += (*iz)*a) % DATA; //当前位乘a再加上低位的进位d//个位放在z[i]中
      d /= DATA; //进位放在d中
    }
        while(d)
        { *iz-- = d % DATA; //个位放在z[i]中
     d /= DATA; //进位放在d中
        }

    //d = 0; //z/=b 的计算，z[i]中每字节存一位数字，从高位向低除。
    for(i0=++iz,ie=z+ARRSIZE; iz<ie; iz++)
    { //上一位除后的余数乘１０加当前位得到被除数
      *iz = (d += *iz) / b; // z[i]中放商
      d = (d % b)*DATA; // d中放余数
    }

        while(*i0==0)i0++;//在z[i]中找到首位非０数//若z[i]中每位都为０则结束
        //if(i0>iz)break;
       
        c=0; //x+=z的计算，从低位向高位计算
    for(ix=x+ARRSIZE-1,ie=iz=z+ARRSIZE-1; i0<=iz; iz--,ix--)
    { *ix=(c+=*iz+*ix)%DATA; //加后和放在c中//和的个位放在当前位置上
      c/=DATA; //c中为进位数１或０
    }
        while((*ix+=c)>=DATA)
        {        c=*ix/DATA;*ix%=DATA;ix--;}

    a++;b+=2;
  }while(i0<=ie);

   //计算完毕，准备显示输出
   ss.Format("改良算法\n用时%g秒",(float)(::GetTickCount()-time1)/1000);
   SetDlgItemText(IDC_BUTTON3,ss);
   ss.Format("(改良)%d",a-1);
   SetDlgItemText(IDC_TIMES,ss);

   delete[] z; int j;
   char *cp=new char[DISPCNT+50],*cc=cp;
   memset(cp,0x30,DISPCNT+50);

   *cc+=*x;cc++;*cc='.';//cc++;
   for(ix=x+1,ie=x+(DISPCNT+PART-1)/PART;ix<ie;ix++)
   {   j=PART;
           do
           {        *(cc+j--)+=*ix%10;
           }while(*ix/=10);
           cc+=PART;
   }
        cp[DISPCNT]=0;
        SetDlgItemText(IDC_OUTPUT,cp);
        delete[] x;delete[] cp;
}
　　后面有我的测试EXE文件，有兴趣的去看看。

[培训]内核驱动高级班，冲击BAT一流互联网大厂工作，每周日13:00-18:00直播授课

上传的附件：

• 计算pi值.rar （50.08kb，5次下载）