# Fast method to calculating signature of a matrix

For calculation of the topological properties of a hamiltonian, sometimes we need the signature of that matrix. This means we only need number of positive eigenvalues. One simple way is to first calculate the eigenvalues and then find difference of number of positive and negative eigenvalues. for example,

H = RandomReal[{-1, 1}, {100, 100}] // (# + #\[ConjugateTranspose]) &;
Total[If[# > 0, 1, -1] & /@ Eigenvalues[H]]


I am interested to know i) is there any function in Mathematica calculate signature ii) is there any fast algorithm within Mathematica to do so fast?

• You mean Sign? This function is Listable i.e. Sign[{-3, 4, 0}] works BTW. Jun 30, 2020 at 3:44
• Total@Sign@Eigenvalues[H] Jun 30, 2020 at 4:49
• The problem is not in sign. The problem is how to get rids of the eigenvalue calculation. Jun 30, 2020 at 5:37
• @Rasoul-Ghadimi First of all, there's no function in Mathematica that directly computes the signature of a matrix. Furthermore, all the definitions for signature I can find lead to eigenvalues. Unless there're definitions that avoid eigenvalues (e.g. ask on Mathematics SE), I don't think there's any faster way to do the calculation. Jun 30, 2020 at 6:52
• Jun 30, 2020 at 17:19

SeedRandom;
H = RandomReal[{-1, 1}, {100, 100}] // (# + #\[ConjugateTranspose]) &;
Short[H, 3]

{{1.269557961,-0.6813805661,<<96>>,-0.3068229219,-0.9273213366},<<98>>,{-<<19>>,<<99>>}}


# Comparison

Eigenvalues is the fastest:

Total@Sign@Eigenvalues[H] // RepeatedTiming

{0.00044, 2}


Descartes' Sign Rule on the characteristic polynomial:

Total[1 - Ratios[Sign@CoefficientList[
CharacteristicPolynomial[H, x], x]
]] - MatrixRank[H] // RepeatedTiming

Total[RealAbs@Differences@Sign@CoefficientList[
CharacteristicPolynomial[H, x], x]] -
MatrixRank[H] // RepeatedTiming

2 (Length[Split@Sign@CoefficientList[
CharacteristicPolynomial[H, x], x]
] - 1) - MatrixRank[H] // RepeatedTiming

{0.0030, 2}
{0.0031, 2}
{0.0032, 2}


$$LDL^\top$$ decomposition (idea, code from @J.M.'stechnicaldifficulties):

Total@Sign@LDLT[H][] // RepeatedTiming

{0.0036, 2}


## Benchmark

Needs["GeneralUtilities"];
BenchmarkPlot[{sig1Eigen, sig2Poly1Ratio, sig2Poly2Diff,
sig2Poly3Split, sig3LDLT},
n \[Function] StatisticsLibraryVectorToSymmetricMatrix[
#[[n + 1 ;;]], #[[;; n]], n
] &@RandomReal[{-1, 1}, Binomial[n + 1, 2]],
"IncludeFits" -> True, TimeConstraint -> 100]


Definitions are the same as above. (Here real symmetric matrices are generated with an undocumented function seen here, which is a bit faster.) Result: It seems that all of these methods are in $$\mathcal{O}(n^3)$$ (Ratios one should be the same), but Eigenvalues one has the smallest coefficient.

## Note: Bug?

One strange thing is that CountRoots doesn't give correct answers here, nor does Reduce or Solve. Is this a bug?

CountRoots[CharacteristicPolynomial[H, x], {x, 0, \[Infinity]}]

21 (* Should be 51 *)

Reduce[CharacteristicPolynomial[H, x] == 0 && x > 0, x] // Length (* or Solve *)

15 (* Should be 51 *)

Solve[{CharacteristicPolynomial[H, x] == 0, x > 0}, x,
Complexes] // Length

51 (* Correct *)

Reduce[CharacteristicPolynomial[H, x] == 0, x] // Length (* or Solve *)

100 (* Correct *)


Otherwise, 2 CountRoots[CharacteristicPolynomial[H, x], {x, 0, \[Infinity]}] - MatrixRank[H] can be used.

### Update

This must be a bug -- Some real-valued roots are considered by CountRoots to have "invisible" imaginary parts, with some of them "$$| Im(x) | \gt 1$$"!

CountRoots[
CharacteristicPolynomial[H, x], {x, -I, 100 + I}] // AbsoluteTiming

{44.9898317, 33} (* Still incorrect *)

CountRoots[
CharacteristicPolynomial[H, x], {x, -5 I,
100 + 5 I}] // AbsoluteTiming

{127.967137, 51} (* Correct, but very slow *)


Solve has slighter peoblem. I think it's due to the machine precision:

Solve[CharacteristicPolynomial[H, x] == 0 && Re[x] > 0, x] //
Length // AbsoluteTiming

{0.0206228, 51} (* Correct *)


### Update 2

According to @MichaelE2, this is because the machine precision is not enough for the large coefficients and high degree of the characteristic polynomial.

• It would be interesting to know how the methods scale with the matrix size $n$. What method is the asymptotic winner for $n\rightarrow\infty$? Jun 30, 2020 at 21:25
• @yarchik, I've updated the answer. Jul 1, 2020 at 5:54
• Nice analysis. I interpret it as the rank-determination is the bottleneck with $\mathcal{O}(n^\omega)$. Jul 1, 2020 at 6:54
• @yarchik Yes, and also the characteristic polynomial, maybe with $\mathcal{O}\left( n^\omega \log n \right)$. Jul 1, 2020 at 8:02
• Probably not a bug. Machine-precision, high-degree polynomials are known to be perfidious. Examine Plot[Evaluate@CharacteristicPolynomial[H, x], {x, 0, 10}]`. Jul 1, 2020 at 21:52