4
$\begingroup$

I compute the eigenvectors and eigenvalues of the same matrix in its dense and sparse version. The eigenvalues are the same to great precision, but the eigenvectors are not: when I compute the eigenvalue of the solution eigenvector, I am losing precision around the 10th decimal place. This has serious consequences for the next computation with that eigenvector (computation of entropy for those who know about quantum information physics), it gives a totally different result when computing with dense or sparse matrices.

dim = 10;

Do[
   If[n == 1,
   ox = {{0., 1}, {1, 0}};
   oz = {{1, 0.}, {0, -1}};
   id = IdentityMatrix[2],

   ox = SparseArray[ox];
   oz = SparseArray[oz];
   id = SparseArray[id];];



   \[Lambda] = 0.1; 
   hx = 0.; hz = 0.; h = 0.;
   Do[
      If[i == 1, a = ox; c = oz, a = id; c = id];
      hxpart = a; hzpart = c;

      Do[
         If[j == i || j == i + 1, b = ox, b = id;];
         hxpart = KroneckerProduct[hxpart, b];
         If[j == i, d = oz, d = id;];
         hzpart = KroneckerProduct[hzpart, d];
      , {j, 2, dim}];
      If[i == dim, hxpart = 0];
      hx = hx + hxpart; hz = hz + hzpart;
   , {i, 1, dim}];

   h = hx + \[Lambda] hz;
   nordre = Ordering[Eigensystem[h, 2][[1]], 1][[1]];
   gsreala = Chop[Eigensystem[h, 2][[2, nordre]]];
   Print[NumberForm[Eigenvalues[h, 4], 15]];
   Print[NumberForm[Conjugate[gsreala].h.gsreala, 15]];

   (*Entropy computation*)

   gsm = Table[gsreala[[2^(dim/2)*i - (2^(dim/2) - 1) + j - 1]], {i, 1, 
   2^(dim/2)}, {j, 1, 2^(dim/2)}];

   singvalues = SingularValueList[gsm];
   dimsing = Length[singvalues];
   eigenvalues = Table[singvalues[[i]]^2, {i, 1, dimsing}];
   epsilon = 10^-10;
   H[x_] := If[Norm[x] > epsilon, -x*Log[2, x], 0.];
   entreal = N[Sum[H[eigenvalues[[i]]], {i, 1, dimsing}]];
   Print[entreal];

, {n, 1, 2}];

These are the results I obtain:

Dense matrix:

Four smallest eigenvalues: {9.03002193787516,-9.03002193787516,-9.03002193767716,9.03002193767714}

Eigenvalue of the smallest eigenvector: -9.03002193787517

Entropy: 1.00001

Sparse matrix:

Four smallest eigenvalues:{-9.03002193787516,9.03002193787516,9.03002193767717,-9.03002193767716}

Eigenvalue of the smallest eigenvector: -9.03002193781702

Entropy: 0.26277

So I am having a loss of precision on the computation of the eigenvector, that ultimately produces a totally different result in the entropy.

Improve this question
$\endgroup$
2
  • $\begingroup$ If you're going to be totally paranoid about this, I would suggest attempting to evaluate to arbitrary precision; that means avoiding machine precision expressions like 0. in your code, and using the second argument of N[], judiciously, with perhaps slightly more precision than what you need. $\endgroup$
    – J. M.'s missing motivation ♦
    Commented Jul 8, 2016 at 14:04
  • $\begingroup$ Thanks, this solved the problem! Putting a precision inside N[ ] of 3 or higher gives me the correct result. $\endgroup$
    – user3154647
    Commented Jul 8, 2016 at 15:03

1 Answer 1

Reset to default
4
$\begingroup$

Specify

Eigensystem[h, 2, Method -> "Direct"]

for both cases. That sparse array may go to a different solver (iterative) than the dense method. You can find more options in the documentation of Eigensystem under the options section.

Improve this answer
$\endgroup$
6
  • $\begingroup$ If memory serves, Arnoldi is indeed performed by default when the matrix is sparse and only the first few or last few eigenpairs are asked for. I would not recommend the use of "Direct" in general. $\endgroup$
    – J. M.'s missing motivation ♦
    Commented Jul 8, 2016 at 13:59
  • $\begingroup$ Thanks user21 for the answer, your method does give the correct result, but it slows the sparse matrix computation to a velocity comparable to the dense matrix computation. What I would like (if possible) is to obtain a good result without losing velocity :) . The problem comes only with the 10th decimal, but as you can see then the entropy computation gives a totally different result which is not acceptable. $\endgroup$
    – user3154647
    Commented Jul 8, 2016 at 14:04
  • $\begingroup$ @user, that's why I said it's not recommended in general; the dense case uses QR by default, and here the setting "Direct" forces the use of QR on your sparse matrix. $\endgroup$
    – J. M.'s missing motivation ♦
    Commented Jul 8, 2016 at 14:10
  • $\begingroup$ @user3154647, there is no general (iterative) solution for this. Try the methods in the documentation, FEAST in a band (if your matrix is Hermitian), set "Arnoldi" with a tolerance. Did you have a look at the documentation? $\endgroup$
    – user21
    Commented Jul 8, 2016 at 14:22
  • 2
    $\begingroup$ I think you have pairs of almost equal eigenvalues. In this case, you may be able to find an eigenspace accurately, but orientation within the eigenspace is arbitrary. Consequently, the two eigenvector solutions may both be correct. $\endgroup$
    – mikado
    Commented Jul 8, 2016 at 19:33

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.