JuliaAlgebra / StaticPolynomials.jl

Fast evaluation of multivariate polynomials
https://juliaalgebra.github.io/StaticPolynomials.jl/latest/
Other
16 stars 4 forks source link

Unexplained Allocation #8

Closed cortner closed 6 years ago

cortner commented 6 years ago

Consider the following, still relatively simple polynomial:

using StaticArrays, BenchmarkTools
using DynamicPolynomials: @polyvar
import StaticPolynomials

@polyvar x1 x2 x3 x4 x5 x6 x7 x8 x9 x10

P = StaticPolynomials.system([   x1 + x2 + x3 + x4, x5 + x6 + x7 + x8 + x9 + x10,
    x1^2 + x2^2 + x3^2 + x4^2,
    x1*x5 + x1*x6 + x1*x7 + x2*x5 + x2*x8 + x2*x9 + x3*x6 + x3*x8 + x3*x10 + x4*x7 + x4*x9 + x4*x10,
    x5^2 + x6^2 + x7^2 + x8^2 + x9^2 + x10^2,
    x5*x6 + x5*x7 + x5*x8 + x5*x9 + x6*x7 + x6*x8 + x6*x10 + x7*x9 + x7*x10 + x8*x9 + x8*x10 + x9*x10,
    x1^3 + x2^3 + x3^3 + x4^3,
    x1^2*x5 + x1^2*x6 + x1^2*x7 + x2^2*x5 + x2^2*x8 + x2^2*x9 + x3^2*x6 + x3^2*x8 + x3^2*x10 + x4^2*x7 + x4^2*x9 + x4^2*x10,
    x1*x2*x5 + x1*x3*x6 + x1*x4*x7 + x2*x3*x8 + x2*x4*x9 + x3*x4*x10,
    x1*x5^2 + x1*x6^2 + x1*x7^2 + x2*x5^2 + x2*x8^2 + x2*x9^2 + x3*x6^2 + x3*x8^2 + x3*x10^2 + x4*x7^2 + x4*x9^2 + x4*x10^2,
    x1*x5*x6 + x1*x5*x7 + x1*x6*x7 + x2*x5*x8 + x2*x5*x9 + x2*x8*x9 + x3*x6*x8 + x3*x6*x10 + x3*x8*x10 + x4*x7*x9 + x4*x7*x10 + x4*x9*x10,
    x5^3 + x6^3 + x7^3 + x8^3 + x9^3 + x10^3,
    (x5^2*x6 + x5^2*x7 + x5^2*x8 + x5^2*x9 + x5*x6^2 + x5*x7^2 + x5*x8^2 +
      x5*x9^2 + x6^2*x7 + x6^2*x8 + x6^2*x10 + x6*x7^2 + x6*x8^2 + x6*x10^2 +
      x7^2*x9 + x7^2*x10 + x7*x9^2 + x7*x10^2 + x8^2*x9 + x8^2*x10 + x8*x9^2 +
      x8*x10^2 + x9^2*x10 + x9*x10^2),
    x5*x6*x7 + x5*x8*x9 + x6*x8*x10 + x7*x9*x10,
    x1^4 + x2^4 + x3^4 + x4^4 ])

xx = @SVector rand(10)
@btime StaticPolynomials.evaluate($P, $xx) #   384.780 ns (17 allocations: 496 bytes)
@btime StaticPolynomials.jacobian($P, $xx) #   300.358 ns (0 allocations: 0 bytes)

Three things are weird here:

saschatimme commented 6 years ago

The reason is an annoying internal tuple limit in combination with the SVector constructor.

SVector(a1,a2,...,a15)  # this allocates
SVector((a1,a2,...,a15))  # this does not allocate

I already run into this at some point earlier but it seems that it slip through in evaluate. This is now fixed on master.

cortner commented 6 years ago

thanks! That also explains a problem I just had with allocating SVectors that are longer than 15 elements :/

saschatimme commented 6 years ago

Yeah, there are really some pitfalls one needs to be aware of.

cortner commented 6 years ago

thanks for the quick fix!