First of all 1e-10
is simply 1*10^-10
which is 0.0000000001
, so numericaly very close to 0
(as it is a result of square rooting, so the actual error in the computation is of row of magnitude 1e-20
). Is it "too big"? Well, library is written in fortran, and is focused on speed, so it is quite acceptable. If you analyze the exact code, you will find out how it is computed:
# fields, Tools for spatial data
# Copyright 2004-2011, Institute for Mathematics Applied Geosciences
# University Corporation for Atmospheric Research
# Licensed under the GPL -- www.gpl.org/licenses/gpl.html
"rdist" <- function(x1, x2) {
if (!is.matrix(x1))
x1 <- as.matrix(x1)
if (missing(x2))
x2 <- x1
if (!is.matrix(x2))
x2 <- as.matrix(x2)
d <- ncol(x1)
n1 <- nrow(x1)
n2 <- nrow(x2)
par <- c(1/2, 0)
temp <- .Fortran("radbas", nd = as.integer(d), x1 = as.double(x1),
n1 = as.integer(n1), x2 = as.double(x2), n2 = as.integer(n2),
par = as.double(par), k = as.double(rep(0, n1 * n2)))$k
return(matrix(temp, ncol = n2, nrow = n1))
}
And the exact answer is hidden in the fortran files (in radfun.f
called from radbas.f
), where you can find the line
if( dtemp.lt.1e-20) dtemp =1e-20
which treats small (even 0) values as 1e-20
, which after taking square root results in 1e-10
. It seems that the motivation was to speed up the process by using logarithm of the value (as a result, square rooting is simply dividing by 2), which of course is not defined for 0
.