Floating-point glitches with the math module. Bug? Or am I missing something?

[Alex Martelli]

Chris Metzler <cmetzler at speakeasy.snip-me.net> wrote:

> On Mon, 20 Sep 2004 19:03:57 -0500, Jeff Epler wrote:
> >
> > Incidentally, the second part of your expression goes to 0 here, with
> > the given values:
> 
> Yeah, I noticed this later (it's sorta implicit in the content of my
> other, followup post).
> 
> This kinda sucks because I need to do some precision surveying calcs
> -- computing lat/lons for a variety of points separated by known
> distances and angles from a point of known lat/lon. If I do everything
> in spherical trig, I run into this issue (that whole dividing-by-earth
> -radius thing).   But if I do things using planar approximations, for
> most of the points I get lat/lons that are too far off from where they
> should be.

I wonder if IBM's Accurate Portable Math Library might help you: they
guarantee precision to the ULP for any function they offer, and of
course you could talk to it from Python via ctypes, pyrex, swig, ...
(haven't tried it out myself, though).

APMathlib's download page is at:
<http://www-124.ibm.com/developerworks/downloads/?group_id=28>

If IEEE 754's precision just ain't enough for you, welcome to the
wonderful but slow world of multi-precision arithmetic!-)

GMP (which gmpy wraps) doesn't do transcendentals.  MPFR, I believe,
does, but I know of no Python interface for it.  Still, check out their
page, it's full of useful links -- <http://www.mpfr.org/> -- worst case,
you can use gmpy and program the transcendentals you need as per
the various well-known series, etc... but you're likely to be able to
find something more direct if you look around a little!


Alex