Остановить аппроксимацию в сложном делении в 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 ответа:
Внутри 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
. Вы, вероятно, хотите удалить чек на равенство.