Алгоритм вычисления логарифма base-10 в Python

Question

Я попытался создать программу для вычисления логарифма по основанию 10 на основе алгоритма на основе рядов Тейлора, описанного в "Руководстве по вычислению математических функций" (я нашел онлайн-копию через библиотеку моего университета).

Аналогичный алгоритм приведен для другого вопроса о Stackru, для которого я не могу сейчас найти ссылку.

10.3.2. Вычисление логарифмов в десятичной основе

Для десятичной основы логарифм base-10 является естественным выбором, и разложение аргумента на экспоненту и дробь дает нам десятичное представление: x = (−1)^s × f × 10^n, либо f = 0 в точности, или f находится в [1/10, 1).

Если f ≤√1/10, установите f = 10 × f и n = n - 1, так что теперь f находится в интервале (√1/10,√10]. Затем введите замену переменной, ряд Тейлора расширение и полиномиальное представление этого расширения:

z = (f - 1)/(f + 1),

f = (1 + z)/(1 - z),

D = 2 log10 (e)

= 2 / log (10)

log10 (f) = D × (z + z3 / 3 + z5 / 5 + z7 / 7 + z9 / 9 + z11 / 11 + · · ·)

≈ D × z + z3Q(z2), полиномиальная подгонка включает D в Q(z2).

Для f в (√1/10,√10] мы имеем z примерно [-0,5195,+0,5195]. Более широкий диапазон z требует более длинных полиномов по сравнению с двоичным случаем, а также делает поправочный член z3Q(z2) относительно больше. Его величина не превышает |0,35z|, поэтому он обеспечивает лишь одну дополнительную десятичную цифру точности вместо двух. Точное вычисление z проще, чем в двоичном случае: просто установите z = fl(fl(f- 12) -12)/ П (F + 1).

Для этого я написал эту программу на Python:

def log10(x):

n = 0.0 #Start exponent of base 10

while (x >= 1.0):
    x = x/10.0
    n+=1


# if x <= sqrt(1/10)
if(x<=0.316227766016838):
    x = x*10.0
    n = n-1

#Produce a change of variable
z = (x-1.0)/(x+1.0)
D = 4.60517018598809 #2*log10(e)

sum = z
for k in range(3,111,2):
    sum+=(z**k)/k

return D*n*sum

Я сравнил результаты с math.log10 функция, и результаты не такие, как ожидалось. Моя самая большая проблема при отладке - это понимание алгоритма и почему он работает.

Answer 1

Вот мой исходный код после предложенных исправлений (изменен оператор возврата на D*sum+n фиксированное значение Dи изменил if(x<=0.316227766016838) в while(x<=0.316227766016838), Я добавил немного if заявления для обработки исключительных случаев.

Приведенный ниже код хорошо работает в пределах моей целевой точности 6 цифр (я тестировал его с очень маленьким вводом, большим вводом).

def log10(x):

    # Handle exceptional cases
    if (x == 1):
        return 0
    if (x == 0):
        return float('-Inf')
    if (x < 0):
        return float('nan')

    n = 0 #Start exponent of base 10

    while (x >= 1.0):
        x = x/10.0
        n+=1

    # if x <= sqrt(1/10)
    while(x<=0.316227766016838):
        x = x*10.0
        n = n-1

    #Produce a change of variable
    z = (x-1.0)/(x+1.0)
    D = 0.868588964 #2*log10(e)

    #Taylor series
    sum = z
    for k in range(3,23,2):
        sum+=(z**k)/k

    return D*sum+n