Preallocating Interpolation object

[tradler]

I am solving a dynamic programing value function using backward induction, and would like to be able to update an Interpolations.jl object (or something similar) which holds an approximation of the value function in each period.

I first preallocate an array which holds the value function, which is modified in each period without any additional memory allocation. My hope is to do the same for the approximation. The reason I am hoping to do this is that in the actual code the approximation for the next-period’s value function is used to compute the current-period value functions, so I am hoping to precompute the approximation only once.

A MWE is below. The goal is for the update_inter! function to not allocate any memory. Any advice in the definition of the interpolation object in the structure would also be greatly appreciated, just because it is a bit clunky.

using BenchmarkTools, Interpolations
struct objects
    vgrid::Vector{Float64}
    agrid::Vector{Float64}
    vfunc::Array{Float64,3}
    interp::Vector{Interpolations.GriddedInterpolation{Float64,2,Matrix{Float64},Gridded{Linear{Throw{OnGrid}}},Tuple{Vector{Float64},Vector{Float64}}}}

    function objects()

        stdV = 2.5
        vpts = 5
        vgrid = collect(range(-stdV, stdV, vpts))
        apts = 10
        agrid = collect(range(0, 1000.0, apts))
        TT = 20
        vfunc = zeros(vpts, apts, TT)
        interp = [Interpolations.interpolate((vgrid, agrid), vfunc[:, :, tt], Gridded(Linear())) for tt = 1:TT]

       return new(vgrid, agrid, vfunc, interp)
    end
end

function update_vfunc!(obj::objects, tt::Integer)
    for vi = eachindex(obj.vgrid), ai = eachindex(obj.agrid)
        obj.vfunc[vi, ai, tt] = obj.vgrid[vi] * obj.agrid[ai] * tt
    end
end

function update_interp!(obj::objects, tt::Integer)
    obj.interp[tt] = Interpolations.interpolate((obj.vgrid, obj.agrid), view(obj.vfunc, :, :, tt), Gridded(Linear()))
end

function solve(obj::objects)
    for tt = 1:obj.tt
        update_vfunc!(obj, tt)
        update_interp!(obj, tt)
    end
end

 obj = objects()
 
 @btime update_vfunc!($obj, 20)
 @btime update_interp!($obj, 20)

update_vfunc! does not allocate any memory, but update_interp! does.

[gdalle]

Hi there!
As far as the struct is concerned, seems fine to me. The important thing is to not have any abstract types lying around, and it looks like you did it right:

julia> isconcretetype(objects)
true

[gdalle]

To find where the allocations are coming from, I profiled your function in VSCode, running it several times to capture enough samples:

 @profview for k in 1:100000; update_interp!(obj, 20); end

Here’s what the result looks like. Every yellow tile is allocating memory.

Unfortunately, it looks like most of it is the fault of Interpolations.jl, so I’m not sure what your margin of improvement is here

[tradler]

Thanks! Yes, I think that is right. I see that in the package there is an interpolate! function which turns the input into an interpolation object, I think without allocating memory, but that is not quite what I want. Hopefully there is a workaround or alterative package.

I currently am using a function I wrote which does a linear approximation each time I need to interpolate something. It might just be easier to keep using that.