An Empty Do

Back when I was first learning Lisp by reading Paul Graham's Ansi Common Lisp, Graham mentioned that sometimes you can do useful work with a DO loop having an empty body. I thought that was pretty neat but I've never come across a case where it was the natural solution. That is until now. As part of a small project I'm working on, I needed to evaluate a polynomial. The efficient way to do that is to use Horner's rule.

If you're a mathematician, you can think of Horner's rule as being a consequence of synthetic division or, if you're a programmer, as an application of strength reduction. Taking the programmer's point of view, we have

a0 + a1x + … + an-1xn-1 + anxn = a0 + x(a1 + x(a2 + … + x(an-1 + anx)….)

Here's Horner's rule in Common Lisp:

(defun horner (a x)
  "Evaluate the polynomial with coefficients in the array a at x."
  (if (= (length a) 1)
      (aref a 0)
      (do* ((k (1- (length a)) (1- k))
            (val (* (aref a k) x) (* (+ val (aref a k)) x)))
           ((<= k 1) (+ val (aref a 0))))))

Notice that the body of the DO is empty. That's a pretty nice example; the implementation is natural and just what you'd write even if you'd never heard of the technique.

After I wrote that, I thought that you could do the same thing with a for loop in C. It turns out, though, that the semantics are just different enough that it doesn't work—or at least I couldn't make it work. The best I could do is

double horner ( double x, double *a, int deg )
    double val;
    int k;

    val = a[ deg ] * x;
    for ( k = deg - 1; k > 0; k-- )
        val = ( val + a[ k ] ) * x;
    return val + a[ 0 ];

You can see that they do the same thing in the same way but Lisp is wonderfully concise. Or, I suppose, horribly obscure depending on your point of view.

Update: Fixed formatting in polynomial strength reduction.

This entry was posted in Programming and tagged , . Bookmark the permalink.
  • R

    Surprisingly I found ECMAScript version is very "Lispy" - still I'm missing a do macro 8P -

    function horner (x, p) {
    var val = 0;
    for (var k = p.length; k > 0; --k,
    val += ("calculate", k === 0? p[k] : p[k] * x));
    return val;

    var pol = [1, 2]; // 2x + 1

    horner(2, pol); // 5