import random
from math import sqrt
def f(x):
return sqrt(1-x**2)
n = 0
N = 10000
for i in range(N):
x = random.uniform(-1.,1.)
y = random.uniform(0.,1.)
if y < f(x):
n += 1
res = n/N * 4
print(res)
3.1348
alternative solution
import random
from math import sqrt
n = 0
N = 10000
for i in range(N):
x = random.uniform(0.,1.)
y = random.uniform(0.,1.)
if (x**2 + y**2) <1 :
n += 1
res = n/N * 4
print(res)
3.14
The SciPy library contains several packages to perform specialized scientific calculations:
scipy.special
)scipy.integrate
)scipy.optimize
)scipy.interpolate
)scipy.fftpack
)scipy.signal
)scipy.linalg
)scipy.sparse.csgraph
)scipy.spatial
)scipy.stats
)scipy.ndimage
)scipy.io
)It is the foundation of python scientific stack.
The basic building block is the numpy.array
data structure. It can be used as a python list of numbers, but it is a specialized efficient way of manipulating numbers in python.
import numpy as np
a = np.array([1, 2, 3, 4], dtype=float)
a
array([1., 2., 3., 4.])
a = range(1000)
%timeit [ i**2 for i in a]
313 µs ± 10.2 µs per loop (mean ± std. dev. of 7 runs, 1000 loops each)
b = np.arange(1000)
%timeit b**2
1.16 µs ± 20.8 ns per loop (mean ± std. dev. of 7 runs, 1000000 loops each)
c = np.array([[1,2],[3,4]])
c
array([[1, 2], [3, 4]])
c.ndim
2
c.shape
(2, 2)
c = np.arange(27)
c.reshape((3,3,3))
array([[[ 0, 1, 2], [ 3, 4, 5], [ 6, 7, 8]], [[ 9, 10, 11], [12, 13, 14], [15, 16, 17]], [[18, 19, 20], [21, 22, 23], [24, 25, 26]]])
np.zeros((2,2))
array([[0., 0.], [0., 0.]])
np.ones((2,1))
array([[1.], [1.]])
a = np.arange(27).reshape((3,3,3))
np.ones_like(a)
array([[[1, 1, 1], [1, 1, 1], [1, 1, 1]], [[1, 1, 1], [1, 1, 1], [1, 1, 1]], [[1, 1, 1], [1, 1, 1], [1, 1, 1]]])
np.eye(3)
array([[1., 0., 0.], [0., 1., 0.], [0., 0., 1.]])
a = np.arange(10)
a
array([0, 1, 2, 3, 4, 5, 6, 7, 8, 9])
a[0]
0
a[-1]
9
a[0:3]
array([0, 1, 2])
a[::2]
array([0, 2, 4, 6, 8])
a = a.reshape(5,2)
a
array([[0, 1], [2, 3], [4, 5], [6, 7], [8, 9]])
a[3,1]
7
a[2,:]
array([4, 5])
A slice or reshape is a view, simply a re-organization of the same data in memory, thus changing one element changes the same element in all views
a = np.arange(9)
b = a.reshape((3,3))
np.shares_memory(a,b)
True
a[3] = -1
b
array([[ 0, 1, 2], [-1, 4, 5], [ 6, 7, 8]])
np.shares_memory(a,b)
True
b = a.copy()
np.shares_memory(a,b)
False
A typical operation done in your daily physics data analysis is to extract from an array the values that match a condition. Consider an array of the energies of particles, and assume you want to use only the energies above a given threshold. Boolean masking comes at a rescue
ene = np.random.exponential(size=10, scale=10.)
ene
array([1.84417793e+01, 1.01995057e-01, 1.15973376e+01, 3.91436666e-04, 1.13769252e+01, 8.49667565e+00, 4.97730541e+00, 8.38086765e+00, 1.65839490e+01, 5.95578226e+00])
mask = ene > 2
mask
array([ True, False, True, False, True, True, True, True, True, True])
ene[mask]
array([18.44177932, 11.59733762, 11.37692521, 8.49667565, 4.97730541, 8.38086765, 16.58394899, 5.95578226])
ene[ene<2]
array([0.10199506, 0.00039144])
ene[ene<2] = 0
ene
array([18.44177932, 0. , 11.59733762, 0. , 11.37692521, 8.49667565, 4.97730541, 8.38086765, 16.58394899, 5.95578226])
Similarly to boolean masks, it is possible to access and modify values directly to an array using a list of indexes
status = np.random.randint(low=0,high=10,size=10)
status
array([2, 2, 2, 6, 2, 0, 2, 4, 8, 5])
status[[0, 3, 5]]
array([2, 6, 0])
status[[0, 3, 5]] = -1
status
array([-1, 2, 2, -1, 2, -1, 2, 4, 8, 5])
a = np.arange(4)
a
array([0, 1, 2, 3])
a+1
array([1, 2, 3, 4])
10**a
array([ 1, 10, 100, 1000])
np.sin(a)
array([0. , 0.84147098, 0.90929743, 0.14112001])
matrix multiplication
a = np.arange(1,5).reshape(2,2)
b = np.arange(5,9).reshape(2,2)
print(a)
print(b)
print("result: \n",a@b)
[[1 2] [3 4]] [[5 6] [7 8]] result: [[19 22] [43 50]]
a = np.random.randint(low=0,high=10,size=4)
a
array([2, 4, 0, 9])
np.sum(a)
15
np.max(a), np.min(a)
(9, 0)
np.argmax(a), np.argmin(a)
(3, 2)
np.mean(a), np.median(a), np.std(a)
(3.75, 3.0, 3.344772040064913)
a = a.reshape(2,2)
a
array([[2, 4], [0, 9]])
np.sum(a,axis=1)
array([6, 9])
m1 = a>3
m1
array([[False, True], [False, True]])
np.all(m1)
False
np.any(m1)
True
do you remember that in python function arguments are passed by copy?
def f(x):
x+=1
y = 1
f(y)
print(y)
1
this is NOT true for numpy arrays!
def f(x):
x+=1
y = np.array([1])
f(y)
print(y)
[2]
numpy arrays are passed by reference, pay attention!
write another code to compute pi using numpy (avoid to use a for loop)
initialize a pseudo-random number generator with:
rng = np.random.default_rng(1234)
example of usage:
rng.uniform(low=0., high=1., size=10)
array([0.97669977, 0.38019574, 0.92324623, 0.26169242, 0.31909706, 0.11809123, 0.24176629, 0.31853393, 0.96407925, 0.2636498 ])