圆周率的几种C语言编程代码

1、计算的方法一、蒙特卡罗法（Monte Carlo）蒙特卡洛方法是一种以概率统计理论为指导的一类非常重要的数值计算方法。蒙特卡罗是摩纳哥公国的一座城市（赌城）。Buffon投针就是蒙特卡洛方法的思想，但是Buffon是蒙特卡洛方法的起源。蒙特卡洛方法同样可以通过随机通过产生随机数的方式来解决计算问题。比如估算和计算不规则图形的面积。图1 蒙特卡罗方法模拟这种方法是一种利用计算机随机数的功能基于“随机数”的算法，通过计算落在单位圆内的点与落在正方形内的点的比值求PI。 图2 蒙特卡罗方法由于图形的对称性，我们靠考虑该图的四分之一部分。 假定一点能够均匀地扔到一个正方形中，计算落入其中的点个数。通过

2、计数其中落入内切圆的点的个数； 如果一共投入N个点，其中有M个落入圆中，则只要点均匀，假定圆周的半径为R，则： 式（1）该方法得到的要得到的精度与投入点的个数有关，一般个数较大时精度比较高。/Monte Carlo Methoddouble rand_pi(int n) int numInCircle = 0;double x, y;double pi;for (int i = 0; i < n; i+)x = rand()*1.0/RAND_MAX;y = rand()*1.0/RAND_MAX;if (x * x + y * y < 1)numInCircle+;pi = (4

3、.0 * numInCircle) / n;return pi;图3 蒙特卡罗示例代码写出来你的代码实际测试结果（10 Marks）表1 蒙特卡罗方法实际测试结果次数10102103104105106107108109计算值运行时间（ms）误差二、数学公式（级数）由数学公式：  式（2）当n取1000时就很接近真实值了n=1000时，pi=3.1406380562059946。/Numerical Series/Sigma(1/n2) = pi/6double math_Pi(int n) int numInCircle = 0;double sum = 0;double pi;fo

4、r (int i = 1; i <= n; i+)sum += 1.0 / (double)i*i);pi = sqrt(sum * 6);return pi;图4 级数方法示例代码写出来你的代码实际测试结果（10 Marks）表2 级数方法实际测试结果次数10102103104105106107108109计算值运行时间（ms）误差三、划分网格计算图5 蒙特卡罗示例代码将图片分为n*n个小方形，统计落在圆内的个数占所有方形的比例。/Grid Countdouble grid_Pi(int n) int i;double sum = 0;for (i = 0; i < n; i+)

5、sum += (int)sqrt(n*(double)n - i*(double)i);return (4.0 * sum) / n / n;图6 网格划分示例代码写出来你的代码实际测试结果（10 Marks）表3 网格划分方法实际测试结果次数10102103104105106107108109计算值运行时间（ms）误差Appendix 时间统计方法一：clock统计在待测试程序的开始部分和结束部分分别加入clock统计，两个之间的差值就是该部分代码运行时间。需要注意的是Windows下和Linux下头文件不一样。clock_tstart;clock_tend;/get current clo

6、ck countstart = clock();/Insert your code here/get current clock countend = clock();printf("The time duration is %d ms.n", end - start);方法二：系统API精确统计时间Windows 和Linux下均有精确计时方法。在Windows下可以用系统提供的API函数 QueryPerformanceFrequency 和 QueryPerformanceCounter 来进行高精度的计时，现在的机器基本上都提供这种高精

7、度的计时啦，所以不用担心。利用该函数可以精确的计时到微妙级别。QueryPerformanceFrequency() 可以得到CPU的时钟频率。QueryPerformanceCounter() 可以通过两次的差值来得到CPU的时钟周期差值。这样就可以计算出两次差值之间所花费的时间，当然这里计算出来的时间是秒，一般都换算为毫秒或者微秒来表示。Windows#ifdef _WIN32LARGE_INTEGER num;long long start, end, freq, time1, time2;#endif/开始计时#ifdef _WIN32QueryPerformanceFrequency

8、(&num);freq = num.QuadPart;QueryPerformanceCounter(&num);start = num.QuadPart;#endif/Insert your code here/插入待测时间代码/ Sleep(1000);/结束计时时间#ifdef _WIN32QueryPerformanceCounter(&num);end = num.QuadPart;time1 = (end - start) * 1000 / freq;printf("nDuration time=%d msn", (end - start

9、) * 1000 / freq);#endif gettimeofday 精确到微秒级别 计时函数 int gettimeofday(structtimeval *tv,struct timezone *tz);说明：gettimeofday将时间保存在结构tv之中，strut timeval gives the number of seconds and microsecondssince the Epoch，(00:00:00 UTC, January 1, 1970), measured inmicroseconds.  tz(即timezone，

10、很少用到)一般我们使用null来代替。保存时间的结构体strut timeval long tv_sec;long tv_usec;注：需要头文件#include <sys/time.h>Linux#ifdef linuxstruct timeval ta, tb;double mseca, msecb;#endif/开始计时#ifdef linuxgettimeofday(&ta, NULL);mseca = ta.tv_sec*1000000.0 + ta.tv_usec;#endif/插入运行代码/Insert your code here/ Sleep(1000);/结