antiloquax
Posts: 406
Joined: Sun Nov 20, 2011 11:37 am
Contact: Website

displaying long need help with floating point numbers

Thu Oct 04, 2012 7:36 pm

Hi,
I have written a little program that works out Pi using atan. I can get it to give me about 50 decimal places, but I really don't understand the "Decimal" module. If anyone can give me some tips on how to get a more precise result, I'd appreciate it.
Here's my code:

Code: Select all

# pi.py - arctan(1) * 4 = pi
from math import *
from decimal import *

err = 0.0000000000000001

def arctan(n, err):
    """Uses Gregory's formula for calcuating atan."""
    temp = n
    atan = 0
    i = 3
    while (abs(atan - n) > err):
        atan = n
        n = n - (pow(temp, i)/i) + ((pow(temp, i + 2)) / (i + 2))
        i += 4
    return n

# set up the start of the fibonacci series
a = 1
b = 2
atan = 0
temp = 5

# this calculates pi/4 using Euler's formula and fibonacci numbers
while (abs(temp - atan) > err):
    temp = atan
    atan += arctan(1/b, err)
    a = b + a
    b = b + a
atan += arctan(1/b, err)
pi = atan * 4
print(Decimal.from_float(pi))
The output I get is:

Code: Select all

3.141592653589793560087173318606801331043243408203125
thanks
mark

BlackJack
Posts: 288
Joined: Sat Aug 04, 2012 8:28 am
Contact: Website

Re: displaying long need help with floating point numbers

Fri Oct 05, 2012 6:59 pm

@antiloquax: You are doing your calculations with `float` objects which are implemented with C's `double` type. Converting values of that type into a `Decimal` object *after* the calculations were done can't magically enhance the precision of the `float` calculations.

Code: Select all

while not self.asleep():
    sheep += 1

antiloquax
Posts: 406
Joined: Sun Nov 20, 2011 11:37 am
Contact: Website

Re: displaying long need help with floating point numbers

Tue Oct 09, 2012 4:52 pm

Thanks BlackJack,
I have had another go at doing this using Decimal. I must have missed something important, because it doesn't make any difference.
This is what I came up with:

Code: Select all

# pi.py - arctan(1) * 4 = pi
from math import *
from decimal import *
getcontext().prec = 50 

def arctan(n):
    """Uses Gregory's formula for calculating atan."""
    temp = n 
    atan = Decimal(0)
    i = Decimal(3)
    while (atan !=  n):
        pow1 = (Decimal.from_float(pow(temp, i)) / i)
        pow2 = (Decimal.from_float(pow(temp, i + 2))) / (i + 2)
        atan = n 
        n = n - pow1 + pow2  
        i += 4
    return n

def euler(a, b):
    """Using Euler's formula."""
    euler = Decimal(0)
    temp = Decimal(5)
    while (temp != euler): 
        temp = euler
        euler += arctan(Decimal(1) / b)
        a = b + a
        b = b + a
    return euler

pi = euler(Decimal(1), Decimal(2)) * Decimal(4)
print(pi)
     
    
I am considering using something like Python's mpmath or bigfloat.
My main aim in this program was just to demonstrate the algorithm, but it would be nice to do a 100 or so places!

BlackJack
Posts: 288
Joined: Sat Aug 04, 2012 8:28 am
Contact: Website

Re: displaying long need help with floating point numbers

Wed Oct 10, 2012 11:53 am

@antiloquax: 1. Don't use `float` values → replace `pow()` by using the power operator (``**``) to keep the precision of `Decimal` in all steps.

2. ``getcontext().prec = 50`` sets the precision to 50. If you want 100 then guess what you have to do. ;-)

Most `Decimal()` calls can be replaced with just integers, because just like a calculation with `float` and `int` results in `float`, a calculation with `Decimal` and `int` results in a `Decimal` value.

The braces around the ``while`` conditions are not necessary.

Code: Select all

# pi.py - arctan(1) * 4 = pi
from decimal import Decimal, localcontext


def arctan(n):
    """Uses Gregory's formula for calculating atan."""
    result = n
    previous = 0
    i = 3
    while previous != result:
        previous = result
        result += (n**(i + 2) / (i + 2)) - (n**i / i)
        i += 4
    return result


def euler(a, b):
    """Using Euler's formula."""
    result = 0
    previous = 5
    while previous != result:
        previous = result
        result += arctan(1 / b)
        a += b
        b += a
    return result


def main():
    with localcontext() as context:
        context.prec = 100
        pi = euler(Decimal(1), Decimal(2)) * 4
        print(pi)


if __name__ == '__main__':
    main()

Code: Select all

while not self.asleep():
    sheep += 1

antiloquax
Posts: 406
Joined: Sun Nov 20, 2011 11:37 am
Contact: Website

Re: displaying long need help with floating point numbers

Wed Oct 10, 2012 12:07 pm

Thanks BlackJack -
I think it was the pow bit that was the main stumbling block
Yeah - I changed it to 50, just because I was trying out a few changes and wanted each run to be a bit quicker.
Thanks for your patience with me! :oops:

EDIT: Yep that worked. I am v. happy with this as I am hoping to do a lesson / part of lesson in the Maths dept. on links between Maths and CSc. I also realised that I only needed one Decimal() - for the value passed to euler as "b". That made all the other values Decimal, as you pointed out.

:mrgreen:

BlackJack
Posts: 288
Joined: Sat Aug 04, 2012 8:28 am
Contact: Website

Re: displaying long need help with floating point numbers

Wed Oct 10, 2012 6:05 pm

I have factored out the Fibonnaci number generation and the „sum the values until the result doesn't change anymore” that was part of both `arctan()` and `euler()`:

Code: Select all

from __future__ import division, print_function
from decimal import Decimal, localcontext
from itertools import count, islice


def iter_fibonacci_numbers():
    a, b = 0, 1
    while True:
        yield a
        a, b = b, a + b


def sum_until_no_change(numbers, start=0):
    result, previous = start, None
    for number in numbers:
        result += number
        if previous == result:
            break
        previous = result
    return result


def arctan(n):
    """Uses Gregory's formula for calculating atan."""
    return sum_until_no_change(
        (
            (n**(i + 2) / (i + 2)) - (n**i / i)
            for i in (k * 4 + 3 for k in count())
        ),
        n
    )


def euler(fibonacci_numbers):
    """Using Euler's formula."""
    return sum_until_no_change(
        arctan(1 / i) for i in islice(fibonacci_numbers, 3, None, 2)
    )


def main():
    with localcontext() as context:
        context.prec = 100
        fibonacci_numbers = (Decimal(i) for i in iter_fibonacci_numbers())
        pi = euler(fibonacci_numbers) * 4
        print(pi)


if __name__ == '__main__':
    main()
The `__future__` import is for Python 2.x, so it works with Python 2 and 3.

Code: Select all

while not self.asleep():
    sheep += 1

antiloquax
Posts: 406
Joined: Sun Nov 20, 2011 11:37 am
Contact: Website

Re: displaying long need help with floating point numbers

Sat Oct 27, 2012 4:26 pm

Hi blackjack, sorry I have not replied before now. Thanks for these improvements. :D

Return to “Python”