9

I wanted to generate a random real (I guess rational) number.

To do this I wanted to use runif(1, min = m, max = M) and my thoughts were to set m; M (absolute) as large as possible in order to make the interval as large as possible. Which brings me to my question:

M <- .Machine$double.xmax
m <- -M
runif(1, m, M)
## which returns
[1] Inf

Why is it not returning a number? Is the chosen interval simply too large?

PS

> .Machine$double.xmax
[1] 1.797693e+308
niko
  • 5,253
  • 1
  • 12
  • 32
  • 6
    Interesting question. A workaround would be `M * runif(1, -1, 1)`. – mt1022 Mar 18 '18 at 11:09
  • 7
    Here is the source of `runif`: https://github.com/wch/r-source/blob/af7f52f70101960861e5d995d3a4bec010bc89e6/src/nmath/runif.c . It looks that `M` is the max double, wihle `M - (-M)` will be large than max double so `b - a` will be `Inf`, which might explain why `a + (b - a) * u` is `Inf`. – mt1022 Mar 18 '18 at 11:19

1 Answers1

1

As hinted by mt1022 the reason is in runif C source:

double runif(double a, double b)
{
    if (!R_FINITE(a) || !R_FINITE(b) || b < a)  ML_ERR_return_NAN;

    if (a == b)
    return a;
    else {
    double u;
    /* This is true of all builtin generators, but protect against
       user-supplied ones */
    do {u = unif_rand();} while (u <= 0 || u >= 1);
    return a + (b - a) * u;
    }
}

In the return argument you can see the formula a + (b - a) * u which transoform uniformly [0, 1] generated random value in the user supplied interval [a, b]. In your case it will be -M + (M + M) * u. So M + M in case it 1.79E308 + 1.79E308 generates Inf. I.e. finite + Inf * finite = Inf:

M + (M - m) * runif(1, 0, 1)
# Inf
Artem
  • 3,304
  • 3
  • 18
  • 41