Benchmark solve_triangular from scipy.linalg

bench_triangular_solve.py

import numpy as np
from scipy import linalg
from datetime import datetime

import gc

mu_sec = 1e-6 # number of seconds in one microseconds

solve_time =[]
solve_triangular_time =[]

dims = np.arange(100, 2020, 20)

for i in dims:
    print "Iteration %s of %s" % (i/20, 100)
    A = np.tril(np.random.randn(i, i)) + 100 * np.eye(i) # to avoid singularities
    Y = np.random.randn(i)

    gc.collect()
    start = datetime.now()
    sol1 = linalg.solve_triangular(A, Y, lower=True)
    delta = datetime.now() - start
    solve_triangular_time.append(delta.seconds + delta.microseconds * mu_sec)
    
    gc.collect()
    start = datetime.now()
    sol2 = linalg.solve(A, Y)
    delta = datetime.now() - start
    solve_time.append(delta.seconds + delta.microseconds * mu_sec)

    # check for consistency
    assert np.linalg.norm(sol1 - sol2) < 0.01

import pylab as pl

pl.plot(dims, solve_triangular_time, label='solve_triangular')
pl.plot(dims, solve_time, label='solve')
pl.xlabel('number of dimensions')
pl.ylabel('seconds')
pl.legend()
pl.axis('tight')
pl.show()