Introduction

Hello nice to meet. This time, I started writing with Qiita. Nice to meet you. Well, I don't know what the first article is, but I decided to write an article about how to implement the Monte Carlo method in Python.

What is the Monte Carlo method?

What kind of method is the Monte Carlo method in the first place? As many of you may know, I will explain it once again. The Monte Carlo method in numerical analysis is often used as a method for finding the probability approximately. You can also find the area surrounded by a closed surface in two dimensions. This time, I will try to find the pi by the Monte Carlo method.

How to find pi by Monte Carlo method

First, suppose you have something like the figure below. The blue line represents the circumference of a circle with a radius of 1. Randomly place points here and count the number of points ($ n ) and the number of points ( p $) placed in the circumference. Then, the following formula holds from the Monte Carlo method.

\frac{r^2\pi}{4S}=\frac{p}{n}

$ S $ in the formula is the area of the square. Since the ratio of $ p $ and $ n $ is the same as the area ratio, the above formula holds. The formula transformation of $ \ pi $ from the above formula is as follows.

\pi=\frac{4Sp}{nr^2}

Of course, the probability is not absolute, so we will only find an approximate value.

Implemented in Python

Let's calculate with a personal computer. As the title says, I implemented it in Python. Below is the source code.

`pi_by_monte.py`


import random

#Variable setting
in_pi = 0

#Number of times to execute
n = 100000

for i in range(n):
    x = random.uniform(0, 1)
    y = random.uniform(0, 1)

    z = x**2 + y**2

    if z <= 1:
        in_pi += 1

pi_by_monte = in_pi / n * 4

print(pi_by_monte)

You could implement it with a surprisingly short code. In the source code, it is executed 100000 times, but in reality, it needs to be executed more.

result

Now let's see the result of the execution.

#Results will vary from time to time
3.14885

For the time being, we succeeded in deriving up to 3.14. But after that it doesn't fit. As mentioned above, the probability is not absolute, so it is just an approximation. I increased the number of executions to about 100 million and executed it.

#Results will vary from time to time
3.14180064

Now it fits up to 3.141. Theoretically, if you increase the number of executions as it is, it will be closer to the circumference ratio, but it will take ** time ** very much. If it was about 100,000 times, it would come out immediately, but if it was done 100 million times, it took more than 3 minutes. It is obvious that the time will increase steadily if we continue to increase it to 1 billion and 10 billion. I think there is a way to make it faster, but this is the limit for this program. By the way, I ran it 10000 times and used matplotlib to actually plot the points. (Figure below) I plotted the inside of the circumference in red and the outside points in blue. The black line is the circumference. Ideally, all of this will be filled, but it will consume a lot of memory and the time to process will be ridiculous.

Summary

This time, I used Python to find the pi by the Monte Carlo method. It's easy to find, but I think it's said that there are some drawbacks in accuracy and processing speed. However, in reality, it will be "What is the Monte Carlo method?" And "Monte Carlo Rosge". Please try it once. If you have any comments or questions, please leave a comment.

Try implementing the Monte Carlo method in Python