# +----------------------------------------------------+
# |                  Lekce 6 - demo 1                  |
# +----------------------------------------------------+
#   MONTE CARLO METHODS 
#   
#  PROBLEM: calculate pi with random throwing into unit circle
#  THIS VERSION: ve srovnani s prednaskou, jen na ctvrtkruhu

import numpy as np
#import scipy.linalg as linalg
import matplotlib.pyplot as plt

n=0    #number of points
npi=0  #number of points inside

#STORING VALUES:
k=0    
kk=10
nk=19
nn=np.zeros(nk)
pp=np.zeros(nk)

ubound=10*2**(nk-1)

while(n<ubound+1):
     n=n+1
     x=np.random.rand(2)
     if (np.linalg.norm(x)<1.0):
         npi=npi+1
     if (n==kk):
         print(n,4*npi/n)
         nn[k]=n
         pp[k]=4*npi/n
         kk=2*kk
         k=k+1
     

#%matplotlib inline
plt.plot(np.log10(nn),np.log10(abs(pp-np.pi)),'bo')
plt.plot(np.log10(nn),np.log10(3/np.sqrt(nn)),'b')
plt.show()