This may be a really dumb question, but here goes: is there any algorithm to compute the signature of a quadratic form (or a symmetric matrix, if you prefer) more efficient (asymptotically or otherwise) than actually computing the eigenvalues?
-
1$\begingroup$ count sign changes in the determinants of the principal minors? (This fails if one of them vanishes, but you can first compute the rank, then do a random coordinate change so you don't get any zeros before that.) $\endgroup$– Noam D. ElkiesCommented Apr 27, 2017 at 0:34
-
$\begingroup$ @NoamD.Elkies Yes, true enough, but I am not sure if this is more efficient than computing the eigenvalues (in fact, I am pretty sure it is not, though I have been wrong before :)) $\endgroup$– Igor RivinCommented Apr 27, 2017 at 0:51
-
$\begingroup$ I don't know. Might depend on whether your matrix has floating-point or exact entries. You might also consider Dodgson condensation, which is inefficient for a single determinant but possibly competitive for this question because it gives you all contiguous minors along the way. $\endgroup$– Noam D. ElkiesCommented Apr 27, 2017 at 0:57
-
1$\begingroup$ If the entries are floating-point, determinants are faster, provided that you compute them all at the same time with a $LU$ factorization (or better $LDL^T$, as suggested in ShakeBaby's answer) instead of one by one. In fact, forget about the determinants, just considering the signs of the diagonal entries in the factorization is enough. $\endgroup$– Federico PoloniCommented Apr 27, 2017 at 6:10
-
$\begingroup$ Note that if the entries are floating-point, in addition to speed another issue is stability. The eigendecomposition (on a symmetric matrix) is probably going to be more stable in cases where you have "almost zero" entries. $\endgroup$– Federico PoloniCommented Apr 27, 2017 at 6:32
3 Answers
Gauss reduction gives you the answer. It writes, quite fast, the quadratic form $q$ as a sum $$\sum_ja_j\ell_j(x)^2$$ where the $\ell_j$'s are independent linear forms. The number of squares gives you the rank of $q$. The signs of the coefficients $a_j$ gives you the signature. I teach that in my undergraduate course in Algebra.
Remark that you cannot calculate the eigenvalues, at least in close form, because this is computing the roots of quite a general polynomial, and this is impossible in dimension $\ge5$.
-
$\begingroup$ Denis, could you please comment on my answer? Your method is what I know from books on integral quadratic forms, and I associate it with the name Hermite. Really just repeated completing the square(s) until all terms are exhausted, and not difficult. Then in 2015, on MSE, I came across the method in my answer, and started asking about it math.stackexchange.com/questions/1388421/… My feeling is that they are really the same method, minor difference in finding either $P$ or $P^{-1}$ first in $P^T AP=D$ $\endgroup$ Commented Apr 27, 2017 at 16:52
-
1$\begingroup$ @Will. There is one difficulty, in that you can at some point have to reduce a quadratic form with no square at all. But then the form represents zero and you just look for a hyperbolic plane. I don't know the terminology Hermite (though he was a French mathematician). Perhaps the method has different names in different countries. $\endgroup$ Commented Apr 27, 2017 at 20:42
-
$\begingroup$ Denis, I think you are right. I am looking at my quadratic forms books, it seems Hermite reduction is far more specific thing. I saw your sum as part of that and incorrectly put his name on this general method. $\endgroup$ Commented Apr 28, 2017 at 0:45
Thursday morning. Noam has suggested that this is equivalent to performing Gram-Schmidt without normalization. That would explain why I could not find any explicit point where anyone wrote "Here is a way to reverse Hermite's type of method." I'm going to try some 2 by 2 and 3 by 3 examples, see if I understand.
Here are explicit example(s) as links, in this first one the form is indefinite;
Igor, there is an easy algorithm that creates $P^T A P = D$ with $D$ diagonal, $\det P = 1$ and the elements of $P$ in the same field as that needed for $A.$
I can describe the way I do it. Let $A_0 = A,$ then $A_{j+1} = P_j^T A_j P_j,$ where $P_j$ is one of three types:
(I) the identity matrix, except for the value $t$ at position $i,j$ in the upper triangle
(II) the identity matrix, except $p_{ii} = 0,$ $p_{jj} = 0,$ $p_{ij} = 1,$ $p_{ji} = -1,$
(III) the identity matrix, except for the fixed value $1$ at a position in the lower triangle.
Oh, after doing several of these, I realized that a bunch of type (I) matrices with the extra off diagonal elements in the same row can be combined into one matrix, as such matrices commute with each other.
I had never seen it before, I asked about it at http://math.stackexchange.com/questions/1388421/reference-for-linear-algebra-books-that-teach-reverse-hermite-method-for-symmetr
In the most favorable cases, $P$ is also upper triangular. Not guaranteed.
-
1$\begingroup$ Very cool! This is clearly superior over integer or rationals, a little harder to tell with floating point! $\endgroup$ Commented Apr 27, 2017 at 2:12
-
1$\begingroup$ Is this basically Gram-Schmidt? (Gram-Schmidt does find an equivalent quadratic form that's diagonal.) $\endgroup$ Commented Apr 27, 2017 at 4:18
-
$\begingroup$ @Noam I don't think so. We don't ever find out the eigenvalues this way. I think it is the same as Denis Serre's answer, which he calls Gauss reduction. That way, you gradually produce $D$ and the rows of some $Q$ with $Q^T DQ = A.$ This method, we gradually produce $P$ and get $P^T AP = D.$ So, the real difference is just $P = Q^{-1},$ unless different choices led to different diagonal $D$ $\endgroup$ Commented Apr 27, 2017 at 16:45
-
1$\begingroup$ Gram-Schmidt doesn't find eigenvalues either. $\endgroup$ Commented Apr 27, 2017 at 17:08
-
1$\begingroup$ G-S doesn't require square roots either unless you insist that your basis be orthonormal, not just orthogonal. $\endgroup$ Commented Apr 27, 2017 at 17:35
There's a large literature on "inertia revealing factorizations" for real or complex matrices. Usually one uses MA57 of the HSL implementation to compute the signature, which does an sparse $LDL^T$ factorization: http://www.hsl.rl.ac.uk/catalogue/ma57.html