Matrix inverse as objective

Specific Domains Optimization (Mathematical)
,
1

Hi everyone,

I am trying to solve a semi-definite program with an objective that contains a matrix inverse.

More precisely, let a \geq 0, and A be a positive definite matrix of size d \times d. \mathbf{1} below is a vector of d ones. I want to solve the following program

\max_{X}\quad [a + \mathbf{1}'(A +X)^{-1} \mathbf{1}]^{-1} \qquad \text{s.t.} \qquad \text{$X$ is a $d\times d$ positive semidefinite matrix}

Is there a way to solve this using the existing functions in Julia and JuMP?

Maybe I have missed it, but I tried looking into the Standard form · MathOptInterface functions and could not find something helpful for my problem.

Any thoughts are appreciated. Thank you!

1 Like
2

Matrix inverse + PSD means there’s usually a trick with Schur’s lemma: 6 Semidefinite optimization — MOSEK Modeling Cookbook 3.3.1

I assume you have some other constraints, so that your objective is always >= 0.

Then it’s equivalent to minimizing the term without the outer inverse.

a is constant so it can be dropped.

You want \min_{X\succeq 0}: 1^\top(A + X)^{-1} 1.

Rewrite to epigraph form: \min_{X\succeq 0,C}: C such that C \succeq 1^\top(A + X)^{-1} 1.

By Schur’s lemma \min_{X\succeq 0,C}: C such that \begin{bmatrix}C & 1^\top\\ 1 &(A+X)\end{bmatrix} \succeq 0.

In JuMP:

julia> using JuMP, SCS

julia> A = [3 -2; -2 4]
2×2 Matrix{Int64}:
  3  -2
 -2   4

julia> d = size(A, 1)
2

julia> model = Model(SCS.Optimizer)
A JuMP Model
├ solver: SCS
├ objective_sense: FEASIBILITY_SENSE
├ num_variables: 0
├ num_constraints: 0
└ Names registered in the model: none

julia> @variable(model, X[1:d, 1:d], PSD)
2×2 Symmetric{VariableRef, Matrix{VariableRef}}:
 X[1,1]  X[1,2]
 X[1,2]  X[2,2]

julia> @constraint(model, LinearAlgebra.tr(X) == 1)  # Or something?
X[1,1] + X[2,2] = 1

julia> @variable(model, C)
C

julia> @objective(model, Min, C)
C

julia> @constraint(model, [C ones(d)'; ones(d) (A + X)] >= 0, PSDCone())
[C  1           1
 1  X[1,1] + 3  X[1,2] - 2
 1  X[1,2] - 2  X[2,2] + 4] ∈ PSDCone()

julia> optimize!(model)
------------------------------------------------------------------
	       SCS v3.2.7 - Splitting Conic Solver
	(c) Brendan O'Donoghue, Stanford University, 2012
------------------------------------------------------------------
problem:  variables n: 4, constraints m: 10
cones: 	  z: primal zero / dual free vars: 1
	  s: psd vars: 9, ssize: 2
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): 9, nnz(P): 0
------------------------------------------------------------------
 iter | pri res | dua res |   gap   |   obj   |  scale  | time (s)
------------------------------------------------------------------
     0| 1.74e+01  1.08e+00  1.90e+01 -7.86e+00  1.00e-01  3.13e-04 
    50| 6.41e-05  1.28e-05  7.62e-05  8.12e-01  1.00e-01  8.29e-04 
------------------------------------------------------------------
status:  solved
timings: total: 8.41e-04s = setup: 1.17e-04s + solve: 7.23e-04s
	 lin-sys: 2.50e-05s, cones: 2.79e-04s, accel: 5.37e-06s
------------------------------------------------------------------
objective = 0.812454
------------------------------------------------------------------

julia> value(C)
0.8124158521710856

julia> ones(d)' * inv(A + value(X)) * ones(d)
0.812469395065637
8 Likes
3

Thank you!!!

1 Like