R help - Feb 2022 - Why does the print method fail for very small numbers?

Dear expeRts,

I'm familiar with IEEE 754. Is there an easy way to explain why even
just printing of small numbers fails?

1e-317 # 1e-317 => fine
1e-318 # 9.999987e-319 => gets tricky; seems to call print() =>
as.character() => format() => paste()
1e-318 == 9.999987e-319 # TRUE
2.48e-324 # prints 4.940656e-324 for me
2.48e-324 == 4.940656e-324 # TRUE
## Relative error as a plot
rel_error <- function(x)
    plot(abs((as.numeric(as.character(x)) - x) / x), type = "l",
         ylab = "Relative error between x and
as.numeric(as.character(x))")
rel_error(seq(0.001, 0.001 + .Machine$double.xmin, length.out = 1001)) # fine
rel_error(seq(0, .Machine$double.xmin, length.out = 1001)) # printing breaks
down

Of course, [0,.Machine$double.xmin] is somewhat of a strange set of numbers to
consider,
and I expect things like "==" to be easily fooled there, but already
the print method (?)

Thanks & cheers,
Marius

sessionInfo()
R version 4.1.2 (2021-11-01)
Platform: x86_64-apple-darwin21.2.0 (64-bit)
Running under: macOS Monterey 12.1
...

Go to https://cs.lmu.edu/~ray/demos/ieee754.html

and enter

2.48e-324

into the encoder box. You get

0000000000000001

which means you are well beyond normalized representation and are approximating
this value with a total of 1 bit of precision.

Do you really want to fault the print method for having difficulty returning
your value to you after this much rounding?

On February 17, 2022 2:40:09 PM PST, Marius Hofert <marius.hofert at
uwaterloo.ca> wrote:
>Dear expeRts,
>
>I'm familiar with IEEE 754. Is there an easy way to explain why even
>just printing of small numbers fails?
>
>1e-317 # 1e-317 => fine
>1e-318 # 9.999987e-319 => gets tricky; seems to call print() =>
as.character() => format() => paste()
>1e-318 == 9.999987e-319 # TRUE
>2.48e-324 # prints 4.940656e-324 for me
>2.48e-324 == 4.940656e-324 # TRUE
>## Relative error as a plot
>rel_error <- function(x)
>    plot(abs((as.numeric(as.character(x)) - x) / x), type = "l",
>         ylab = "Relative error between x and
as.numeric(as.character(x))")
>rel_error(seq(0.001, 0.001 + .Machine$double.xmin, length.out = 1001)) #
fine
>rel_error(seq(0, .Machine$double.xmin, length.out = 1001)) # printing breaks
down
>
>Of course, [0,.Machine$double.xmin] is somewhat of a strange set of numbers
to consider,
>and I expect things like "==" to be easily fooled there, but
already the print method (?)
>
>Thanks & cheers,
>Marius
>
>sessionInfo()
>R version 4.1.2 (2021-11-01)
>Platform: x86_64-apple-darwin21.2.0 (64-bit)
>Running under: macOS Monterey 12.1
>...
>
>
>______________________________________________
>R-help at r-project.org mailing list -- To UNSUBSCRIBE and more, see
>https://stat.ethz.ch/mailman/listinfo/r-help
>PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
>and provide commented, minimal, self-contained, reproducible code.
-- 
Sent from my phone. Please excuse my brevity.

Values this small are denormalized numbers and only have a few binary
digits of accuracy.  The smallest positive value an IEEE double
precision number can have is 2^-1074.  The second smallest IEEE double
precision number is 2^-1073, and the third smallest is 2^-1073 +
2^-1074.

The numbers 2.48e-324 and 4.940656e-324 are both read in as 2^-1074:
> 2.48e-324 == 2^-1074
[1] TRUE
> 4.940656e-324 == 2^-1074
[1] TRUE

You can see the binary representations with

example(numToBits) ## for bitC
bitC(2^-1074)
bitC(2^-1073)
bitC(2^-1073 + 2^-1074)

Look up 'denormalized numbers' if you want to learn more.

Best,

luke

On Thu, 17 Feb 2022, Marius Hofert wrote:
> Dear expeRts,
>
> I'm familiar with IEEE 754. Is there an easy way to explain why even
> just printing of small numbers fails?
>
> 1e-317 # 1e-317 => fine
> 1e-318 # 9.999987e-319 => gets tricky; seems to call print() =>
as.character() => format() => paste()
> 1e-318 == 9.999987e-319 # TRUE
> 2.48e-324 # prints 4.940656e-324 for me
> 2.48e-324 == 4.940656e-324 # TRUE
> ## Relative error as a plot
> rel_error <- function(x)
>    plot(abs((as.numeric(as.character(x)) - x) / x), type = "l",
>         ylab = "Relative error between x and
as.numeric(as.character(x))")
> rel_error(seq(0.001, 0.001 + .Machine$double.xmin, length.out = 1001)) #
fine
> rel_error(seq(0, .Machine$double.xmin, length.out = 1001)) # printing
breaks down
>
> Of course, [0,.Machine$double.xmin] is somewhat of a strange set of numbers
to consider,
> and I expect things like "==" to be easily fooled there, but
already the print method (?)
>
> Thanks & cheers,
> Marius
>
> sessionInfo()
> R version 4.1.2 (2021-11-01)
> Platform: x86_64-apple-darwin21.2.0 (64-bit)
> Running under: macOS Monterey 12.1
> ...
>
>
> ______________________________________________
> R-help at r-project.org mailing list -- To UNSUBSCRIBE and more, see
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide
http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>
-- 
Luke Tierney
Ralph E. Wareham Professor of Mathematical Sciences
University of Iowa                  Phone:             319-335-3386
Department of Statistics and        Fax:               319-335-3017
    Actuarial Science
241 Schaeffer Hall                  email:   luke-tierney at uiowa.edu
Iowa City, IA 52242                 WWW:  http://www.stat.uiowa.edu

It should also be noted that format(x, digits = 17) instead of
as.character(x) won't lose any accuracy.

On Thu, Feb 17, 2022, 17:41 Marius Hofert <marius.hofert at uwaterloo.ca>
wrote:
> Dear expeRts,
>
> I'm familiar with IEEE 754. Is there an easy way to explain why even
> just printing of small numbers fails?
>
> 1e-317 # 1e-317 => fine
> 1e-318 # 9.999987e-319 => gets tricky; seems to call print() =>
> as.character() => format() => paste()
> 1e-318 == 9.999987e-319 # TRUE
> 2.48e-324 # prints 4.940656e-324 for me
> 2.48e-324 == 4.940656e-324 # TRUE
> ## Relative error as a plot
> rel_error <- function(x)
>     plot(abs((as.numeric(as.character(x)) - x) / x), type = "l",
>          ylab = "Relative error between x and
as.numeric(as.character(x))")
> rel_error(seq(0.001, 0.001 + .Machine$double.xmin, length.out = 1001)) #
> fine
> rel_error(seq(0, .Machine$double.xmin, length.out = 1001)) # printing
> breaks down
>
> Of course, [0,.Machine$double.xmin] is somewhat of a strange set of
> numbers to consider,
> and I expect things like "==" to be easily fooled there, but
already the
> print method (?)
>
> Thanks & cheers,
> Marius
>
> sessionInfo()
> R version 4.1.2 (2021-11-01)
> Platform: x86_64-apple-darwin21.2.0 (64-bit)
> Running under: macOS Monterey 12.1
> ...
>
>
> ______________________________________________
> R-help at r-project.org mailing list -- To UNSUBSCRIBE and more, see
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide
> http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>
	[[alternative HTML version deleted]]

You say that you already understand IEEE floating point
arithmetic, but you have placed the problem at the wrong
end.  The problem is at the INPUT end.  The smallest
number that can be represented in IEEE double format
is 4.940656458412465d-324
but you asked for 2.48e-324.
The two representable numbers closest to that are
+0.0 and 4.94...e-324.  Of those, 4.94...e-324 is
a little bit closer.  So 4.94...e-324 was THE BEST
ANSWER YOU COULD POSSIBLY GET.

The number having been rounded up to the nearest
representable number, that representable number was
printed pretty darned well.  So:
- you have NOT found an error in printing
- you have NOT found an error in comparison
- you HAVE found that rounding in the range of
  denormalised numbers is pretty harsh.

One more time: the print function DID NOT FAIL
for tiny numbers.  It printed what you gave it
and printed it well.  The problem was at the
INPUT end.

On Fri, 18 Feb 2022 at 11:41, Marius Hofert <marius.hofert at
uwaterloo.ca>
wrote:
> Dear expeRts,
>
> I'm familiar with IEEE 754. Is there an easy way to explain why even
> just printing of small numbers fails?
>
> 1e-317 # 1e-317 => fine
> 1e-318 # 9.999987e-319 => gets tricky; seems to call print() =>
> as.character() => format() => paste()
> 1e-318 == 9.999987e-319 # TRUE
> 2.48e-324 # prints 4.940656e-324 for me
> 2.48e-324 == 4.940656e-324 # TRUE
> ## Relative error as a plot
> rel_error <- function(x)
>     plot(abs((as.numeric(as.character(x)) - x) / x), type = "l",
>          ylab = "Relative error between x and
as.numeric(as.character(x))")
> rel_error(seq(0.001, 0.001 + .Machine$double.xmin, length.out = 1001)) #
> fine
> rel_error(seq(0, .Machine$double.xmin, length.out = 1001)) # printing
> breaks down
>
> Of course, [0,.Machine$double.xmin] is somewhat of a strange set of
> numbers to consider,
> and I expect things like "==" to be easily fooled there, but
already the
> print method (?)
>
> Thanks & cheers,
> Marius
>
> sessionInfo()
> R version 4.1.2 (2021-11-01)
> Platform: x86_64-apple-darwin21.2.0 (64-bit)
> Running under: macOS Monterey 12.1
> ...
>
>
> ______________________________________________
> R-help at r-project.org mailing list -- To UNSUBSCRIBE and more, see
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide
> http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>
	[[alternative HTML version deleted]]