Rmpfr可以做到使用mpfr_set_str字符串轉換...
val <- mpfr("1e309")
## 1 'mpfr' number of precision 17 bits
## [1] 9.999997e308
# set a precision (assume base 10)...
est_prec <- function(e) floor(e/log10(2)) + 1
val <- mpfr("1e309", est_prec(309))
## 1 'mpfr' number of precision 1027 bits
## [1]1000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000
.mpfr2bigz(val)
## Big Integer ('bigz') :
## [1] 1000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000
# extract exponent from a scientific notation string
get_exp <- function(sci) as.numeric(gsub("^.*e",'', sci))
# Put it together
sci2bigz <- function(str) {
.mpfr2bigz(mpfr(str, est_prec(get_exp(str))))
}
val <- sci2bigz(paste0(format(Const("pi", 1027)), "e309"))
identical(val, .mpfr2bigz(Const("pi",1027)*mpfr(10,1027)^309))
## [1] TRUE
## Big Integer ('bigz') :
## [1] 3141592653589793238462643383279502884197169399375105820974944592307816406286208998628034825342117067982148086513282306647093844609550582231725359408128481117450284102701938521105559644622948954930381964428810975665933446128475648233786783165271201909145648566923460348610454326648213393607260249141273724587004
至於爲什麼上存儲的數大於.Machine$double.xmax
,上浮點編碼在IEEE規範的文檔時,R常見問題和維基百科進入所有的行話,但我發現它有助於只定義條款(使用?'.Machine'
)...
double.xmax
(最大標準化浮點數)=
(1 - double.neg.eps) * double.base^double.max.exp
其中
double.neg.eps
(一個小的正浮點數x使得1 - X = 1)!= double.base^double.neg.ulp.digits
其中
double.neg.ulp.digits
=最大負整數,使得1 - double.base^i != 1
和
double.max.exp
= double.base溢出的最小正功率和
double.base
(浮點數表示的基數)= 2(二進制數)。
考慮有限浮點數可以與另一個區別開來嗎?在IEEE規範告訴我們,對於一個binary64數11位都將用於指數,所以我們的2^(11-1)-1=1023
最大的指數,但我們希望的是溢出最大指數所以double.max.exp
是1024
# Maximum number of representations
# double.base^double.max.exp
base <- mpfr(2, 2048)
max.exp <- mpfr(1024, 2048)
# This is where the big part of the 1.79... comes from
base^max.exp
## 1 'mpfr' number of precision 2048 bits
## [1] 179769313486231590772930519078902473361797697894230657273430081157732675805500963132708477322407536021120113879871393357658789768814416622492847430639474124377767893424865485276302219601246094119453082952085005768838150682342462881473913110540827237163350510684586298239947245938479716304835356329624224137216
# Smallest definitive unit.
# Find the largest negative integer...
neg.ulp.digits <- -64; while((1 - 2^neg.ulp.digits) == 1)
neg.ulp.digits <<- neg.ulp.digits + 1
neg.ulp.digits
## [1] -53
# It makes a real small number...
neg.eps <- base^neg.ulp.digits
neg.eps
## 1 'mpfr' number of precision 2048 bits
## [1] 1.11022302462515654042363166809082031250000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000e-16
# Largest difinitive floating point number less than 1
# times the number of representations
xmax <- (1-neg.eps) * base^max.exp
xmax
## 1 'mpfr' number of precision 2048 bits
## [1] 179769313486231570814527423731704356798070567525844996598917476803157260780028538760589558632766878171540458953514382464234321326889464182768467546703537516986049910576551282076245490090389328944075868508455133942304583236903222948165808559332123348274797826204144723168738177180919299881250404026184124858368
identical(asNumeric(xmax), .Machine$double.xmax)
## [1] TRUE
您可能會發現閱讀有用的文檔,特別是'?as.bigz'的Note部分。 – joran 2013-02-11 18:22:43
謝謝喬蘭。我錯過了最後一行。令人煩惱的是,我現在不能使用科學記數法! – 2013-02-11 18:42:31
是的,但你可以做'as.bigz(10)^ 309'。事實上,你可以這樣做:''%e%'< - function(x,y)as.bigz(x)* 10^as.bigz(y); 1%e%309' – 2013-02-11 19:18:15