support linear algebra on ArrayPartition

Author: ChrisRackauckas

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"

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

src/array_partition.jl

@@ -293,3 +293,9 @@ common_number(a, b) =
     (b == 0 ? a :
      (a == b ? a :
       throw(DimensionMismatch("number of partitions must be equal"))))
+
+## Linear Algebra
+
+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))

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

test/upstream.jl

@@ -0,0 +1,29 @@
+using OrdinaryDiffEq, NLsolve, RecursiveArrayTools, Test
+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 u0 .* u0' .* false isa Matrix
+tspan = (0.0,100.0)
+prob = ODEProblem(lorenz,u0,tspan)
+sol = solve(prob,Tsit5());
+sol = solve(prob,AutoTsit5(Rosenbrock23(autodiff=false)));
+
+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]))