MonteCarlo_Pi_numpy_matplotlib.py

# Heitetään pisteitä satunnaisesti neliölle
# Osa pisteistä osuu neliön sisällä olevalle ympyrälle
# Punaiset pisteet ympyrän sisällä, siniset ulkona
# Lasketaan tuloksista Piin arvo
# Juhani Kaukoranta 31.12.2019
import numpy as np
import matplotlib.pyplot as plt
import time

def monte_carlo_pi(r,n):
    # pisteiden x- ja y-koordinaatit satunnaisvektoreina
    # ja x- ja ykoordinaattien neliöt
    x= 2*r*np.random.random_sample(n)
    y= 2*r*np.random.random_sample(n)
    x2 = (x - r)**2
    y2 = (y - r)**2
    r2 = r**2
    # circle-ehdon täyttävät pisteet ympyrän sisällä
    circle = x2+ y2 <= r2  # ympyrän sisällä olevien ehto
    nocircle = x2 + y2 > r2 # ympyrän ulkopuolella
    accum = np.sum(circle)  # ympyrän sisällä olevien pisteiden määrä
    print("Pii = ",4*accum/n)
    xout = x[nocircle]
    yout = y[nocircle]
    xcircle = x[circle]
    ycircle = y[circle]
    # plottaus, aluksi ympyrä
    alfa = np.linspace(0, 2*np.pi, 100)
    x1 = r+r*np.cos(alfa)
    y1 = r+r*np.sin(alfa)
    fig,ax = plt.subplots(1)
    ax.plot(x1, y1)  # piirtää ympyrän
    ax.set_aspect(1)  # asettaa näytön oikein
    # plottaus pisteet ympyrän sisällä ja ulkona
    ax.plot(xcircle, ycircle,'r.')  # sisällä punaiset pisteet
    ax.plot(xout,yout,'b.')  # ulkona siniset pisteet
    plt.xlim(0,2*r)
    plt.ylim(0,2*r)
    plt.xlabel("Satunnaispisteet ja ympyrä")
    fig.tight_layout()
    plt.show()
    # accum sisältää ympyräneljännekseen osuvien pisteiden lukumäärän
    # Piin arvon laskeminen
    return

pisteita = int(input("Kuinka monta pistettä arvotaan "))
radius = float(input("Anna ympyrän säde "))
time0 = time.perf_counter()
print("Pii = ",monte_carlo_pi(radius,pisteita))
time1 = time.perf_counter()
print("Aikaa kului ",time1-time0," sekuntia laskentaan, piirtämiseen ja kuvan katsomiseen")