Come aumentare le prestazioni per la stima di `Pi`in Python

Question

Ho scritto il seguente codice in Python, per stimare il valore di Pi. Si chiama metodo Monte Carlo. Ovviamente aumentando il numero di campioni il codice diventa più lento e presumo che la parte più lenta del codice sia nella parte di campionamento. Come posso renderlo più veloce?

from __future__ import division 
import numpy as np 

a = 1 
n = 1000000 

s1 = np.random.uniform(0,a,n) 
s2 = np.random.uniform(0,a,n) 

ii=0 
jj=0 

for item in range(n): 
    if ((s1[item])**2 + (s2[item])**2) < 1: 
     ii = ii + 1 

print float(ii*4/(n)) 

Suggerite altri (presumibilmente più veloci) codici?

Answer 1

Il collo di bottiglia qui è in realtà il ciclo for. I loop Python for sono relativamente lenti, quindi se hai bisogno di iterare oltre un milione di elementi, puoi guadagnare molta velocità evitandoli del tutto. In questo caso, è abbastanza facile. Invece di questo:

for item in range(n): 
    if ((s1[item])**2 + (s2[item])**2) < 1: 
     ii = ii + 1 

fare questo:

ii = ((s1 ** 2 + s2 ** 2) < 1).sum() 

Questo funziona perché numpy ha supporto incorporato per le operazioni di matrice ottimizzate. Il ciclo si verifica in c anziché in python, quindi è molto più veloce. Ho fatto un test rapido in modo da poter vedere la differenza:

>>> def estimate_pi_loop(x, y): 
...  total = 0 
...  for i in xrange(len(x)): 
...   if x[i] ** 2 + y[i] ** 2 < 1: 
...    total += 1 
...  return total * 4.0/len(x) 
... 
>>> def estimate_pi_numpy(x, y): 
...  return ((x ** 2 + y ** 2) < 1).sum() 
... 
>>> %timeit estimate_pi_loop(x, y) 
1 loops, best of 3: 3.33 s per loop 
>>> %timeit estimate_pi_numpy(x, y) 
100 loops, best of 3: 10.4 ms per loop 

Ecco alcuni esempi dei tipi di operazioni che sono possibili, solo così si ha un senso di come funziona.

quadratura una matrice:

>>> a = numpy.arange(5) 
>>> a ** 2 
array([ 0, 1, 4, 9, 16]) 

array Aggiunta:

>>> a + a 
array([0, 2, 4, 6, 8]) 

array confronto:

>>> a > 2 
array([False, False, False, True, True], dtype=bool) 

sommatori valori:

>>> (a > 2).sum() 
2 

Come probabilmente ti rendi conto, ci sono modi più veloci per stimare Pi, ma ammetterò che ho sempre ammirato la semplicità e l'efficacia di questo metodo.

Answer 2

Sono stati assegnati array numpy, quindi è consigliabile utilizzarli a proprio vantaggio.

for item in range(n): 
    if ((s1[item])**2 + (s2[item])**2) < 1: 
     ii = ii + 1 

diventa

s1sqr = s1*s1 
s2sqr = s2*s2 
s_sum = s1sqr + s2sqr 
s_sum_bool = s_sum < 1 
ii = s_sum_bool.sum() 
print float(ii*4/(n)) 

dove si effettua la quadratura del array, li somma, controllando se la somma è inferiore a 1 e poi sommando i valori booleani (false = 0, true = 1) per ottenere il numero totale che soddisfa i criteri.

Answer 3

ho upvoted risposta senderle s', ma nel caso in cui non si vuole modificare il codice troppo:

numba è una libreria progettata per questo scopo.

basta definire l'algoritmo in funzione, e aggiungere il @jit decoratore:

from __future__ import division 
import numpy as np 
from numba import jit 

a = 1 
n = 1000000 

s1 = np.random.uniform(0,a,n) 
s2 = np.random.uniform(0,a,n) 

@jit 
def estimate_pi(s1, s2): 
    ii = 0 
    for item in range(n): 
     if ((s1[item])**2 + (s2[item])**2) < 1: 
      ii = ii + 1 
    return float(ii*4/(n)) 

print estimate_pi(s1, s2) 

Sul mio portatile, guadagna circa un aumento di velocità 20x per n = 100000000, e un aumento di velocità 3x per n = 1000000.