Extended floating point precision in R with Rmpfr

I learnt from a recent post on John Cook’s excellent blog that it’s really easy to do extended floating point computations in R using the Rmpfr package. Rmpfr is R’s wrapper around the C library MPFR, which stands for “Multiple Precision Floating-point Reliable”.

The main function that users will interact with is the mpfr function: it converts numeric values into (typically) high-precision numbers, which can then be used for computation. The function’s first argument is the numeric value(s) to be converted, and the second argument, precBits, represents the maximal precision to be used in numbers of bits. For example, precBits = 53 corresponds to double precision.

In his blog post, Cook gives an example of computing \pi to 100 decimal places by multiplying the arctangent of 1 by 4 (recall that \tan (\pi / 4) = 1, so \text{arctan}(1) = \pi / 4):

4 * atan(mpfr(1, 333))
# 1 'mpfr' number of precision  333   bits 
# [1] 3.14159265358979323846264338327950288419716939937510582097494459230781640628620899862803482534211706807

Why does he set the precision to 333 bits? This link suggests that with b bits, we get d = \log_{10}(2^b) \approx 0.3010b decimal digits of precision. (Reality for floating point numbers is not quite as straightforward as that: see this for a discussion. But for our purposes, this approximation will do.) Hence, to get 100 decimal places, we need around b = 100 / 0.3010 \approx 332.2 bits, so he rounds it up to 333 bits.

The first argument to mpfr can be a vector as well:

mpfr(1:10, 5)
# 10 'mpfr' numbers of precision  5   bits 
# [1]  1  2  3  4  5  6  7  8  9 10

As the next code snippet shows, R does NOT consider the output of a call to mpfr a numeric variable.

x <- sin(mpfr(1, 100))
x
# 1 'mpfr' number of precision  100   bits 
# [1] 0.84147098480789650665250232163005
is.numeric(x)
# [1] FALSE

We can use the asNumeric function to convert it to a numeric:

y <- asNumeric(x)
y
# [1] 0.841471
is.numeric(y)
# [1] TRUE

Can we use the more familiar as.numeric instead? According to the function’s documentation, as.numeric coerces to both “numeric” and to a vector, whereas asNumeric() should keep dim (and other) attributes. We can see this through a small example:

x <- mpfr(matrix(1:4, nrow = 2), 10)
x
# 'mpfrMatrix' of dim(.) =  (2, 2) of precision  10   bits 
# [,1]   [,2]  
# [1,] 1.0000 3.0000
# [2,] 2.0000 4.0000
asNumeric(x)
# [,1] [,2]
# [1,]    1    3
# [2,]    2    4
as.numeric(x)
# [1] 1 2 3 4

Leave a Reply

Fill in your details below or click an icon to log in:

WordPress.com Logo

You are commenting using your WordPress.com account. Log Out /  Change )

Google photo

You are commenting using your Google account. Log Out /  Change )

Twitter picture

You are commenting using your Twitter account. Log Out /  Change )

Facebook photo

You are commenting using your Facebook account. Log Out /  Change )

Connecting to %s