# A Tale of Two Trapezoids

In our discussion of flow control we wrote a python function to compute a trapezoidal rule approximation to the integral function we specified over some interval we also specified. That approximation was given by $${\int}_{a}^{b}\text{fct}\left(x\right)dx=\frac{h}{2}\left[\text{fct}\left(a\right)+\text{fct}\left(b\right)+2\sum _{i=1}^{n-1}\text{fct}\left({x}_{i}\right)\right]=\frac{h}{2}\left[-\text{fct}\left(a\right)-\text{fct}\left(b\right)+2\sum _{i=0}^{n}\text{fct}\left({x}_{i}\right)\right],$$ where $a={x}_{0}$ and $b={x}_{n}$. The function we used to evaluate this approximation was

def trapezoid(fct,a,b,n): if n<0 or n==0: print "Error: n must be positive" return False h = (double(b)-double(a))/double(n) total = 0.0 for i in arange(n+1): total += fct(a+i*h) total *= 2.0 total -= fct(a)+fct(b) return total*h/2.0

It seems interesting to find out how fast this function is. We can
do that using the `clock` function from the `time` module.
This function gives the current time to some level of precision. On our
computer it gave the time to the nearest hundredth of a second. We
should be able to find the time before we start our trapezoid function,
then run the trapezoid function, then find the time again. The difference
between the two times is roughly how long our function took to run.
We can write a little function to do this for us.

>>> import myfunctions as mf >>> import time >>> def ttrap(fct): ... t1=time.clock() ... fct(mf.f,0,1,1000000) ... t2=time.clock()-t1 ... print t2 ... >>> ttrap(mf.trapezoid) 11.58

Our function is called `ttrap`, and it takes one argument:
the name of the function it is going to time. It finds the current
time and puts it in a variable called `t1`. Then it runs
the function you supplied as an argument. After that it finds
the time again, subtracts `t1` to find the elapsed time, and
prints that. We find that running our trapezoid function with
a million subintervals ($n=1000000$)
takes eleven and a half seconds.

You might wonder why we were so cagey about the name of the function
we wanted to time - why not just type `trapezoid` instead of
`fct` and then supplying `trapezoid` as an argument.
You might also wonder whether we could write a faster function. That,
of course is the point. We can write a function that runs faster,
and we'll give it a different name so we can use this `ttrap`
function to time it.

The thing that takes all the time in the `trapezoid` function
is the "for" construction. A million times it has to evaluate
`f(x)` for some scalar `x`. We could write this instead
to do the calculations on a numpy array. Maybe that would go faster.
Our new function looks like this.

def arraytrapezoid(fct,a,b,n): if n<0 or n==0: print "Error: n must be positive" return False h = (double(b)-double(a))/double(n) x = arange(a,b+h,h) y = 2.0*fct(x) y[0] = y[0]/2.0 y[n] = y[n]/2.0 return sum(y)*h/2.0

Everything is the same up through the definition of
`h`. At that point, instead of starting a "for"
construction, we define a vector `x` of points of the
partition, in such a way that
${x}_{i}=a+ih$
for each
$i=\text{0, 1, ...,}n$.
After that, we have only to evaluate the function for
each point of that array, which we do simply by feeding the
entire array to `fct` as an argument. The properties
of arrays take care of the elementwise multiplication that
is required. After that we just subtract the extra function
value at $a$ and $b$,
and return the result.

Evidently this version of the code needs a great deal more memory,
since it has to store the entire vector `x` as well as
the function evaluation over that: `y`. However, this
version of the code is spectacularly faster than our earlier
version. Remember that our other version took over eleven
seconds to run.

>>> ttrap(mf.arraytrapezoid) 0.05

Why is this more than twenty times faster? It is the difference
between *interpreted* code and *compiled* code. These
are essentially the two ways to run computer programs. Interpreters
look at the text commands you type, parse their meaning, translate
them into machine code, and then execute them. Compilers translate
an entire file of code into machine language, which then only needs
to be loaded and run full speed. Interpreted code thus runs much more
slowly than that compiled, as it were pretranslated code.

In our case, the `trapezoid` function relies on that
"for" construction. Thus each time it executes the code, the
python interpreter must translate that command to add to the total,
then execute it. By contrast, in the `arraytrapezoid` function,
everything happens through the numpy `array` construction.
This, like many things in numpy, is compiled. Thus in this case, instead
of one million interpreted lines, there is only one, and then the
computation happens in compiled code, so it is much faster. In general, because
the numpy package contains much compiled code, the more we rely on it
and not native python, the faster our programs will run.

Assignment 7 is posted.

A solution to the test is
available.