From cdb13fc677349ab8e23c1e6f081817a48841a65b Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Beno=C3=AEt=20Legat?= Date: Wed, 3 Jul 2024 16:01:20 +0200 Subject: [PATCH 1/2] More substitution --- src/arithmetic.jl | 10 ++++++++++ test/runtests.jl | 3 +++ 2 files changed, 13 insertions(+) diff --git a/src/arithmetic.jl b/src/arithmetic.jl index e806d82..46e7134 100644 --- a/src/arithmetic.jl +++ b/src/arithmetic.jl @@ -99,6 +99,16 @@ function (p::_AE)(args::MP.AbstractSubstitution...) return MP.substitute(MP.Eval(), p, args...) end +function (p::_AE)(x::NTuple{N,<:Number}) where {N} + return (MP.polynomial(p))(x) +end + +function (p::_AE)(x::AbstractVector{<:Number}) + return (MP.polynomial(p))(x) +end + +(p::_AE)(x::Number...) = (MP.polynomial(p))(x...) + function MP.differentiate(p::_AE, args...) return MP.differentiate(MP.polynomial(p), args...) end diff --git a/test/runtests.jl b/test/runtests.jl index 763c2c0..067c35f 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -102,6 +102,9 @@ function api_test(B::Type{<:MB.AbstractMonomialIndexed}, degree) @test typeof(MB.algebra_element(sum(x))) == MA.promote_operation(MB.algebra_element, typeof(sum(x))) @test const_alg_el(x => ones(length(x))) == const_mono + @test const_alg_el(ones(length(x))) == const_mono + @test const_alg_el(ones(length(x))...) == const_mono + @test const_alg_el(tuple(ones(length(x))...)) == const_mono @test subs(const_alg_el, x => ones(length(x))) == const_mono @test differentiate(const_alg_el, x) == differentiate(const_mono, x) end From fffee1c2870b91123a3eb3747404ef40e8452de7 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Beno=C3=AEt=20Legat?= Date: Wed, 3 Jul 2024 19:47:56 +0200 Subject: [PATCH 2/2] Not specific to Float64 --- src/scaled.jl | 22 +++++++++++----------- 1 file changed, 11 insertions(+), 11 deletions(-) diff --git a/src/scaled.jl b/src/scaled.jl index b5d6c7a..1eb89aa 100644 --- a/src/scaled.jl +++ b/src/scaled.jl @@ -51,13 +51,13 @@ end _float(::Type{T}) where {T<:Number} = float(T) # Could be for instance `MathOptInterface.ScalarAffineFunction{Float64}` # which does not implement `float` -_float(::Type{F}) where {F} = F +_float(::Type{F}) where {F} = MA.promote_operation(+, F, F) _promote_coef(::Type{T}, ::Type{ScaledMonomial}) where {T} = _float(T) -function MP.polynomial(f::Function, basis::SubBasis{ScaledMonomial}) +function MP.polynomial(f::Function, basis::SubBasis{ScaledMonomial}, ::Type{T}) where {T} return MP.polynomial( - i -> scaling(basis.monomials[i]) * f(i), + i -> scaling(T, basis.monomials[i]) * f(i), basis.monomials, ) end @@ -69,12 +69,12 @@ function Base.promote_rule( return SubBasis{Monomial,M,V} end -function scaling(m::MP.AbstractMonomial) - return √( +function scaling(::Type{T}, m::MP.AbstractMonomial) where {T} + return √(T( factorial(MP.degree(m)) / prod(factorial, MP.exponents(m); init = 1), - ) + )) end -unscale_coef(t::MP.AbstractTerm) = MP.coefficient(t) / scaling(MP.monomial(t)) +unscale_coef(t::MP.AbstractTerm) = MP.coefficient(t) / scaling(MP.coefficient_type(t), MP.monomial(t)) function SA.coeffs( t::MP.AbstractTermLike, ::FullBasis{ScaledMonomial}, @@ -84,10 +84,10 @@ function SA.coeffs( return MP.term(mono * MP.coefficient(t), mono) end function MP.coefficients(p, ::FullBasis{ScaledMonomial}) - return unscale_coef.(MP.terms(p)) + return unscale_coef.(MP.coefficient_type(p), MP.terms(p)) end function MP.coefficients(p, basis::SubBasis{ScaledMonomial}) - return MP.coefficients(p, basis.monomials) ./ scaling.(MP.monomials(p)) + return MP.coefficients(p, basis.monomials) ./ scaling.(MP.coefficient_type(p), MP.monomials(p)) end function SA.coeffs( p::MP.AbstractPolynomialLike, @@ -109,7 +109,7 @@ function SA.coeffs!( SA.unsafe_push!( res, target[Polynomial{Monomial}(mono)], - v * scaling(mono), + v * scaling(eltype(res), mono), ) end MA.operate!(SA.canonical, res) @@ -128,7 +128,7 @@ function SA.coeffs!( SA.unsafe_push!( res, target[Polynomial{ScaledMonomial}(mono)], - v / scaling(mono), + v / scaling(eltype(T), mono), ) end MA.operate!(SA.canonical, res)