SIMDPolynomials

---

SIMDPolynomials.jl uses bit-packed monomials so that most of the operations on multivariate monomials can be done in a few CPU instructions. Currently, it is only optimized for relatively small polynomials. Contributions, especially optimizations for large polynomials, are welcome!

Examples:

julia> using SIMDPolynomials, BenchmarkTools

julia> x, y, z, t = [PackedMonomial{4,7}(i) for i in 0:3]; # PackedMonomial with maximum of 4 variables and 7 bits of exponents.

julia> p = x * y + 3 * (z * t)
x₀x₁ + 3x₂x₃

julia> q = (p + 1) * p
x₀²x₁² + 6x₀x₁x₂x₃ + 9x₂²x₃² + x₀x₁ + 3x₂x₃

julia> @btime gcd($p, $q)
  4.019 μs (94 allocations: 8.06 KiB)
x₀x₁ + 3x₂x₃

julia> begin
           c1 = 10*(x * z + x)
           c2 = 2*(x^2 + z)
           c3 = 2*(2 - z  )
           c4 = 20*(x * z^2)
           e1 = 0
           e2 = 5
           e3 = 7
           e4 = 10
           p = c1 * y^e1 + c2 * y^e2 + c3 * y^e3 + c4 * y^e4
           q = prod(i->p + i, 0:3)
       end;

julia> @btime for i in 0:3
           gcd($p + i, $q)
       end
  350.086 μs (1159 allocations: 588.06 KiB)

The same micro-benchmark using AbstractAlgebra:

julia> using AbstractAlgebra, BenchmarkTools

julia> R, (x, y, z, t) = PolynomialRing(AbstractAlgebra.Integers{Int}(), [:x, :y, :z, :t], ordering=:deglex);

julia> p = x * y + 3 * (z * t)
x*y + 3*z*t

julia> q = (p + 1) * p
x^2*y^2 + 6*x*y*z*t + 9*z^2*t^2 + x*y + 3*z*t

julia> @btime gcd($p, $q) # SIMDPolynomials.jl is 30x faster
  119.795 μs (1320 allocations: 89.17 KiB)
x*y + 3*z*t

julia> begin
           c1 = 10*(x * z + x)
           c2 = 2*(x^2 + z)
           c3 = 2*(2 - z  )
           c4 = 20*(x * z^2)
           e1 = 0
           e2 = 5
           e3 = 7
           e4 = 10
           p = c1 * y^e1 + c2 * y^e2 + c3 * y^e3 + c4 * y^e4
           q = prod(i->p + i, 0:3)
       end;

julia> @btime for i in 0:3 # SIMDPolynomials.jl is 14x faster
           gcd($p + i, $q)
       end
  4.934 ms (32235 allocations: 3.43 MiB)

The same micro-benchmark using DynamicPolynomials:

julia> using DynamicPolynomials, BenchmarkTools

julia> @polyvar x y z t;

julia> p = x * y + 3 * (z * t)
xy + 3zt

julia> q = (p + 1) * p
x²y² + 6xyzt + 9z²t² + xy + 3zt

julia> @btime gcd($p, $q)  # SIMDPolynomials.jl is 65x faster
  264.561 μs (4962 allocations: 298.19 KiB)
xy + 3zt

julia> begin
           c1 = 10*(x * z + x)
           c2 = 2*(x^2 + z)
           c3 = 2*(2 - z  )
           c4 = 20*(x * z^2)
           e1 = 0
           e2 = 5
           e3 = 7
           e4 = 10
           p = c1 * y^e1 + c2 * y^e2 + c3 * y^e3 + c4 * y^e4
           q = prod(i->p + i, 0:3)
       end;

julia> @btime for i in 0:3 # SIMDPolynomials.jl is 82x faster
           gcd($p + i, $q)
       end
  28.943 ms (529642 allocations: 31.20 MiB)