math module for Decimals

[Steven D'Aprano]

On Sat, 27 Dec 2008 16:02:51 -0800, jerry.carl.mi wrote:

> Hi,
> 
> I have been looking for a Python module with math functions that would
> both eat and spit Decimals. The standard math module eats Decimals
> allright but spits floats... herefore exp(sin(Decimal())) produces exp
> () of a float :-(


You can write your own wrappers. Assuming you always want to return 
decimals instead of floats:

from decimal import Decimal
import math

def exp(x):
    if isinstance(x, Decimal):
        return x.exp()
    f = math.exp(x)
    return make_decimal(f)    

def sin(x):
    f = math.sin(x)
    return make_decimal(f)    

You can get fancier if you like:

def cos(x):
    try:
        return x.cos()
    except AttributeError:
        f = math.cos(x)
        return make_decimal(f)

That version will automatically use Decimal.cos() in the future, assuming 
Decimals grow a cos() method. If they don't, it will continue to work 
(although a tad slower than the sin() above.


Converting floats to decimals is tricky, because you have to decide how 
much precision you care about. You can keep all the precision of floats. 
Depending on your needs, this is either a feature or a gotcha:

>>> (1.0/3).as_integer_ratio()  # expecting (1, 3)
(6004799503160661L, 18014398509481984L)


Anyway, here are at least three ways to write make_decimal:


# Python 2.6 only
def make_decimal(f):  # exact, no loss of precision
    # however, will fail for f = NAN or INF
    a, b = f.as_integer_ratio()
    return Decimal(a)/Decimal(b)

def make_decimal(f):
    return Decimal("%r" % f)

def make_decimal(f, precision=16):
    # choose how many decimal places you want to keep
    return Decimal.from_float(f, precision)


You choose which suits you best.



> So far, i found:
> 
> -AJDecimalMathAdditions (http://www.ajgs.com/programming/
> PythonForDownload.html)
> -decimalfuncs (http://pypi.python.org/pypi/decimalfuncs/1.4) -and dmath
> (http://pypi.python.org/pypi/dmath/0.9)
> 
> I tried using the AJDecimalMathAdditions, but ran into issues like dSin
> (1E4) would take forever to calculate and would result in sin() > 1 ...
> If i understand it correctly, the program is using maclaurin series to
> calculate the value and since it does not chop off all the multiples of
> 2*pi, the maclaurin approximation becomes useless when its too far from
> x=0.


You might try something like this:

pi = Decimal("%r" % math.pi)
twopi = 2*pi
pi2 = pi/2


def sin(x):
    #  *** UNTESTED ***
    if x < 0:
        return -sin(-x)
    elif x > twopi:
        x = remainder(x, twopi)
        return sin(x)
    elif x > pi:
        return -sin(x - pi)
    elif x > pi2:
        return sin(pi - x)
    else:
        return AJDecimalMathAdditions.dSin(x)


You can treat cos and tan similarly but not identically.

However, you might find that it quicker and simpler to just do your maths 
calculations as floats, then convert to decimals, as in my earlier 
examples.



> I also ran into some issues with the decimalfuncs, but i have never
> tried the dmath thing.
> 
> Now, i can go and start modifying these modules to behave more like the
> standard math module but since i am neither mathematician or programer,
> i really hesitate.

And so you should. Getting numerical functions correct is damn hard.

However, if all you need is "close enough", then don't be scared -- the 
sort of tricks I've shown are probably "close enough", although you may 
want a disclaimer about not using them for running nuclear reactors or 
flying aircraft. The basic tools Python supplies are pretty good, and you 
can probably get a 95% solution with them without being an expert in 
computational mathematics.

 
> So my questions are:
> 
> (1) what do folks use when they need to calculate something like exp
> (sin(Decimal())) or even more complex things? Any recommendations? Or am
> I completely missing something?
> 
> (2) Is there any plan to provide a standard python module that would do
> that? (Python 4.0? ;-)


Guido is very reluctant to put numeric libraries in the standard library, 
precisely because it is so hard to get it right. Don't expect Python to 
become Mathematica any time soon!



-- 
Steven