Остановить аппроксимацию в сложном делении в Python


Я писал некоторый код, чтобы перечислить гауссовы целочисленные делители рациональных целых чисел в Python. (Относящийся к задаче 153 проекта Эйлера)

Я, кажется, достиг некоторой проблемы с некоторыми числами, и я думаю, что это связано с Python, аппроксимирующим деление комплексных чисел.

Вот мой код для функции:

def IsGaussian(z):
    #returns True if the complex number is a Gaussian integer
    return complex(int(z.real), int(z.imag)) == z

def Divisors(n):
    divisors = []

    #Firstly, append the rational integer divisors
    for x in range(1, int(n / 2 + 1)):
        if n % x == 0:
            divisors.append(x)

    #Secondly, two for loops are used to append the complex Guassian integer divisors
    for x in range(1, int(n / 2 + 1)):
        for y in range(1, int(n / 2 + 1)):
            if IsGaussian(n / complex(x, y)) == n:
                divisors.append(complex(x, y))
                divisors.append(complex(x, -y))

    divisors.append(n)

    return divisors

Когда я выполняю Divisors(29), я получаю [1, 29], но при этом упускаю четыре других делителя, один из которых (5 + 2j), который можно ясно видеть, чтобы разделите на 29.

При запуске 29 / complex(5, 2), Python дает (5 - 2.0000000000000004j)

Этот результат неверен, как и должно быть (5 - 2j). Есть ли способ как-то обойти приближение питона? И почему эта проблема не поднялась для многих других рациональных целых чисел меньше 100?

Заранее спасибо за вашу помощь.

2 2

2 ответа:

Внутри CPython использует пару поплавков двойной точности для комплексных чисел. Поведение численных решений в целом слишком сложно, чтобы суммировать здесь, но некоторая ошибка неизбежна в численных расчетах.

Например:

>>>print(.3/3)
0.09999999999999999
Таким образом, при тестировании решений такого рода часто правильно использовать приближенное равенство, а не фактическое равенство.

Функцияisclose в модуле cmath доступна именно по этой причине.

>>>print(.3/3 == .1)
False
>>>print(isclose(.3/3, .1))
True

Этот вид вопросов является областью численного анализа; это может быть полезным тегом для дальнейших вопросов на эту тему.

Обратите внимание, что считается "питоническим", если идентификаторы функций находятся в snake_case.

from cmath import isclose
def is_gaussian(z):
    #returns True if the complex number is a Gaussian integer
    rounded = complex(round(z.real), round(z.imag))
    return isclose(rounded, z)

Вы можете определить Эпсилон, используя round для округления до нужного числа десятичных знаков / точности (например, 10):

def IsGaussian(z, prec=10):
    # returns True if the complex number is a Gaussian integer
    # rounds the input number to the `prec` number of digits
    z = complex(round(z.real,prec), round(z.imag,prec))
    return complex(int(z.real), int(z.imag)) == z

У вашего кода есть еще одна проблема:

if IsGaussian(n / complex(x, y)) == n:
Это даст результаты только для n = 0 или n = 1. Вы, вероятно, хотите удалить чек на равенство.