[Tutor] greater precision?
Peter Otten
__peter__ at web.de
Mon Oct 29 11:43:33 CET 2012
John Collins wrote:
> Hi,
> OS: XP
> Py: 2.7.2
>
> I am running some old scripts, not written by me, that produce
> substantial calculated values data as decimals, 12 significant
> figures. AFAIK, the keys calls are;
[On modern hardware] Python uses IEEE 754 double-precision
<http://en.wikipedia.org/wiki/Double_precision> internally which gives 15 to
17 digits. But of course errors may accumulate.
> def crosspoint(xa1,ya1,xa2,ya2,xb1,yb1,xb2,yb2):
> "Give the intersection point of the (possibly extrapolated) lines\n"\
> "segments (xa1,ya1)-(xa2,ya2) and (xb1,yb1)-(xb2,yb2)."
> dxa = xa2-xa1
> dya = ya2-ya1
> dxb = xb2-xb1
> dyb = yb2-yb1
> if dya * dxb == dxa * dyb:
> return None
> if dxa == 0:
> return (xa1, (xa1 - xb1) * dyb / dxb + yb1)
> if dxb == 0:
> return (xb1, (xb1 - xa1) * dya / dxa + ya1)
> if dya == 0:
> return ((ya1 - yb1) * dxb / dyb + xb1, ya1)
> if dyb == 0:
> return ((yb1 - ya1) * dxa / dya + xa1, yb1)
>
> det = dya * dxb - dyb * dxa
> xtop = dxb * dxa * (yb1-ya1) + dya * dxb * xa1 - dyb * dxa * xb1
> ytop = dya * dyb * (xa1-xb1) + dxb * dya * yb1 - dxa * dyb * ya1
>
> return (xtop / det, ytop / det)
>
>
> I am not sure what == means,
Equality.
>>> 1 == 0
False
>>> 1 == 1
True
Because of rounding errors this is a dangerous operation for floating point
numbers:
>>> print 1.0 == .1 + .1 + .1 + .1 + .1 + .1 + .1 + .1 + .1 + .1
False
> nor if any 'special' maths functions are
> called from crosspoint.py (to me, it appears not?), so, as I understand
> it, the
>
> from math import
>
> line *is* the 'math engine' if you will, and is the source of the 12
> sig fig limit, yes?
It is a toolbox rather than an engine.
math is basically a wrapper around the C library, and these functions all
use C's double and [most of the time] your hardware's floating point unit.
> One other odd thing that occurs, is when the scripts are run, a weird,
> small, non ASCII file called crosspoint.pyc is created? Is this the
> interpreter 'compiling' crosspoint.py ( either before, during, or
> after?) the 2nd script calls it?
Yes, Python compiles the source to bytecode before it executes it, and for
all modules except the main script this bytecode is cached in a .pyc file.
> Sorry to bother, but if I can squeeze out *at least* 15 sig figs, (30
> or more would be better!) I'd be a happy camper!
As you correctly observed crospoints() uses only +-*/ and ==, so you can
feed it every type that supports these operations (this is called "duck
typing"). For example:
>>> from fractions import Fraction
>>> from crosspoint import crosspoint
>>> crosspoint(Fraction(0, 1), Fraction(1, 3), Fraction(1, 1), Fraction(1,
3), Fraction(5, 7), Fraction(0, 1), Fraction(5, 7), Fraction(1, 1))
(Fraction(5, 7), Fraction(1, 3))
>>> p = crosspoint(Fraction(0, 1), Fraction(1, 3), Fraction(1, 1),
Fraction(1, 3), Fraction(5, 7), Fraction(0, 1), Fraction(5, 7), Fraction(1,
1))
>>> print "%s, %s" % p
5/7, 1/3
Yay, infinite precision ;)
Of more practical relevance may be something like gmpy
<http://code.google.com/p/gmpy/>
More information about the Tutor
mailing list