Last digit of floating point numbers truncated in R -> Julia connection

[chriselrod]

Hi,

This seems to round the inputs.

library(XRJulia)
Sys.setenv(JULIA_BIN = "/home/celrod/Documents/julia-static-llvm/usr/bin/julia",
           LD_LIBRARY_PATH="/usr/local/lib64");
findJulia()

id <- JuliaFunction("identity")

id(pi) - pi
# -3.108624e-15

id(2*pi - id(pi)) - pi # was hoping to compensate
# 7.105427e-15

id(exp(1)) - exp(1)
# 4.884981e-15

id(1/3) - 1/3
# -3.330669e-16

nf <- JuliaFunction("nextfloat")
nf(pi) - pi
# -2.664535e-15
nf(nf(pi)) - pi
# -2.664535e-15
nf(nf(nf(pi))) - pi
# -2.664535e-15

I got this on two computers, both running Ubuntu 16.04 and Julia 0.6.2.

EDIT: To say something positive, I noticed this issue because of this test problem:

assess_approx <- function (ints, target){
  (ints[1] / ints[2]) %>% (function(approx) c(approx, approx-target))
}
juliaEval("function check_approx!(best_inds, i, j, target, best)
            diff = abs(i / j - target)
            @inbounds if diff < best
              best_inds[1] = Int(i)
              best_inds[2] = j
            return diff
            end
            best
          end")
juliaEval("function slow_approx(n::Int, target=pi)
            best = Inf
            best_inds = [0,0]
            for j in 1:n
              i_mid = j*target
              best = check_approx!(best_inds, floor(i_mid), j, target, best)
              best = check_approx!(best_inds, ceil(i_mid), j, target, best)
            end
            best_inds
          end")

sa_jl <- JuliaFunction("slow_approx")
slow_approx_jl <- function(x) juliaGet(sa_jl(as.integer(x)))
slow_approx_jl(1e3)
system.time(exactjl <- slow_approx_jl(1e9))
#   user  system elapsed 
#  0.006   0.000   7.239
assess_approx(exactjl,pi)
# 3.141593 0.000000
slow_approx_jl2 <- function(x, y = pi) juliaGet(sa_jl(as.integer(x),y))
slow_approx_jl2(1e3)
system.time(closejl <- slow_approx_jl2(1e9))
#   user  system elapsed 
#  0.002   0.004   7.363
assess_approx(closejl,pi)
# 3.141593e+00 -3.108624e-15
assess_approx(closejl,id(pi))
# 3.141593 0.000000

versus

library(Rcpp); library(magrittr)
"#include <Rcpp.h>
#include <math.h>
#include <limits>
#include <tuple>
using namespace Rcpp;
typedef std::tuple<double, double, int> best_tup;
  best_tup check_approx(double i, int j, double target, best_tup best){
  double diff = std::abs(i / j - target);
  return diff < std::get<0>(best) ? std::make_tuple(diff, i, j) : best;
}
// [[Rcpp::export]]
IntegerVector brutish_tup(int n, double target) {
  best_tup best = std::make_tuple(std::numeric_limits<double>::infinity(),0.0,0);
  for (int j=1; j <= n; j++ ){
    double i_mid = j*target;
    best = check_approx( std::floor(i_mid), j, target, best);
    best = check_approx( std::ceil(i_mid), j, target, best);
  }
  return IntegerVector::create((int) std::get<1>(best), std::get<2>(best));
}" %>% sourceCpp(code = .) #requires c++11
brutish_tup(500, pi)
## [1] 355 113
system.time(exact <- brutish_tup(1e9, pi))
#   user  system elapsed 
# 14.904   0.011  14.916
assess_approx(exact, pi)
# 3.141593 0.000000

On my other computer the best (for C++) relative time I got was about 17 vs 5 seconds. I don't know C++ well, so perhaps I'm doing something very wrong that is killing performance -- when I tried a version that looks more like the Julia code, C++ took closer to 50 seconds. Theoretically, C++ should not be slower, so I'll look into that more. I just bring this up to say, XRJulia seems to be working well -- in this case it's looking much better than Rcpp.

EDIT: https://github.com/johnmchambers/XRJulia/issues/12#issuecomment-374330083 "have to format and then parse data every time." It looks like the problem here is that it truncates the last digit.

julia> pi |> Float64
3.141592653589793

julia> 3.14159265358979 - pi
-3.1086244689504383e-15

julia> ℯ |> Float64
2.718281828459045

julia> 2.71828182845904 - ℯ
-4.884981308350689e-15

julia> exp(1)
2.718281828459045

julia> 1/3
0.3333333333333333

julia> 0.333333333333333 - 1/3
-3.3306690738754696e-16