diff --git a/.travis.yml b/.travis.yml index 2e762b8f..084ebb99 100644 --- a/.travis.yml +++ b/.travis.yml @@ -3,11 +3,11 @@ language: julia os: - linux julia: - - 1.0 - - nightly -matrix: - allow_failures: - - julia: nightly + - 1 +# - nightly +#matrix: +# allow_failures: +# - julia: nightly notifications: email: false # uncomment the following lines to override the default test script diff --git a/Project.toml b/Project.toml index 382dedf5..0e33e5ac 100644 --- a/Project.toml +++ b/Project.toml @@ -5,6 +5,9 @@ version = "2.2.0" [deps] ArrayInterface = "4fba245c-0d91-5ea0-9b3e-6abc04ee57a9" +LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e" +NLsolve = "2774e3e8-f4cf-5e23-947b-6d7e65073b56" +OrdinaryDiffEq = "1dea7af3-3e70-54e6-95c3-0bf5283fa5ed" RecipesBase = "3cdcf5f2-1ef4-517c-9805-6587b60abb01" Requires = "ae029012-a4dd-5104-9daa-d747884805df" StaticArrays = "90137ffa-7385-5640-81b9-e52037218182" @@ -12,12 +15,12 @@ Statistics = "10745b16-79ce-11e8-11f9-7d13ad32a3b2" ZygoteRules = "700de1a5-db45-46bc-99cf-38207098b444" [compat] -ArrayInterface = "1.2, 2.0" +ArrayInterface = "2.7" RecipesBase = "0.7, 0.8, 1.0" Requires = "0.5, 1.0" StaticArrays = "0.10, 0.11, 0.12" ZygoteRules = "0.2" -julia = "1" +julia = "1.3" [extras] Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40" diff --git a/src/RecursiveArrayTools.jl b/src/RecursiveArrayTools.jl index 664c1171..e3ca5ee5 100644 --- a/src/RecursiveArrayTools.jl +++ b/src/RecursiveArrayTools.jl @@ -3,7 +3,7 @@ __precompile__() module RecursiveArrayTools using Requires, RecipesBase, StaticArrays, Statistics, - ArrayInterface, ZygoteRules + ArrayInterface, ZygoteRules, LinearAlgebra abstract type AbstractVectorOfArray{T, N, A} <: AbstractArray{T, N} end abstract type AbstractDiffEqArray{T, N, A} <: AbstractVectorOfArray{T, N, A} end diff --git a/src/array_partition.jl b/src/array_partition.jl index 63864ffa..e4ec24d3 100644 --- a/src/array_partition.jl +++ b/src/array_partition.jl @@ -233,6 +233,8 @@ ArrayPartitionStyle(::Val{N}) where N = ArrayPartitionStyle{Broadcast.DefaultArr function Broadcast.BroadcastStyle(::ArrayPartitionStyle{AStyle}, ::ArrayPartitionStyle{BStyle}) where {AStyle, BStyle} ArrayPartitionStyle(Broadcast.BroadcastStyle(AStyle(), BStyle())) end +Broadcast.BroadcastStyle(::ArrayPartitionStyle, ::Broadcast.DefaultArrayStyle{0}) = Broadcast.DefaultArrayStyle{1}() +Broadcast.BroadcastStyle(::ArrayPartitionStyle, ::Broadcast.DefaultArrayStyle{N}) where N = Broadcast.DefaultArrayStyle{N}() combine_styles(args::Tuple{}) = Broadcast.DefaultArrayStyle{0}() combine_styles(args::Tuple{Any}) = Broadcast.result_style(Broadcast.BroadcastStyle(args[1])) @@ -252,7 +254,7 @@ end ArrayPartition(f, N) end -@inline function Base.copyto!(dest::ArrayPartition, bc::Broadcast.Broadcasted) +@inline function Base.copyto!(dest::ArrayPartition, bc::Broadcast.Broadcasted{ArrayPartitionStyle{Style}}) where Style N = npartitions(dest, bc) for i in 1:N copyto!(dest.x[i], unpack(bc, i)) @@ -293,3 +295,10 @@ common_number(a, b) = (b == 0 ? a : (a == b ? a : throw(DimensionMismatch("number of partitions must be equal")))) + +## Linear Algebra + +ArrayInterface.zeromatrix(A::ArrayPartition) = ArrayInterface.zeromatrix(reduce(vcat,vec.(A.x))) +LinearAlgebra.ldiv!(A::LinearAlgebra.LU,b::ArrayPartition) = ldiv!(A,Array(b)) +LinearAlgebra.ldiv!(A::LinearAlgebra.QR,b::ArrayPartition) = ldiv!(A,Array(b)) +LinearAlgebra.ldiv!(A::LinearAlgebra.SVD,b::ArrayPartition) = ldiv!(A,Array(b)) diff --git a/test/partitions_test.jl b/test/partitions_test.jl index 0f9f3ff3..3f05c534 100644 --- a/test/partitions_test.jl +++ b/test/partitions_test.jl @@ -26,7 +26,21 @@ a = 5 @. p = p*p2 K = p.*p2 -@test_broken p.*rand(10) +x = rand(10) +y = p.*x +@test y[1:5] == p.x[1] .* x[1:5] +@test y[6:10] == p.x[2] .* x[6:10] +y = p.*x' +for i in 1:10 + @test y[1:5,i] == p.x[1] .* x[i] + @test y[6:10,i] == p.x[2] .* x[i] +end +y = p .* p' +@test y[1:5,1:5] == p.x[1] .* p.x[1]' +@test y[6:10,6:10] == p.x[2] .* p.x[2]' +@test y[1:5,6:10] == p.x[1] .* p.x[2]' +@test y[6:10,1:5] == p.x[2] .* p.x[1]' + b = rand(10) c = rand(10) copyto!(b,p) @@ -94,7 +108,7 @@ S = [ ] for sizes in S - x = ArrayPartition( randn.(sizes[1]) ) + x = ArrayPartition( randn.(sizes[1]) ) y = ArrayPartition( zeros.(sizes[2]) ) y_array = zeros(length(x)) copyto!(y,x) #testing Base.copyto!(dest::ArrayPartition,A::ArrayPartition) @@ -102,5 +116,3 @@ for sizes in S @test all([x[i] == y[i] for i in eachindex(x)]) @test all([x[i] == y_array[i] for i in eachindex(x)]) end - - diff --git a/test/runtests.jl b/test/runtests.jl index a8d6f996..46ceadaf 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -7,4 +7,5 @@ using Test @time @testset "VecOfArr Indexing Tests" begin include("basic_indexing.jl") end @time @testset "VecOfArr Interface Tests" begin include("interface_tests.jl") end @time @testset "StaticArrays Tests" begin include("copy_static_array_test.jl") end + @time @testset "Upstream Tests" begin include("upstream.jl") end end diff --git a/test/upstream.jl b/test/upstream.jl new file mode 100644 index 00000000..88f2a2e4 --- /dev/null +++ b/test/upstream.jl @@ -0,0 +1,30 @@ +using OrdinaryDiffEq, NLsolve, RecursiveArrayTools, Test, ArrayInterface +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 +u0 = ArrayPartition([1.0,0.0],[0.0]) +@test ArrayInterface.zeromatrix(u0) isa Matrix +tspan = (0.0,100.0) +prob = ODEProblem(lorenz,u0,tspan) +sol = solve(prob,Tsit5()) +sol = solve(prob,AutoTsit5(Rosenbrock23(autodiff=false))) +sol = solve(prob,AutoTsit5(Rosenbrock23())) + +function f!(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 = f!(F, u0) + +# To show the NLsolve error that results with ArrayPartitions: +nlsolve(f!, ArrayPartition([0.1; 1.2], [0.1; 1.2]))