Unexpected Allocations in Julia 0.5

General Usage
1

I am working on an ODE/DAE solver in Julia (GitHub - JuliaGNI/GeometricIntegrators.jl: Geometric Numerical Integration in Julia), where I encountered a problem with unexpected allocations when calling a function which is part of a type, e.g.,

type Data{T}
    t::T          # time t
    x::Vector{T}  # solution x
    y::Vector{T}  # vector field y=f(t,x)
    f::Function   # function f
end

called by

function f(t, x, fx)
    fx[1] = x[2]
    fx[2] = sin(x[1])
    nothing
end

t = 1.
x = [1., 2.]
y = zeros(x) 
d = Data(t, x, y, f)

@time d.f(d.t, d.x, d.y)

which (after warmup) leads to

  0.000002 seconds (1 allocation: 16 bytes)

The problem is the t parameter. If fixed, e.g.,

@time d.f(1., d.x, d.y)

the allocation disappears and the result becomes

  0.000001 seconds

Similarly, if the t parameter is removed altogether or wrapped into an array, the allocation disappears as well, e.g.,

function f(x, fx)
    fx[1] = x[2]
    fx[2] = sin(x[1])
    nothing
end

type Data{T}
    t::T          # time t
    x::Vector{T}  # solution x
    y::Vector{T}  # vector field y=f(t,x)
    f::Function   # function f
end

t = 1.
x = [1., 2.]
y = zeros(x) 
d = Data(t, x, y, f)

@time d.f(d.x, d.y)

which leads to

  0.000000 seconds

or

function f(t, x, fx)
    fx[1] = x[2]
    fx[2] = sin(x[1])
    nothing
end

type Data{T}
    t::Vector{T}
    x::Vector{T}
    y::Vector{T}
    f::Function
end

t = [1.]
x = [1., 2.]
y = zeros(x) 
d = Data(t, x, y, f)

@time d.f(d.t, d.x, d.y)

which leads to

  0.000001 seconds

So Julia seems to have troubles figuring out the type of the scalar. This problem persists even if I fix the type in the function header, i.e.,

function f{T}(t::T, x::Vector{T}, y::Vector{T})
   ...
end

or

function f(t::Float64, x::Vector{Float64}, y::Vector{Float64})
   ...
end

Below is the full “minimal” example,

# define function to integrate
function f(t, x, fx)
    fx[1] = x[2]
    fx[2] = sin(x[1])
    nothing
end

# define datatype which holds time t, solution x,
# vector field y=f(t,x), and function f
type Data{T}
    t::T
    x::Vector{T}
    y::Vector{T}
    f::Function
end

# function for timing f with data from data type
function call_f(d::Data)
    d.f(d.t, d.x, d.y)
    
    @time d.f(d.t, d.x, d.y)
    @time d.f(1.0, d.x, d.y)
end

# create data type and call f with wrapped data
function time_f1()
    t = 1.
    x = [1., 2.]
    y = zeros(x) 
    my_data = Data(t, x, y, f)
    call_f(my_data)
end

# function for timing f without wrapping data
function time_f2()
    t = 1.
    x = [1., 2.]
    y = zeros(x)

    f(t, x, y)
    
    @time f(t, x, y)
    @time f(1., x, y)
end

# time f, calling it with and without data type
time_f1()
time_f2()

Calling this script leads to

  0.000002 seconds (1 allocation: 16 bytes)
  0.000001 seconds
  0.000000 seconds
  0.000000 seconds
2

What happens if you

const t = 1.
3

const will not work in local scope.

Your issue is not from t but f. The call on d.f cannot be inferred so it’s a dynamic dispatch. t needs to be boxed as a result of that. Allocation is only needed if it was not boxed already.

If you want to optimize for performance, you should remove the dynamic dispatch which should remove the allocation too although the allocation you see should not cause any noticeable performance issue.

4

Instead:

type Data{T,F}
    t::T          # time t
    x::Vector{T}  # solution x
    y::Vector{T}  # vector field y=f(t,x)
    f::F   # function f
end

Each function has its own type, thus to make the field concrete, you need to give it a type parameter.

2 Likes
5

@yuyichao, @mauro3
Thanks for making this clear. Now I see the problem. @mauro3’s fix solved the issue. Regarding performance, for simple functions like in my example, the difference is usually about a factor of ~1.5, but can be up to ~2. For more complicated functions I expect it to be smaller, though.