Author: Gareth McCaughan
Date: 17:40:01 12/11/01
Go up one level in this thread
On December 11, 2001 at 17:49:24, Dann Corbit wrote:
> However, I broadly prefer the template implementation. Take a gander
> at the STL complex template sometime, and see the incredible power of
> the template idea. I can create a template that uses a big number class
> as easily as a float, double, or long double complex type.
>
> Here are some templates I wrote which are freely available for any use:
> A Kahan adder (see http://docs.sun.com/htmlcoll/coll.648.2/iso-8859-1/NUMCOMPGD/ncg_goldberg.html
> under "Theorem 8 Kahan Summation Formula"):
> ftp://cap.connx.com/pub/tournament_software/Kahan.Hpp
Just for laughs, here's a rather different view of how to use
generic programming and generative programming to do Kahan adders
efficiently. Take a deep breath...
(defun kahan-adder-flet-clause (sum-var func-name type carry-var)
(let ((y-var (gensym))
(temp-var (gensym))
(addend-var (gensym)))
`(,func-name (,addend-var)
(let* ((,y-var (- ,addend-var ,carry-var))
(,temp-var (+ ,sum-var ,y-var)))
,@(when type
`((declare (type ,type ,y-var ,temp-var))))
(setf ,carry-var (- ,temp-var ,sum-var ,y-var)
,sum-var ,temp-var)))))
(defmacro with-kahan-adders (adders &body body)
(let ((carries (loop for x in adders collect (gensym))))
`(let ,(loop for (var func type) in adders for c in carries
nconc (let ((zero (if type `(coerce 0 ',type) '0)))
`((,var ,zero) (,c ,zero))))
,@(loop for (var func type) in adders for c in carries
nconc (when type
`((declare (type ,type ,var ,c)))))
(flet ,(loop for (var func type) in adders for c in carries
collect (kahan-adder-flet-clause var func type c))
(declare (inline ,@(mapcar #'second adders)))
,@body))))
This little (and admittedly ugly) snippet of Common Lisp
lets you say things like this:
(with-kahan-adders ((x add-to-x double-float)
(y add-to-y double-float))
(loop for i from 1 to 1000 do
(let ((delta (/ (coerce i 'double-float))))
(add-to-x delta) <--- \ These are magic!
(add-to-y (- 1 delta)))) <--- /
(list x y))
and the nice thing is that the sum and carry for each of
those adders will just be an ordinary local variable, and
there will be no function call overhead for any of the
additions. The temporaries used for the Kahan summing
will be ordinary local variables too. This may be a good
thing or a bad thing, depending on what sort of numbers
you're working with. For the common case of ordinary
floating-point numbers, it's good.
In other words, that code above turns -- at compile time --
into something equivalent to the following C++:
{
double x,x_carry;
double y,y_carry;
for (int i=0; i<=1000; ++i) {
double delta = 1./(double)i;
// Next bit is the translation of (add-to-x delta).
{
double diff = delta - x_carry;
double temp = x + diff;
x_carry = (temp-x) - diff;
x = temp;
}
// Next bit is the translation of (add-to-y (- 1 delta)).
{
double diff = (1-delta) - y_carry;
double temp = y + diff;
y_carry = (temp-y) - diff;
y = temp;
}
}
return make_list_2(x, y);
}
I compared the Lisp code with Dann's C++ code for the task of
adding the reciprocals of the first 10^7 positive integers,
as double-precision reals, 100 times.
The C++ code (compiled using g++ 2.95.3, flags -O9)
took 24.9 seconds elapsed, 22.4 seconds CPU.
The Common Lisp code (compiled using CMU CL 18c, speed 3, safety 0)
took 26.4 seconds elapsed, 24.3 seconds CPU.
I wouldn't want to claim that the CL implementation of Kahan
adders, above, is a thing of beauty, but macro definitions
don't have to be (any more than template definitions do).
Much of the oddness you probably feel when reading it is
just unfamiliarity, but defining macros is never going to
be really pretty.
The two bits of code aren't strictly comparable, by the way.
I've omitted the ability to add two Kahan adders together.
In an application where that's important, an approach more
like Dann's (making the Kahan adders real objects, not just
syntactic phantoms) would be appropriate. You can do that
in CL too, and it would actually look nicer than the code
above, but it would be a bit slower.
--
g
This page took 0 seconds to execute
Last modified: Thu, 15 Apr 21 08:11:13 -0700
Current Computer Chess Club Forums at Talkchess. This site by Sean Mintz.