Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Fix monomials for LexOrder #159

Merged
merged 4 commits into from
May 2, 2024
Merged

Fix monomials for LexOrder #159

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
4 commits
Select commit Hold shift + click to select a range
8ca687c
Fix monomials for LexOrder
blegat Apr 30, 2024
06f7e1b
Fixes
blegat Apr 30, 2024
e0f9594
Fix
blegat May 2, 2024
30d73a5
Fix
blegat May 2, 2024
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
103 changes: 79 additions & 24 deletions src/monomial_vector.jl
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -12,9 +12,6 @@ struct MonomialVector{V,M} <: AbstractVector{Monomial{V,M}}
@assert !iscomm(V) || issorted(vars, rev = true)
@assert all(z -> length(z) == length(vars), Z)

_isless = let M = M
(a, b) -> MP.compare(a, b, M) < 0
end
return new{V,M}(vars, Z)
end
end
Expand Down Expand Up @@ -125,68 +122,126 @@ function _error_for_negative_degree(deg)
end
end

function fillZfordeg!(Z, n, deg, ::Type{Commutative}, filter::Function, ::Int)
const _Lex = Union{MP.LexOrder,MP.InverseLexOrder}

_last_lex_index(n, ::Type{MP.LexOrder}) = n
_prev_lex_index(i, ::Type{MP.LexOrder}) = i - 1
_not_first_indices(n, ::Type{MP.LexOrder}) = n:-1:2
_last_lex_index(_, ::Type{MP.InverseLexOrder}) = 1
_prev_lex_index(i, ::Type{MP.InverseLexOrder}) = i + 1
_not_first_indices(n, ::Type{MP.InverseLexOrder}) = 1:(n-1)

function _fill_exponents!(Z, n, degs, ::Type{Commutative}, M::Type{<:_Lex}, filter::Function)
_error_for_negative_degree.(degs)
maxdeg = maximum(degs, init = 0)
I = _not_first_indices(n, M)
z = zeros(Int, n)
while true
deg = sum(z)
if deg in degs && filter(z)
push!(Z, z)
z = copy(z)
end
if deg == maxdeg
i = findfirst(i -> !iszero(z[i]), I)
if isnothing(i)
break
end
j = I[i]
z[j] = 0
z[_prev_lex_index(j, M)] += 1
else
z[_last_lex_index(n, M)] += 1
end
end
end

function _fill_exponents!(Z, n, deg, ::Type{Commutative}, M::Type{<:_Lex}, filter::Function, ::Int)
_error_for_negative_degree(deg)
I = _not_first_indices(n, M)
z = zeros(Int, n)
z[end] = deg
z[_last_lex_index(n, M)] = deg
while true
if filter(z)
push!(Z, z)
z = copy(z)
end
if z[1] == deg
i = findfirst(i -> !iszero(z[i]), I)
if isnothing(i)
break
end
i = findfirst(i -> !iszero(z[i]), n:-1:2)
j = (n:-1:2)[i]
j = I[i]
p = z[j]
z[j] = 0
z[end] = p - 1
z[j-1] += 1
z[_last_lex_index(n, M)] = p - 1
z[_prev_lex_index(j, M)] += 1
end
end
function fillZrec!(Z, z, i, n, deg, filter::Function)

function _fill_noncomm_exponents_rec!(Z, z, i, n, deg, ::Type{MP.LexOrder}, filter::Function)
if deg == 0
if filter(z)
push!(Z, copy(z))
end
else
for i in i:i+n-1
z[i] += 1
fillZrec!(Z, z, i, n, deg - 1, filter)
_fill_noncomm_exponents_rec!(Z, z, i, n, deg - 1, LexOrder, filter)
z[i] -= 1
end
end
end
function fillZfordeg!(

function _fill_exponents!(
Z,
n,
deg,
::Type{NonCommutative},
::Type{MP.LexOrder},
filter::Function,
maxdeg::Int,
)
_error_for_negative_degree(deg)
_error_for_negative_degree(maxdeg)
z = zeros(Int, maxdeg * n - maxdeg + 1)
start = length(Z) + 1
fillZrec!(Z, z, 1, n, deg, filter)
_fill_noncomm_exponents_rec!(Z, z, 1, n, deg, MP.LexOrder, filter)
return reverse!(view(Z, start:length(Z)))
end
# List exponents in decreasing Graded Lexicographic Order
function getZfordegs(

function _fill_exponents!(Z, n, deg, ::Type{V}, ::Type{MP.Reverse{M}}, args...) where {V,M}
prev = lastindex(Z)
_fill_exponents!(Z, n, deg, V, M, args...)
reverse!(view(Z, (prev + 1):lastindex(Z)))
return
end

function _fill_exponents!(
Z::Vector{Vector{Int}},
n,
degs::AbstractVector{Int},
::Type{V},
::Type{M},
::Type{MP.Graded{M}},
filter::Function,
) where {V,M}
Z = Vector{Vector{Int}}()
# For non-commutative, lower degree need to create a vector of exponent as large as for the highest degree
maxdeg = isempty(degs) ? 0 : maximum(degs)
maxdeg = maximum(degs, init = 0)
for deg in sort(degs)
fillZfordeg!(Z, n, deg, V, filter, maxdeg)
_fill_exponents!(Z, n, deg, V, M, filter, maxdeg)
end
return
end

# List exponents in decreasing Graded Lexicographic Order
function _all_exponents(
n,
degs::AbstractVector{Int},
::Type{V},
::Type{M},
filter::Function,
) where {V,M}
Z = Vector{Vector{Int}}()
_fill_exponents!(Z, n, degs, V, M, filter)
_isless = let M = M
(a, b) -> MP.compare(a, b, M) < 0
end
Expand All @@ -202,7 +257,7 @@ function MonomialVector(
vars = unique!(sort(vars, rev = true))
return MonomialVector(
vars,
getZfordegs(
_all_exponents(
length(vars),
degs,
Commutative,
Expand All @@ -222,7 +277,7 @@ function MonomialVector(
filter::Function = x -> true,
) where {M}
vars = unique!(sort(vars, rev = true))
Z = getZfordegs(
Z = _all_exponents(
length(vars),
degs,
NonCommutative,
Expand All @@ -248,11 +303,11 @@ function MP.monomials(vars::Tuple{Vararg{Variable}}, args...)
end

#function MP.monomials(vars::TupOrVec{Variable{true}}, degs::AbstractVector{Int}, filter::Function = x->true)
# Z = getZfordegs(length(vars), degs, true, z -> filter(Monomial(vars, z)))
# Z = _all_exponents(length(vars), degs, true, z -> filter(Monomial(vars, z)))
# [Monomial{true}(vars, z) for z in Z]
#end
#function MP.monomials(vars::TupOrVec{<:Variable{<:NonCommutative}}, degs::AbstractVector{Int}, filter::Function = x->true)
# Z = getZfordegs(length(vars), degs, false, z -> filter(Monomial(vars, z)))
# Z = _all_exponents(length(vars), degs, false, z -> filter(Monomial(vars, z)))
# v = isempty(Z) ? vars : getvarsforlength(vars, length(first(Z)))
# [Monomial(v, z) for z in Z]
#end
Expand Down
62 changes: 62 additions & 0 deletions test/comp.jl
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -44,3 +44,65 @@ end
@test ordering(x[1]) == order
@test issorted(monomials(x[1], 0:2))
end

function _test_less(a, b)
@test a < b
@test b > a
end

function _test_monomials(vars, degs, exp)
# Without `collect`, `exp` is promoted to a `MonomialVector`
# which sorts it so it doesn't test the order
@test collect(monomials(vars, degs)) == exp
end

@testset "LexOrder" begin
@polyvar x y monomial_order = LexOrder
_test_less(y, y^2)
_test_less(x^0, y)
_test_less(y^2, x)
_test_less(x * y^2, x^2)
_test_monomials([x, y], 3, [y^3, y^2 * x, y * x^2, x^3])
_test_monomials([x, y], 1:2, [y, y^2, x, x * y, x^2])
_test_monomials([x, y], [0, 1, 3], [1, y, y^3, x, y^2 * x, y * x^2, x^3])
end

@testset "InverseLexOrder" begin
@polyvar x y monomial_order = InverseLexOrder
_test_less(y, y^2)
_test_less(x^0, y)
_test_less(x, y^2)
_test_less(x^2, x * y^2)
_test_monomials([x, y], 3, [x^3, x^2 * y, x * y^2, y^3])
end

@testset "Reverse{LexOrder}" begin
@polyvar x y monomial_order = Reverse{LexOrder}
_test_less(y^2, y)
_test_less(y, x^0)
_test_less(x, y^2)
_test_less(x^2, x * y^2)
_test_monomials([x, y], 3, [x^3, x^2 * y, x * y^2, y^3])
end

@testset "Reverse{InverseLexOrder}" begin
@polyvar x y monomial_order = Reverse{InverseLexOrder}
_test_less(y^2, y)
_test_less(y, x^0)
_test_less(y^2, x)
_test_less(x * y^2, x^2)
_test_monomials([x, y], 3, [y^3, y^2 * x, y * x^2, x^3])
_test_monomials([x, y], 1:2, [y^2, x * y, y, x^2, x])
_test_monomials([x, y], [0, 1, 3], [y^3, x*y^2, x^2*y, y, x^3, x, 1])
end

@testset "Graded{Reverse{InverseLexOrder}}" begin
@polyvar x y monomial_order = Graded{Reverse{InverseLexOrder}}
_test_less(y, y^2)
_test_less(x^0, y)
_test_less(x, y^2)
_test_less(x^2, x * y^2)
_test_monomials([x, y], 3, [y^3, y^2 * x, y * x^2, x^3])
_test_monomials([x, y], 1:2, [y, x, y^2, x * y, x^2])
_test_monomials([x, y], [0, 1, 3], [1, y, x, y^3, x*y^2, x^2*y, x^3])
end
Loading