Merge pull request #105 from JuliaAlgebra/st/arrays

Better support arrays of polynomials.

Author: saschatimme

src/operators.jl

@@ -27,6 +27,17 @@ for (op, fun) in [(:+, :plusconstant), (:-, :minusconstant), (:*, :multconstant)
     @eval Base.$op(p::APL, α) = $fun(p, α)
     @eval Base.$op(α, p::APL) = $fun(α, p)
 end
+
+# Special case AbstractArrays of APLs
+# We add these instead of relying on the broadcasting API since the above method definitions are very wide.
+# In particular, there is support for Matrices as coefficents. In order to avoid issues like #104 we therefore
+# explicitly define this (instead of implictly getting unexpected results).
+for op in [:+, :-, :*]
+    @eval Base.$op(p::APL, A::AbstractArray{<:APL}) = map(f -> $op(p, f), A)
+    @eval Base.$op(A::AbstractArray{<:APL}, p::APL) = map(f -> $op(f, p), A)
+end
+Base.:/(A::AbstractArray{<:APL}, p::APL) = map(f -> f / p, A)
+
 Base.isapprox(p::APL, α; kwargs...) = isapprox(promote(p, α)...; kwargs...)
 Base.isapprox(α, p::APL; kwargs...) = isapprox(promote(p, α)...; kwargs...)
 

src/substitution.jl

@@ -75,6 +75,10 @@ substitute(st::AST, q::RationalPoly, s::Substitutions) = substitute(st, q.num, s
 # subs(x, x=>x+y, y=>2) would call substitute(Subs(), x+y, y=>2)
 #substitute(st::AST, p::Union{APL, RationalPoly}, s::AbstractSubstitution...) = substitute(st, p, s)
 
+# Substitute Arrays
+function substitute(st::AST, A::AbstractArray{<:APL}, s::Substitutions)
+    map(p -> substitute(st, p, s), A)
+end
 ## Everything else
 substitute(::AST, x, s::Substitutions) = x
 # subs(x, x=>1, y=>2) would call substitute(Subs(), 1, y=>2)

test/algebra.jl

@@ -65,4 +65,14 @@
 
     @test iszero((x+x-2*x) * (x * (x^2 + y^2)))
     @test iszero((0*x) * (x*y * (x^2 + y^2)))
+
+    @testset "Scalar - Array" begin
+        Mod.@polyvar x y
+        @test x + [x^2+y y; y x*y] == [x+x^2+y x+y; x+y x+x*y]
+        @test [x^2+y y; y x*y] + x == [x^2+x+y y+x; y+x x*y+x]
+
+        @test x * [1+y y] == [x*(1+y) x*y]
+        @test [1+y y] * x == [(1+y)*x y*x]
+        @test [1+y y] / x == [(1+y)/x y/x]
+    end
 end

test/substitution.jl

@@ -63,4 +63,8 @@ import Test: @inferred
     @test p(variables(p) => valsv) == 13
     @test p(varst => valst) == 13
     @test p(varsv => valsv) == 13
+
+
+    Mod.@polyvar x y
+    @test subs([x^2 + y, x + y], x => y) == [y^2+y, 2y]
 end