前提
巨大な値の計算をPythonで行ってみました。
実現したいこと
ideoneのPyPy 2.7.13のdecimalで
(2^((2^(519-1))-1))*(2-(2^-10889035741470030830827987437816582766072))
を計算してみましたが、途中で誤差が出て正しい値が出力されているのか分かりません。
正しい値を知りたいです。先頭100桁とE+ほにゃららで構いません。
ちなみに桁数が10^156ありますので、Python3.10.8やGoogleColabのdecimalでは桁数の桁数が18桁までなので計算できません。
ideoneのPyPy 2では桁数が10^10^5程度までは扱えるようです。(が、一部で誤差が出る。そのためにPython3では桁数制限がかかっていると思われる。)
よろしくお願いいたします。
8.438485139005238033716204113978142132776783207915429197719643649741250839391892405849510296078232762E+258313751232903212140244172706732768962283773495639777019746650509389263410185892774365178070439343478843471144414568958296922314778842870932034756007552127
で正しいかどうかを知りたいわけです。
発生している問題・エラーメッセージ
途中の計算結果:
0E-1000000098
該当のソースコード
Python
1#(2^((2^(w-1))-1))*(2-(2^-t)) 2#2^128倍精度 3#(2^((2^(519-1))-1))*(2-(2^-10889035741470030830827987437816582766072)) 4from decimal import * 5getcontext().prec = 100 6from decimal import localcontext 7with localcontext() as ctx: 8 print(ctx.Emax) # デフォルトの指数上限は999999なので、今回の計算だと足りない 9 ctx.Emax = 10**10**5 # 指数上限を増やす 10 #指数部 11 w = 519 12 #仮数部 13 t = 10889035741470030830827987437816582766072 14 a = ctx.power(2,(w - 1)) 15 b = ctx.power(2,(a - 1)) 16 print(b) 17 c = t * -1 18 d = ctx.power(2,c) 19 print(d) 20 e = (b * 2 - b * d) 21 #最大値 22 print("最大値") 23 print(e) 24 #必要に応じて有効桁を調整する 25 format_str = "{:.1e}".format(e) 26 print(format_str) 27 f = e.log10() 28 #桁数-1 29# print("桁数-1") 30# print(f.quantize(Decimal('0'), rounding=ROUND_FLOOR)) 31 #桁数 32# print("桁数") 33# print(f.quantize(Decimal('0'), rounding=ROUND_CEILING)) 34 g = f.log10() 35 #桁数の桁数-1 36# print("桁数の桁数-1") 37# print(g.quantize(Decimal('0'), rounding=ROUND_FLOOR)) 38 #ウルフラムアルファ確認用 39 print("ウルフラムアルファ確認用") 40 print(g) 41 h = ctx.power(2,t) 42 i = h.log10() 43 #有効桁数 44 print("有効桁数") 45 print(i.quantize(Decimal('0'), rounding=ROUND_FLOOR)) 46
Python+mpmath+gmpy2
1#2^128倍精度 2#(2^((2^(519-1))-1))*(2-(2^-10889035741470030830827987437816582766072)) 3from mpmath import * 4mp.dps = 1000000 5(mpf(2)**((mpf(2)**(mpf(519)-mpf(1)))-mpf(1)))*(mpf(2)-(mpf(2)**mpf(-10889035741470030830827987437816582766072))) 6
試したこと
Python3.10.8(VSCode)やGoogleColabではdecimalでは桁数の桁数が18桁までなので計算できません。
mpmath+gmpy2でも計算してみましたが、dpsが100000まで程度では値の先頭の方、桁数ともに誤差が出て、正しい値がはっきりしません。dpsを1000000にしたら結果がコンソールからはみ出て値の先頭の方が確認できませんでした。また桁数の値がPyPy 2の場合と相違しました。
Python+mpmath+gmpy2でdps = 100の場合
Python+mpmath+gmpy2
1>>> from mpmath import * 2>>> mp.dps = 100 3>>> (mpf(2)**((mpf(2)**(mpf(519)-mpf(1)))-mpf(1)))*(mpf(2)-(mpf(2)**mpf(-10889035741470030830827987437816582766072))) 4mpf('1.738741545572851906007851776954032208812769161318139539364439351610202064301416872947694837078743209829e+258313751232903212140244172706732768962283773495639777019746650509389263410185892774365178070439343470463409842392677007886631633129570799399839688363235724') 5
Python+mpmath+gmpy2でdps = 1000の場合
Python+mpmath+gmpy2
1>>> from mpmath import * 2>>> mp.dps = 1000 3>>> (mpf(2)**((mpf(2)**(mpf(519)-mpf(1)))-mpf(1)))*(mpf(2)-(mpf(2)**mpf(-10889035741470030830827987437816582766072))) 4mpf('8.693707727864259530039258884770161044063845806590697696822196758051010321507084364738474185393716049143428269592664416299506985821607020224451757384081772736694085185953261476092256949585533792382065398255366123346945632455365114496387996812763405070547543636242943582509351695964148881539829278098776930760208512040241459131050857534909186967359333680690876095696574392978402776982242910476888722508873652625966431429199319419162400479389670802250530127697191701685799046457340098916039874145467184794888167675164383011930352737537848617171247730837116311290657456022018847615163555388577340388773985857825851948069000990304471913408925992764573057043682223975592937342873908916124950471446485565735968453329290029508092694774450347618896698260229705612494938176332564765641437441682011291579724298194210916645951676005349723798753644015878111182523064775510486058470694029884469024542601898715230728898371787332915920789625287283111060329377477628692626381070127590574256433161812199884160657947907722e+258313751232903212140244172706732768962283773495639777019746650509389263410185892774365178070439343470463409842392677007886631633129570799399839688363235723') 5
回答2件
あなたの回答
tips
プレビュー
バッドをするには、ログインかつ
こちらの条件を満たす必要があります。
2022/11/11 06:48
2022/11/11 13:07