Logo Questions Linux Laravel Mysql Ubuntu Git Menu
 

Double precision (64-bit) representation of numeric value in R (sign, exponent, significand)

R FAQ states that:

The only numbers that can be represented exactly in R’s numeric type are integers and fractions whose denominator is a power of 2. All other numbers are internally rounded to (typically) 53 binary digits accuracy.

R uses IEEE 754 double-precision floating-point numbers which is

  • 1 bit for sign
  • 11 bits for exponent
  • 52 bits for mantissa (or significand)

which sums up to 64-bits.

For the numeric number 0.1, R represents

sprintf("%.60f", 0.1)
[1] "0.100000000000000005551115123125782702118158340454101562500000"

Double (IEEE754 Double precision 64-bit) gives us this binary representation for 0.1 :

00111111 10111001 10011001 10011001 10011001 10011001 10011001 10011010

How we can get this representation in R and how does it relate to the output given by sprintf in our example?

like image 885
989 Avatar asked May 07 '18 15:05

989


People also ask

What is significand and exponent?

Generally speaking, the significand can be thought of as an integer, as a fraction, or as some other fixed-point form by choosing an appropriate exponent offset (that is, an appropriate bias). What's more, a decimal or subnormal binary significand may also contain leading zeros.

What is the number of exponent bits for double precision number?

It is commonly known simply as double. The IEEE 754 standard specifies a binary64 as having: Sign bit: 1 bit. Exponent: 11 bits.

What is the significand in binary?

sign. fraction (aka significand) — the valuable digits (the meaning, the payload) of the number. exponent — controls how far and in which direction to move the decimal point in the fraction.

Is 64 bit double precision?

The XDR standard defines the encoding for the double-precision floating-point data type as a double. The length of a double is 64 bits or 8 bytes. Doubles are encoded using the IEEE standard for normalized double-precision floating-point numbers.


3 Answers

The answer to the question raised by @chux in the comments is "yes"; R supports the %a format:

sprintf("%a", 0.1)
#> [1] "0x1.999999999999ap-4"

If you want to access the underlying bit pattern, you will have to reinterpret the double as a 64bit integer. For this task one can use C++ via Rcpp:

Rcpp::cppFunction('void print_hex(double x) {
    uint64_t y;
    static_assert(sizeof x == sizeof y, "Size does not match!");
    std::memcpy(&y, &x, sizeof y);
    Rcpp::Rcout << std::hex << y << std::endl;
}', plugins = "cpp11", includes = "#include <cstdint>")
print_hex(0.1)
#> 3fb999999999999a

This hexadecimal representation is identical to your binary representation. How does one get to the decimal representation?

  • The first bit is zero, hence the sign is positive
  • The exponent is 0x3fb, i.e. 1019 in decimal. Given the exponent bias this corresponds to an actual exponent of -4.
  • The mantissa is 0x1999999999999a × 2^-52 including the implicit 1, i.e. 2^−52 × 7,205,759,403,792,794.
  • In total this gives 2^−56 × 7,205,759,403,792,794:

    sprintf("%.60f", 2^-56 * 7205759403792794)
    #> [1] "0.100000000000000005551115123125782702118158340454101562500000"
    
like image 102
Ralf Stubner Avatar answered Oct 19 '22 23:10

Ralf Stubner


Take for example 0.3 into account. Run in R console

> sprintf("%a", 0.3)
[1] "0x1.3333333333333p-2"

Mantissa or Significand

The hex representation 3333333333333 to binary would give us the mantissa (or significand) part. That is

0011001100110011001100110011001100110011001100110011

Exponent

The exponent part (11 bits) should be the offset from 2^(11-1) - 1 = 1023 so as the trailing 3 is p-2 (in the output given by sprintf) we have

-2 + 1023 = 1021

and its binary representation fixed in 11 bits is

01111111101

Sign

As for the sign bit, its 0 for positive and 1 otherwise

Double Precision Representation

So the complete representation is

0 | 01111111101 | 0011001100110011001100110011001100110011001100110011

Another example:

> sprintf("%a", -2.94)
[1] "-0x1.7851eb851eb85p+1"

# Mantissa or Significand
(7851eb851eb85) # base 16 
(0111100001010001111010111000010100011110101110000101) # base 2

# Exponent
1 + 1023 = 1024 # base 10
10000000000 # base 2

# So the complete representation is
1 | 10000000000 | 0111100001010001111010111000010100011110101110000101
like image 20
989 Avatar answered Oct 20 '22 00:10

989


From decimal to normalized double precion:

library(BMS)

from10toNdp <- function(my10baseNumber) {
out <- list()

# Handle special cases (0, Inf, -Inf)
if (my10baseNumber %in% c(0,Inf,-Inf)) {
if (my10baseNumber==0)    { out <- "0000000000000000000000000000000000000000000000000000000000000000" }
if (my10baseNumber==Inf)  { out <- "0111111111110000000000000000000000000000000000000000000000000000" }
if (my10baseNumber==-Inf) { out <- "1111111111110000000000000000000000000000000000000000000000000000" }
} else {

signBit <- 0 # assign initial value

from10to2 <- function(deciNumber) {
  binaryVector <- rep(0, 1 + floor(log(deciNumber, 2)))
  while (deciNumber >= 2) {
    theExpo <- floor(log(deciNumber, 2))
    binaryVector[1 + theExpo] <- 1
    deciNumber <- deciNumber - 2^theExpo  }
  binaryVector[1] <- deciNumber %% 2
  paste(rev(binaryVector), collapse = "")}

#Sign bit
if (my10baseNumber<0) { signBit <- 1 
} else { signBit <- 0 }

# Biased Exponent
BiasedExponent <- strsplit(from10to2(as.numeric(substr(sprintf("%a", my10baseNumber), which(strsplit( sprintf("%a", my10baseNumber), "")[[1]]=="p")+1, length( strsplit( sprintf("%a", my10baseNumber), "")[[1]]))) + 1023), "")[[1]] 
BiasedExponent <- paste(BiasedExponent, collapse='')
if (nchar(BiasedExponent)<11) {BiasedExponent <-  paste(c(  rep(0,11-nchar(BiasedExponent)), BiasedExponent),collapse='')    }

# Significand
significand <- BMS::hex2bin(substr( sprintf("%a", my10baseNumber) , which(strsplit( sprintf("%a", my10baseNumber), "")[[1]]=="x")+3, which(strsplit( sprintf("%a", my10baseNumber), "")[[1]]=="p")-1))

significand <- paste(significand, collapse='')
if (nchar(significand)<52) {significand <-  paste(c( significand,rep(0,52-nchar(significand))),collapse='')    }

out <- paste(c(signBit, BiasedExponent, significand), collapse='')
}

out
}

Hence,

from10toNdp(0.1)
# "0011111110111001100110011001100110011001100110011001100110011010"
like image 24
Erdogan CEVHER Avatar answered Oct 20 '22 00:10

Erdogan CEVHER