Project Euler #645 -- speed up Monte-Carlo simulation in Python I am trying to solve Q645. While the logic used for my code seems to be appropriate, the code itself is way too slow for the large number required in this question. May I ask for suggestion(s) to improve my code's performance?
The question is as in the link: https://projecteuler.net/problem=645
My Python code is as follow:
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 &lt; iters:
        n_emperor = Exp(D)
        yield n_emperor
        iter += 1

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

And my logic is as follow:
For the day_list inside the Exp(D) function, where D is the number of days in a year, zeros mean no holiday, and ones mean holiday. Initially the day_list is all zeros since there is no holiday to begin with.
The rules of defining a random day (d) as a holiday is as follow:

At the beginning of the reign of the current Emperor, his birthday is declared a holiday from that year onwards.

If both the day before and after a day d are holidays, then d also becomes a holiday.


I then subsequently implement the rules stated for the question, to gradually add holidays (ones) into the day_list. After num_emperor number of emperors, all the days (d) in day_list will become 1, i.e. all days will become holiday. This is the point to quit the while_loop in Exp(D) function and count the number of emperors required. To get the average number of emperors required for all the days to become holidays (avg_n_emperor), I then apply the monte-carlo method.
For my current code, the time takes is as follow:
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 which the time takes increase approx. linearly with the iters.
However,
avg_n_emperor = my_mean(monte_carlo(iters=1000,D=365)) #about 68 seconds

already takes about 68 seconds, and the question is asking for D=10000. Not to mention that the iters required for the answer to be accurate within 4 digits after the decimal points (as required by the question) would be much larger than 1000000 too...
Any suggestion(s) to speed up my code would be appreciated! :)
Welcome to Code Review. Nice implementation, easy to read and understand.
Optimization
There are some &quot;expensive&quot; operations that can be simplified. Below I commented the relevant parts:
def Exp(D):
    # the method &quot;all&quot; 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

By \$O(D)\$ I mean that in the worst case such operation will iterate D times, where D is the number of days.
The condition in the while loop can be simplified by checking if the number of holidays is &lt; days:
def Exp(D):
    holidays = 0
    while holidays &lt; D:
        # increment holidays 
    return num_emperor

The second optimization is to avoid the inner loops. Once the new birthday has been calculated, it's enough to &quot;look around&quot; that specific day:
def Exp(D):
    # ..
    while holidays &lt; 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

The % operator prevents to overflow the array, so you don't need to catch exceptions.
Running the average:
avg_n_emperor = my_mean(monte_carlo(iters=1000,D=365))
# Output: 1173.786
# Running time: around 2 seconds

Regarding the style, @Peilonrayz already provided an excellent review.