我想知道是否有這樣做的演算法。
考慮sin(60.9^100)例如，
如何60.9^100 - 2pi *N 在2 pi范圍內制作。
uj5u.com熱心網友回復：
這里有兩種半自動方式。他們仍然需要手動配置，具體取決于輸入。也許可以完全自動化，但至少它們提供了某種方式（以及檢查進一步方法的結果）。
使用 Python 的decimal模塊及其pi配方，顯然 60.9^100 % 2pi 大約是 0.4826（然后可以給sin）。以 180 到 290 位精度計算的結果（代碼在末尾）：
180 0.52113386128181643243087541425218797675113893601959695702254815952665...
190 0.48262316221535366629016856570286348468626721388098119253199769818223...
200 0.48262316221828443267196371207773451732899712100881145938907295835606...
210 0.48262316221828443267246775208563277802202330286500415343966588161647...
220 0.48262316221828443267246775208566344687793590859019274697600998645752...
230 0.48262316221828443267246775208566344687793479998648411772237563268709...
240 0.48262316221828443267246775208566344687793479998648362494580989872984...
250 0.48262316221828443267246775208566344687793479998648362494580991864728...
260 0.48262316221828443267246775208566344687793479998648362494580991864728...
270 0.48262316221828443267246775208566344687793479998648362494580991864728...
280 0.48262316221828443267246775208566344687793479998648362494580991864728...
290 0.48262316221828443267246775208566344687793479998648362494580991864728...
Wolfram Alpha失敗了，計算“零”。但對于指數 30，它仍然顯示有效結果，我們匹配：
WolframAlpha: 6.0148312092022347033088447399833343520115646793565705028401966310...
Mine: 6.01483120920223470330884473998333435201156467935657050284019663107410...
另一種方法，使用從某個站點復制的 pi 的前 1001 位數字，并使用整數直到最后，給出 0.48262316221828444（在線嘗試！）：
a, b = 609, 100
pi = 31415926535897932384626433832795028841971693993751058209749445923078164062862089986280348253421170679821480865132823066470938446095505822317253594081284811174502841027019385211055596446229489549303819644288109756659334461284756482337867831652712019091456485669234603486104543266482133936072602491412737245870066063155881748815209209628292540917153643678925903600113305305488204665213841469519415116094330572703657595919530921861173819326117931051185480744623799627495673518857527248912279381830119491298336733624406566430860213949463952247371907021798609437027705392171762931767523846748184676694051320005681271452635608277857713427577896091736371787214684409012249534301465495853710507922796892589235420199561121290219608640344181598136297747713099605187072113499999983729780499510597317328160963185950244594553469083026425223082533446850352619311881710100031378387528865875332083814206171776691473035982534904287554687311595628638823537875937519577818577805321712268066130019278766111959092164201989
print(a**b * 10**(1000-b) % (2*pi) / 10**1000)
這對大整數進行操作，“放大”10 ^ 1000，直到最終的縮小除法給出一個浮點數。
第三種方式，使用 Python 的fraction，同樣得到 0.48262316221828444（在線試用！）：
from fractions import Fraction
a, b = Fraction('60.9'), 100
pi = Fraction('3.1415926535897932384626433832795028841971693993751058209749445923078164062862089986280348253421170679821480865132823066470938446095505822317253594081284811174502841027019385211055596446229489549303819644288109756659334461284756482337867831652712019091456485669234603486104543266482133936072602491412737245870066063155881748815209209628292540917153643678925903600113305305488204665213841469519415116094330572703657595919530921861173819326117931051185480744623799627495673518857527248912279381830119491298336733624406566430860213949463952247371907021798609437027705392171762931767523846748184676694051320005681271452635608277857713427577896091736371787214684409012249534301465495853710507922796892589235420199561121290219608640344181598136297747713099605187072113499999983729780499510597317328160963185950244594553469083026425223082533446850352619311881710100031378387528865875332083814206171776691473035982534904287554687311595628638823537875937519577818577805321712268066130019278766111959092164201989')
print(float(a**b % (2*pi)))
使用十進制的代碼（在線嘗試！）：
from decimal import *
a, b = '60.9', 100
def pi():
"""Compute Pi to the current precision.
>>> print(pi())
3.141592653589793238462643383
"""
getcontext().prec = 2 # extra digits for intermediate steps
three = Decimal(3) # substitute "three=3.0" for regular floats
lasts, t, s, n, na, d, da = 0, three, 3, 1, 0, 0, 24
while s != lasts:
lasts = s
n, na = n na, na 8
d, da = d da, da 32
t = (t * n) / d
s = t
getcontext().prec -= 2
return s # unary plus applies the new precision
for prec in range(100, 300, 10):
setcontext(Context(prec=prec, Emax=MAX_EMAX, Emin=MIN_EMIN))
x = Decimal(a) ** b
try:
print(prec, str(x % (2*pi()))[:70] '...')
except:
pass
uj5u.com熱心網友回復：
使用任意精度浮點數和600小數位進行除法精度會導致：
a = 60.89999999999999999999999999999999999999999999999999999999999999999999999999999999999999999999999999999999999999999999999999999999999999999999999999999999999999999999999999999999999996385120202345673823
b = 100
pi2 = 6.28318530717958647692528676655900576839433879875021164194988918461563281257241799725606965068423413596429617302656461329418768921910116446345071881625696223490056820540387704221111928924589790986076393
a^b = 28955379630372944287405172428351524098232652637160696453610571103951984391040712391370032514984003111097712417114477560530145388673931838194276269928567117043004481189909961875857.75649514478978553782864231122006591641737915748789056749200619756952612910126893035405495135890665595449781429534911495188619709183484163296347386211939610629876033802394440562845818548327528930373773
a^b mod pi2 = 0.48090443944555460103493387524923864059291388194583637380891303988368325373623569292130225575175455300500465326465846619858936543877561381898189967500531981888897878580530179993696579483550649680773777
這將作為double計算的參考：
pow(60.9,100) = 2.895537963037287827000000000000000000000e 178
正如您所看到的，舍入誤差的大小~10^160遠遠大于所需的模數，導致結果完全偏離......
通過平方和模（modpow）使用冪來保持子結果很小，就像 Blackgaurd 的回答一樣：
//---------------------------------------------------------------------------
double mod(double a,double p)
{
a/=p;
a-=floor(a);
a*=p;
return a;
}
//---------------------------------------------------------------------------
double modpow(double a,unsigned int b,double p)
{
int i;
double d;
a=mod(a,p);
d=mod(1.0,p);
for (i=0;i<32;i )
{
d=mod(d*d,p);
if (DWORD(b&0x80000000))
{
d=mod(d*a,p);
}
b<<=1;
}
return mod(d,p);
}
//---------------------------------------------------------------------------
看起來稍微好一點，但仍然很差，因為在每次迭代中我們丟失了尾數的精度位：
modpow(60.9,100,2.0*M_PI) = 0.000000000004526729
然而，正如 Kelly Bundy 指出的那樣，浮點數的模不能在pow迭代內部或之前使用，因為它破壞了數學（與整數不同），使得這種方法無法使用。
因此，計算它的唯一方法（除非有一些很好的數學恒等式可用于此）是使用更大的位寬數（定點或浮點數）。的結果60.9^100有 178 個十進制數字，因此您至少需要使用：
178/log10(2) = 591.3 bits
在上面的計算中，我為尾數選擇了 600 位（通過平方而不在迭代中使用 mod）。請注意，您需要PI在相似的位寬中保持恒定才能使這項作業。
所以你的問題的答案是：不，我看不出有辦法計算這個double，你必須使用大數字。
這取決于您要在...中實作此功能的語言
[Edit1]這里是根據不同精度位寬計算的東西：
600 bits: a^b mod p = 0.48090443944543787059955763581639220817922063926527443008616966832526276298682799066248021870151726073022025792461359895510004764246892419475030025711547915571989225146813845497084456179528624932006095
1000 bits: a^b mod p = 0.48262316221828443267246775208566344687793479998648362494580991864728550767313153551910251537428629715050683896930935124670316833678465659261329197846231478482857003648065425433322518882397051713086293
10000 bits: a^b mod p = 0.48262316221828443267246775208566344687793479998648362494580991864728550767313153551910251537428629715050683896930935124670316834694077218601842581288592371259021844778135739430850760935223120750400005
uj5u.com熱心網友回復：
完整的 Python 實作將我的第一個答案的半自動decimal方法與 Spektre 對如何計算所需精度的觀察相結合。
我們的輸入sin是a b % 2π，它位于[0, 2π)，所以它的整數部分有一位。除了那個數字之外，我想要 15 個小數位，因為 Python 的float精度約為 16 位。因此，我使用比 a b的整數部分多 30 位的數字來獲得 30 個小數位（大約 15 位就足夠了，但再多一些也無妨，而且我不想嘗試準確計算需要多少位）。
from math import log10, sin
from decimal import setcontext, Context, getcontext, Decimal
def bigsin(a, b):
prec = int(log10(a)*b) 30
setcontext(Context(prec=prec, Emax=prec 10))
a = Decimal(str(a))
return sin(float(a**b % (2*pi())))
def pi():
"""Compute Pi to the current precision.
>>> print(pi())
3.141592653589793238462643383
"""
getcontext().prec = 2 # extra digits for intermediate steps
three = Decimal(3) # substitute "three=3.0" for regular floats
lasts, t, s, n, na, d, da = 0, three, 3, 1, 0, 0, 24
while s != lasts:
lasts = s
n, na = n na, na 8
d, da = d da, da 32
t = (t * n) / d
s = t
getcontext().prec -= 2
return s # unary plus applies the new precision
print(bigsin(60.9, 100))
輸出（在線嘗試！）：
0.4641043156966329
請注意，這與我的第一個答案不同，因為在這里我計算了正弦值，而不僅僅是模冪。
60.9另請注意，雖然我以a開頭float（不完全代表 60.9），但str將其轉換為string 60.9，因此Decimalthen確實代表 60.9。
uj5u.com熱心網友回復：
以下是@Marat 在 C 中的評論示例：
double dmod(double x, double y) {
return x - (int)(x/y) * y;
}
double bin_pow(double a, double exp, double mod) {
if (!exp) return 1.0;
double ret = 1;
while (exp > 0) {
if (exp & 1) ret = dmod(ret * a, mod);
a = dmod(a * a, mod);
exp >>= 1;
}
return ret;
}
轉載請註明出處，本文鏈接：https://www.uj5u.com/gongcheng/486308.html