1-dimensional equations are converted to scalar equations

[alina-artcibasova]

As discussed here, @dlekshmi and myself are having trouble running 1D example from the CRAN ODE vignette.

The original code with function definition through R

f <- function(u,p,t) {
  return(1.01*u)
}

u0 = 1/2
tspan <- list(0.0,1.0)
sol = diffeqr::ode.solve(f,u0,tspan)

gives the following error

Error: Error happens in Julia. MethodError: Cannot convert an object of type OrdinaryDiffEq.ODECompositeSolution{Float64,1,Array{Float64,1},Nothing,Nothing,Array{Float64,1},Array{Array{Float64,1},1},ODEProblem{Float64,Tuple{Float64,Float64},false,Nothing,ODEFunction{false,RCall.var"#11#12"{RObject{ClosSxp}},UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},Base.Iterators.Pairs{Union{},Union{},Tuple{},NamedTuple{,Tuple{}}},DiffEqBase.StandardODEProblem},CompositeAlgorithm{Tuple{Tsit5,Rosenbrock23{0,false,DefaultLinSolve,DataType}},AutoSwitch{Tsit5,Rosenbrock23{0,false,DefaultLinSolve,DataType},Rational{Int64},Int64}},OrdinaryDiffEq.CompositeInterpolationData{ODEFunction{false,RCall.var"#11#12"{RObject{ClosSxp}},UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing},Array{Float64,1},Array{Float64,1},Array{Array{Float64,1},1},OrdinaryDiffEq.CompositeCache{Tuple{OrdinaryDif

And a modified version with function definition through JuliaCall

f <- JuliaCall::julia_eval("
function f(du,u,p,t)
  du = 1.01*u
end")

tspan=list(0.0,1.0)
u0=c(0.5)
sol = diffeqr::ode.solve('f',u0,tspan)

gives

Error: Error happens in Julia. MethodError: no method matching similar(::Float64) Closest candidates are: similar(!Matched::Sundials.NVector) at C:\Users\alinaa.julia\packages\Sundials\SKP8f\src\nvector_wrapper.jl:69 similar(!Matched::Array{T,1}) where T at array.jl:375 similar(!Matched::Array{T,2}) where T at array.jl:376 ... Stacktrace: [1] alg_cache(::Tsit5, ::Float64, ::Float64, ::Type{T} where T, ::Type{T} where T, ::Type{T} where T, ::Float64, ::Float64, ::ODEFunction{true,typeof(f),UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing}, ::Float64, ::Float64, ::Float64, ::Nothing, ::Bool, ::Val{true}) at C:\Users\alinaa.julia\packages\OrdinaryDiffEq\lGGkK\src\caches\low_order_rk_caches.jl:356 [2] (::OrdinaryDiffEq.var"#172#173"{Float64,Float64,DataType,DataType,DataType,Float64,Float64,ODEFunction{true,typeof(f),UniformScaling{Bool},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,N

Issue can be circumvented by turning the function into a system of equations. I put together a more detailed R notebook with my session info/system of equations hack, and previously suggested fixes.

[ChrisRackauckas]

f <- function(u,p,t) {
  return(1.01*u)
}

u0 = 1/2
tspan <- list(0.0,1.0)
sol = diffeqr::ode.solve(f,u0,tspan)

works with the latest versions, but

f <- JuliaCall::julia_eval("
function f(du,u,p,t)
  du .= 1.01.*u
end")

does not because in R there is no difference between scalars and arrays, so the R->Julia connection converts arrays of length 1 into scalars, meaning on the Julia side it should have scalar semantics instead of array semantics. That similar error is precisely because it's trying to create a cache array that's similar to a scalar, when it should've been an array. I'm not sure how fixable this is given that R doesn't distinguish between these two semantics. Thankfully, this is just something that will occur in the 1D case.