2016-07-27 1 views
0

Следующие примеры кода моя проблема, которая не возникает между 10 мощностью 10 и 10 мощностью 11, но для примера приведена в коде и над ним.Python - Как избежать расхождения базы y ** log base y of x, in gmpy2

Я не вижу, где в моем коде я неправильно обрабатываю исходное значение. Может быть, я просто пропустил что-то простое.

Я должен быть уверен, что смогу восстановить x от log x для различных баз. Вместо того, чтобы полагаться на библиотечную функцию, такую ​​как gmpy2, есть ли обратный алгоритм анти-журнала, который гарантирует, что, например, 2**log2(x) он даст x.

Я могу видеть, как напрямую разработать журнал, но не как вернуться, например, серия Тейлора требует много терминов ... How can I write a power function myself? и @ dan04 ответ. Далее следует код.

from gmpy2 import gcd, floor, next_prime, is_prime  
from gmpy2 import factorial, sqrt, exp, log,log2,log10,exp2,exp10  
from gmpy2 import mpz, mpq, mpfr, mpc, f_mod, c_mod,lgamma  
from time import clock  
import random  
from decimal import getcontext 
x=getcontext().prec=1000 #also tried 56, 28 
print(getcontext()) 

def rint():#check accuracy of exp(log(x)) 
    e=exp(1) 
    l2=log(2) 
    l10=log(10) 
    #x=random.randint(10**20,10**21) --replaced with an actual value on next line 
    x=481945878080003762113 
    # logs to different bases 
    x2=log2(x) 
    x10=log10(x) 
    xe=log(x) 
    # logs back to base e 
    x2e=xe/l2 
    x10e=xe/l10 
    # 
    e2=round(2**x2) 
    e10=round(10**x10) 
    ex=round(e**xe) 
    # 
    ex2e=round(2**x2e) 
    ex10e=round(10**x10e) 
    error=5*x-(e2+e10+ex+ex2e+ex10e) 
    print(x,"error sum",error) 
    #print(x,x2,x10,xe) 
    #print(x2e,x10e) 
    print(e2,e10,ex) 
    print(ex2e,ex10e) 
rint() 
+0

Во-первых, я не знаю, как вы получили «random.randint» для работы с этими ограничениями. Это не на моей машине ... Но в любом случае я подозреваю, что это связано с арифметической точностью float. – Aguy

+0

Я думаю, что вы правы, и в этом случае мой вопрос: как мне обойти это, как в python, так и в gmpy2 или в какой-то другой библиотеке? – ocopa

+0

Я могу получить анти-журнал по методу, описанному на https://en.wikipedia.org/wiki/Exponentiation_by_squaring, поэтому для журнала 2 это вопрос не перехода через естественный журнал, а непосредственно вычисление базы журнала 2. – ocopa

ответ

1

Примечание: Я поддерживаю библиотеку gmpy2.

В вашем примере вы используете getcontext() из модуля decimal. Вы не меняете точность, используемую gmpy2. Поскольку точность по умолчанию gmpy2 составляет 53 бита, а для вашего значения x требуется 69 бит, ожидается, что у вас есть ошибка.

Вот скорректированная версия вашего примера, которая иллюстрирует изменение накопленной ошибки при увеличении точности.

import gmpy2 

def rint(n): 
    gmpy2.get_context().precision = n 
    # check accuracy of exp(log(x)) 
    e = gmpy2.exp(1) 
    l2 = gmpy2.log(2) 
    l10 = gmpy2.log(10) 
    x = 481945878080003762113 
    # logs to different bases 
    x2 = gmpy2.log2(x) 
    x10 = gmpy2.log10(x) 
    xe = gmpy2.log(x) 
    # logs back to base e 
    x2e = xe/l2 
    x10e = xe/l10 
    # 
    e2 = round(2**x2) 
    e10 = round(10**x10) 
    ex = round(e**xe) 
    # 
    ex2e = round(2**x2e) 
    ex10e = round(10**x10e) 
    error = 5 * x - (e2 + e10 + ex + ex2e + ex10e) 
    print("precision", n, "value", x, "error sum", error) 

for n in range(65, 81): 
    rint(n) 

И вот результаты.

precision 65 value 481945878080003762113 error sum 1061 
precision 66 value 481945878080003762113 error sum 525 
precision 67 value 481945878080003762113 error sum -219 
precision 68 value 481945878080003762113 error sum 181 
precision 69 value 481945878080003762113 error sum -79 
precision 70 value 481945878080003762113 error sum 50 
precision 71 value 481945878080003762113 error sum -15 
precision 72 value 481945878080003762113 error sum -14 
precision 73 value 481945878080003762113 error sum 0 
precision 74 value 481945878080003762113 error sum -2 
precision 75 value 481945878080003762113 error sum 1 
precision 76 value 481945878080003762113 error sum 0 
precision 77 value 481945878080003762113 error sum 0 
precision 78 value 481945878080003762113 error sum 0 
precision 79 value 481945878080003762113 error sum 0 
precision 80 value 481945878080003762113 error sum 0 
+0

спасибо за это – ocopa

1

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

from decimal import Decimal, getcontext 

getcontext().prec = 1000 

# Just a different method to get the random number: 
x = Decimal(round(10**20 * (1 + 9 * random.random()))) 

x10 = Decimal.log10(x) 
e10 = 10**x10 

e10 - x 
#outputs: Decimal('5.2E-978') 

Для различных баз вы можете захотеть использовать логарифмическую формулу:

x2 = Decimal.log10(x)/Decimal.log10(Decimal('2')) 
e2 = 2**x2 

e2 - x 
#outputs: Decimal('3.9E-978') 
+0

OK , спасибо, это работает для Decimal, log10, но что работает вообще для «разных баз», т. е. разных баз, любая целочисленная база, отличная от естественного log, например 2, для восстановления x? – ocopa

+0

@ocopa - см. Мое правление относительно разных оснований. Просто используйте log10 (x)/log10 (base). – Aguy

+0

Да, я знал это, но я пытался выяснить, почему я получил неправильный ответ от gmpy2, и каков был базовый код. Кажется, gmpy2 использует базу данных журнала за сценой, или, может быть, есть проблема с поплавком, как мы думали.я буду отмечать это как ответ после того, как я напишу пару полезных ссылок, в ближайшее время. – ocopa

0

Aguy решила мою проблему, как признано. Я не принимал во внимание необходимость более 15 цифр точности. Этот ответ на другой вопрос охватывает эту почву. gmpy2 log2 not accurate after 16 digits