`"mpfr"`

is the class of **M**ultiple **P**recision
**F**loatingpoint numbers with **R**eliable arithmetic.

For the high-level user, `"mpfr"`

objects should behave
as standard **R**'s `numeric`

*vectors*. They would just
print differently and use the prespecified (typically high) precision
instead of the double precision of ‘traditional’ **R** numbers
(with `class(.) == "numeric"`

and
`typeof(.) == "double"`

).

`hypot(x,y)`

computes the hypothenuse length *z* in a rectangular
triangle with “leg” side lengths *x* and *y*, i.e.,

*z = hypot(x,y) = sqrt(x^2 + y^2),*

in a numerically stable way.

1 |

`x,y` |
an object of class |

`rnd.mode` |
a 1-letter string specifying how |

Objects are typically created by `mpfr(<number>, precBits)`

.

Internally, `"mpfr"`

objects just contain standard **R**
`list`

s where each list element is of class
`"mpfr1"`

, representing *one* MPFR number, in a structure
with four slots, very much parallelizing the C `struc`

in the
`mpfr`

C library to which the Rmpfr package interfaces.

An object of class `"mpfr1"`

has slots

`prec`

:`"integer"`

specifying the maxmimal precision in**bits**.`exp`

:`"integer"`

specifying the base-**2**exponent of the number.`sign`

:`"integer"`

, typically`-1`

or`1`

, specifying the sign (i.e.`sign(.)`

) of the number.`d`

:an

`"integer"`

vector (of 32-bit “limbs”) which corresponds to the full mantissa of the number.

- abs
`signature(x = "mpfr")`

: ...- atan2
`signature(y = "mpfr", x = "ANY")`

, and- atan2
`signature(x = "ANY", y = "mpfr")`

: compute the arc-tangent of two arguments:`atan2(y, x)`

returns the angle between the x-axis and the vector from the origin to*(x, y)*, i.e., for positive arguments`atan2(y, x) == atan(y/x)`

.- lbeta
`signature(a = "ANY", b = "mpfrArray")`

, is*log(abs(B(a,b)))*where*B(a,b)*is the Beta function,`beta(a,b)`

.- beta
`signature(a = "mpfr", b = "ANY")`

,- beta
`signature(a = "mpfr", b = "mpfr")`

, ..., etc: Compute the beta function*B(a,b)*, using high precision, building on internal`gamma`

or`lgamma`

. See the help for**R**'s base function`beta`

for more. Currently, there,*a,b >= 0*is required. Here, we provide (non-`NaN`

) for all numeric`a, b`

.When either

*a*,*b*, or*a+b*is a negative*integer*,*Γ(.)*has a pole there and is undefined (`NaN`

). However the Beta function can be defined there as “limit”, in some cases. Following other software such as SAGE, Maple or Mathematica, we provide finite values in these cases. However, note that these are not proper limits (two-dimensional in*(a,b)*), but useful for some applications. E.g.,*B(a,b)*is defined as zero when*a+b*is a negative integer, but neither*a*nor*b*is. Further, if*a > b > 0*are integers,*B(-a,b)= B(b,-a)*can be seen as*(-1)^b * B(a-b+1,b)*.- dim<-
`signature(x = "mpfr")`

: Setting a dimension`dim`

on an`"mpfr"`

object makes it into an object of class`"mpfrArray"`

or (more specifically)`"mpfrMatrix"`

for a length-2 dimension, see their help page; note that`t(x)`

(below) is a special case of this.- Ops
`signature(e1 = "mpfr", e2 = "ANY")`

: ...- Ops
`signature(e1 = "ANY", e2 = "mpfr")`

: ...- Arith
`signature(e1 = "mpfr", e2 = "missing")`

: ...- Arith
`signature(e1 = "mpfr", e2 = "mpfr")`

: ...- Arith
`signature(e1 = "mpfr", e2 = "integer")`

: ...- Arith
`signature(e1 = "mpfr", e2 = "numeric")`

: ...- Arith
`signature(e1 = "integer", e2 = "mpfr")`

: ...- Arith
`signature(e1 = "numeric", e2 = "mpfr")`

: ...- Compare
`signature(e1 = "mpfr", e2 = "mpfr")`

: ...- Compare
`signature(e1 = "mpfr", e2 = "integer")`

: ...- Compare
`signature(e1 = "mpfr", e2 = "numeric")`

: ...- Compare
`signature(e1 = "integer", e2 = "mpfr")`

: ...- Compare
`signature(e1 = "numeric", e2 = "mpfr")`

: ...- Logic
`signature(e1 = "mpfr", e2 = "mpfr")`

: ...- Summary
`signature(x = "mpfr")`

: The S4`Summary`

group functions,`max`

,`min`

,`range`

,`prod`

,`sum`

,`any`

, and`all`

are all defined for MPFR numbers.`mean(x, trim)`

for non-0`trim`

works analogously to`mean.default`

.- median
`signature(x = "mpfr")`

: works via

- quantile
`signature(x = "mpfr")`

: a simple wrapper of the`quantile.default`

method from stats.- Math
`signature(x = "mpfr")`

: All the S4`Math`

group functions are defined, using multiple precision (MPFR) arithmetic, from`getGroupMembers("Math")`

, these are (in alphabetical order):`abs`

,`sign`

,`sqrt`

,`ceiling`

,`floor`

,`trunc`

,`cummax`

,`cummin`

,`cumprod`

,`cumsum`

,`exp`

,`expm1`

,`log`

,`log10`

,`log2`

,`log1p`

,`cos`

,`cosh`

,`sin`

,`sinh`

,`tan`

,`tanh`

,`acos`

,`acosh`

,`asin`

,`asinh`

,`atan`

,`atanh`

,`gamma`

,`lgamma`

,`digamma`

, and`trigamma`

.Currently,

`trigamma`

is not provided by the MPFR library and hence not yet implemented.

Further, the`cum*()`

methods are*not yet*implemented.- factorial
`signature(x = "mpfr")`

: this will`round`

the result when`x`

is integer valued. Note however that`factorialMpfr(n)`

for integer`n`

is slightly more efficient, using the MPFR function mpfr_fac_ui.- Math2
`signature(x = "mpfr")`

:`round(x, digits)`

and`signif(x, digits)`

methods. Note that these do not change the formal precision (`'prec'`

slot), and you may often want to apply`roundMpfr()`

in addition or preference.- as.numeric
`signature(x = "mpfr")`

: ...- as.vector
`signature(x = "mpfrArray")`

: as for standard`array`

s, this “drops” the`dim`

(and`dimnames`

), i.e., transforms`x`

into an ‘MPFR’ number vector, i.e., class`mpfr`

.- [[
`signature(x = "mpfr", i = "ANY")`

, and- [
`signature(x = "mpfr", i = "ANY", j = "missing", drop = "missing")`

: subsetting aka “indexing” happens as for numeric vectors.- format
`signature(x = "mpfr")`

, further arguments`digits = NULL, scientific = NA`

, etc: returns`character`

vector of same length as`x`

; when`digits`

is`NULL`

, with*enough*digits to recreate`x`

accurately. For details, see`formatMpfr`

.- is.finite
`signature(x = "mpfr")`

: ...- is.infinite
`signature(x = "mpfr")`

: ...- is.na
`signature(x = "mpfr")`

: ...- is.nan
`signature(x = "mpfr")`

: ...- log
`signature(x = "mpfr")`

: ...- show
`signature(object = "mpfr")`

: ...- sign
`signature(x = "mpfr")`

: ...- Re, Im
`signature(z = "mpfr")`

: simply return`z`

or`0`

(as`"mpfr"`

numbers of correct precision), as mpfr numbers are ‘real’ numbers.- Arg, Mod, Conj
`signature(z = "mpfr")`

: these are trivial for our ‘real’ mpfr numbers, but defined to work correctly when used in**R**code that also allows complex number input.- all.equal
`signature(target = "mpfr", current = "mpfr")`

,- all.equal
`signature(target = "mpfr", current = "ANY")`

, and- all.equal
`signature(target = "ANY", current = "mpfr")`

: methods for numerical (approximate) equality,`all.equal`

of multiple precision numbers. Note that the default`tolerance`

(argument) is taken to correspond to the (smaller of the two) precisions when both main arguments are of class`"mpfr"`

, and hence can be considerably less than double precision machine epsilon`.Machine$double.eps`

.- coerce
`signature(from = "numeric", to = "mpfr")`

:`as(., "mpfr")`

coercion methods are available for`character`

strings,`numeric`

,`integer`

,`logical`

, and even`raw`

. Note however, that`mpfr(., precBits, base)`

is more flexible.- coerce
`signature(from = "mpfr", to = "bigz")`

: coerces to biginteger, see`bigz`

in package gmp.- coerce
`signature(from = "mpfr", to = "numeric")`

: ...- coerce
`signature(from = "mpfr", to = "character")`

: ...- unique
`signature(x = "mpfr")`

: and- duplicated
`signature(x = "mpfr")`

: just work as with numbers.- t
`signature(x = "mpfr")`

: makes`x`

into an*n x 1*`mpfrMatrix`

.- which.min
`signature(x = "mpfr")`

: gives the index of the first minimum, see`which.min`

.- which.max
`signature(x = "mpfr")`

: gives the index of the first maximum, see`which.max`

.

Many more methods (“functions”) automagically work for
`"mpfr"`

number vectors (and matrices, see the
`mpfrMatrix`

class doc),
notably
`sort`

, `order`

, `quantile`

,
`rank`

.

Martin Maechler

The `"mpfrMatrix"`

class, which extends the
`"mpfr"`

one.

`roundMpfr`

to *change* precision of an `"mpfr"`

object which is typically desirable *instead* of or in addition
to `signif()`

or `round()`

;
`is.whole()`

etc.

Special mathematical functions such as some Bessel ones, e.g., `jn`

;
further, `zeta(.)`

*(= ζ(.))*, `Ei()`

etc.
`Bernoulli`

numbers and the Pochhammer function
`pochMpfr`

.

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 | ```
## 30 digit precision
str(x <- mpfr(c(2:3, pi), prec = 30 * log2(10)))
x^2
x[1] / x[2] # 0.66666... ~ 30 digits
## indexing - as with numeric vectors
stopifnot(identical(x[2], x[[2]]),
## indexing "outside" gives NA (well: "mpfr-NaN" for now):
is.na(x[5]),
## whereas "[[" cannot index outside:
is(try(x[[5]]), "try-error"),
## and only select *one* element:
is(try(x[[2:3]]), "try-error"))
## factorial() & lfactorial would work automagically via [l]gamma(),
## but factorial() additionally has an "mpfr" method which rounds
f200 <- factorial(mpfr(200, prec = 1500)) # need high prec.!
f200
as.numeric(log2(f200))# 1245.38 -- need precBits >~ 1246 for full precision
##--> see factorialMpfr() for more such computations.
##--- "Underflow" **much** later -- exponents have 30(+1) bits themselves:
mpfr.min.exp2 <- - (2^30 + 1)
two <- mpfr(2, 55)
stopifnot(two ^ mpfr.min.exp2 == 0)
## whereas
two ^ (mpfr.min.exp2 * (1 - 1e-15))
## 2.38256490488795107e-323228497 ["typically"]
##--- "Assert" that {sort}, {order}, {quantile}, {rank}, all work :
p <- mpfr(rpois(32, lambda=500), precBits=128)^10
np <- as.numeric(log(p))
stopifnot(all(diff(sort(p)) >= 0),
identical(order(p), order(np)),
identical(rank (p), rank (np)),
all.equal(sapply(1:9, function(Typ) quantile(np, type=Typ, names=FALSE)),
sapply(lapply(1:9, function(Typ) quantile( p, type=Typ, names=FALSE)),
function(x) as.numeric(log(x))),
tol = 1e-3),# quantiles: interpolated in orig. <--> log scale
TRUE)
m0 <- mpfr(numeric(), 99)
xy <- expand.grid(x = -2:2, y = -2:2) ; x <- xy[,"x"] ; y <- xy[,"y"]
a2. <- atan2(y,x)
stopifnot(identical(which.min(m0), integer(0)),
identical(which.max(m0), integer(0)),
all.equal(a2., atan2(as(y,"mpfr"), x)),
max(m0) == mpfr(-Inf, 53), # (53 is not a feature, but ok)
min(m0) == mpfr(+Inf, 53),
sum(m0) == 0, prod(m0) == 1)
``` |

Questions? Problems? Suggestions? Tweet to @rdrrHQ or email at ian@mutexlabs.com.

All documentation is copyright its authors; we didn't write any of that.