Get positive eigenvalues of density matrix

I would like to compute the eigenvector of the normalised Laplacian of a graph, which I think has properties of density matrix, that is trace = 1 and eigenvalues between 0 and 1. However when I compute the eigenvalues some are negative. Is there a problem with my code?

My initial adjacency matrix is:

A1 =  
     0     1     1     1     0     0     0     0
     1     0     1     1     0     0     0     0
     1     1     0     1     0     0     0     0
     1     1     1     0     0     0     0     0
     0     0     0     0     0     1     1     1
     0     0     0     0     1     0     1     1
     0     0     0     0     1     1     0     1
     0     0     0     0     1     1     1     0

My code then is:

[n,n] = size(A1);
D = zeros(n,1);
D(:,1) = sum(A1);
Norm1 = sum(D(:, 1));
L1 = (D(:, 1) - A1)/Norm1;

[V,D] = eig(L1)

The resulting eigenvalues are:

     -0.1250
      0.0417
      0.0417
      0.0417
      0.0417
      0.0417
      0.0417
      0.8750

The first one is negative, but since L1 is a density matrix I think there is a problem with my code.

[James Tursa]

How was this matrix built? Where did the values come from? What is this matrix used for downstream in your code?

[Walter Roberson]

I just tried with symbolic toolbox. You are not just hitting round-off error. If those entries are to be interpreted as 1/8 and 1/12 then the first eigenvalue is -1/8; the next is 7/8; the rest are 1/24

[David Goodmanson]

the eigenvector for the negative eigenvalue is interesting: const x [1 1 1 1 -1 -1 -1 -1]', which is related to the fact that L has the block form

[A B
 B A]

[RM]

Thank you for the feedback! I have edited the question so that it is more clear. Perhaps I am using matrix algebra wrong and not getting a density matrix?

[David Goodmanson]

Hi RM, I believe the problem is that the initial vector D is correct, but when it comes to defining the Laplacian matrix, D is supposed to be not a vector but a diagonal matrix with the elements of D on the diagonal, zero everywhere else. The diag command can create the one from the other. Defining L1 as

D = sum(A1);      % sum down columns of A1 to make a row vector
Dmat = diag(D);   % diagonal matrix
L1 = Dmat - A1;

does give a matrix whose eigenvalues are all positive or zero. Normalization is another question, since the Wikipedia definition, anyway, appears to differ from what you have.

Incidentally, the expression

D(:, 1) - A1

subtracts a matrix from a column vector. Before Matlab 2016b this would have given an error message. Now, because of the explicit expansion feature it does not (as long as the number of rows of D matches up with the number of rows in A1). It's a good feature but can lead to a whole new territory of unintended behavior.

[RM]

Perfect! Thank you very much. This solves my problem.

Answer 1

[John D'Errico]

You cannot compute eigenvalues for a matrix subject to a constraint. I'm sorry, but that simply does not make mathematical sense.

In fact, even if I replace the elements of your matrix with 1/12, the eigenvalues still have one negative eigenvalue.

So lets see if it is just a precision thing? I've turned your matrix into a symbolic one, assuming that you really intended 1/12 when you show 0.0833.

Ms =
[  1/8, 1/12, 1/12, 1/12,  1/8,  1/8,  1/8,  1/8]
[ 1/12,  1/8, 1/12, 1/12,  1/8,  1/8,  1/8,  1/8]
[ 1/12, 1/12,  1/8, 1/12,  1/8,  1/8,  1/8,  1/8]
[ 1/12, 1/12, 1/12,  1/8,  1/8,  1/8,  1/8,  1/8]
[  1/8,  1/8,  1/8,  1/8,  1/8, 1/12, 1/12, 1/12]
[  1/8,  1/8,  1/8,  1/8, 1/12,  1/8, 1/12, 1/12]
[  1/8,  1/8,  1/8,  1/8, 1/12, 1/12,  1/8, 1/12]
[  1/8,  1/8,  1/8,  1/8, 1/12, 1/12, 1/12,  1/8]

eig(Ms)
ans =
 -1/8
 1/24
 1/24
 1/24
 1/24
 1/24
 1/24
  7/8

Lets see if this makes sense. I'll compute the eigenvector for that eigenvalue.

[V,D] = eig(M);
V(:,1)
ans =
      0.35355
      0.35355
      0.35355
      0.35355
     -0.35355
     -0.35355
     -0.35355
     -0.35355

M*V(:,1)
ans =
    -0.044194
    -0.044194
    -0.044194
    -0.044194
     0.044194
     0.044194
     0.044194
     0.044194

So, in fact, that negative eigenvalue makes complete sense. If I multiply your matrix by the vector V(:,1), then I get -1/8*V(:,1) back. That is what an eigenvalue/eigenvector pair means.

I would suggest you are in error when you claim that ALL of the eigenvalues of this matrix must be positive.

[RM]

Thank you for your answer. I wrote the question wrong I think. I may be computing the Laplacian of the adjacency matrix wrong and therefore having negative eigenvalues.

[John D'Errico]

That could be. This matrix clearly has a negative eigenvalue.