Skip to content

Commit fd6020e

Browse files
authored
adding draft of monte carlo chapter (#115)
* adding draft of monte carlo chapter. * removing unnecessary tag language tag * changing images to gif. * fixing typo * fixing typo. * fixing a few mroe smallscale typos
1 parent 19db513 commit fd6020e

File tree

8 files changed

+128
-0
lines changed

8 files changed

+128
-0
lines changed

SUMMARY.md

Lines changed: 1 addition & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -28,6 +28,7 @@
2828
* [Tree Traversal](chapters/tree_traversal/tree_traversal.md)
2929
* [Euclidean Algorithm](chapters/euclidean_algorithm/euclidean.md)
3030
* [Multiplication](chapters/multiplication/multiplication.md)
31+
* [Monte Carlo](chapters/monte_carlo/monte_carlo.md)
3132
* [Matrix Methods](chapters/matrix_methods/matrix_methods.md)
3233
* [Gaussian Elimination](chapters/matrix_methods/gaussian_elimination/gaussian_elimination.md)
3334
* [Thomas Algorithm](chapters/matrix_methods/thomas/thomas.md)
Lines changed: 27 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,27 @@
1+
# function to determine whether an x, y point is in the unit circle
2+
function in_circle(x_pos::Float64, y_pos::Float64, radius::Float64)
3+
if (x_pos^2 + y_pos^2 < radius^2)
4+
return true
5+
else
6+
return false
7+
end
8+
end
9+
10+
# function to integrate a unit circle to find pi via monte_carlo
11+
function monte_carlo(n::Int64, radius::Float64)
12+
13+
pi_count = 0
14+
for i = 1:n
15+
point_x = rand()
16+
point_y = rand()
17+
18+
if (in_circle(point_x, point_y, radius))
19+
pi_count += 1
20+
end
21+
end
22+
23+
pi_estimate = 4*pi_count/(n*radius^2)
24+
println("Percent error is: ", signif(100*(pi - pi_estimate)/pi, 3), " %")
25+
end
26+
27+
monte_carlo(10000000, 0.5)

chapters/monte_carlo/monte_carlo.md

Lines changed: 100 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,100 @@
1+
# Monte Carlo Integration
2+
3+
Monte Carlo methods were some of the first methods I ever used for research, and when I learned about them, they seemed like some sort of magic.
4+
Their premise is simple: random numbers can be used to integrate arbitrary shapes embedded into other objects.
5+
Nowadays, "monte carlo" has become a bit of a catch-all term for methods that use random numbers to produce real results, but it all started as a straightforward method to integrate objects.
6+
No matter how you slice it, the idea seems a bit crazy at first.
7+
After all, random numbers are random.
8+
How could they possibly be used to find non-random values?
9+
10+
Well, imagine you have a square.
11+
The area of the square is simple, $$\text{Area}_{\text{square}} = \text{length} \times \text{width}$$.
12+
Since it's a square, the $$\text{length}$$ and $$\text{width}$$ are the same, so the formula is technically just $$\text{Area}_{\text{square}} = \text{length}^2$$.
13+
If we embed a circle into the square with a radius $$r = \text{length}$$ (shown below), then it's area is $$\text{Area}_{\text{circle}}=\pi r^2$$.
14+
For simplicity, we can also say that $$\text{Area}_{\text{square}}=4r^2$$.
15+
16+
<p align="center">
17+
<img src="res/square_circle.png" width="300"/>
18+
</p>
19+
20+
Now, let's say we want to find the area of the circle without an equation.
21+
As we said before, it's embedded in the square, so we should be able to find some ratio of the area of the square to the area of the circle:
22+
23+
$$
24+
\text{Ratio} = \frac{\text{Area}_{\text{circle}}}{\text{Area}_{\text{square}}}
25+
$$
26+
27+
This means,
28+
29+
$$
30+
\text{Area}_{\text{circle}} = \text{Area}_{\text{square}}\times\text{Ratio} = 4r^2 \times \text{ratio}
31+
$$
32+
33+
So, if we can find the $$\text{Ratio}$$ and we know $$r$$, we should be able to easily find the $$\text{Area}_{\text{circle}}$$.
34+
The question is, "How do we easily find the $$\text{Ratio}$$?"
35+
Well, one way is with *random sampling*.
36+
We basically just pick a bunch of points randomly in the square, and
37+
each point is tested to see whether it's in the circle or not:
38+
39+
{% method %}
40+
{% sample lang="jl" %}
41+
[import:2-8, lang:"julia"](code/julia/monte_carlo.jl)
42+
{% endmethod %}
43+
44+
If it's in the circle, we increase an internal count by one, and in the end,
45+
46+
$$
47+
\text{Ratio} = \frac{\text{count in circle}}{\text{total number of points used}}
48+
$$
49+
50+
If we use a small number of points, this will only give us a rough approximation, but as we start adding more and more points, the approximation becomes much, much better (as shown below)!
51+
52+
<p align="center">
53+
<img src="res/monte_carlo.gif" width="400"/>
54+
</p>
55+
56+
The true power of monte carlo comes from the fact that it can be used to integrate literally any object that can be embedded into the square.
57+
As long as you can write some function to tell whether the provided point is inside the shape you want (like `in_circle()` in this case), you can use monte carlo integration!
58+
This is obviously an incredibly powerful tool and has been used time and time again for many different areas of physics and engineering.
59+
I can gaurantee that we will see similar methods crop up all over the place in the future!
60+
61+
# Example Code
62+
Monte carlo methods are famous for their simplicity.
63+
It doesn't take too many lines to get something simple going.
64+
Here, we are just integrating a circle, like we described above; however, there is a small twist and trick.
65+
Instead of calculating the area of the circle, we are instead trying to find the value of $$\pi$$, and
66+
rather than integrating the entire circle, we are only integrating the upper right quadrant of the circle from $$0 < x,y < 1$$.
67+
This saves a bit of computation time, but also requires us to multiply our output by $$4$$.
68+
69+
That's all there is to it!
70+
Feel free to submit your version via pull request, and thanks for reading!
71+
72+
{% method %}
73+
{% sample lang="jl" %}
74+
### Julia
75+
[import, lang:"julia"](code/julia/monte_carlo.jl)
76+
{% endmethod %}
77+
78+
79+
<script>
80+
MathJax.Hub.Queue(["Typeset",MathJax.Hub]);
81+
</script>
82+
$$
83+
\newcommand{\d}{\mathrm{d}}
84+
\newcommand{\bff}{\boldsymbol{f}}
85+
\newcommand{\bfg}{\boldsymbol{g}}
86+
\newcommand{\bfp}{\boldsymbol{p}}
87+
\newcommand{\bfq}{\boldsymbol{q}}
88+
\newcommand{\bfx}{\boldsymbol{x}}
89+
\newcommand{\bfu}{\boldsymbol{u}}
90+
\newcommand{\bfv}{\boldsymbol{v}}
91+
\newcommand{\bfA}{\boldsymbol{A}}
92+
\newcommand{\bfB}{\boldsymbol{B}}
93+
\newcommand{\bfC}{\boldsymbol{C}}
94+
\newcommand{\bfM}{\boldsymbol{M}}
95+
\newcommand{\bfJ}{\boldsymbol{J}}
96+
\newcommand{\bfR}{\boldsymbol{R}}
97+
\newcommand{\bfT}{\boldsymbol{T}}
98+
\newcommand{\bfomega}{\boldsymbol{\omega}}
99+
\newcommand{\bftau}{\boldsymbol{\tau}}
100+
$$

chapters/monte_carlo/res/13311.png

177 KB
Loading

chapters/monte_carlo/res/195583.png

154 KB
Loading

chapters/monte_carlo/res/31.png

39.2 KB
Loading
1.6 MB
Loading
6.76 KB
Loading

0 commit comments

Comments
 (0)