求对圆周率代码进行优化,该如何处理
求对圆周率代码进行优化
设计多种解法计算圆周率π，并进行方法比较。例如可以考虑下述方法：（1）随机数法，思路是取一个边长为1的单位正方形，在其中做它的内切圆，再向正方形内扔点，点落在圆内则计数，落在圆外不计数。扔到5000个点后停止，用落入圆内的点数的4倍除以总的扔的点数，就得到π的一个近似值。（2）用我国古代数学家祖冲之的方法，即用圆内接正多边形逼近。可以从圆内接正六边形出发，迭代计算园内接正12、24…… 边形的边长。（3）采用级数：π/2 = 1 + 1/3 + (1*2)/(3*5) +(1*2*3)/(3*5*7) + …… +(1*2*…n)/(3*5*…(2n+1)) =1+1/3*(1+2/5*(1+……+(n-1)/(2n-1)*(1+n/(2n+1)......)。说明解法（3）最有利于高精度计算，例如你用这一解法编写的的程序应该能够求出圆周率π精确到小数点后100位，求出 圆周率π = 3. 1415926535 8979323846 2643383279 5028841971 6939937510 5820974944 5923078164 0628620899 8628034825 3421170677。
提示：（2）可以推导得出若圆内接正i边形边长的一半为xi , 则正2i边形边长的一半是： x(2i)={根号下（2-2根号下（1-xi*xi））}/2
#include<stdio.h>
#include<time.h>
#include<stdlib.h>
#include<math.h>
#define N 5000
void main()
{
void jishu();
void diedai();
void suiji();
suiji();
diedai();
jishu();
}
void suiji()
{
double x,y;
int a=0,b=0;
srand(time(0));
while(a++<=N)//投5000次点
{
x=(double)rand()/RAND_MAX;//产生0~1之间的浮点数
y=(double)rand()/RAND_MAX;//产生0~1之间的浮点数
if((x-0.5)*(x-0.5)+(y-0.5)*(y-0.5)<=0.25)//判断所产生的点是否在圆内
b++;//汇总落在圆内的点数
}
printf("用随机数法求得pi=%lf\n",4.0*b/N);
}
void diedai()
{
double e=0.1,b=0.5,d;
long int i;//正多边形边数
for(i=6;;i=i*2)
{
d=1.0-sqrt(1.0-b*b);
b=0.5*sqrt(b*b+d*d);
if(2*i*b-i*e<1e-15)
break;
e=b;
}
printf("祖冲之多边形迭代所求得的pi值为：%.15lf\n",2*i*b);
}
void jishu()
{
float s;
int x,b,c,i,j,n,e,a[4000];
for(s=0,n=1;n<=4000;n++)//累加确定计算的项数
{
s=s+log((2*n+1)/n);
if(s/log(10)>100) break;
}
for(i=0;i<=101;i++) a[0]=0;
for(x=1,b=n;b>=1;b--)//按公式分布计算
{
c=2*b+1;
for(i=0;i<=100;i++)
{
a[i]=x/c; x=(x%c)*10+a[i+1];
a[101]=x/c;
}
for(j=0,i=101;i>=0;i--)
{
a[i]=a[i]*b+j;
j=a[i]/10;
a[i]=a[i]%10;
}
a[0]=a[0]+1;x=a[0];
}
for(j=0,i=101;i>=0;i--)
{
a[i]=a[i]*2+j;j=a[i]/10;a[i]=a[i]%10;
}
printf("用级数方法求得pi=%d.",a[0]);//输出结果
for(e=10,i=1;i<=100;i++)
{
printf("%d",a[i]);e++;
if(e%30==0) printf("\n ");
}
printf("\n");
}
------解决方案--------------------
设计多种解法计算圆周率π，并进行方法比较。例如可以考虑下述方法：（1）随机数法，思路是取一个边长为1的单位正方形，在其中做它的内切圆，再向正方形内扔点，点落在圆内则计数，落在圆外不计数。扔到5000个点后停止，用落入圆内的点数的4倍除以总的扔的点数，就得到π的一个近似值。（2）用我国古代数学家祖冲之的方法，即用圆内接正多边形逼近。可以从圆内接正六边形出发，迭代计算园内接正12、24…… 边形的边长。（3）采用级数：π/2 = 1 + 1/3 + (1*2)/(3*5) +(1*2*3)/(3*5*7) + …… +(1*2*…n)/(3*5*…(2n+1)) =1+1/3*(1+2/5*(1+……+(n-1)/(2n-1)*(1+n/(2n+1)......)。说明解法（3）最有利于高精度计算，例如你用这一解法编写的的程序应该能够求出圆周率π精确到小数点后100位，求出 圆周率π = 3. 1415926535 8979323846 2643383279 5028841971 6939937510 5820974944 5923078164 0628620899 8628034825 3421170677。
提示：（2）可以推导得出若圆内接正i边形边长的一半为xi , 则正2i边形边长的一半是： x(2i)={根号下（2-2根号下（1-xi*xi））}/2
#include<stdio.h>
#include<time.h>
#include<stdlib.h>
#include<math.h>
#define N 5000
void main()
{
void jishu();
void diedai();
void suiji();
suiji();
diedai();
jishu();
}
void suiji()
{
double x,y;
int a=0,b=0;
srand(time(0));
while(a++<=N)//投5000次点
{
x=(double)rand()/RAND_MAX;//产生0~1之间的浮点数
y=(double)rand()/RAND_MAX;//产生0~1之间的浮点数
if((x-0.5)*(x-0.5)+(y-0.5)*(y-0.5)<=0.25)//判断所产生的点是否在圆内
b++;//汇总落在圆内的点数
}
printf("用随机数法求得pi=%lf\n",4.0*b/N);
}
void diedai()
{
double e=0.1,b=0.5,d;
long int i;//正多边形边数
for(i=6;;i=i*2)
{
d=1.0-sqrt(1.0-b*b);
b=0.5*sqrt(b*b+d*d);
if(2*i*b-i*e<1e-15)
break;
e=b;
}
printf("祖冲之多边形迭代所求得的pi值为：%.15lf\n",2*i*b);
}
void jishu()
{
float s;
int x,b,c,i,j,n,e,a[4000];
for(s=0,n=1;n<=4000;n++)//累加确定计算的项数
{
s=s+log((2*n+1)/n);
if(s/log(10)>100) break;
}
for(i=0;i<=101;i++) a[0]=0;
for(x=1,b=n;b>=1;b--)//按公式分布计算
{
c=2*b+1;
for(i=0;i<=100;i++)
{
a[i]=x/c; x=(x%c)*10+a[i+1];
a[101]=x/c;
}
for(j=0,i=101;i>=0;i--)
{
a[i]=a[i]*b+j;
j=a[i]/10;
a[i]=a[i]%10;
}
a[0]=a[0]+1;x=a[0];
}
for(j=0,i=101;i>=0;i--)
{
a[i]=a[i]*2+j;j=a[i]/10;a[i]=a[i]%10;
}
printf("用级数方法求得pi=%d.",a[0]);//输出结果
for(e=10,i=1;i<=100;i++)
{
printf("%d",a[i]);e++;
if(e%30==0) printf("\n ");
}
printf("\n");
}
------解决方案--------------------
- C/C++ code
#include <stdio.h>
long a=10000;
long b;
long c=2800;
long d;
long e;
long f[2801];
long g;
int main() {
for(;b-c;) f[b++]=a/5;
// while (1) {
// if (0==b-c) break;
// f[b]=a/5;
// b++;
// }
//f[0 - 2800] = 10000/5
for(;d=0,g=c*2;c-=14,printf("%.4d",e+d/a),e=d%a)
for(b=c;d+=f[b]*a,f[b]=d%--g,d/=g--,--b;d*=b);
// while (1) {
// d=0;
// g=c*2;
// if (0==g) break;
// b=c;
// while (1) {
// d+=f[b]*a;
// f[b]=d%--g;
// d/=g--;
// --b;
// if (0==b) break;
// d*=b;
// }
// c-=14;
// printf("%.4d",e+d/a);
// e=d%a;
// }
return 0;
}