Is there a better way to do run time analysis than this?

Question

I currently have 2 different functions with options to vectorise them: acc_rej_sine(max_iter, algorithm=None) and analytical_sine(max_iter, algorithm=None) for which I'm trying to compare their run time against the number of iterations computed. i.e. compare all 4 methods; 2 looped, 2 vectorised. Essentially my code goes something like this:

def analytical_sine(max_iter, algorithm=None):
    if algorithm is None:
        count = 0
        analytical_hist = []
        for i in range(max_iter):
            count += 1
            progress = round((count/max_iter)*100)
            sys.stdout.write('\r' + str(progress) + '%')
            uni_dist = np.random.uniform(0, 1)
            arccos = np.arccos(1 - 2*uni_dist)
            analytical_hist.append(arccos)
    elif algorithm is "vectorise":
        analytical_hist = np.arccos(1 - 2*np.random.uniform(0, 1, max_iter))

    return analytical_hist

def acc_rej_sine(max_iter, algorithm=None):
    x = np.random.uniform(0, np.pi, max_iter)
    y = np.random.rand(max_iter)
    if algorithm is None:
        accepted_x = []
        j = count = 0
        for i in range(max_iter):
            count += 1
            progress = round((count/max_iter)*100)
            sys.stdout.write('\r' + str(progress) + '%')
            if y[i] <= np.sin(x[i]):
                accepted_x.append(x[i])
                j += 1
    elif algorithm is "vectorise":
        accepted_x = np.extract((y <= np.sin(x)), x)

    return accepted_x

def runtime(func, runs, max_iter, algorithm=None):
    time = []
    for i in range(runs):
        start = timer()
        func()
        end = timer()
        time.append((end-start))
    error = np.std(time)
    time = sum(time)/runs

    return time, error

def time_analysis():
    time1, time2, time3, time4 = [], [], [], []
    error1, error2, error3, error4 = [], [], [], []
    for i in np.arange(1, 8, 1):
        max_iter = 10**i
        time, error = runtime(analytical_sine, 5, int(max_iter))
        time1.append(time)
        error1.append(error)
        time, error = runtime(analytical_sine, 5, int(max_iter), "vectorise")
        time2.append(time)
        error2.append(error)
        time, error = runtime(acc_rej_sine, 5, int(max_iter))
        time3.append(time)
        error3.append(error)
        time, error = runtime(acc_rej_sine(max_iter), 5, int(max_iter), "vectorise")
        time4.append(time)
        error4.append(error)
    return [time1, time2, time3, time4], [error1, error2, error3, error4]

# to run the code I would probably do something like this
time, error = time_analysis()
#then if I wanna plot number of iterations vs run time with errors I would probably do something along the lines of
plt.plot(max_iter, time) # max_iter would be a list of [10**i for i in np.arange(1, 8, 1)]
plt.errorbar(error)

So the idea is my runtime() function would allow me to put in any of the 4 functions that I want to compare(which currently still isn't working yet and I can't/haven't figured out why), and run it for 5 times, work out the mean run time and standard deviation as the error. In return my time_analysis() function would run runtime() for different functions for different max_iter(max iterations) which goes like [10, 1E2, 1E3, 1E4, 1E5, 1E6, 1E7] so I can plot max iterations against run time. However, this whole method seems quite cumbersome and inelegant, as my time_analysis() requires me to repeatedly work out time and error and append it to a list. Is there a better way of timing this?(Also my runtime() doesn't work yet lol because my algorithm argument seems to be making something not callable)

Answer 1

Welcome to SE SciComp. First of all, I would suggest using Jupyter so that you have access to IPython and its nice timing magics (see %time and %timeit magic function). These magics take into account that you need to run the code a couple of times to get reliable measurements. If you need even more in-depth comparison, you also need to take into account memory allocation and caching for which you can use the %prun magic.

Second, it's a bit unfair to compare timings between a vectorized function (in your case a NumPy call, which uses C as a backend) and a routine where you also perform I/O (i.e. printing your progress bar on screen). I can understand the comparison between the vanilla non-vectorized routine and the NumPy routine, but you shouldn't penalize the vanilla routine by having print outs to the screen. You're comparing grapes to apples...

Third, try to avoid the if-then-else clause in the routines. Write separate routines for each case (vectorized vs non-vectorized) that you want to time/compare.

Fourth, in your non-vectorized algorithm you append a value to a list in each iteration. Since you use NumPy, why not pre-allocate a vector of size max_iter?

See code below (note that it has to be run inside an IPython environment, i.e. a Jupyter notebook or Qtconsole or plain IPython console). Is this want you want to achieve?

import numpy as np

def analytical_sine(max_iter):
    analytical_hist = np.zeros(max_iter)
    for i in range(max_iter):
        uni_dist = np.random.uniform(0, 1)
        analytical_hist[i] = np.arccos(1 - 2*uni_dist)
    return analytical_hist

def analytical_sine_vec(max_iter):
    return np.arccos(1 - 2*np.random.uniform(0, 1, max_iter))

print("Non-vectorized")
for i in range(10):
    print(f"Max_iter = {100*(i+1)}")
    %timeit res = analytical_sine(100*(i+1))
print("Vectorized")
for i in range(10):
    print(f"Max_iter = {100*(i+1)}")
    %timeit res = analytical_sine_vec(100*(i+1))