JuMP Efficient Quadratic Objective Formulation

Specific Domains Optimization (Mathematical)
1

Hi everyone,

I (eventually) want to run an MIQP where the objective is a simple sum of squares say sum((Y-XB)^2). Without even getting to the constraints or the integer part of the problem, I’m running into issues simply formulating the objective for large X and Y (say X~500x1000 and Y~500x1). Specifically, I get an “OutofMemoryError()”, note that I’m not solving anything yet just setting up the model. I have a simple function snippet below I’ve been using to see how big I can make X before things break. The use of “dot()” may be silly/naive to more experienced julia users, it was just working faster then writing out the whole thing as a sum. So far I’ve made it up to n=500 and p=900 before the error appears (applying @time to this shows the output: “64.652337 seconds (1.81 M allocations: 18.305 GB, 19.55% gc time)”)
I’m working on a cluster and have to request memory for jobs, so theoretically I could request more though that puts me further down the queue. So if there are glaring inefficiencies here or nice tricks to use, any input would be great!

function test(n,p)
X=randn(n,p)
Y=randn(n)
m=Model(solver=CplexSolver())
@variable(m, B[1:p])
W = Y-X*B
@objective(m,Min,dot(W,W))
end

1 Like
2

Try re-writting this way:

function test2(n,p)
    X=randn(n,p)'
    Y=randn(n)
    m=Model(solver=CplexSolver())
    @variable(m, B[1:p])
    @variable(m, w[1:n])
    @constraint(m, ctr[i in 1:n], w[i]==Y[i]-sum(X[j,i]*B[j] for j in 1:p))
    @objective(m,Min,sum(w[i]*w[i] for i in 1:n))
    nothing    
end

p=n=1000 takes 1s on my notebook

2 Likes
3

Thanks joaquimg! Works great, just curious but is there any particular reason to transpose X and go across rows in sum(X[j,i]*B[j] rather then columns? Are longer matrices faster to process/store then wider ones?

4

That was not the key feature in the rewrite, but inner loop should go column-wise in julia.
Actually the I should have written: X=randn(p,n).

5

Got it, and, if I’m understanding correctly, it seems like the key features would be defining placeholder variables and working with looped operations rather then straight matrix/vector multiplication.