undefined | Better HN

story

0 pointsChrisRackauckas3y ago0 comments

Ahh then it's all based on the speed of the objective function and the gradient calculation code. Do you have an example to look at for that? Or did you try modelingtoolkitizing and running an auto-simplified form?

0 comments

celrod3y ago

I've noticed that R often defaults to much higher tolerances than Julia, even when it's wrappers to the same C library, like cubature.

R cubature[0]: 1e-5 Cubature.jl [1]: 1e-8

The difference for NLopt in R vs Julia is smaller. `NLopt.DEFAULT_OPTIONS`[2] in Julia shows `1e-7` for `ftol_rel`, `xtol_rel`, and `constrtol_abs`, while in R `xtol_rel` is `1e-6` and the others are `0.0`[3]. So, the options at least aren't the same with nlopt. Anyway, I always recommend confirming that you're comparing the same settings.

And of course, in Julia, you'll probably want to `JET.report_opt` your function and fix and glaring performance issues.

NLopt seems like it may be a bit of an exception, but I noticed this is pretty common pattern elsewhere, uniroot[4] being another example, with eps()^(1/4) default tolerance, far higher than Julia root solvers will use.

[0] https://cran.r-project.org/web/packages/cubature/cubature.pd... [1] https://github.com/JuliaMath/Cubature.jl [2] https://github.com/JuliaOpt/NLopt.jl/blob/6ade25740362895bbf... [3] https://cran.r-project.org/web/packages/nloptr/nloptr.pdf [4] https://www.rdocumentation.org/packages/stats/versions/3.6.2...

BrandonS1133y ago

We take of that by setting all the NLOPT options to be the same across calling languages. We get pretty much the same number of calls in julia, r and c.

BrandonS1133y ago

Yes, well no gradients. It a simple parametric function (6 parameters) applied to a vector of lengths 10 to 20. All vectorizes in the language of R. All powers, logs, exponentials, ratios, and sums. Large number of local maxima, and places where function cannot be evaluated. So probabilistic optimizer is very important.

Come to think of it, for this sort of calculations, R and Julia should take the same time.

DNF23y ago

Of course, this is piques one's curiosity. It might be, if the function is simple enough, that there is little advantage to Julia here.

But if you are combining multiple operations on a vector, there could be opportunities for Julia, in-place operations, fusing, simd. Maybe even StaticArrays.

Any chance of sharing that little piece of code?

BrandonS1133y ago

sure, for parameters par (what is optimized for), data vector x (typical length from 10 to 20), constants n and n2, a typical function is

  if((1 - par[3]^2)<0 ) return(100)
  if(par[1] + par[2] \* par[5] \* sqrt(1 - par[3]^2)<0) 
  return(500 )
  tmp <- (x - par[4])^2 + par[5]^2
  tmp2 <- par[1] + par[2] \* (par[3] \* (x - par[4]) + 
  sqrt(tmp))
  result <- ((sqrt(tmp2) - Ce) / n2)^2
  result <- sqrt(sum(n \* result) / n)
  if(is.na(result)) return(1111)
  return(result)

DNF23y ago

I'm not that good at reading R. But if the Julia code is similar, then this code is type unstable, sometimes returning an Int, sometimes a Float. That harms performance.

Generally, it looks like a function where Julia could have a significant performance advantage.

1 more reply

j / k navigate · click thread line to collapse

0 comments

celrod3y ago

I've noticed that R often defaults to much higher tolerances than Julia, even when it's wrappers to the same C library, like cubature.

R cubature[0]: 1e-5 Cubature.jl [1]: 1e-8

And of course, in Julia, you'll probably want to `JET.report_opt` your function and fix and glaring performance issues.

BrandonS1133y ago

We take of that by setting all the NLOPT options to be the same across calling languages. We get pretty much the same number of calls in julia, r and c.

BrandonS1133y ago

Come to think of it, for this sort of calculations, R and Julia should take the same time.

DNF23y ago

Of course, this is piques one's curiosity. It might be, if the function is simple enough, that there is little advantage to Julia here.

But if you are combining multiple operations on a vector, there could be opportunities for Julia, in-place operations, fusing, simd. Maybe even StaticArrays.

Any chance of sharing that little piece of code?

BrandonS1133y ago

sure, for parameters par (what is optimized for), data vector x (typical length from 10 to 20), constants n and n2, a typical function is

  if((1 - par[3]^2)<0 ) return(100)
  if(par[1] + par[2] \* par[5] \* sqrt(1 - par[3]^2)<0) 
  return(500 )
  tmp <- (x - par[4])^2 + par[5]^2
  tmp2 <- par[1] + par[2] \* (par[3] \* (x - par[4]) + 
  sqrt(tmp))
  result <- ((sqrt(tmp2) - Ce) / n2)^2
  result <- sqrt(sum(n \* result) / n)
  if(is.na(result)) return(1111)
  return(result)

DNF23y ago

I'm not that good at reading R. But if the Julia code is similar, then this code is type unstable, sometimes returning an Int, sometimes a Float. That harms performance.

Generally, it looks like a function where Julia could have a significant performance advantage.

1 more reply

j / k navigate · click thread line to collapse