Defining gradient

[rajatsanya...@gmail.com]

I am confused how to define a gradient of the following form-

(the equation has been written in latex format)

\begin{equation}
\nabla_X f(X) = \begin{bmatrix}
0 & & \lambda_1\max\{0,2X_{1,2}\} & & 0 & & \cdots & & 0 \\
0 & & 0 & & \lambda_2\max\{0,2X_{2,3}\} & & \cdots & & 0 \\
\vdots & & \vdots & & & & & & \vdots \\
0 & & 0 & & 0 & & \cdots & & \lambda_n\max\{0,2X_{n-1,n}\} \\
\end{bmatrix}
\end{equation}

Should I directly define it or is there any other smarter way to do it?

[BM]

Hello, 

I think, you could go ahead and directly define it in the function problem.egrad. Let us know how it goes. How big is "n" by the way? 

Regards,

BM

[rajatsanya...@gmail.com]

It is in the order of thousands.

[BM]

Okay... The for loop would be a way forward, but may not be suitable for those orders of magnitude loops.

What is X_{1,2} and what is the structure of X?

It seems to me that \nabla_X f(X) is mostly a sparse matrix and only the off diagonal elements have the max operations.... Is my understanding correct? If that is so, then the 'diag' function of Matlab would be helpful. 

Follow the following four steps. The text in blue are Matlab commands.

1) diag(X, 1) would collect all the first off diagonal elements of X. I assume that you would have the lambda vector with [lambda_1; lambda_2; ...; lambda_{n-1}] (@Rajat: verify the length of the lambda vector. I believe, the last element should be lambda_{n-1} and not lambda_n.)

The multiplication of lambda and max operation is defined by lambda.*max([diag(X,1), zeros(n-1,1)], [], 2). Lets call this vector b.

2) Lets create the indices of the off diagonal entries. (To be done only once!)

I = [1 : 1 : n-1];

J = [2: 1 : n];

IJ =  sub2ind([n n], I, J); % linear indices.

3) Create nabla_X = zeros(n,n); (Edit: To be done only once!)

4) Replace the off diagonal elements of nabla_X with b. This can be done by nabla_X(IJ) = b;

Let us know if this works out for you.

Regards,

BM

[Nicolas Boumal]

For reference, this is the LaTeX bit in the original post:

[Nicolas Boumal]

I think you want to use spdiags (look it up in Matlab): it helps create sparse matrices directly in Matlab, by specifying the main diagonal, and secondary diagonals as well. Should do the trick I suppose.

With which manifold do you use this?

Best,

Nicolas

[rajatsanya...@gmail.com]

Yes, there is a slip. Yes, it should be \lambda_{n-1} and the gradient should be-

\begin{equation}
\nabla_X f(X) = \begin{bmatrix}
0 & & \lambda_1\max\{0,2X_{1,2}\} & & 0 & & \cdots & & 0 \\
0 & & 0 & & \lambda_2\max\{0,2X_{2,3}\} & & \cdots & & 0 \\
\vdots & & \vdots & & & & & & \vdots \\

0 & & 0 & & 0 & & \cdots & & \lambda_{n-1} \max\{0,2X_{n-1,n}\} \\
0 & & 0 & & 0 & & \cdots & & 0
\end{bmatrix}
\end{equation}

Thanks for your input, and yes the gradient matrix is sparse. Only (n-1) non-zero elements.

[rajatsanya...@gmail.com]

Sir, I am using this with spectrahedronfactory.

[rajatsanya...@gmail.com]

I am getting the following error-

----------------------------------x---------------------x-------------------
Error using trace (line 13)
Matrix must be square.

Error in Prog>@(X)trace(B*X) (line 11)
problem.cost = @(X) trace(B*X);

Error in getCost (line 38)
cost = problem.cost(x);

Error in getCostGrad (line 56)
cost = getCost(problem, x, storedb, key);

Error in trustregions (line 388)
[fx, fgradx] = getCostGrad(problem, x, storedb, key);

Error in Prog (line 14)
[u0, xcost, info, options] = trustregions(problem);
------------------------x-------------------x-----------------------------

Do you have any idea why I am getting this error? By the way, I checked my input matrix dimension, and it is a square symmetric matrix.

The code is as follows-

---------------------------x-----------------x-----------------------------
function [u0,xcost,info,options] = Prog(B)
% Input:
% B - Symmetric matrix
%% Initialize few parameters
n = size(B,1);
%% MANOPT
% Create the problem structure.
manifold = spectrahedronfactory(n,1);
problem.M = manifold;
% Define the problem cost function and its Euclidean gradient.
problem.cost = @(X) trace(B*X);
problem.egrad = @(X) B;
% Solve.
[u0, xcost, info, options] = trustregions(problem);
end
---------------------------x------------x------------------------------

[Nicolas Boumal]

Hello,

Points in the spectrahedron are indeed square psd matrices, but they are represented using tall-and-skinny matrices (so, rectangular, with more rows than columns), such that X = YY'.

So your cost function, gradient, etc. are supposed to receive Y as input, not X. You cost would thus be trace(B*Y*Y'), or, more efficiently written, trace(Y'*B*Y) (which can be further simplified, but this is particularly readable). The gradient is 2*B*Y.

See this example:

http://www.manopt.org/reference/examples/maxcut.html

and look in the code for the maxcut_fixedrank function.

Hope this helps,

Nicolas

[rajatsanya...@gmail.com]

Thank you again sir.