Implementation of expl()
brde at optusnet.com.au
Mon Dec 10 02:16:33 PST 2007
On Mon, 10 Dec 2007, David Schultz wrote:
> On Mon, Dec 10, 2007, Bruce Evans wrote:
>> Better algorithms exist, like someone named das@ used for exp2(), but
>> would be harder to debug and maintain.
> I'm not worried about maintenance, since I don't expect God to add
> any major new features to e^x any time soon.
But every different precision requires a significantly different
implementation (mainly for magic number, but table-driven methods
give a huge number of those) since we don't have algorithms or
tools to generate all the magic numbers from the precision parameter.
> Writing exp2() took a
> lot of reading papers and tinkering, and it's a pain in the neck
> to generate the constants and figure out the resulting error for
> each interval. I seem to recall trying the same tricks for expl(),
> but there were problems with rescaling in base e instead of base 2
> without losing accuracy. In any case, I don't have the kind of time
> needed to fix all of that stuff now.
The problems are well understood :-). I didn't understand them when
you started working on exp2() but do now. Essentially, they are due
to the strange phenomenon that linear functions are usually harder to
approximate (to within 1 ULP) than transcendental functions expanded
in a Taylor series about 0. The linear scaling step with a not-exactly
representable factor like log2(e) ends up being the most difficult point
(unless it is done inaccurately).
> I don't really care how expl() is implemented; anything that works
> is better than anything that doesn't. If someone comes up with a
> better, cleverer scheme later, that's great, but we've been
> talking about a lot of this stuff forever and there's still
> nothing in the tree.
I only agree if "working" includes no new bugs and good accuracy at little
We still don't support 64-bit precision on i386. Fixing that has higher
More information about the freebsd-standards