Basis for unstable manifold of a matrix - Mathematica Stack Exchange most recent 30 from mathematica.stackexchange.com 2019-11-14T23:39:21Z https://mathematica.stackexchange.com/feeds/question/179645 https://creativecommons.org/licenses/by-sa/4.0/rdf https://mathematica.stackexchange.com/q/179645 3 Basis for unstable manifold of a matrix erfink https://mathematica.stackexchange.com/users/43617 2018-08-07T18:40:41Z 2018-08-07T23:13:22Z <p>Given a square matrix <code>A</code>, how can I generate a basis for the <a href="http://mathworld.wolfram.com/GeneralizedEigenvector.html" rel="nofollow noreferrer"><em>generalized eigenspace</em></a> corresponding to all eigenvectors $\lambda_i$ such that $\vert \lambda_i \vert &gt; 1$? I.e., if $A$ is $n \times n$ and acts on $\mathbb{R}^n$, then $\mathbb{R}^n = E^u \oplus E^1 \oplus E^s$, corresponding to the action of eigenvalues of modulus greater than / less than / equal to 1. </p> <p>If <code>A</code> is diagonalizable, this can easily be accomplished with the <code>Eigensystem</code> command:</p> <pre><code>{evals, evecs}=Eigensystem[A]; Take[evecs, Length@Select[evals, Abs[#] &gt; 1 &amp;]] </code></pre> <p>This approach relies on the <code>Eigenvalue</code> command returning the eigenvalues (and hence eigenvectors) in decreasing order of their moduli. Also note that MMA returns the eigenvalues with repetition according to their algebraic multiplicity, and that the eigenvectors are returned as "row vectors" (technically, a list of vectors).</p> <p>My problem arises when <code>A</code> is <strong>not</strong> diagonalizable. The natural function to look at would be <code>JordanDecomposition</code>, which retuns the generalized eigenvectors (N.B. as columns in a matrix!) and the Jordan Canonical Form. However, the ordering is a bit...abstract. As best as I've been able to tell from taking the <code>JordanDecomposition</code> of random matrices, the ordering of the Jordan blocks seems to be some combination of the order in which roots are returned by the <code>Root</code> command for irreducible factors of the characteristic polynomial, the ordinary <code>Greater</code> sorting, the size of the Jordan block, and whether it corresponds to a complex eigenvalue. </p> <p>The kludge that I've been able to come up with is the following:</p> <pre><code>{s, j} = JordanDecomposition[A]; bigevals=Union@Select[Eigenvalues[A], Abs[#] &gt; 1 &amp;]; cols= Flatten[ Position[Diagonal[j], #]&amp; /@ bigevals]; s[[All, cols]] </code></pre> <p>This seems to return the appropriate generalized eigenvectors. </p> <h3>The Questions</h3> <ol> <li>Surely there's a better way? The above kludge seems...kludgy. Perhaps something like <code>NullSpace[MatrixPower[A-lambda IdentityMatrix[First@Dimensions@A], n]]</code> upon detecting that an eigenvalue <code>lambda</code> has algebraic multiplicity <code>n</code>, i.e., going directly to the definition?</li> <li>For my application, the matrix <code>A</code> will only ever be real valued. Thus all complex eigenvalues occur as complex conjugates (perhaps with multiplicity!). I would prefer to have a real basis. For example, if <code>A={{0,2},{-2,0}</code> then the eigenvalues are $\lambda_{\pm} = \pm 2 i$. Then $\left\{ \begin{bmatrix} 1 \\ 0 \end{bmatrix} , \begin{bmatrix} 0 \\ 1 \end{bmatrix} \right\}$ is a basis for the unstable manifold of $A$. <code>Nullspace[MatrixPower[A,2]+4 IdentityMatrix]</code> gives the correct result here. Note that this question is a tertiary concern at best; there are other ways of dealing with this in "post-production."</li> </ol> https://mathematica.stackexchange.com/questions/179645/-/179662#179662 5 Answer by Hector for Basis for unstable manifold of a matrix Hector https://mathematica.stackexchange.com/users/8803 2018-08-07T23:13:22Z 2018-08-07T23:13:22Z <p>Remember that the diagonal elements of the Jordan normal form are the eigenvalues. There is no need to compute the eigenvalues and the decomposition. Instead, the following code creates a picker list from the diagonal elements.</p> <pre><code>A = 1/8 {{27, 48, 81}, {-6, 0, 0}, {1, 0, 3}}; {s, j} = JordanDecomposition[A]; Pick[Transpose[s], Map[Abs[N[#]] &gt; 1 &amp;, Transpose[j, {1, 1}]]] </code></pre> <p>From there, you might want to use <code>Orthogonalize</code> …</p>