使用蒙特卡洛估算Pi的值

蒙特卡洛估计
蒙特卡洛方法是一类广泛的计算算法，它们依赖于重复的随机采样来获得数值结果。蒙特卡洛算法入门的基本示例之一是Pi的估计。

Pi的估计
这个想法是模拟2D平面中的随机(x，y)点，其域为边1单位的平方。想象一个在相同区域内的，直径相同，并刻在正方形上的圆。然后，我们计算在圆内撒谎的点数与生成的点总数之比。请参考下图：

生成随机点，其中只有很少一部分位于假想圆之外

我们知道正方形的面积是1单位平方，而圆形的面积是$\pi \ast (\frac{1}{2})^{2} = \frac{\pi}{4}$ 。现在，对于大量生成的点，

$\frac{\textrm{area of the circle}}{\textrm{area of the square}} = \frac{\textrm{no. of points generated inside the circle}}{\textrm{total no. of points generated or no. of points generated inside the square}}$那是， $\pi = 4 \ast \frac{\textrm{no. of points generated inside the circle}}{\textrm{no. of points generated inside the square}}$

该算法的优点在于，我们不需要任何图形或模拟即可显示生成的点。我们只需生成随机(x，y)对，然后检查是否$x^{2} + y^{2} \leqslant 1$ 。如果是，我们增加出现在圆内的点数。在像蒙特卡洛这样的随机和仿真算法中，迭代次数越多，结果越准确。因此，标题是“估计Pi的值”而不是“计算Pi的值”。下面是该方法的算法：

算法
1.将circle_points，square_points和间隔初始化为0。
2.生成随机点x。
3.生成随机点y。
4.计算d = x * x + y * y。
5.如果d <= 1，则增加circle_points。
6.递增square_points。
7.递增间隔。
8.如果增量 9.计算pi = 4 *(circle_points / square_points)。
10.终止。

该代码不等待通过stdin进行的任何输入，因为可以根据所需的迭代次数更改宏INTERVAL。迭代次数是INTERVAL的平方。另外，我已经使用getch()暂停了前10次迭代的屏幕，并且每次迭代的输出都以下面给出的格式显示。您可以根据需要更改或删除它们。

x y circle_points square_points - pi 

例子：

INTERVAL = 5
Output : Final Estimation of Pi = 2.56

INTERVAL = 10
Output : Final Estimation of Pi = 3.24

INTERVAL = 100
Output : Final Estimation of Pi = 3.0916

C++

/* C++ program for estimation of Pi using Monte
   Carlo Simulation */
#include 
  
// Defines precision for x and y values. More the
// interval, more the number of significant digits
#define INTERVAL 10000
using namespace std;
  
int main()
{
    int interval, i;
    double rand_x, rand_y, origin_dist, pi;
    int circle_points = 0, square_points = 0;
  
    // Initializing rand()
    srand(time(NULL));
  
    // Total Random numbers generated = possible x
    // values * possible y values
    for (i = 0; i < (INTERVAL * INTERVAL); i++) {
  
        // Randomly generated x and y values
        rand_x = double(rand() % (INTERVAL + 1)) / INTERVAL;
        rand_y = double(rand() % (INTERVAL + 1)) / INTERVAL;
  
        // Distance between (x, y) from the origin
        origin_dist = rand_x * rand_x + rand_y * rand_y;
  
        // Checking if (x, y) lies inside the define
        // circle with R=1
        if (origin_dist <= 1)
            circle_points++;
  
        // Total number of points generated
        square_points++;
  
        // estimated pi after this iteration
        pi = double(4 * circle_points) / square_points;
  
        // For visual understanding (Optional)
        cout << rand_x << " " << rand_y << " " << circle_points
             << " " << square_points << " - " << pi << endl << endl;
  
        // Pausing estimation for first 10 values (Optional)
        if (i < 20)
            getchar();
    }
  
    // Final Estimated Value
    cout << "\nFinal Estimation of Pi = " << pi;
  
    return 0;
}

Python

import random
  
INTERVAL= 1000
  
circle_points= 0
square_points= 0
  
# Total Random numbers generated= possible x
# values* possible y values
for i in range(INTERVAL**2):
  
    # Randomly generated x and y values from a
    # uniform distribution
    # Rannge of x and y values is -1 to 1
    rand_x= random.uniform(-1, 1)
    rand_y= random.uniform(-1, 1)
  
    # Distance between (x, y) from the origin
    origin_dist= rand_x**2 + rand_y**2
  
    # Checking if (x, y) lies inside the circle
    if origin_dist<= 1:
        circle_points+= 1
  
    square_points+= 1
  
    # Estimating value of pi,
    # pi= 4*(no. of points generated inside the 
    # circle)/ (no. of points generated inside the square)
    pi = 4* circle_points/ square_points
  
##    print(rand_x, rand_y, circle_points, square_points, "-", pi)
##    print("\n")
  
print("Final Estimation of Pi=", pi)

输出：

Final Estimation of Pi = 3.16116