JuMP with user-defined non linear functions

Specific Domains Optimization (Mathematical)
1

Hi everyone!

I am new to Julia and I am trying to solve a highly non linear optimization problem. I am encountering a problem when registering the objective function and very little idea about what could be going on. Heres the piece of code:

function S0(kappa,age)
    Add_0 = zeros(Real,size(age,1),1);
    Add_0 = kappa[1] .+ kappa[2].*(age .== 2) .+ kappa[3].*(age .== 3) .+ kappa[4].*(age .== 4) .+
    kappa[5].*(age.== 5) .+ kappa[6].*(age .== 6);
    return Add_0;
end

where Age is a Nx1 array and kappa is a 6x1 array.

function S_t(x,InitialAdd,a)

Add = zeros(Real,size(a,1),size(a,2));
   for t in 1:size(a,2)
         if t == 1
            Add[:,1] .= (1-x[2]) .* InitialAdd;
         else
            Tmax = t-1;
            aux = Array{Real}(undef,Tmax,1);
         for tt = 1:Tmax
            aux[tt,1] = (1-x[2]).^(Tmax-tt);
         end
         Add[:,t] .= ((1-x[2])^t) .* InitialAdd .+ x[1] .* (a[:,1:t-1]*aux[:,1]);
      end
   end
return Add;
end

where InitialAdd is the output of S0(kappa,age), x is a 2x1 array with parameters and a is a NxT vector so the output of this functions is a NxT array.
The previous two are two user defined functions that I know use in the optimization routine as follows


m = Model(Ipopt.Optimizer)
@variable(m, x[1:8])

function Stinner2(x...)
    aux1 = Array{Real}(undef,N,T)
    aux = S_t(x[1:2],S0(x[3:8],Age),a);
    for i in 1:N
        for t in 1:T
            aux1[i,t] = aux[i,t]^2;
        end
    end
    return sum(aux1[i,j] for i = 1:N, j = 1:T);
end
JuMP.register(m, :Stinner2, 8, Stinner2; autodiff=true)

@NLobjective(m, Min, Stinner2(x[1],x[2],x[3],x[4],x[5],x[6],x[7],x[8]))


JuMP.optimize!(m)

The error that appears is

Unable to register the function :Stinner2 because it does not support differentiation via ForwardDiff.

Weird thing is that I sometimes run it and it works and others it does not. There are some suggestions in the documentation like avoiding vector operations (which I do have in the two previously defined functions. However, I have tried to build them using loops and scalars and the same problem still persists), creating Real temporary arrays but it seems that none of them are working. If I input some numbers into Stinner2 I do get an answer but still the program fails to register it. All the variables (N,T, etc.) that appear in the code are previously defined so that seems not to be the issue either.

Any help would be greatly appreciated,

Best!

2

I don’t have enough information to reproduce this. What is a and Age in Stinner2?

2 Likes
3

You can tidy up the functions to remove a lot of the temporary matrices. This should point you in the right direction. I haven’t tested because I don’t know Age or a, so there might be some typos etc.

S0(kappa,age) = kappa[1] .+ sum(kappa[i] .* (age .== i) for i in 2:6)

function S_t(x::Vector{T}, InitialAdd, a) where {T}
    output = zero(T)
    for t in 1:size(a, 2)
        if t == 1
            for i in 1:size(a, 1)
                output += ((1 - x[2]) * InitialAdd[i])^2
            end
        else
            Tmax = t - 1
            aux = T[(1 - x[2]) .^ (Tmax - tt) for tt in 1:Tmax]
            for i in 1:size(a, 1)
                y = (1 - x[2])^t * InitialAdd[i] + x[1] * a[i, 1:t-1] * aux
                output += y^2
            end
        end
    end
    return output
end

# TODO: what is Age, a?
Stinner2(x...) = S_t(x[1:2], S0(x[3:8], Age), a)

m = Model(Ipopt.Optimizer)
@variable(m, x[1:8])
JuMP.register(m, :Stinner2, 8, Stinner2; autodiff=true)
@NLobjective(m, Min, Stinner2(x...))
optimize!(m)
2 Likes
4

Thanks a lot for the swift response! With a few modifications it worked. Age and a were Nx1 and NxT arrays, I had to adjust some of the code to incorporate this and it worked. thanks!

2 Likes