A Jacobi Algorithm in Phase Space: Diagonalizing (skew-) Hamiltonian and Symplectic Matrices with Dirac-Majorana Matrices
aa r X i v : . [ m a t h - ph ] A ug A Jacobi Algorithm in Phase Space: Diagonalizing (skew-) Hamiltonian andSymplectic Matrices with Dirac-Majorana Matrices
C. Baumgarten
Paul Scherrer Institute, Switzerland ∗ (Dated: September 1, 2020)Jacobi’s method is a well-known algorithm in linear algebra to diagonalize symmetric matricesby successive elementary rotations. We report here about the generalization of these elementaryrotations towards canonical transformations acting in Hamiltonian phase spaces. This generalizationallows to use Jacobi’s method in order to compute eigenvalues and eigenvectors of Hamiltonian (andskew-Hamiltonian) matrices with either purely real or purely imaginary eigenvalues by successiveelementary “decoupling” transformations. The real importance of Einstein’s workwas that he introduced Lorentz transforma-tions as something fundamental in physics –P.A.M. Dirac [1]
I. INTRODUCTION
The problem of eigenvector and eigenvalue computa-tion of (skew-) symmetric, Hamiltonian and symplecticmatrices received considerable attention in the past .Here we describe a method that is entirely base on pureHamiltonian (symplectic) notions. We shall develop thereal Clifford algebra Cl (3 ,
1) from algebraic Hamiltoniansymmetries and demonstrate that the group of linearcanonical transformations in classical phase space anda generalized Lorentz group are isomorphic.Jacobi’s Method is a well known numerical methodthat allows to diagonalize symmetric matrices . Themethod is most easily explained when we start with somearbitrary real symmetric 2 × A : A = (cid:18) a a a a (cid:19) (1)This matrix can be diagonalized by an orthogonal Matrix R T = R − of the form R = (cid:18) c − ss c (cid:19) (2)where c = cos θ und s = sin θ , so that ˜A = R A R − ˜ a = a c + a s − a s c ˜ a = a ( c − s ) + ( a − a ) s c ˜ a = a c + a s + 2 a s c (3)where ˜ a simplifies to:˜ a = a C + a − a S (4) ∗ [email protected] See Refs. [2–8] and references therein. See for instance Sec. 9.2 in Ref. [9]. where S = 2 s c = sin 2 θ and C = c − s = cos 2 θ . Thetransformed matrix ˜A is diagonal, if ˜ a = 0 and hence,if the angle of rotation θ is:0 = a C + a − a Sθ = 12 arctan ( 2 a a − a ) (5)This method allows to diagonalize real symmetric n × n -matrices by successive rotations, if the orthogonal rota-tion matrix R kl is a unit matrix except for the k -th and l -th rows and column. This means that r ii = 1 for i / ∈ l, k and r ll = r kk = c . The only non-vanishing off-diagonalelements are − r lk = r kl = s : R = . . . c − s . . . s c . . . (6)This rotation then allows to zero any selected non-diagonal element a kl = 0. The Jacobi method demandsto chose the order of the successive rotations such thatalways the dominant non-diagonal element a kl is selectedfor the next rotation – until convergence. The successiverotations are collected in a matrix RR = R N , R N − , . . . , R (7)so that ˜A = Diag( α , . . . , α n )= R A R T A = R T ˜A R . (8)Since real symmetric matrices have real eigenvalues, all α i are real. II. PHASE SPACE AND HAMILTONIANMATRICES
The conventional Jacobi method operates on matricesin n -dimensional Euklidean space in the sense that or-thogonal transformations do not affect the n -dimensionalEuklidean norm. In contrast to those geometrical spaces( R n ) phase spaces are dynamical spaces and the mostimportant transformations in phase spaces are canonicaltransformations, which preserve Hamilton’s equations ofmotion.Given a coordinate ψ = ( q , p , . . . , q n , p n ) T in a clas-sical phase space of n degrees of freedom, we define aHamiltonian function H ( ψ ) by the quadratic form H ( ψ ) = 12 ψ T A ψ (9)with a 2 n × n -dimensional real symmetric matrix A .Hamilton’s equations of motion are obtained by the con-straint ˙ H ( ψ ) = n X k =1 ∂ H ∂q k ˙ q k + ∂ H ∂p k ˙ p k = 0 (10)with the general solution˙ q k = ∂ H ∂p k ˙ p k = − ∂ H ∂q k (11)This can be written in vectorial form as˙ H ( ψ ) = ∇ ψ H ( ψ ) · ˙ ψ = 0˙ ψ = γ ∇ ψ H ( ψ )= γ A ψ (12)where the constraint that H ( ψ ) = const leads to:˙ H ( ψ ) = ψ T A γ A ψ = 0 (13)which is true whenever A γ A is skew-symmetric, i.e.whenever γ is skew-symmetric. Typically the so-called symplectic unit matrix (SUM) γ is chosen to either havethe form γ = Diag( η, η, . . . , η ) = n × n ⊗ η (14)with the 2 × η = (cid:18) − (cid:19) (15) It suffices to refer to the quadratic form of a more general Taylorseries of H ( ψ ) since we are only interested in linear transforma-tions. or γ = η ⊗ n × n = (cid:18) n × n − n × n (cid:19) . (16)The former case corresponds to an alternating or-der of canonical coordinates and momenta ( ψ =( q , p , q , p , . . . ) T ), in the latter case one uses an or-dering where first all coordinates appear and then allmomenta: ψ = ( q , q , . . . , q n , p , p , . . . , p n ) T . The dif-ference is mostly notational, but the algorithm that weaim to describe, favours alternate ordering.The evolution of the system is, for A = const given bythe matrix exponential: ψ ( t ) = exp ( H τ ) ψ (0) . (17)The matrix H that appears in the exponent H = γ A (18)is called a Hamiltonian matrix and it obeys the relation H T = γ H γ (19)The simplest case is certainly that of a diagonal matrix A , and the unit matrix is the simplest diagonal matrix.Hence, if we consider this simplest case, then A = ω ,where ω is the (eigen-) frequency. Then Eq. 17 results in ψ ( t ) = exp ( γ ω τ ) ψ (0) = X k ( γ ω τ ) k k ! ψ (0) . (20)which can be splitted in the even and odd powers:exp ( γ ω τ ) = X k ( γ ω τ ) k k ! + X k ( γ ω τ ) k +1 (2 k + 1)!= X k ( − ) k ( ω τ ) k k ! + γ X k ( − ) k ( ω τ ) k +1 (2 k + 1)!= cos ( ω τ ) + γ sin ( ωτ ) . (21)In this derivation we made use of nothing but the factthat γ = − . Hence, the matrix exponential of sym-metric (Hamiltonian) matrices γ a with γ a = 1 yields thecorresponding result:exp ( γ a ω τ ) = X k ( γ a ω τ ) k k ! + X k ( γ a ω τ ) k +1 (2 k + 1)!= X k ( ω τ ) k k ! + γ a X k ( ω τ ) k +1 (2 k + 1)!= cosh ( ω τ ) + γ a sinh ( ωτ ) . (22)It is well-known (and easily proven) in linear Hamilto-nian theory that the matrix exponentials of Hamiltonianmatrices are symplectic and the matrix logarithm of anysymplectic matrix is Hamiltonian [10]. A symplectic ma-trix M represents a linear canonical transformation ofthe phase space coordinates and fulfills the following def-inition: M γ M T = γ , (23)and the corresponding transformation is˜ ψ = M ψ ˙˜ ψ = M H M − ˜ ψ (24)Note that the matrix γ is both, Hamiltonian and sym-plectic.From Eq. 23 one obtains: M − = − γ M T γ , (25)so that the transformed Hamiltonian ˜H is given by ˜H = M H M − . (26)and must be Hamiltonian again, i.e.: ˜H T = ( M H M − ) T = ( M − ) T H T M T = − ( γ M T γ ) T H T M T = − γ T M γ T γ H γ M T = γ ( M H M − ) γ = γ ˜H γ (27)where Eq. 23 has been used. This means that symplec-tic transformations are structure preserving: The formand structure of all matrices and specifically the form ofHamilton’s equations of motion is preserved under sym-plectic similarity transformations. III. THE REAL PAULI ALGEBRA
The canonical pair is the smallest meaningful elementin phase space and accordingly the smallest Hamiltonianmatrix has size 2 × H = (cid:18) h h h h (cid:19) (28)The trace of the product of two matrices holdsTr( A T B ) = Tr( A B T ) (29)Hence, if A is symmetric and B skew-symmetric, thenthe left side of Eq. 29 givesTr( A T B ) = Tr( A B ) (30)and the right sideTr(
A B T ) = − Tr(
A B ) (31) which can only be true, if the trace of such a productvanishes: The trace of the product of some symmetric A and some skew-symmetric matrix B is zero:Tr( A B ) = 0 (32)Hence all Hamiltonian matrices have a vanishing trace,i.e. h + h = 0, so that a Hamiltonian 2 × H = (cid:18) h h h − h (cid:19) (33)Furthermore it follows from Eq. 19 that symmetries playan essential role in linear Hamiltonian theory. If a generalHamiltonian matrix is splitted into it’s symmetric H s andskew-symmetric H a parts H s = 12 ( H + H T )= γ H γ − γ H ) H a = 12 ( H − H T )= − γ H γ + γ H ) , (34)then it is evident that symmetries and (anti-) commuta-tion properties are related: the part of H which com-mutes with γ , is skew-symmetric and the part of H which anti-commutes with γ , is symmetric. Hence it isrequired to distinguish symmetric from skew-symmetriccomponents: H = h (cid:18) − (cid:19) + h (cid:18) (cid:19) + h (cid:18) − (cid:19) (35)In consequence the parameters h , h and h representdynamical structures with specific well-defined symmetryproperties.Since any 2 n × n Hamiltonian matrix may depend on2 n (2 n + 1) / × Cl (1 , γ = η = (cid:18) − (cid:19) η = (cid:18) (cid:19) η = η η = (cid:18) − (cid:19) (36)The Pauli matrices mutually anti-commute and square to ± . In order to make this system of matrices η i complete(i.e. a group), we define η = × . IV. THE REAL DIRAC ALGEBRA
The original Jacobi algorithm picks out two “spatialcoordinates” of R n and “rotates them” by the use of 2 × × c and s are replaced by 2 × × × × ˜H = R H R − = (cid:18) G × K × (cid:19) (37)Let us first remark, before going into the details, thatHamiltonian matrices have in general the following prop-erty: If λ is an eigenvalue of some Hamiltonian matrix H , in the general case being complex, then − λ as well asthe complex conjugate values ± ¯ λ are also eigenvalues of H [10]. Since 2 × × H has complex eigenvalues off axis (i.e. neitherpurely real nor purely imaginary), then symplectic block-diagonalization is impossible since similarity transforma-tions, regardless whether they are symplectic or not, pre-serve the eigenvalues.Furthermore, let us remark that 2 × H = ( h η + h η + h η ) = − ( h − h − h ) × (38)Hence the eigenvalues λ of H × are readily obtained by λ = ± p − h + h + h . Eq. 38 provides evidence thatthe Hamiltonian matrix of an oscillatory degree of free-dom, normalized to the frequency, squares to − × ± not only have rather simplematrix exponentials but also unique symmetry proper-ties, it is nearby to use a Clifford algebraic parameteriza-tion of the required 4 × The computation of matrix exponentials is in general signifi-cantly more involved [13].
The Kronecker product of two matrices A = { a ij } and B = { b kl } is given by: C = A ⊗ B = (cid:18) a B a B a B a B (cid:19) = a b a b a b a b a b a b a b a b a b a b a b a b a b a b a b a b , (39)i.e. the Kronecker product is a method to systematicallywrite down all possible products between all elements of A and B , respectively. The rules of Kronecker matrixproducts are [11, 14]:( A ⊗ B ) T = A T ⊗ B T A ⊗ ( B + C ) = A ⊗ B + A ⊗ C ( A ⊗ B ) ( C ⊗ D ) = A C ⊗ B D
Tr( A ⊗ B ) = Tr( A ) Tr( B )( A ⊗ B ) − = A − ⊗ B − (40)If some set η k of matrices represents a Clifford algebra,as in case of the real Pauli matrices, then it is straight-forward to verify that the set of all Kronecker products η i ⊗ η j is again a representation of some Clifford algebra .Hence the required set of sixteen 4 × γ k with k ∈ [0 , . . . ,
15] where the (skew-symmetric) symplectic unitmatrix (SUM) is γ and the unit matrix is γ = × .In the following we shall investigate whether one can finda logical order for the remaining matrices.As well known, the elements of Clifford algebrasare generated (from products) of their mutually anti-commuting basis elements. A Clifford algebra of dimen-sion N has N = p + q basis elements, called “vectors”, p of which square to (i.e. have signature 1) and q squareto − (i.e. have signature − N which is the product of all basis vectors and iscalled “pseudo-scalar”. Real matrices of a given dimen-sion 2 m × m always represent some Clifford algebra, butcan only represent CAs of maximal possible dimension N = p + q = 2 m , if p − q = 0 , , , (41)which, since N = p + q is even in our case, reduces to p − q = 0 , . (42) More generally one finds (without proof) that all real matrix rep-resentations of Clifford algebras can be generated from (multiple)Kronecker products of the real Pauli matrices.
This is usually called Bott’s periodicity [15–17].The N generators (or vectors) can be used to obtainnew elements by multiplication since products of two (ormore) different generators γ i γ j are unique elements, dif-ferent from the unit element and different from each fac-tor. It follows from combinatorics that there are (cid:16) Nk (cid:17) k -vectors, so that one has X k (cid:18) Nk (cid:19) = 2 N (43)elements in total. This means that the structure of Clif-ford algebras is closely related to Pascal’s triangle.The question then is whether and how the Cliffordalgebraic structure might support the objective to de-velop a Jacobi method for phase space. In order to clar-ify this structure, let us first remark that all elementsof any (real rep of a) Clifford algebra either commuteor anti-commute with any other element, that any ele-ment is either purely symmetric or skew-symmetric andthat therefore any element is either Hamiltonian or skew-Hamiltonian.If we denote Hamiltonian elements by S k and skew-Hamiltonian elements by C k , then it is straightforwardto show that S S − S S C C − C C C S + S CS n +1 ⇒ Hamiltonian S S + S S C C + C C C S − S CS n C n ⇒ skew − Hamiltonian (44)Hence the commutator of Hamiltonian elements is againHamiltonian.From Eq. 34 we know that any Clifford element thatcommutes with γ , is skew-symmetric (with square − )and any component of H that anti-commutes with γ ,is symmetric (with square to + ). By definition all ele-ments of a Clifford basis mutually anti-commute. Hence,if γ is a basis element, then all other basis elements must anti-commute with γ . It follows that all other basis el-ements are either Hamiltonian and symmetric or skew-Hamiltonian and skew-symmetric.A purely Hamiltonian basis is hence only possible forClifford algebras of type Cl ( N − ,
1) so that Bott’s pe-riodicity (Eq. 42 results in :( N − − N − , . (45) From a physical point of view, this might be an essential insight:It follows that Hamiltonian notions applied to real Clifford alge-bras require vector spaces in which the metric is of the Minkowskitype.
Type Symbol Symmetry Hamiltonian1-vector γ - + ~γ + +2-vector γ ~γ + + γ γ ~γ - +3-vector γ γ = γ γ γ - - γ ~γ + -4-vector γ = γ γ γ γ - -(pseudo-scalar)skalar γ = + -TABLE I. Structure of the symplectic Clifford algebra Cl (3 , This includes the real Dirac algebra with N = 4. Hencea pure Hamiltonian basis exists and generates a uniquerepresentation of Cl (3 , (cid:0) (cid:1) = 6 bi-vectors, which are,according to Eq. 44, all Hamiltonian.It is straightforward to verify that products of two anti-commuting symmetric elements are skew-symmetric:( γ γ ) T = γ T γ T = γ γ = − γ γ (46)while products of γ and γ k ( k = 1 , ,
3) are symmetric bi-vectors. It follows that, when a SUM γ and three moremutually anti-commuting symmetric basis elements γ , γ and γ are selected, the structure of the HamiltonianClifford algebra is completely determined.Since symmetric 4 × : 1) There are 4 basis elements,the SUM γ and 3 symmetric elements γ , γ , γ . 2)There are 3 symmetric bi-vectors: γ = γ γ , γ = γ γ ,and γ = γ γ . 3) There are 3 skew-symmetric bi-vectors: γ = γ γ , γ = γ γ , and γ = γ γ . Notethat γ commutes with the skew-symmetric bi-vectorsand anti-commutes with the symmetric bi-vectors.The remaining elements of the Clifford algebra are four3-vectors and the pseudo-skalar 4-vector γ = γ γ γ γ ,which are all skew-Hamiltonian. The skew-symmetric bi-vectors can also be written as γ = γ γ γ γ = γ γ γ γ = γ γ γ . (47)It is sort of nearby to use a vectorial notation for thetriple ~γ = ( γ , γ , γ ) T , so that we have the “vectors”( γ , ~γ ) and the bi-vectors γ ~γ and γ γ ~γ , respectively.The question then is, whether and how the ten symplec- To our knowledge, there is no such stringent logic to be found incomplex representations of Clifford algebras. tic generators γ , . . . , γ that have been identified, whenregarded as generators of symplectic transformations, al-low to diagonalize real Hamiltonian 4 × H andwhether this provides the means for a symplectic Jacobialgorithm.Let us therefore return to Eq. 22 and investigate theeffect of transformations of the form R i ( τ ) = exp ( γ i τ /
2) (48) where i ∈ [0 , . . . , R i γ j R − i = ( c + γ i s ) γ j ( c − γ i s )= c γ j − γ i γ j γ i s + ( γ i γ j − γ j γ i ) c s ) (49)where c = cos ( τ /
2) and s = sin ( τ /
2) for γ i = − and c = cosh ( τ /
2) and s = sinh ( τ /
2) for γ i = . If γ i and γ j commute, this results in: R i γ j R − i = γ j , (50)and if γ i and γ j anti-commute, one obtains: R i γ j R − i = γ j ( c + γ i s ) + γ i γ j c s . (51)In case of a skew-symmetric γ i , the final result is : R i γ j R − i = γ j cos ( τ ) + γ i γ j sin ( τ ) . (52)and in case of a symmetric γ i : R i γ j R − i = γ j cosh ( τ ) + γ i γ j sinh ( τ ) . (53)Hence skew-symmetric γ i generate rotation-like transfor-mations and symmetric γ i generate boost-like symplectictransformations.Note that we obtained the complete structure of thealgebra from only two requirements, namely that γ isa basis element of the Clifford algebra and that all ele-ments of the basis should be Hamiltonian as defined byEq. 19. These two (nearby) conditions inevitably lead tothe unique structure of Cl (3 ,
1) as listed in Tab. I.
V. THE ELECTROMAGNETIC EQUIVALENCE
Since any real-valued 4 × M can be writtenas M = X k =0 m k γ k , (54) Note that the inverse transformation is given by the negativeargument: R − i = R ( − τ ) = exp ( − γ i τ/ The half-angle arguments has been chosen in order to obtainfull-angle arguments here. any real-valued
Hamiltonian × H can accord-ingly be written as: H = X k =0 h k γ k (55)Since the unit element γ = × is the only γ -matrixwith non-vanishing trace, the trace of M is 4 m . Thecoefficients m k therefore can be obtained from the traceof the product of γ Tk and M m k = 14 Tr( γ Tk M ) , (56)since γ Tk γ k = for all k ∈ [0 , . . . , H are given by h k = 14 Tr( γ Tk H ) . (57)The structure of Cl (3 ,
1) is, beyond the sign of the metric,closely related to Cl (1 ,
3) which is used in the usual pre-sentation of the Dirac equation. Indeed the “real Diracmatrices” are, up to multiplication with the unit imagi-nary, identical to the Majorana matrices [18]. We shallnot discuss here whether (and which) physical insightsthis might eventually provide, but the transformationproperties of the ten relevant parameters provide a strictformal correspondence to quantumelectrodynamics andtherefore it is convenient to use this isomorphism by thefollowing notation: h = E ( h , h , h ) T = ~P ( h , h , h ) T = ~E ( h , h , h ) T = ~B (58)so that any 4 × P and the bi-vectorcomponents F as follows: P = E γ + ~P · ~γ = E γ + P x γ + P y γ + P z γ F = γ ~E · ~γ + γ γ ~B · ~γ = E x γ + E y γ + E z γ + B x γ + B y γ + B z γ H = P + F (59)With the basis given in App. A one obtains the follow-ing explicit form for the vectors: P = − P z E − P x P y −E − P x P z P y P y − P z E + P x P y −E + P x P z (60) For the bi-vectors one obtains: F = − E x E z + B y E y − B z B x E z − B y E x − B x − E y − B z E y + B z B x E x E z − B y − B x − E y + B z E z + B y − E x (61) γ γ γ γ ˜ γ = γ ˜ γ = c γ .. − s γ − s γ + s γ ˜ γ = c γ .. − s γ + s γ − s γ ˜ γ = c γ .. − s γ − s γ + s γ ˜ γ = c γ .. + s γ − s γ γ ˜ γ = c γ .. + s γ + s γ − s γ ˜ γ = c γ .. + s γ − s γ + s γ ˜ γ = c γ .. − s γ + s γ ˜ γ = c γ .. + s γ − s γ ˜ γ = c γ .. − s γ + s γ TABLE II. Table of symplectic rotations of two degrees offreedom. a indicates the rows and b the column: γ ′ a =exp ( γ b τ / γ a exp ( − γ b τ / γ a and γ b anticommute, thenthe result is γ ′ a = c γ a + s γ a γ b where c and s are the sine-and cosine-function of τ . γ γ γ γ γ γ ˜ γ = c γ .. + s γ + s γ + s γ − s γ − s γ − s γ ˜ γ = c γ .. + s γ − s γ − s γ ˜ γ = c γ .. − s γ + s γ − s γ ˜ γ = c γ .. + s γ − s γ − s γ ˜ γ = c γ .. + s γ + s γ − s γ ˜ γ = c γ .. + s γ − s γ + s γ ˜ γ = c γ .. + s γ + s γ − s γ ˜ γ = c γ .. − s γ + s γ − s γ + s γ ˜ γ = c γ .. + s γ − s γ + s γ − s γ ˜ γ = c γ .. − s γ + s γ − s γ + s γ TABLE III. Table of symplectic boosts with symmetric gen-erator of two degrees of freedom. a indicates the rows and b the column: ˜ γ a = exp ( γ b τ / γ a exp ( − γ b τ / γ a and γ b anticommute, then the result is ˜ γ a = c γ a + s γ a γ b where c and s are the hyperbolic sine- and cosine-function of τ . If γ a and γ b commute, then γ ′ a = γ a . The explicite computation of the effects of symplec-tic similarity transformations are given in Tab. II andTab. III. The calculations verify that the gyroscopicHamiltonian terms associated here with ~B generate ro-tations which are isomorph to spatial rotation of vec-tors. Likewise the symmetric bi-vector terms ~E generateboosts that are isomorphic to Lorentz-Boosts. This iso-morphism is indeed helpful to grasp the geometric con-tent which is needed to obtain block-diagonalization. VI. BLOCK-DIAGONALIZATION
The analogy with the conventional Jacobi algorithm isachieved by a block-diagonalization of H . It is achieved ˜ ε r ˜ ε g ˜ ε b γ ε r c + ε g s ε g c − ε r s ε b c + ~P − ~E s γ ε r C − ( ~b ) x S ε g ε b C − ( ~r ) x Sγ ε r C − ( ~b ) y S ε g ε b C − ( ~r ) y Sγ ε r C − ( ~b ) z S ε g ε b C − ( ~r ) z Sγ ε r ε g C + ( ~b ) x S ε b C + ( ~g ) x Sγ ε r ε g C + ( ~b ) y S ε b C + ( ~g ) y Sγ ε r ε g C + ( ~b ) z S ε b C + ( ~g ) z S TABLE IV. Table of “scalar products” after symplectic trans-formations generated by symplectic boosts and γ , respec-tively (left column). See Eqns. 63, 64 and 65. when the transformed matrix ˜H = M H M − (62)has block-diagonal form, which means that ˜ P y = ˜ E y =˜ B x = ˜ B z = 0 (see Eq. 60 and Eq. 61). In other words,the “vectors” ~P and ~E must be orthogonal to ~B . Let’sdefine the following auxiliary “scalar products” ε r = ~E · ~Bε g = ~B · ~Pε b = ~E · ~P (63)The first objective that has to be achieved, is therefore ε r = ε g = 0. The second objective is to rotate the “sys-tem” such that ~B = B y . For this we introduce the fol-lowing auxiliary “vectors”: ~r ≡ E ~P + ~B × ~E~g ≡ E ~E + ~P × ~B~b ≡ E ~B + ~E × ~P , (64)The “scalar products” of Eq. 64 are invariant under spa-tial rotations. Hence we need to analyze the transforma-tion behavior of these scalar products only under sym-plectic boosts and under the (yet unnamed) rotation gen-erated by γ . We introduce the following abbreviationsfor a better readability c = cos ( τ ) s = sin ( τ ) c = cos (2 τ ) s = sin (2 τ ) C = cosh ( τ ) S = sinh ( τ ) C = cosh (2 τ ) S = sinh (2 τ ) (65)The results are summarized in Tab. IV. The inspectionof this tabel reveals that ε r is invariant under the actionof R , , while ε g is invariant under the action of R , , .This findings suggests the following possible strategies:Firstly, we can use R to either make ε r = 0, then (oneof) R , , to make ε g = 0 or we could first use R tomake ε g = 0, followed by R , , to make ε r = 0. In bothcases, it is possible to use spatial rotations ( R , , ) toalign ~b along one axis, preferably the y -axis since in thechosen matrix representation, b x = b z = 0 will eventuallylead to the required “vector orientation” of ~B = B y , i.e.provide B x = B z = 0.Hence we can identify four steps that are needed toachieve block-diagonalization. After each step, the co-efficients h k , the auxiliary parameters ε i and the “vec-tors” ~r, ~g,~b have to be re-evaluated. The required stepsare: 1) Rotate using R (see Eq. 48) and angle τ = − arctan ( ε r /ε g ). This result in ˜ ε r = 0. 2) Rotate ~b us-ing R with angle τ = − arctan ( b x /b y ) such that ˜ b x = 0.3) Rotate ~b using R with angle τ = arctan ( b z /b y ) suchthat ˜ b z = 0. 4) The next step is to apply a boost with R and τ = − artanh( ε g / b y ). This last step requiresthat | ε g /b y | < ˜H is thengiven by: ˜H = R R R R H R − R − R − R − (66)Alternatively we could also proceed as follows: 1) Ro-tate using R (see Eq. 48) and angle τ = arctan ( ε g /ε r ).This should make ˜ ε g = 0. 2) Rotate ~b using R withangle τ = − arctan ( b x /b y ) such that ˜ b x = 0. 3) Ro-tate ~b using R with angle τ = arctan ( b z /b y ) such that˜ b z = 0. 4) The last step is to apply a boost with R and τ = artanh( ε r / b y ). Again, this is only possible for | ε r /b y | < ˜H = R R R R H R − R − R − R − (67)In order to block-diagonalize 2 n × n Hamiltonian ma-trices, the algorithm follows the same logic as the usualJacobi algorithm and requires, before each step, to selectthe “dominant” off-diagonal 2 × × × n ,i.e. with O ( n ) steps [19]. We do not provide a rigorousproof of convergence. VII. DIAGONALIZATION
After successfull block-diagonalization, the problem isreduced to that of Hamiltonian 2 × H × = h η + h η + h η ˜H × = exp ( η τ / H × exp ( − η τ /
2) (68) with τ = − arctan ( h /h ) one obtains: ˜H × = h + ˜ h − h + ˜ h ! ˜ h = h p h /h ) (69)In the next step, we transform using the generator η : ˜H × → exp ( η τ / H × exp ( − η τ / h − ˜ h ! τ = log p h − ˜ h p h + ˜ h ! ˜ h = ω = q h − h − h (70)This is still not a diagonal matrix, but it contains allrequired information and it is all that can be done by theuse of symplectic transformations. The form of ˜H × , asgiven in the first line of Eq. 70, is the normal form of aone-dimensional oscillator.Since the eigenvalues of ˜H × are directly given as ± i ω , it is possible – though not required – to take thenext step and diagonalize the matrix. We provide thisstep only for completeness. The (almost trivial) matrix V of eigenvectors and the final transformation to diago-nal form are: V = 1 √ i − i ! V − = 1 √ − i i ! V ˜H × V − = i ω − i ω ! (71)This last step is not a symplectic transformation.Nonetheless it is an interesting step as it reveals the last(or first) step that transforms between a “classical” realand a seemingly “non-classical” (complex) form of Hamil-tonian dynamics, which consists in replacement˜ ψ = V ψ = 1 √ i − i ! qp ! = 1 √ p + i qp − i q ! (72)The equations of motion then read, after multiplicationwith the unit imaginary and ~ : i ~ ˙˜ ψ = − ~ ω ~ ω ! ˜ ψ , (73) VIII. EIGENVALUES
The eigenvalues of H × can also be obtained by a dif-ferent method. It is well known in linear algebra that thetrace of a matrix is the sum of it’s eigenvalues, the trace ofthe square of a matrix is the sum of squares of it’s eigen-values. Since the sum of the eigenvalues of a Hamiltonianmatrix (and it’s odd powers) vanishes, the trace of theeven powers do not vanish. Recall that those Hamilto-nian matrices that allow for block-diagonalization, haveeither purely real or purely imaginary eigenvalues. Let ± λ and ± Λ be the four eigenvalues, then:Tr( H ) = 2 λ + 2 Λ Tr( H ) = 2 λ + 2 Λ (74)If we introduce the definitions K = Tr( H ) / λ + Λ K = Tr( H ) / − K / λ + Λ − λ + Λ + 2 λ Λ
16= ( λ − Λ )
16 (75)so that one obtains: K + 2 p K = λ + Λ λ − Λ λ K − p K = λ + Λ − λ − Λ (76)and hence the eigenvalues can be determined from K and K alone: λ = ± q K + 2 p K Λ = ± q K − p K (77)Inserting Eq. 60 and Eq. 61 yields [21]: K = −E − ~B + ~P + ~E K = ( E ~B + ~E × ~P ) − ( ~E · ~B ) − ( ~P · ~B ) = ~b − ε r − ε g (78)There are two cases of special interest, namely if H = P (only vector components) K = −E + ~P K = 0 , (79)and H = F (only bi-vector components), K = ~E − ~B K = − ( ~E · ~B ) . (80)Recall that symplectic eigenvalues are invariants . IX. SYMPLECTIC MATRICES
According to linear Hamiltonian theory, every sym-plectic matrix M is the matrix exponentials of someHamiltonian matrix H [10]: M = exp ( H τ )= + H τ + 12 H τ + . . . . (81)Since any analytic function of block-diagonal matricesis again a block-diagonal matrix it is evident that M is blockdiagonal whenever H is block- diagonal. FromEq. 44 we know that all odd terms of the Taylor seriesEq. 81 are Hamiltonian and all even terms are skew-Hamiltonian. Hence, we can use the brute force methodand remove all skew-Hamiltonian terms from M : ˜M = 12 ( M + γ M T γ ) . (82)Whatever similarity transformation is used to block-diagonalize ˜M , it will automatically block-diagonalize M as well.In other words: If a Hamiltonian matrix can be writtenas H = V D V − , (83)then the matrix exponential of H is: M = exp ( H ) = V exp ( D ) V − , (84)so that an eigenvalue λ of H is replaced by matrix expo-nentiation with e λ in M . X. SKEW-HAMILTONIAN MATRICES
Skew-Hamiltonian matrices C , when expressed withthe real Dirac algebra, can be written as (see Tab. I): C = X k =10 c k γ k . (85)The four components ( c , c , c , c ) are sometimescalled pseudo-vector or “axial” vector and they transformaccordingly. The pseudo-scalar c and the scalar c areinvariants. For instance H / c = ~P · ~B ( c , c , c ) T = E ~B + ~E × ~Pc = ~E · ~Bc = ( −E + ~P − ~B + ~E ) / . (86)The reader will have noticed that c = ε g , c = ε r and( c , c , c ) T = ~b have been defined before in Eq. 630and Eq. 64. Hence the simplest way to determine thoseparameters, is to compute H and extract those compo-nents: ε g = Tr( γ T10 H ) / ε r = Tr( γ T14 H ) / b x = Tr( γ T11 H ) / b y = Tr( γ T12 H ) / b z = Tr( γ T13 H ) / c = c = c = c = 0, or, in the notation introduced above that ε r = ε g = b x = b z = 0. This is exactly the strategy of thedescribed algorithm. Hence the developed algorithm is already able to block-diagonalize skew-Hamiltonian ma-trices. XI. APPLICATIONS
Many problems in engineering and physics involveHamiltonian matrices, often of size greater 2 ×
2. Anexample from accelerator physics can, for instance, befound in Ref. [22]. Since most accelerators make use ofa median plane symmetry, the vertical degree of freedomdoes, in linear approximation, not couple to the trans-verse horizontal and the longitudinal degrees of freedom.Then the involved 6 × × × ψ of individualparticles are the dynamical variables and they represent(small) deviations of a particle’s position and momen-tum relative to some reference orbit. But the “observ-ables” that are relevant to describe the collective behav-ior are not individual particle positions but statisticalaverages of particle ensembles as for instance the root-mean-square width of the beam, which is encoded in thematrix of second moments:Σ = h ψψ T i . (88)The evolution of Σ in time is given by˙Σ = h ˙ ψψ T i + h ψ ˙ ψ T i = h H ψψ T i + h ψψ T H T i = H Σ + Σ H T = H Σ + Σ γ H γ . (89)If we multiply the last line with γ T from the right, weobtain: ˙Σ γ T = H (Σ γ T ) − (Σ γ T ) H ˙ S = H S − S H , (90) with the Hamiltonian matrix S = Σ γ T . The last line tellsus that H and S are a so-called “Lax pair”. As Peter Laxhas proven[23], it follows from the validity of Eq. 90 thatthe following expressionTr( S k ) = const (91)is invariant for any integer k .Since the evolution of ψ in time is given by the sym-plectic transformation ψ ( t ) = M ( t ) ψ (0) = exp ( H τ ) ψ (0) , (92)then we obtain for the the matrix of second moments:Σ( t ) = M Σ(0) M T S ( t ) = M S (0) γ M T γ T = M S (0) M − (93)Hence the evolution of the matrix S in time is describedby a symplectic similarity transformation, i.e. the eigen-values of S are conserved (which follows already fromLax’ theorem).A typical problem in accelerator physics requires to de-termine a so-called matched beam: Given that the (sym-plectic) transfer matrix M of some beamline (or accel-erator ring) is known, the problem is to find a stable(“matched”) distribution S for given beam emittances.Matrices that share a system of eigenvectors, commutewith each other. Hence, provided that S and M share asystem of eigenvectors, we obtain from Eq. 93: S ( t ) = S (0) . (94)Hence, in the matched case, the same transformationsthat allow to block-diagonalize (the Hamiltonian part of)the matrix M , also diagonalize S . Provided one has analgorithm to generate single-variate Gaussian distribu-tions, then the described block-diagonalization enablesto generate arbitrary multivariate matched distributionsby application of the reverse transformation [24]. XII. SUMMARY AND OUTLOOK
The described Jacobi method allows to bring (skew-)Hamiltonian and symplectic matrices to block-diagonalform. This enables to compute matrix exponentials andlogarithms, eigenvalues and eigenvectors of Hamiltonianmatrices with eigenvalues on the real and/or the imagi-nary axis.The presented analysis of the structure of Hamilto-nian phase spaces reveils an isomorphism between four-component spinors and points in a classical Hamilto-nian phase space for two degrees of freedom. It is thissame structure that has been described as the “square1root of geometry” and which is also used in quantum-electrodynamics .It is known (though not well-known) since long thatthe symplectic group Sp (4) and the algebraic structureof the Dirac algebra are closely related [26], but few au-thors made use of this remarkable fact. We emphasizethat the construction of both the real Pauli as well asthe real Dirac algebra was derived exclusively from thesymmetry properties of a general classical Hamiltonianphase space. It follows that the restricted Lorentz groupis isomorphic to symplectic transformations of classicalphase space ensembles with generators of the bi-vector-type. Since the presented algorithm is based solely onreal Hamiltonian 4 × ACKNOWLEDGMENTS L A TEX has been used to write this article,Mathematica R (cid:13) has been used for parts of the sym- bolic calculations. Appendix A: The Real Dirac Matrices
In this article we used the following representation: γ = − − γ = − − γ = γ = − − However, any representation allows for an equivalent al-gorithm, firstly due to Pauli’s fundamental theorem ofthe Dirac matrices, but secondly since a different choiceof matrices mostly stems from a different ordering of thecanonical coordinates and momenta in ψ . [1] P.A.M. Dirac P.A.M. Dirac “Why We Believe in the Ein-stein Theory” in: Bruno Gruber and Richard S. Mill-man (Eds.) “SYMMETRIES IN SCIENCE” PLENUMPRESS, NEW YORK AND LONDON (1980), DOI:10.1007/978.1.4684.3833.8.[2] Alan J. Laub and Kenneth Meyer “Canonical Forms forSymplectic and Hamiltonian Matrices”; Cel. Mech. (1974), pp. 213-23.[3] Gene H. Golub and Henk A. van de Vorst “Eigenvaluecomputation in the 20th century”; J. Comp. Appl. Math. (2000), pp. 35-65.[4] H. Fassbender, D.S. Mackey and N. Mackey: “Hamil-ton and Jacobi come full circle: Jacobi algorithms forstructured Hamiltonian eigenproblems”: Lin. Alg. Appl.332-334 (2001) 37-80.[5] M. Kleinsteuber “A sort-Jacobi algorithm for semisimpleLie algebras”; Lin. Alg. Appl. (2009), pp 155-175.[6] V. Mehrmann, C. Schr¨oder and D.S. Watkins “A newblock method for computing the Hamiltonian Schurform”; Lin. Alg. Appl. 431 (2009), pp. 350-368.[7] P. Benner, H. Fassbender and M. Stoll “A HamiltonianKrylov-Schur-type method based on the symplectic Lanc-zos process”; Lin. Alg. Appl. 435 (2011) 578-600.[8] Ru-Ru Ma, Zhi-Gang Jia and Zheng-Jian Bai “Astructure-preserving Jacobi algorithm for quaternion “No one fully understands spinors. Their algebra is formallyunderstood but their general significance is mysterious. In somesense they describe the square root of geometry and, just asunderstanding the square root of -1 took centuries, the samemight be true of spinors” [25]. Hermitian eigenvalue problems”; Comp. Math. Appl. (2018), pp. 809-820.[9] Philipp O.J. Scherer “Computational Physics”; 2nd Ed.Springer (2013).[10] K.R. Meyer, G.R. Hall and D. Offin “Introduction toHamiltonian Dynamical Systems and the N-Body Prob-lem”; 2nd. Ed., Springer, New York, 2000.[11] Karim B. Abadir and Jan R. Magnus; “Matrix Algebra”;Cambridge University Press 2005.[12] Hermann Weyl “The Classical Groups: Their Invariantsand Representations”; Princeton University Press (1939).[13] C. Baumgarten “Analytic Expressions for Ex-ponentials of Specific Hamiltonian Matrices”;https://arxiv.org/abs/1703.02893.[14] Charles F. Van Loan “The ubiquitous Kronecker prod-uct”; J. Comp. Appl. Math. (2000), pp. 85-100.[15] Raoul Bott “The Stable Homotopy of the ClassicalGroups” PNAS No. 10 (1957), pp. 933-935.[16] M. Atiyah, R. Bott and A. Shapiro “Clifford Modules”;Topology Suppl. 1,(1964), pp. 3-38.[17] Susumu Okubo “Real representations of finite Cliffordalgebras”, Part 1 and 2; J. Math. Phys. 32, (1991), pp.1657-1668 and pp. 1669-1673.[18] Ettore Majorana; Il Nuovo Cimento 5 (1937), 171-184.[19] C. Baumgarten “Geometrical method of decou-pling”; Phys. Rev. ST Accel. Beams. 15, 124001(2012). A Mathematica notebook is available here: .[20] C. Baumgarten “Use of Real Dirac Matrices in 2-dimensional Coupled Linear Optics”; Phys. Rev. ST Ac-cel. Beams. 14, 114002 (2011).[21] C. Baumgarten; “Old Game, New Rules: Rethink-ing the Form of Physics”, Symmetry 2016, 8(5), 30;2