Nothing
#### Low level stuff - debugging etc
#### ========= =========
require("Rmpfr")
options(warn = 2)# warning -> error
identical3 <- function(x,y,z) identical(x,y) && identical (y,z)
identical4 <- function(a,b,c,d) identical(a,b) && identical3(b,c,d)
## sane state [when re-source()ing this file]:
.mpfr_erange_set("Emin", -(2^30-1))
.mpfr_erange_set("Emax", +(2^30-1))
###----- _1_ mpfr1 , import, xport etc -----------------------------------------
i8 <- mpfr(-2:5, 32)
x4 <- mpfr(c(NA, NaN, -Inf, Inf), 32); x4 # NA -> NaN as well
stopifnot(identical3(is.na(x4), is.nan(x4), c(T,T,F,F)))
o1 <- as(x4[1], "mpfr1")
stopifnot(is(o1, "mpfr1")) # failed previously
validObject(o1) # ditto (failed on 64-bit only)
stopifnot(
getPrec("0xabc", base=16, doNumeric=FALSE) == 3*4,
getPrec( "abc", base=16, doNumeric=FALSE) == 3*4,
getPrec("0b1001", base=2, doNumeric=FALSE) == 4,
getPrec( "1001", base=2, doNumeric=FALSE) == 4,
identical3(mpfr("0b101", base= 2),
mpfr( "101", base= 2), mpfr(5, precBits = 3))
,
identical3(mpfr("0xabc", base=16),
mpfr( "abc", base=16), mpfr(2748, base=16, precBits = 12))
)
## save initial (Emin, Emax) eranges :
erangesOrig <- .mpfr_erange()
###----- _2_ Debugging, changing MPFR defaults, .. -----------------------------
## NB: Currently mostly *not* documented, not even .mpfr_erange()
stopifnot(Rmpfr:::.mpfr_debug() == 0 # the default level
## Activate debugging level 1:
, Rmpfr:::.mpfr_debug(1) == 0 # the previous level
## and check it :
, Rmpfr:::.mpfr_debug() == 1 # the current level
)
r <- mpfr(7, 100)^-1000
r
## (same as without debugging)
## where as this does print info: -- notably the very large values [3..6]:
.eranges <- function() sapply(.mpfr_erange_kinds, .mpfr_erange, USE.NAMES=FALSE)
## now, mpfr_erange() works with a *vector* of args:
.erange2 <- function() .mpfr_erange(.mpfr_erange_kinds)
## now returning *double* - which loses some precision [ending in '04' instead of '03']:
formatC(.eranges(), format="fg")
stopifnot(identical(.eranges(), .erange2()))
.mpfr_minPrec()
.mpfr_maxPrec()# debug printing shows the long integer (on 64 bit)
## Now, level 2 :
stopifnot(Rmpfr:::.mpfr_debug(2) == 1)
r
## with quite a bit of output
if(FALSE) # on Winbuilder [2019-08-08, both 32 and 64 bit]:
.mpfr_erange_set("Emax", 1073741823)
r2 <- r^100
r2
L <- r^-100000
L3 <- L^3
str(L3, internal=TRUE)
## Class 'mpfr' [package "Rmpfr"] of length 1 and precision 100
## internally @.Data: List of 1
## $ :Formal class 'mpfr1' [package "Rmpfr"] with 4 slots
## .. ..@ prec: int 100
## .. ..@ exp : int [1:2] 842206477 0
## .. ..@ sign: int 1
## .. ..@ d : int [1:4] 268435456 761715680 1492345294 -1000766770
str(L3)
## lots of debugging output, then
## 1.00989692356e+253529412
## ^^~~~~~~~~~ 10 ^ 253'529'412 that is humongous
if(!interactive()) # not seg.faulting, but printing a *huge* line [no longer!]
show(L3)
## segmentation fault -- randomly; 2017-06: no longer see any problem, not even with
if(FALSE) ## well, not really, definitely not interactively for now
if(interactive())
for(i in 1:256) show(L3)
##
## quite platform dependent {valgrind ==> bug? even in mpfr/gmp/.. ?}
str(.mpfr2list(x4))
## slightly nicer ["uniformly not worse"] (still very similar) :
str(x4, internal=TRUE)
x4 ## "similar info" as .mpfr2list(.)
## Increase maximal exponent:
tools:::assertWarning(
.mpfr_erange_set("Emax", 5e18)) # too large {FIXME why only warning and not error ??}
.mpfr_erange("Emax") # is unchanged
if(4e18 < .mpfr_erange("max.emax")) {
.mpfr_erange_set("Emax", 4e18) # now ok:
stopifnot(.mpfr_erange("Emax") == 4e18)
}
## revert to no debugging:
stopifnot(Rmpfr:::.mpfr_debug(0) == 2)
.mpfr_maxPrec()
L / (r2^-1000)# 1.00000....448 (could be more accurate?)
stopifnot(exprs = {
all.equal(L, r2^-1000, tol= 1e-27) # why not more accurate?
all.equal(log(L), -100000 * (-1000) * log(7), tol = 1e-15)
})
## Now, our experimental "transport vehicle":
stopifnot(length(rv <- c(r, r2, L)) == 3)
str(mpfrXport(rv))
str(mpfrXport(mpfr(2, 64)^(-3:3)))
str(mpfrXport(Const("pi")* 2^(-3:3)))
## and a very large one
mil <- mpfr(1025, 111)
str(mm <- mpfrXport(xx <- mil^(2^25)))
stopifnot(all.equal(log2(xx) * 2^-25, log2(mil), tol=1e-15))
## even larger -- strictly needs extended erange:
if(.mpfr_erange("min.emin") <= -2^40) {
.mpfr_erange_set("Emin", - 2^40)
show(xe <- 2^mpfr(-seq(1,70, by=3)*8e8, 64))
## used to print wrongly {because of integer overflow in .mpfr2str()$exp},
## with some exponents large positive
stopifnot(exprs = {
! .mpfr_erange_is_int() # as 'exp's now are double
(ee <- as.numeric(sub(".*e","", formatMpfr(xe)))) < -240e6
(diff(ee) + 722471990) %in% 0:1
})
} else {
cat(sprintf(
"Cannot set 'Emin' to -2^40 (= %g), as .mpfr_erange(\"min.emin\") is larger,
namely %g.\n",
- 2^40, .mpfr_erange("min.emin")))
}
## Bill Dunlap's example (with patch about convert S_alloc bug):
## (precision increases, then decreases)
z <- c(mpfr(1,8)/19, mpfr(1,32)/19, mpfr(1,24)/19)
cbind(fz <- format(z))
stopifnot(identical(fz, rev(format(rev(z)))))
stopifnot(identical(fz, c("0.05273",
"0.052631578947",
"0.0526315793"))) # << smaller prec, again since 2019-08-09
e.xx. <- .mpfr2exp(xx)
e.z. <- .mpfr2exp(z)
## revert to original 'erange' settings (which gives integer 'exp'):
.mpfr_erange_set("Emax", erangesOrig[["Emax"]]) # typically 2^30 - 1 = 1073741823
.mpfr_erange_set("Emin", erangesOrig[["Emin"]])
e.xx <- .mpfr2exp(xx)
e.z <- .mpfr2exp(z)
stopifnot(exprs = {
.mpfr_erange_is_int()
e.xx == e.xx.
e.xx == 335591572
e.z == e.z.
e.z == -4
is.integer(e.xx) # but e.xx. is double
is.integer(e.z)
})
k1 <- mpfr( c(123, 1234, 12345, 123456), precBits=2)
(N1 <- asNumeric(k1))# 128 1024 12288 131072 -- correct
str(sk1 <- .mpfr2str(k1))
str(sk1. <- .mpfr2str(k1, maybe.full=TRUE))
str(sk1.2 <- .mpfr2str(k1, digits=2, base=2))
str(sk1.2F <- .mpfr2str(k1, maybe.full=TRUE, base=2))
stopifnot(exprs = {
identical(sk1 [1:2], list(str = c("13", "10", "12", "13"), exp = 3:6))
identical(sk1.[1:2], list(str = c("128", "1024", "12288", "131072"), exp = 3:6))
identical(sk1.2, list(str = c("10", "10", "11", "10"),
exp = c( 8L, 11L, 14L, 18L),
finite = rep(TRUE, 4), is.0 = rep(FALSE, 4)))
all.equal(sk1.2[2:4], .mpfr_formatinfo(k1), tol=0) # not identical(): int <-> double
identical(formatMpfr(k1, base=2, digits=20, drop0trailing=TRUE),
with(sk1.2, paste0(str, sapply(exp - nchar(str), strrep, x="0"))))
identical(formatMpfr(k1, base=2, digits=2, exponent.plus=FALSE),
c("1.0e7", "1.0e10", "1.1e13", "1.0e17"))
})
## MM: --> need_dig is fine but is not used in the string that is returned !!
(fk1sF <- formatMpfr(k1, scientific=FALSE)) # "the bug" --- now fixed! ==> new "Bug" in new Rmpfr ????
## was "128." "1024." "12288." "131072." , but now obeying internal precision gives
## "1.e2" "1.e3" "1.e4" "1.e5"
(fk1 <- formatMpfr(k1, digits=6))
stopifnot(exprs = {
N1 == as.numeric(fk1)
## FIXME: This should change again "1024"
identical(format(k1, digits=3), c("128.", "1020.", "1.23e+4", "1.31e+5"))
})
##
digs <- setNames(1:6, 1:6)
## Each of these are 4 x 6 matrices
ffix <- sapply(digs, function(d) format(k1, digits = d, scientific = FALSE)) ## *not* good at all ..
## ==> need a maybe.full=TRUE even here
ff <- sapply(digs, function(d) format(k1, digits = d))# sci..fic = NA -- digits=1 failing for '128'
fsci <- sapply(digs, function(d) format(k1, digits = d, scientific = TRUE)) # perfect
stopifnot(exprs = {
length(dd <- dim(ff)) == 2
identical(dd, dim(ffix))
identical(dd, dim(fsci))
all.equal(asNumeric(fsci), asNumeric(ffix) -> dmat, tol=0)
all.equal(asNumeric(ff), asNumeric(ffix), tol=0)
})
rE <- 1 - dmat / asNumeric(k1)
i <- 1:5
summary(fm <- lm(log10(colMeans(abs(rE)))[i] ~ i))
stopifnot(exprs = {
rE[ cbind(FALSE, upper.tri(rE)[,-6]) ] == 0
abs(residuals(fm)) < 0.15
})
## formatting / printing :
tenth <- mpfr(-12:12, 52)/10
cents <- mpfr(-11:11, 64)/100
(kxi <- sort(c(k1, x4, i8, tenth, cents), na.last=FALSE))
mstr <- .mpfr2str (kxi)
mfi <- .mpfr_formatinfo(kxi)
es <- mstr$exp # base 10 ; with '0' when !is.finite or is0
ef <- mfi $exp # base 2 ; "undefined" when !is.finite or is0
j2 <- c("finite", "is.0")
dxi <- cbind(x = asNumeric(kxi), prec = .getPrec(kxi),
as.data.frame(mstr, stringsAsFactors = FALSE))
stopifnot(is.data.frame(dxi), identical(mstr$str, dxi[,"str"]),
identical(mstr[j2], mfi[j2]),
identical(ef, .mpfr2exp(kxi)))
dxi ## 2019-08-09: again *varying* size of 'str' rather than only growing !!
## Show that *order* no longer matters:
n <- length(ixk <- rev(kxi))
dix <- cbind(x = asNumeric(ixk), prec = .getPrec(ixk),
as.data.frame(.mpfr2str(ixk), stringsAsFactors = FALSE))[n:1,]
attr(dix, "row.names") <- .set_row_names(n)
stopifnot(identical(dxi, dix))
## somewhat (but not so much) revealing :
cbind(prec = .getPrec(kxi), kxi = asNumeric(kxi), str = es,
fi.10 = ceiling(ef/log2(10)), str.2 = as.integer(es*log2(10)), fi = ef)
## Bug example from RMH 2018-03-16 :
(x <- mpfr(c(65, 650, 6500, 65000, 650000), precBits=6))
data.frame(fDec = formatDec(x), f = formatMpfr(x))
x. <- as.numeric(xDec <- formatDec(x))
stopifnot(abs(x - x.) <= c(0, 0, 2, 12, 360))
cat("Checking compatibility .mpfr_formatinfo() <--> .mpfr2str(*, base=2) :\n")
for(nm in ls())
if(is(OO <- get(nm), "mpfr")) {
cat(nm,": str(*) :\n"); str(OO); cat("compatibility: ")
I <- .mpfr_formatinfo(OO)
S <- .mpfr2str(OO, base = 2L)
if(identical(I, S[-1]))
cat("[Ok]\n")
else {
if(any(B <- !I$finite)) I$exp[B] <- S$exp[B]
if(any(B <- I $ is.0)) I$exp[B] <- S$exp[B]
if(identical(I, S[-1]))
cat(" after fixup [Ok]\n")
else
stop(".mpfr_formatinfo(*) and .mpfr2str(*, base=2) do not match")
}
}
Any scripts or data that you put into this service are public.
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.