Plotting Timing Experiments
Like any other person I regulary need to run experiments to check how fast or slow a particular algorithm/implementation is for a given problem. The natural choice is to plot the data. This way you at least get some more or less pretty picture out of the tendious experience of having to wait for the experiment to finish. I used to write crappy code to generate these pictures myself and I could not convince myself to remember the appropriate commands for matplotlib and R. Today I sat down and learned the five lines of code necessary to have decent plots for my experiments. I’m putting examples here for no good reason except maybe to show off Sage’s new HNF code which I use as a showcase.
First lets compare how long it takes to compute the Hermite Normal Form for a given random $n \times n$ matrix with (possibly negative) integer entries of size bounded absolute by $2^{16}$.
n = 10
b = 16
st =[]
mt = []
x = [20*i for i in range(n)]
for i in range(n):
A = random_matrix(ZZ,20*i,20*i, x=-2**b, y=2**b)
t = cputime()
E = A.echelon_form()
st.append(cputime(t))
AM = A._magma_()
t = magma.cputime()
EM = AM.EchelonForm()
mt.append(magma.cputime(t))
import pylab
pylab.clf() # clear the figure first
pylab.figure(1)
# plot some data and add a legend
pylab.plot(x,st,label="Sage")
pylab.plot(x,mt,label="Magma")
pylab.legend() # print the legend
pylab.title("HNF for Random Matrices with $%d$-bit Integer Entries: Sage vs. Magma"%b)
pylab.ylabel("execution time $t$") # label the axes
pylab.xlabel("n for n x n matrix")
pylab.savefig('foo.png',dpi=72) # fire!
Now lets use R to see how the runtimes vary for random $160 \times 160$ matrices with (possible negative) integer entries bounded absolute by $2^{10}$.
b = 10
st = []
for i in range(500):
A = random_matrix(ZZ,160,160, x=-2**b, y=2**b)
t = cputime()
E = A.echelon_form()
st.append(cputime(t))
from rpy import r
r.png('histogram.png',width=640,height=480)
r.hist(st,r.seq(1.2,3.7,0.02),main="SAGE HNF Histogram",col="lightblue", prob=True, xlab="seconds")
r.lines(r.density(st,bw=0.05),col="black")
r.rug(st)
r.dev_off()
Neat, isn’t it? Btw. Pygments is also neat, thanks rpw.

