Help tracking down an InexactError()


#1

I’m trying to solve a dynamic programming problem and I am running into an InexactError() that I have not been able to track down. Can you please help me out?

const nᵢ = 2

const δ = 0.8
const σ = 1.5
const ρ = 5.0

const S = 1:3
const γc = 1.2
const πmat = [  0.9 0.05 0.05;
                0.05 0.9 0.05;
                0.05 0.05 0.9 ]
const hi = (0.52 * γc) ^ (1 - σ) / (1 - σ)
const lo = (0.35 * γc) ^ (1 - σ) / (1 - σ)
const w = [     hi lo lo; # is
                lo hi lo ]

function J(m, y, s)
    σ / (1 - σ) * sum( ((1 - δ) * (y[i] + m[i])) ^ (1 / σ) for i in 1:nᵢ ) * m[nᵢ + 1] ^ ((σ - 1) / σ) -
    sum( m[i] * w[i, s] for i in 1:nᵢ ) + m[nᵢ + 1] * γc
end

const grid_m = 0:0.1:1.5
const grid_y′ = 0:0.1:1

using Interpolations
function interp1(S, grid_y′1, grid_y′2, Darr)
    itp = interpolate(Darr, BSpline(Cubic(Line())), OnGrid())
    sitp = scale(itp, S, grid_y′1, grid_y′2)
    return sitp
end

interp(S, grid_y′1, grid_y′2, Darr, s, y1, y2) = interp1(S, grid_y′1, grid_y′2, Darr)[s, y1, y2]

using Optim

function main()
    s = 2
    y1 = 0.0
    y2 = sqrt(abs(1.0 - y1^2.0))

    Darr = [J([1.0, 1.0, 1.0], [y1, y2], s) for s in S, y1 in grid_y′, y2 in grid_y′]

    z = x -> begin
                m1 = x[1]
                m2 = x[2]
                m3 = x[3]
                y′1 = x[4]
                y′2 = sqrt(abs(1.0 - y′1^2.0))
                return J([m1, m2, m3], [y1, y2], s) + δ * interp(S, grid_y′, grid_y′, Darr, s, y′1, y′2)
            end
        twicedif = TwiceDifferentiable(z, [0.0, 0.0, 0.0, 0.0]; autodiff = :forward)
        opt = optimize(twicedif, [1.0, 1.0, 1.0, 1.0], Newton(linesearch = LineSearches.Static()))
end

Running main() gives the error.


#2

I think you are running into the float indices being too large to be represented in Int64, eg like in

julia> trunc(Int64, typemax(Int64) + 1.0)
ERROR: InexactError()
Stacktrace:
 [1] trunc(::Type{Int64}, ::Float64) at ./float.jl:672

It is possible that you are leaving the grid when looking for the optimum. Since B-splines don’t necessarily do well outside the grid, I would just restrict the policy choice to something sensible. Also, note that optimize minimizes, is this what you want?


#3

The error is thrown in getindex and when it happens you have (s, y1, y2) = (2, NaN, NaN) so I guess trying to getindex the return value of interp1 with NaN doesn’t work so well.


#4

Yes, I’m looking for a minimum (unfortunate and inconsistent choice of naming optimize, since also to obtain the results one calls minimize and minimum). Is there way to provide bounds to optimize?


#5

NaN? When does the program go to those values?


#6

Yes, using box minimization.


#7

I printed the values inside interp and that’ what they had when the error was thrown


#8

So it seems that was the problem, now with constraints I obtain the solution with no problem. Thanks!