JuliaAlgebra / DynamicPolynomials.jl

Multivariate polynomials implementation of commutative and non-commutative variables
Other
60 stars 21 forks source link

Multiplication of Polynomial{C1,Polynomial{C1,T}} #69

Open David-Berghaus opened 3 years ago

David-Berghaus commented 3 years ago

Since you have defined the coefficient vectors of polynomials as arbitrary types, I have been playing around with "polynomials over polynomials". Although this type of usage was probably not intended by you, many functionalities of this are actually already supported:

julia> @ncpolyvar x0 x1 y0 y1
(x0, x1, y0, y1)

julia> tmp = Polynomial([y0*y1+y1],[x0*x1])
(y0*y1 + y1)x0x1

julia> tmp+tmp
(2*y0*y1 + 2*y1)x0x1

julia> tmp*tmp
ERROR: InexactError: convert(Term{false,Polynomial{false,Int64}}, y0*y1*y0*y1*x0*x1*x0*x1 + y0*y1^2*x0*x1*x0*x1 + y1*y0*y1*x0*x1*x0*x1 + y1^2*x0*x1*x0*x1)

As you can see however, multiplication of these objects is not possible at this point. I could obviously write a custom function to deal with this concrete scenario, however I wonder if you would be interested to add this functionality to your library :)

Best and thanks in advance, David

blegat commented 3 years ago

This seem to resolve one design decision I had to make and I apparently did wrong ^^ Given a coefficient a and a monomial m, how do we build a term with a and m ? I thought that as it can be done with a * m, there is no need for term(a, m). However, in the case where a is a polynomial, the two are different as a * m would give a polynomial in x0, x1, y0, y1 while term(a, m) gives a term with coefficient a.

Would you like to give it a try ?

David-Berghaus commented 3 years ago

Sorry for the delay, GitHub somehow did not send me a notification email :) So what I currently did is to treat the specific case of T being a polynomial:

Base.:*(t1::Term{C1,T1}, t2::Term{C2,T2}) where {C1, C2, T1<:Union{AbstractPolynomial,AbstractTerm}, T2<:Union{AbstractPolynomial,AbstractTerm}} = Term(t1.α*t2.α,t1.x*t2.x)
Base.:*(p1::Polynomial{C1, T1}, t2::Term{C2, T2}) where {C1, C2, T1<:Union{AbstractPolynomial,AbstractTerm}, T2<:Union{AbstractPolynomial,AbstractTerm}} = Polynomial(p1.a*t2.α,p1.x*t2.x)
Base.:*(t2::Term{C2, T2}, p1::Polynomial{C1, T1}) where {C1, C2, T1<:Union{AbstractPolynomial,AbstractTerm}, T2<:Union{AbstractPolynomial,AbstractTerm}} = Polynomial(t2.α*p1.a,t2.x*p1.x)
function Base.:*(p1::Polynomial{C1, T1}, p2::Polynomial{C2, T2}) where {C1, C2, T1<:Union{AbstractPolynomial,AbstractTerm}, T2<:Union{AbstractPolynomial,AbstractTerm}}
    res = 0
    @inbounds for i = 1:length(p1.a)
        @inbounds for j = 1:length(p2.a)
            res += Term(p1.a[i]*p2.a[j],p1.x[i]*p2.x[j])
        end
    end
    return res
end

This workaround is probably not the most efficient way though...

blegat commented 3 years ago

It should work with https://github.com/JuliaAlgebra/DynamicPolynomials.jl/pull/75 and https://github.com/JuliaAlgebra/MultivariatePolynomials.jl/pull/149