I’m currently involved in a project for a paper I’m writing where I’m trying to find the auto-covariance function of a queue, which involves evaluating modified Bessel functions of the first kind; in essence these are the solutions to a particular form of differential equations which tell you how queues work, and each evaluation takes a relatively large amount of time: they are the limiting step in getting towards the answer I’m after in my code. I’m very close to generating the results I think I need to finish the paper, and was getting very frustrated by how slowly the code ran.
For my application, the Bessel function, implemented in matlab and (with a really weird bug initially) in Octave, as
besseli is ideal. Essentially, I need to evaluate this over and over again; the Bessel function has two parameters, an index, which in my case is an integer, and another parameter, in my case a positive real number, which is fixed for each function evaluation.
The thing I learnt last night is the extreme time saving in doing things in a sensible order- I’ve had quite a lot of coding training, but to be honest, never really come across situations where using it is needed. So essentially my code was doing this for each z,
out+=m*(besseli(m-n,z)-besseli(m+n,z) +...% (more functions of besseli)
This function has many calls to
besseli, each of which is an expensive step. Sometimes we even evaluate the same function twice. By replacing this function with a vectorized call to
besseli, we should reduce the overall calls to
b=besseli((0:(n+N+50)),z); % <-Vectorized version of besseli
out+=m*(b(abs(m-n)+1)-b(abs(m+n)+1) +...%(more functions of b)
And it’s as simple as that; it helps in this case that the function is symmetric (and
besseli(n,z)-besseli(-n,z). The results are that the code runs about four times faster. For those of you that programme more regularly, and as it should have been for me, this is quite obvious. I write this because
- With a little bit of thought, I could have saved myself many hours of computing time, and have finished this paper about a month ago, probably.
- It might help someone else out.
- As a reminder to myself to think carefully about things like this before any serious coding happens.