Project Euler # 645: velocizza la simulazione Monte-Carlo in Python

Question

Sto cercando di risolvere Q645. Mentre la logica utilizzata per il mio codice sembra essere appropriata, il codice stesso è troppo lento per il numero elevato richiesto in questa domanda. Posso chiedere suggerimenti per migliorare le prestazioni del mio codice?

La domanda è come nel link: https://projecteuler.net/problem=645

Il mio codice Python è il seguente:

def Exp(D):
    day_list = [0]*D
    num_emperor = 0
    while all((d == 1 for d in day_list)) == False:
        #the birthday of the emperors are independent and uniformly distributed throughout the D days of the year
        bday = np.random.randint(0,D)
        day_list[bday] = 1
        num_emperor+=1
        #indices of d in day_list where d == 0
        zero_ind = (i for i,v in enumerate(day_list) if v == 0)
        for ind in zero_ind:
            try:
                if day_list[ind-1] and day_list[ind+1] == 1:
                    day_list[ind] = 1
            except IndexError:
                if ind == 0:
                    if day_list[-1] and day_list[1] == 1:
                        day_list[0] = 1
                elif ind == len(day_list)-1:
                    if day_list[len(day_list)-2] and day_list[0] == 1:
                        day_list[len(day_list)-1] = 1
    return num_emperor

def my_mean(values):
    n = 0
    summ = 0.0
    for value in values:
        summ += value
        n += 1
    return summ/n

def monte_carlo(iters, D):
    iter = 0
    n_emperor = 0
    while iter < iters:
        n_emperor = Exp(D)
        yield n_emperor
        iter += 1

avg_n_emperor = my_mean(monte_carlo(iters,D))
print(avg_n_emperor)

E la mia logica è la seguente:

Per il day_list all'interno della funzione Exp (D) , dove D è il numero di giorni in un anno, gli zeri non indicano festività e quelli significano festività. Inizialmente il day_list è tutti zeri poiché non ci sono festività con cui iniziare .

Le regole per definire un giorno casuale ( d ) come festività sono le seguenti:

1. All'inizio del regno dell'attuale Imperatore, il suo compleanno è dichiarato festivo da quell'anno in poi.

2. Se sia il giorno prima che quello dopo un giorno d sono festivi, anche d diventa festivo.

Successivamente ho implementato le regole indicate per la domanda, per aggiungere gradualmente le festività (quelle) nella day_list . Dopo num_emperor numero di imperatori, tutti i giorni ( d ) in day_list diventeranno 1, cioè tutti i giorni diventeranno festivi . Questo è il punto per chiudere il while_loop nella funzione Exp (D) e contare il numero di imperatori richiesti. Per ottenere il numero medio di imperatori necessari affinché tutti i giorni diventino vacanze ( avg_n_emperor ), applico quindi il metodo monte-carlo.

Per il mio codice attuale, il tempo necessario è il seguente:

avg_n_emperor = my_mean(monte_carlo(iters=100000,D=5)) #6-7 seconds

avg_n_emperor = my_mean(monte_carlo(iters=1000000,D=5)) #about 62 seconds

in cui il tempo impiega aumentare di ca. linearmente con gli iter .

Tuttavia,

avg_n_emperor = my_mean(monte_carlo(iters=1000,D=365)) #about 68 seconds

impiega già circa 68 secondi e la domanda richiede D = 10000. Per non parlare del fatto che gli iter richiesti affinché la risposta sia accurata entro 4 cifre dopo i punti decimali (come richiesto dalla domanda) sarebbero anche molto più grandi di 1000000 ...

Qualunque suggerimento per velocizzare il mio codice sarebbe apprezzato! :)

Answer 1

Benvenuto in Code Review. Bella implementazione, facile da leggere e da capire.

Ottimizzazione

Ci sono alcune operazioni "costose" che possono essere semplificate. Di seguito ho commentato le parti rilevanti:

def Exp(D):
    # the method "all" takes O(D)
    while all((d == 1 for d in day_list)) == False:
        # O(D)
        zero_ind = (i for i,v in enumerate(day_list) if v == 0) 
        # O(D)
        for ind in zero_ind:
            # Here there are only O(1) operations
    return num_emperor

Di \$O(D)\$Voglio dire che nel peggiore dei casi tale operazione itererà i Dtempi, dov'è Dil numero di giorni.

La condizione nel ciclo while può essere semplificata controllando se il numero di ferie è <giorni:

def Exp(D):
    holidays = 0
    while holidays < D:
        # increment holidays 
    return num_emperor

La seconda ottimizzazione è evitare i loop interni. Una volta calcolato il nuovo compleanno, è sufficiente "guardarsi intorno" in quel giorno specifico:

def Exp(D):
    # ..
    while holidays < D:
        bday = np.random.randint(0,D)
        # Increment holidays only if birthday is not in a holiday
        if day_list[bday] == 0:
            holidays += 1
        day_list[bday] = 1
        num_emperor+=1

        yesterday = (bday - 1) % D
        day_before_yesterday = (bday - 2) % D
        if day_list[day_before_yesterday] == 1 and day_list[yesterday] == 0:
            day_list[yesterday] = 1
            holidays += 1

        tomorrow = (bday + 1) % D
        day_after_tomorrow = (bday + 2) % D
        if day_list[day_after_tomorrow] == 1 and day_list[tomorrow] == 0:
            day_list[tomorrow] = 1
            holidays += 1
    return num_emperor

L' %operatore impedisce l'overflow dell'array, quindi non è necessario rilevare eccezioni.

Esecuzione della media:

avg_n_emperor = my_mean(monte_carlo(iters=1000,D=365))
# Output: 1173.786
# Running time: around 2 seconds

Per quanto riguarda lo stile, @Peilonrayz ha già fornito un'ottima recensione.

Answer 2

In primo luogo, facciamo in modo che il codice sia un po 'più pulito:

• Puoi usare statistics.meanpiuttosto che creare my_mean.

• Dovresti usare un forciclo piuttosto che un ciclo while monte_carlo.

• Non è necessario eseguire alcuna assegnazione n_emperernella funzione.

• Expe Ddovrebbe essere lower_snake_case. Questo è perché sono funzioni e variabili.

• Dovresti mettere degli spazi intorno a tutti gli operatori.

• Dovrebbe esserci uno spazio dopo le virgole.

• Dovresti avere nomi migliori che day_listpotrebbero essere days, Dpotrebbero anche essere qualcosa come days, summpossono essere total, iterspotrebbero essere amounts.

• Puoi semplicemente usare all(day_list)piuttosto che all((d == 1 for d in day_list)).

• Non utilizzare ==per confrontare i singleton come False. Sarebbe meglio se invece usassi not.

• Questo non controlla se entrambi i valori sono 1, controlla se il primo è vero e il secondo è uno. Ciò significa che se imposti day_list[index - 1]a due sarebbe ancora vero.

  day_list[ind - 1] and day_list[ind + 1] == 1
  

  Per verificare che siano entrambi uguali a quello che dovresti usare:

  day_list[ind - 1] == 1 and day_list[ind + 1] == 1
  

  Qui invece vorrei solo controllare se sono veritiere.

• Non è necessario if ind == 0:che se indfosse 0, allora ind - 1lo sarà -1.

• Puoi semplicemente usare (ind + 1) % len(days)per rimuovere la necessità di elif index == len(days)-1:.

import random
import statistics


def simulate(days_in_year):
    days = [0] * days_in_year
    emperors = 0
    while not all(days):
        days[random.randrange(len(days))] = 1
        emperors += 1
        for index, value in enumerate(days):
            if value:
                continue
            if days[index - 1] and days[(index + 1) % len(days)]:
                days[index] = 1
    return emperors


def monte_carlo(amount, days):
    for _ in range(amount):
        yield simulate(days)


print(statistics.mean(monte_carlo(amount, days)))

Ora che il codice è piccolo e carino, possiamo concentrarci su ciò che causa problemi di prestazioni.

1. Quanto segue anyviene eseguito in \$O(n)\$ora, dove \$n\$è la lunghezza di days. Ciò significa che nel peggiore dei casi verrà eseguito per quanto lunghi siano i giorni ogni volta che lo chiami.

   not all(days)
   

   Possiamo fare di meglio aggiungendo una variabile in tali incrementi ogni volta che cambiamo uno 0 in 1. Possiamo quindi confrontarla con days_in_yearper vedere se l'elenco è pieno. Questo verrà eseguito in \$O(1)\$ tempo determinando un notevole risparmio.

2. Se un nuovo imperatore nasce in una vacanza già esistente, non verranno fatte vacanze extra.

3. Quando nasce un nuovo imperatore non è necessario verificare se ogni zero può essere modificato, basta invece controllarne due. Questo taglierà un altro \$O(n)\$operazione a \$O(1)\$.
   Supponiamo di avere quanto segue come days:

   0123456
   1000010
   

   Se il nuovo compleanno è:

   • 6 - Poiché sia ​​5 che 0 sono già 1, non è possibile effettuare festività aggiuntive.

   • 3 - Poiché 4 è uno 0 e 5 è un 1, 4 può diventare un 1. Poiché 2 è uno 0 ma 1 è uno 0, 3 non può diventare un 1.

     Questo non può propagarsi verso l'esterno.

Answer 3

Realisticamente, la mia recensione dovrebbe essere "Non funzionerà, non otterrai la precisione richiesta con un esperimento del genere. Hai bisogno di un approccio diverso" .

Ma ecco una simulazione del tempo O (D). Invece di generare potenzialmente compleanni già avvenuti più e più volte, mi concentro solo sui nuovi compleanni. Cioè, all'inizio mescolo tutti i compleanni possibili e poi li passo. Ovviamente questo significa che non posso semplicemente farlo emperors += 1. Invece, aggiungo il numero previsto di nuovi imperatori necessari per incontrare un nuovo compleanno.

Con 1000 simulazioni, il mio laptop impiega circa 0,6 secondi per D = 365, 1,8 secondi per D = 1000 o 19 secondi per D = 10000.

from random import sample
from statistics import mean

def Exp(D):
    emperors = 0
    holidays = set()
    for i, day in enumerate(sample(range(D), D)):
        emperors += D / (D - i)
        holidays.add(day)
        if (day + 2) % D in holidays:
            holidays.add((day + 1) % D)
        if (day - 2) % D in holidays:
            holidays.add((day - 1) % D)
        if len(holidays) == D:
            return emperors

print(mean(Exp(365) for _ in range(1000)))

Meh. Ho appena provato anche in questo emperor += 1modo, ci sono voluti circa 1,35, 4,1 e 62 secondi:

from random import randrange
from statistics import mean

def Exp(D):
    emperors = 0
    holidays = set()
    while len(holidays) < D:
        emperors += 1
        day = randrange(D)
        if day not in holidays:
            holidays.add(day)
            if (day + 2) % D in holidays:
                holidays.add((day + 1) % D)
            if (day - 2) % D in holidays:
                holidays.add((day - 1) % D)
    return emperors

print(mean(Exp(365) for _ in range(1000)))