Nuclear norm minimization objective using JuMP

Specific Domains Optimization (Mathematical)
2

Hi @MLL, welcome to JuMP.

First, here’s how you can solve your model with JuMP. (Note that I made up B, because it wasn’t in your example code.)

julia> using JuMP

julia> import Random

julia> import SCS

julia> begin
           rng = Random.MersenneTwister(1234)
           n, r = 20, 2
           Y = randn(rng, n, r) * randn(rng, r, n)
           X = Y .* (rand(rng, n, n) .< 1 / 2)
           mask = X .!= 0
           B = rand(rng, n, n)
           model = Model(SCS.Optimizer)
           @variable(model, X[1:n, 1:n])
           @constraint(model, X[mask] .== B[mask])
           @variable(model, t)
           @constraint(model, [t; vec(X)] in MOI.NormNuclearCone(n, n))
           @objective(model, Min, t)
           optimize!(model)
       end
------------------------------------------------------------------
	       SCS v3.2.6 - Splitting Conic Solver
	(c) Brendan O'Donoghue, Stanford University, 2012
------------------------------------------------------------------
problem:  variables n: 821, constraints m: 1006
cones: 	  z: primal zero / dual free vars: 185
	  l: linear vars: 1
	  s: psd vars: 820, ssize: 1
settings: eps_abs: 1.0e-04, eps_rel: 1.0e-04, eps_infeas: 1.0e-07
	  alpha: 1.50, scale: 1.00e-01, adaptive_scale: 1
	  max_iters: 100000, normalize: 1, rho_x: 1.00e-06
	  acceleration_lookback: 10, acceleration_interval: 10
	  compiled with openmp parallelization enabled
lin-sys:  sparse-direct-amd-qdldl
	  nnz(A): 1046, nnz(P): 0
------------------------------------------------------------------
 iter | pri res | dua res |   gap   |   obj   |  scale  | time (s)
------------------------------------------------------------------
     0| 7.78e+01  1.00e+00  3.78e+02 -1.72e+02  1.00e-01  1.76e-03 
   125| 3.73e-05  1.47e-05  3.48e-07  2.03e+01  6.29e-01  2.06e-02 
------------------------------------------------------------------
status:  solved
timings: total: 2.06e-02s = setup: 4.78e-04s + solve: 2.02e-02s
	 lin-sys: 1.85e-03s, cones: 1.63e-02s, accel: 8.10e-05s
------------------------------------------------------------------
objective = 20.349979
------------------------------------------------------------------

I note that we don’t have this example in our documentation, so it’s not surprising that you didn’t find it! I’ll fix by adding a new section to Modeling with cones · JuMP

Edit: see [docs] add Norm{Nuclear,Spectral}Cone examples by odow · Pull Request #3887 · jump-dev/JuMP.jl · GitHub

Your other attempts didn’t work for a number of reasons:

You need to use op_nuclearnorm in @objective, not the function nuclearnorm.

But operators take splatted scalar arguments, so you cannot pass a matrix X.

The result would be a nonlinear program, so you couldn’t use ProxSDP, you’d need a solver like Ipopt.

All together, it’d be something like this:

using JuMP, Ipopt
model = Model(Ipopt.Optimizer)
@variable(model, X[1:n, 1:n])
nuclearnorm(x...) = sum(LinearAlgebra.svdvals(reshape(collect(x), n, n)))
@operator(model, op_nuclearnorm, n^2, nuclearnorm)
@objective(model, Min, op_nuclearnorm(X...))

except this would still be broken, because we use ForwardDiff.jl to automatically compute the Jacobian of nuclearnorm, and ForwardDiff.jl does not support differentiating through svdvals.

So all in all, it’s just a lot easier to use MOI.NormNuclearCone!

3 Likes