qwr
qwr

Reputation: 10891

Python Decimal type precision error

I have a quite puzzling issue that I suspect has to do with scientific notation and decimal precision. Here is part of my code:

    def atan(x):
        # Calculate arctan(1/x)
        x = Decimal(x)
        current_value = Decimal(0)
        divisor = 1
        x_squared = x * x
        current_term = 1 / x


        while True:
            current_value += current_term

            divisor += 2
            current_term = (-current_term / x_squared) / divisor
            print(current_term)

            # The issue
            if current_term == Decimal(0):
                break

        return current_value

    print(atan(5))

This is based on the formula atan(1/x) = 1/x - 1/(3x^3) + 1/(5x^5) - ...

However, I discovered that current_term, which gets smaller every loop iteration, is going into values like 4E-80000. Since I've set my decimal precision getcontext().prec to 20, current term should not even support these values. I think somehow current_term is not of decimal type but of scientific notation/float type, but python tells me it's still decimal type.

The correct value for arctan(1/5) is about 0.1973955. I get a value of 0.1973545, which is wrong starting at the 5th digit. Even if I manually break the loop the value is still wrong for some reason. Any help fixing this issue is appreciated.

Upvotes: 2

Views: 1070

Answers (1)

Raymond Hettinger
Raymond Hettinger

Reputation: 226171

You code doesn't match the formula. It got a little too tricky with inferring one term from the next ;-) The 1/(5x^5) term isn't a multiple of 1/(3x^3) term.

Here is code that models the formula directly:

from decimal import Decimal

def atan_recip(x):
    # Calculate arctan(1/x)
    x = Decimal(x)

    total = Decimal(0)
    sign = 1
    for i in range(1, 35, 2):
        total += sign / (i * x ** i)
        sign = -sign
        print(total)

atan_recip(5)

The output is what you expected:

0.2
0.1973333333333333333333333333
0.1973973333333333333333333333
0.1973955047619047619047619047
0.1973955616507936507936507936
0.1973955597889754689754689754
0.1973955598519908535908535908
0.1973955598498063202575202575
0.1973955598498834214339908457
0.1973955598498806620234645299
0.1973955598498807618878454823
0.1973955598498807582406246127
0.1973955598498807583748423407
0.1973955598498807583698713137
0.1973955598498807583700564416
0.1973955598498807583700495142
0.1973955598498807583700497745

Upvotes: 3

Related Questions