# How can I compute a basis for a matrix Lie algebra given a finite set of generators?

Given an arbitrary set of (numerical) square complex matrices $\mathcal{A}=\{A_1,A_2,\cdots,A_m\}$, I am interested in computing the real matrix Lie algebra generated by $\mathcal{A}$, call it $\mathcal{L_\mathcal{A}}$. That is, I would like a basis for $$\mathcal{L_\mathcal{A}} = \mathbb{span_R}\{B:B\in\cup_{k=1}^{\infty}\mathcal{C}_k\}$$ where $\mathcal{C}_k$ is defined recursively as $\mathcal{C_1}=\mathcal{A}$, and $\mathcal{C_{k+1}}=\{[X,Y]:X,Y\in\cup_{j=1}^k\mathcal{C_j}\}$ for $k\geq 1$.

This calculation comes up in (quantum) control theory.

Currently I am using a method found here which searches only through repeated Lie brackets (i.e. ones of the form $[A_{j_1},[A_{j_2},[A_{j_3},\cdots[A_{j_{n-1}},A_{j_n}]\cdots]]]$), and is guaranteed to terminate. However I'm interested to know if there are any other (faster) methods. Perhaps using P. Hall bases? Perhaps a recursive algorithm? My default language at the moment is Matlab.

• I am guessing that your original generators are Hermitian. Is this true? If so, I would imagine the first step would be to compare the eigenspaces of the generators, as commutators are only nonzero when the eigenspaces differ. Feb 1 '12 at 18:02
• @JackPoulson Yes, the A's come from Hamiltonians, and so are skew-Hermitian (not Hermitian because they are multiplied by the i in Schroedinger's equation). I'm not sure I understand why this would be a good first step. Wouldn't calculating the commutators and checking to see if they are non-zero be faster than fiddling with eigenspaces? Feb 1 '12 at 18:15
• For a single level of commutators, probably yes. But there is a combinatorial explosion when you start considering several levels of commutators. I do not know of an algorithm, but usually it is a good idea to exploit as much structure as possible. I would carefully think about whether you knew any other properties that relate your generators as well. Feb 1 '12 at 19:00

On an only somewhat related note, if I were implementing this I would worry about the numerical instability of testing linear dependency. Make sure to use a method for testing independence of new matrices that allows for numerical inaccuracy - maybe comparing the norm of $A - p(A)$ to the norm of $A$, where $p$ is the projection onto the space of matrices you've found before.
• @IanHincks: What I meant was to compare $\Vert A - p(A) \Vert$ to $\Vert A \Vert$, but that wasn't based on any particularly deep thoughts. You'll need to experiment. The numerically best criterion may be to view all matrices as $\mathbb{R}^{n^2}$-vectors and do a sparse SVD of the rectangular $n^2 \times k$ matrix obtained by putting them next to each other, then discarding the "vector" added last if the smallest singular value is very small. But that will be very expensive computationally. First see if you really need it - and if so, maybe do a cheap test first. Feb 2 '12 at 20:37