Python - как избежать несоответствия базы y ** логарифма y x в gmpy2

В следующих примерах кода моя проблема не возникает между 10 степенями 10 и 10 степенями 11, но возникает для примера, приведенного в коде и выше.

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

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

Я вижу, как напрямую разработать журнал, но не вижу, как вернуться, например, для серии Тейлора требуется много терминов ... Как мне самому написать степенную функцию? и ответ @ 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()

python logarithm gmpy
person ocopa    schedule 27.07.2016    source источник
comment
Во-первых, я не знаю, как вам random.randint удалось работать с этими ограничениями. На моей машине этого нет ... Но в любом случае я подозреваю, что это связано с точностью арифметических вычислений с плавающей запятой.   -  person Aguy    schedule 27.07.2016
comment
Я думаю, что вы правы, и в этом случае у меня вопрос, как мне обойти это, в python, gmpy2 или какой-либо другой библиотеке?   -  person ocopa    schedule 27.07.2016
comment
Я могу получить анти-журнал с помощью метода, описанного на странице en.wikipedia.org/wiki/Exponentiation_by_squaring , поэтому для журнала 2 это вопрос не использования натурального журнала, а прямого вычисления базы журнала 2.   -  person ocopa    schedule 27.07.2016

Ответы (3)


arrow_upward
1
arrow_downward

Примечание. Я поддерживаю библиотеку 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
person casevh    schedule 29.07.2016

arrow_upward
1
arrow_downward

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

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')
person Aguy    schedule 27.07.2016
comment
Хорошо, спасибо, это работает для Decimal, log10, но что в целом работает для различных оснований, т.е. разных оснований, любой целочисленной базы, кроме натурального логарифма, например 2, для восстановления x? - person ocopa; 28.07.2016
comment
@ocopa - см. мою правку о разных базах. Просто используйте log10 (x) / log10 (base). - person Aguy; 28.07.2016
comment
Да, я знал это, но я пытался понять, почему я получаю неправильный ответ от gmpy2 и каков основной код. Похоже, что gmpy2 использует базу журналов e за сценой, или, возможно, есть проблема с точностью с плавающей запятой, как мы думали. Я отмечу это как ответ после того, как вскоре опубликую пару полезных ссылок. - person ocopa; 28.07.2016
comment
Полезные ссылки на логарифмы en.wikipedia.org/wiki/Binary_logarithm, ответ от Dan04 на stackoverflow.com/questions/2882706/ ответьте log0 в stackoverflow.com/questions/3719631/ - person ocopa; 28.07.2016

arrow_upward
0
arrow_downward

Как признался, Агуй решил мою проблему. Я не учел необходимости более 15 знаков точности. Этот ответ на другой вопрос охватывает это основание. gmpy2 log2 неточен после 16 цифр

person ocopa    schedule 28.07.2016