שמות:

In [6]:
from IPython.display import HTML
import matplotlib.pyplot as plt

<font color = blue>
חישוב שטחים בשיטת מונטה קרלו
</font>

<font color = blue>

בתרגיל זה ניישם שיטה חישובית המכונה מונטה קרלו. בעזרת שיטה זו נחשב בהמשך שטח של אליפסה. אבל בשלב ראשון ניישם את השיטה המתוארת בהמשך לחישוב $\pi$.
</font>

<font color = green>
קראו בעיון את הקטע שבהמשך

In [7]:
HTML('''<iframe src="https://he.wikipedia.org/wiki/%D7%A9%D7%99%D7%98%D7%AA_%D7%9E%D7%95%D7%A0%D7%98%D7%94_%D7%A7%D7%A8%D7%9C%D7%95"
     width=900 height=700"></iframe>''')
Out[7]:

<font color = green>

בהמשך קטע קוד המיישם את האלגוריתם המתואר למעלה למציאת $\pi$. הוסיפו ליד כל סימן סולמית הערה המסבירה את תפקיד משפט הקוד שמשמאלה.
</font>

In [5]:
from random import uniform #הבאת פונקציה שתביא מספר רנדומאלי לא שלם בתחום שנכתב
N=10000#מספר הנסיונות\נקודות
count=0#מספר הקליעות
for i in range(N):
    x=uniform(-1,1)#הגרלה של שיעור ה-איקס בין 1- ל1
    y=uniform(-1,1)#הגרלה של שיעור הוואי בין 1- ל1
    if x**2+y**2<=1:# בדיקה אם הוא נמצא בתוך או על המעגל
        count=count+1
float(count)/N*4
Out[5]:
3.1408

<font color = blue>

חישוב $\pi$ בשיטה המתוארת כאן מכונה חישוב בדגימה ישירה
</font>

<font color = green>

הריצו את הפונקציה תוך כדי הגדלת ערכי N מ-10 ועד ל- 10000. מה המסקנה? הסבירו במה במה תלוי הדיוק בחישוב $\pi$?
</font>

In [7]:
from random import uniform
def pi_Evaluation(N):
    count=0
    for i in range(N):
        x=uniform(0,1)
        y=uniform(0,1)
        if x**2+y**2<=1:
            count=count+1
    return float(count)/N*4
pi_Evaluation(100000)
#הדיוק בחישוב פאי תלוי במספר ההרצות
Out[7]:
3.13716

<font color =green>

מידת הפיזור של תוצאות "ניסוי" מכונה סטיית התקן בסידרת התרגילים שבהמשך תכתבו פונקציות שמטרתן לחשב את סטיית התקן המתקבלת בחישוב $\pi$ בדגימה ישירה. בשלב ראשון עליכם לכתוב פונקציה המריצה עבור אותו N את הפונקציה pi_Evaluation מספר פעמים המועבר אליה כפרמטר. הפונקציה צריכה להחזיר רשימה של כל ערכי ה-$\pi$ שהתקבלו.

השלימו את קטע הקוד שבהמשך.

</div><

In [9]:
def get_trails(k,N):
    pi_list=[pi_Evaluation(N) for i in range(k)]
    return pi_list
get_trails(5,10000)
Out[9]:
[3.1468, 3.1424, 3.13, 3.12, 3.1432]

<font color = green>

שורת הפקודות שלמטה משרטט תרשים שכיחויות של תוצאות של הרצות שונות עבור אותו גודל של N.הריצו את הקוד עבור ערכי N שונים ובחנו כיצד משתנה פיזור התוצאות.

In [12]:
import matplotlib.pyplot as plt
N=5000
trails=10000
o=plt.hist(get_trails(trails,N),bins=30)
Out[12]:

<font color = green>

העזרו בחבילה אותה פיתחתם לחישוב פונקציות סטיסטיות כדי לענות על השאלות שבהמשך.

<font color = green>

בעזרת הפונקציה המחשבת ממוצע חשבו ערך ממוצע של $\pi$ המתקבל לאחר ביצוע של k ניסויים בהם בוחרים N נקודות.

In [14]:
import py_compile
py_compile.compile('my_statistics.py')
In [17]:
from my_statistics import mean
def mean_for_pi(k,N):
    return mean(get_trails(k,N))
mean_for_pi(5,10000)
Out[17]:
3.14032

<font color =green>

הקישור שלמטה הוא קישור לדף ווקיפדיה של המונח סטיית תקן. קראו על המונח . העזרו בחבילה אותה פיתחתם כדי לשרטט גרף של סטיית התקן כתלות ב- N. שנו את N מ- 10 ל-100 במרווחים של 10. קבעו את k ל-1000. השתמשו בפונקציה זו כדי לשרטט גרף של סטיית התקן כתלות ב-N. מה המסקנה?
</font>

In [2]:
from my_statistics import std
import matplotlib.pyplot as plt

get_trails(5,1000)
 10001

<font color = green>

שרטטו גרף של סטיית התקן כתלות ב: $\sqrt{N}\over{N}$. האם קיבלתם בקרוב יחס ישר? אם כן כיצד סטיית התקן תלויה ב-N?

In [0]:
 

<div dir =rtl> <font color = blue>אליפסה היא המקום הגיאומטרי של כל הנקודות שסכום מרחקיהן משתי נקודות קבועות המכונות מוקדי האליפסה הוא קבוע ושווה לקוטר הגדול של האליפסה, 2a, הקוטר הקטן של האליפסה שווה ל-2b. </font>
<font color =green> כתבו פונקציה בשם ellipse_area המקבל שלושה פרמטרים a, b ו- tol המחשבת שטח של אליפסה שחצי הקוטר הגדול שלה הוא a וחצי הקוטר הקטן הוא b בדיוק של tol. השגיאה בחישוב השטח היא מסדר גודל של סטיית התקן.</font></div>

In [45]:
import random
from math import pi
def ellipse_area(a,b,tol):
    rec_area = 2*a*2*b
    con1=0
    con2=0
    ellipse_area1=pi*a*b#חישוב לפי נוסחא
    for i in range(tol):
        x=uniform(-a,a)
        y=uniform(-b,b)
        if ((x*x)/(a*a) + (y*y)/(b*b) <= 1):
            count=count+1
    ellipse_area2=(float(count)/tol)*rec_area#חישוב לפי שיטת מונטה קרלו
    return ellipse_area1, ellipse_area2#השוואה בין השטחים
ellipse_area(6,3,10000)
Out[45]:
(56.548667764616276, 56.455200000000005)
In [0]: