Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Some DifferentialEquations.jl solvers fail when using ArrayPartition #44

Closed
Sleort opened this issue Oct 25, 2018 · 4 comments · Fixed by #88
Closed

Some DifferentialEquations.jl solvers fail when using ArrayPartition #44

Sleort opened this issue Oct 25, 2018 · 4 comments · Fixed by #88

Comments

@Sleort
Copy link

Sleort commented Oct 25, 2018

Note the use of ArrayPartition below

julia> using DifferentialEquations,RecursiveArrayTools

julia> function lorenz(du,u,p,t)
           du[1] = 10.0*(u[2]-u[1])
           du[2] = u[1]*(28.0-u[3]) - u[2]
           du[3] = u[1]*u[2] - (8/3)*u[3]
       end
lorenz (generic function with 1 method)

julia> u0 = ArrayPartition([1.0,0.0],[0.0])
([1.0, 0.0], [0.0])

julia> tspan = (0.0,100.0)
(0.0, 100.0)

julia> prob = ODEProblem(lorenz,u0,tspan)
ODEProblem with uType ArrayPartition{Float64,Tuple{Array{Float64,1},Array{Float64,1}}} and tType Float64. In-place: true
timespan: (0.0, 100.0)
u0: [1.0, 0.0][0.0]

julia> sol = solve(prob,Tsit5());

julia> sol = solve(prob,AutoTsit5(Rosenbrock23()));
ERROR: MethodError: no method matching ldiv!(::LinearAlgebra.LU{Float64,Array{Float64,2}}, ::ArrayPartition{Float64,Tuple{Array{Float64,1},Array{Float64,1}}})
Closest candidates are:
  ldiv!(::LinearAlgebra.Transpose{#s561,#s560} where #s560<:LinearAlgebra.LowerTriangular where #s561, ::AbstractArray{T,1} where T) at /buildworker/worker/package_linux64/build/usr/share/julia/stdlib/v1.0/LinearAlgebra/src/triangular.jl:1217
  ldiv!(::LinearAlgebra.Transpose{#s561,#s560} where #s560<:LinearAlgebra.LowerTriangular where #s561, ::AbstractArray{T,1} where T, ::AbstractArray{T,1} where T) at /buildworker/worker/package_linux64/build/usr/share/julia/stdlib/v1.0/LinearAlgebra/src/triangular.jl:1201
  ldiv!(::LinearAlgebra.Transpose{#s561,#s560} where #s560<:LinearAlgebra.UnitLowerTriangular where #s561, ::AbstractArray{T,1} where T) at /buildworker/worker/package_linux64/build/usr/share/julia/stdlib/v1.0/LinearAlgebra/src/triangular.jl:1235
  ...

(Everything works fine without ArrayPartition, i.e. u0 = [1.0, 0.0, 0.0])

I'm using DifferentialEquations v5.3.1, RecursiveArrayTools v0.18.3, and Julia 1.0.1 on Ubuntu 16.04.

@YingboMa
Copy link
Member

There is no linear algebra defined on ArrayPartition. We should fix that. Thanks for reporting!

@vnvasquez
Copy link

vnvasquez commented Apr 7, 2020

Adding another MWE to this discussion showing the current incompatibility of ArrayPartitions with NLsolve:

using NLsolve
using RecursiveArrayTools
function mymodel(F, vars)
    x = vars.x[1]    
    F.x[1][1] = (x[1]+3)*(x[2]^3-7)+18
    F.x[1][2] = sin(x[2]*exp(x[1])-1)
    y=vars.x[2]
    F.x[2][1] = (y[1]+3)*(y[2]^3-7)+18
    F.x[2][2] = sin(y[2]*exp(y[1])-1)
end

# To show that the function works
F = ArrayPartition([0.0 0.0],[0.0, 0.0])
u0= ArrayPartition([0.1; 1.2], [0.1; 1.2])
result = mymodel(F, u0) 

# To show the NLsolve error that results with ArrayPartitions: 
nlsolve(f!, ArrayPartition([0.1; 1.2], [0.1; 1.2]))

@ChrisRackauckas
Copy link
Member

@YingboMa what do you think about just defining a slow fallback for now? I don't think it would actually be a major performance hit since the factorization is just a normal array factorization, we just need to define ldiv! to convert the ArrayPartition to an Array and then do the quick step?

@ChrisRackauckas
Copy link
Member

These examples a lot more are all fixed. Go ham.

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

Successfully merging a pull request may close this issue.

4 participants