Computing signature - MathOverflow most recent 30 from http://mathoverflow.net 2013-06-18T06:06:06Z http://mathoverflow.net/feeds/question/107593 http://www.creativecommons.org/licenses/by-nc/2.5/rdf http://mathoverflow.net/questions/107593/computing-signature Computing signature Igor Rivin 2012-09-19T16:57:31Z 2012-09-20T16:44:36Z <p>I have a feeling that this might have already been asked, but can't find the question. Anyway, the question is: given a symmetric $n\times n$ matrix, is there a faster way to compute its signature than computing the full eigen-decomposition? (the gt tag is because my application is in computing signatures of four-manifolds, and some others might want to do the same computation...)</p> http://mathoverflow.net/questions/107593/computing-signature/107595#107595 Answer by Terry Loring for Computing signature Terry Loring 2012-09-19T17:21:54Z 2012-09-19T22:17:05Z <p>This seems to be an $O(n^3)$ problem, so if your matrix is smaller than 100 by 100 you might as well use an eigensolver. </p> <p>For bigger matrices you can use Newton's method (see Higham's papers or book) to get a real almost orthogonal matrix and take the trace and round off. To get near the orthogonal in the polar decomposition iteratively replace $A$ by the average of $A$ with the transpose of the inverse of $A$.</p> <p>That is, let <code>$A_0 = A$</code> and <code>$A_{n+1} = \left( \frac{1}{2}A + \frac{1}{2}\left( A^{\mathrm{T}} \right)^{-1} \right)$</code> and then the signature is <code>$\lim_{n \rightarrow \infty} \mathrm{Trace} \left( A_N \right )$</code></p> <p>If you have a priori information on the norm of $A$ and the size of the spectral gap at zero, you can figure how many iterations to do easily. If not many, this is very fast.</p> <p>Nicholas J. Higham and Pythagoras Papadimitriou, A parallel algorithm for computing the polar decomposition, Elsevier, Parallel Computing, 20, (1994).</p> http://mathoverflow.net/questions/107593/computing-signature/107607#107607 Answer by Ryan Budney for Computing signature Ryan Budney 2012-09-19T18:57:25Z 2012-09-20T06:15:45Z <p>In the software <a href="http://regina.sourceforge.net/" rel="nofollow">Regina</a> what I do is apply <a href="http://en.wikipedia.org/wiki/Descartes%27_rule_of_signs" rel="nofollow">Descartes Rule of Signs</a> on the characteristic polynomial. The algorithm for computing the determinant we use has complexity of the order of $n^4$, but it seems to work quite effectively for the signature computations we've been doing so far. </p> <p>The primary reason we use these techniques is they're fairly easy to implement. Also, we're worried about rigour so we want to avoid floating-point data types that could have round-off errors. We use arbitrary precision integers for the computations, via the GNU Multi-Precision libraries. </p> http://mathoverflow.net/questions/107593/computing-signature/107608#107608 Answer by Ian Agol for Computing signature Ian Agol 2012-09-19T19:02:17Z 2012-09-20T16:44:36Z <p>As Fernando Muro points out in the comments, <a href="http://en.wikipedia.org/wiki/Sylvester_law" rel="nofollow">Sylvester's law of inertia</a> is probably the easiest way to determine the signature. You diagonalize the symmetric matrix by the <a href="http://en.wikipedia.org/wiki/Gram-schmidt" rel="nofollow">Gram-Schmidt process</a>. This is essentially as easy as Gaussian elimination. </p> <p>A symmetric matrix is equivalent to a bilinear form, with basis <code>$\{v_1, \ldots, v_n\}$</code>, and matrix entries inner products $\langle v_i,v_j\rangle, i,j\leq n$. If $\langle v_1,v_1\rangle \neq 0$, then replace $v_i$ with $v_i-v_1 \langle v_1,v_i\rangle/\langle v_1,v_1\rangle$, and proceed by induction (this has the effect of performing the same row and column operations, so as to keep the matrix symmetric). The only variation on the standard Gram-Schmidt is what to do when $\langle v_1,v_1\rangle=0$, which clearly means the form is indefinite. If $\langle v_1, v_i\rangle=0$ for all $i\leq n$, then you can proceed by induction. Otherwise, by permuting indices, you may assume that $\langle v_1,v_2\rangle \neq 0$. Replace $v_1$ with $v_2+tv_1$, such that $\langle v_2,v_2\rangle+2t\langle v_1,v_2\rangle\neq 0$, and proceed as with the usual Gram-Schmidt. Once you have diagonalized in this way, Sylvester's law implies that the signature of the diagonal matrix is the same as that of the original one. </p> <p><strong>Note added:</strong> Since you claim you're interested in computing the signatures of 4-manifolds, there's some special points one might be able to take advantage of. Certainly then the form will be integral unimodular, so presumably your symmetric matrix has integer entries and determinant $\pm 1$. If the 4-manifolds are smooth, then by Donaldson's theorem, if the intersection form is (positive) definite, then it is integrally equivalent to the diagonal form <code>$1\oplus 1\oplus \cdots\oplus 1=1^n$</code>. Otherwise, if the form is indefinite, then it is either odd and integrally equivalent to a diagonal form $1^k\oplus (-1)^{n-k}$, or it is even and there is a unique equivalence class for a given $n$ and signature. One can then perform the diagonalization process over $\mathbb{Z}$ using some version of the Euclidean algorithm, except in the even case (in which case there is a standard tridiagonal integral form). Actually, in the even case one could take a sum $\oplus \pm 1$ to make the form odd, and then diagonalize over the integers (a similar trick taking a sum with a hyperbolic form makes things integrally diagonalizable in the general case if you are also considering non-smooth manifolds which might have non-diagonalizable forms). See Milnor-Husemoller or <a href="http://books.google.com/books?id=upYwZ6cQumoC&amp;lpg=PP1&amp;pg=PP1#v=onepage&amp;q&amp;f=false" rel="nofollow">Conway-Sloane</a> for more details. </p> http://mathoverflow.net/questions/107593/computing-signature/107619#107619 Answer by Denis Serre for Computing signature Denis Serre 2012-09-19T20:25:08Z 2012-09-19T20:25:08Z <p>If $n$ is not too big, compute the characteristic polynomial $P$ by the method of Leverrier, modified by Preparata &amp; Sarwate; this takes $O(n^{3.5})$ operations, see the 2nd edition of my book <em>Matrices</em> (GTM <strong>216</strong>, Springer-Verlag). Then apply the Routh-Hurwitz algorithm to find the number of ots positive roots. This second step, and variants, are very well described by Gantmacher in the 2nd volume of his book on matrices.</p> <p>Notice that this takes finitely many operations. Actually the whole book of Gantmacher is motivated by stability issues in dynamical systems. Thus its main purpose is to check whether a given matrix has its spectrum in a given half of the complex plane.</p> http://mathoverflow.net/questions/107593/computing-signature/107622#107622 Answer by suVRit for Computing signature suVRit 2012-09-19T21:24:33Z 2012-09-19T21:24:33Z <p>Here is a textbook solution that runs in time $O(n^3)$ (assuming that Gaussian Elimination can be done in that time, which might not be immediate unless one takes steps to rule out intermediate coefficient growth). The first part of the statement below is merely LU decomposition. The last sentence is what you're looking for.</p> <blockquote> <p><strong>Theorem</strong> (<a href="http://books.google.de/books?id=FtC9rHW7ORMC&amp;lpg=PA193&amp;ots=nDB3pdYYTz&amp;dq=compute%2520inertia%2520of%2520matrix&amp;pg=PA193#v=onepage&amp;q=compute%2520inertia%2520of%2520matrix&amp;f=false" rel="nofollow">Stewart</a>) Let $A$ be symmetric and suppose that the leading principal submatrices of $A$ are nonsingular. Then, there exists a unit lower triangular matrix $L$ and upper triangular matrix $U$ such that $A=LU$. Moreover, if $D_U=\text{diag}(u_{11},\ldots,u_{nn})$, then $\text{Inertia}(D_U)=\text{Inertia}(A)$.</p> </blockquote> <p>Unfortunately, if one pivots (which in practice one must), the above rule breaks. In this case, I would recommend first tridiagonalizing the matrix using orthgonal transformations---this will cost $O(n^3)$. Then, you can compute the eigenvalues of the tridiagonal matrix in time $O(n^2)$ (pseudo-polynomial time, though), and then you'll have the inertia.</p> http://mathoverflow.net/questions/107593/computing-signature/107675#107675 Answer by Bazin for Computing signature Bazin 2012-09-20T12:14:27Z 2012-09-20T12:14:27Z <p>May I attract your attention on the following formula, where $B$ is a $n\times n$ non-singular symmetric real-valued matrix $$\int_{\mathbb R^n}e^{-2i\pi x\cdot \xi} e^{i\pi Bx\cdot x}dx=e^{\frac{i\pi}{4}\text{signature } B}\vert\det B\vert^{-1/2}e^{-i\pi B^{-1}\xi\cdot \xi},$$ so that $e^{\frac{i\pi}{4}\text{signature } B}=\vert\det B\vert^{1/2} \int_{\mathbb R^n} e^{i\pi Bx\cdot x}dx.$ Of course the Fresnel-type integral on the rhs of last formula is not absolutely converging but, as an oscillatory integral, will become so after integration by parts.</p>