A note on trigonometric identities involving non-commuting matrices

TL;DR

This paper introduces an algorithm for approximating trigonometric functions of sums of non-commuting matrices, useful in quantum mechanics, involving nested commutators and converging within the Zassenhaus formula domain.

Contribution

It presents a novel iterative method to compute trigonometric functions of non-commuting matrices with convergence guarantees.

Findings

Algorithm converges within the Zassenhaus domain

Expressions involve nested commutators

Applicable to perturbative quantum mechanics

Abstract

An algorithm is presented for generating successive approximations to trigonometric functions of sums of non-commuting matrices. The resulting expressions involve nested commutators of the respective matrices. The procedure is shown to converge in the convergent domain of the Zassenhaus formula and can be useful in the perturbative treatment of quantum mechanical problems, where exponentials of sums of non-commuting skew-Hermitian matrices frequently appear.

Equations

$$\frac{d^{2}y}{d{t}^2}+A^{2}y=0,y(0)=y_0,y^{\prime }(0)=y_0^{\prime },$$

$$y(t)=\cos (tA)y_0+A^{-1}\sin (tA)y_0^{\prime }.$$

$$\cos (A)$$

$$\sin (A)$$

$${\mathrm{e}}^{iA}=\cos (A)+i\sin (A),$$

$$\cos (A)=\frac{{\mathrm{e}}^{iA}+{\mathrm{e}}^{-iA}}{2},\sin (A)=\frac{{\mathrm{e}}^{iA}-{\mathrm{e}}^{-iA}}{2i}$$

$${\cos }^2(A)+{\sin }^2(A)=I.$$

$$\cos (2X)=2{\cos }^2(X)-I,$$

$$\cos ((A+B)t)$$

$$\sin ((A+B)t)$$

$$A=\pi\left(\begin{array}[]{cr}0&\ \alpha\\ -1/\alpha&0\end{array}\right),\qquad\quad B=\pi\left(\begin{array}[]{cc}0&\ (10+4\sqrt{6})\alpha\\ (-10+4\sqrt{6})/\alpha&0\end{array}\right).$$

$$Z=\log (e^A{\mathrm{e}}^B)=A+B+\sum _{m=2}^{\infty }{Z}_m(A,B)$$

$$m=1:Z_1$$

$$m=2:Z_2$$

$$m=3:Z_3$$

$$m=4:Z_4$$

$$m=5:Z_5$$

$$\cos (A+B)$$

$$\sin (A+B)$$

$${\mathrm{e}}^{X+Y}={\mathrm{e}}^X{\mathrm{e}}^Y\prod _{n=2}^{\infty }{\mathrm{e}}^{C_n(X,Y)}={\mathrm{e}}^X{\mathrm{e}}^Y{\mathrm{e}}^{C_2(X,Y)}{\mathrm{e}}^{C_3(X,Y)}\cdots {\mathrm{e}}^{C_k(X,Y)}\cdots ,$$

$$C_2(X,Y)=-\frac{1}{2}[X,Y]$$

$$C_3(X,Y)=\frac{1}{3}[Y,[X,Y]]+\frac{1}{6}[X,[X,Y]]$$

$$C_4(X,Y)=-\frac{1}{24}[X,[X,[X,Y]]]-\frac{1}{8}[Y,[X,[X,Y]]]-\frac{1}{8}[Y,[Y,[X,Y]]].$$

$$\begin{array}[]{l}\mbox{Define}\;f_{1,k}\;\mbox{by}\\ \;\;\;\displaystyle f_{1,k}=\sum_{j=1}^{k}\frac{(-1)^{k}}{j!(k-j)!}{\rm ad}_{Y}^{k-j}{\rm ad}_{X}^{j}Y,\\ \displaystyle C_{2}=\frac{1}{2}\,f_{1,1},\\ \mbox{Define}\;f_{n,k}\quad n\geq 2,\;k\geq n\;\;\mbox{by}\\ \;\;\;\displaystyle f_{n,k}=\sum_{j=0}^{[k/n]-1}\frac{(-1)^{j}}{j!}{\rm ad}_{C_{n}}^{j}f_{n-1,k-nj},\\ \displaystyle C_{n}=\frac{1}{n}f_{[(n-1)/2],n-1}\quad n\geq 3.\end{array}$$

$${\mathrm{ad}}_{A}B=[A,B],{\mathrm{ad}}_A^{j}B=[A,{\mathrm{ad}}_A^{j-1}B],{\mathrm{ad}}_A^{0}B=B.$$

$$\underset{n\to \infty }{\lim }{\mathrm{e}}^X{\mathrm{e}}^Y{\mathrm{e}}^{C_2}\cdots {\mathrm{e}}^{C_n}={\mathrm{e}}^{X+Y}$$

$$\Vert f_{n,k}\Vert \le d_{n,k},\Vert C_n\Vert \le {\delta }_n=\frac{1}{n}{d}_{[(n-1)/2],n-1}$$

$${\mathrm{e}}^{X+Y}\approx {\mathrm{e}}^X{\mathrm{e}}^Y{e}^{C_2(X,Y)}{\mathrm{e}}^{C_3(X,Y)}\cdots {\mathrm{e}}^{C_n(X,Y)}.$$

$${\mathrm{e}}^{i(X+Y)}$$

$${\mathrm{e}}^{-i(X+Y)}$$

$${\hat{C}}_n=i^n{C}_n,{\tilde{C}}_n=(-i{)}^n{C}_n,n\ge 2$$

$${\hat{C}}_{2k}={\tilde{C}}_{2k}=(-1{)}^k{C}_{2k},$$

$${\hat{C}}_{2k+1}=-{\tilde{C}}_{2k+1}=(-1{)}^{k}i{C}_{2k+1},k\ge 1.$$

$$\cos(X+Y)=\frac{1}{2}\big{(}{\rm e}^{i(X+Y)}+{\rm e}^{-i(X+Y)}\big{)},\qquad\sin(X+Y)=\frac{1}{2i}\big{(}{\rm e}^{i(X+Y)}-{\rm e}^{-i(X+Y)}\big{)},$$

$$z_{1,1}\equiv {\mathrm{e}}^{iX}{\mathrm{e}}^{iY}=\cos (X)\cos (Y)-\sin (X)\sin (Y)$$

$$\displaystyle\quad+i\big{(}\cos(X)\sin(Y)+\sin(X)\cos(Y)\big{)}$$

Full text

1 Introduction

Trigonometric matrix functions appear naturally when solving systems of second order differential equations

[TABLE]

whose solution is expressed by

[TABLE]

for all $n\times n$ matrices $A$ higham08fom . When $A$ is singular, (2) is interpreted by expanding the matrix cosine and sine functions as power series of $A$:

[TABLE]

Equation (1) arises in finite element semidiscretizations of the wave equation, whereas similar equations with a non-vanishing right-hand side of the form $g(t,y(t),y^{\prime}(t))$ appear in highly oscillatory problems, control theory, etc.

In this case one has also the matrix analogue of Euler’s formula,

[TABLE]

so that

[TABLE]

and

[TABLE]

Different algorithms exist in the literature for the practical computation of the matrix cosine and sine (see e.g. almohy15naf ; higham08fom and references therein). Several of them make use of the double angle formula,

[TABLE]

to construct an approximation $Y$ to $\cos(A)$ by first considering a matrix $X=2^{-s}A$ with small norm and then approximating $\cos(X)$ by a function $r(X)$ (a truncated Taylor series, a Padé approximant, etc.). $Y$ is then determined by applying formula (5) $s$ times.

Identity (5) is a special case of the addition formulae

[TABLE]

which hold if and only if $AB=BA$ (higham08fom, , p. 287). This is not necessary the case, however, when $t=1$, as the following pair of matrices illustrate frechet52lsn :

[TABLE]

Although $AB\neq BA$ for all $\alpha\neq 0$, a straightforward calculation shows that, indeed, equations (1) with $t=1$ are still valid here. For general matrices $A$ and $B$, however, one cannot expect them to hold unless their commutator $[A,B]\equiv AB-BA$ vanishes. This property is of course related through eq. (4) with the celebrated Baker–Campbell–Hausdorff (BCH) formula bonfiglioli12tin . Roughly speaking, ${\rm e}^{A}\,{\rm e}^{B}={\rm e}^{A+B+C}$, where the additional term $C$ is due to the non-commutativity of $A$ and $B$. More in detail, the BCH theorem establishes that ${\rm e}^{A}\,{\rm e}^{B}={\rm e}^{Z}$, with

[TABLE]

and $Z_{m}(A,B)$ is a linear combination (with rational coefficients) of nested commutators involving $m$ operators $A$ and $B$. The first terms read explicitly

[TABLE]

An efficient algorithm for generating explicit expressions of $Z_{m}(A,B)$ up to an arbitrary $m$ in terms of independent commutators is presented in casas09aea . At this point it is natural to raise the following question: since formulae (1) do not hold in general for non-commutative matrices, is it still possible to express $\cos(A+B)$ in terms of the cosine and sine of $A$ and $B$ for general matrices when $[A,B]\neq 0$? And if the answer is in the affirmative, can this be done in a systematic (and hopefully efficient) way?

It is the purpose of this note to develop an algorithm that generalizes identities (1) to non-commuting operators, thus providing successive approximations to $\cos(A+B)$ and $\sin(A+B)$ involving $n$-nested commutators of $A$ and $B$ for $n=1,2,\ldots$. As an illustration, if $A$ and $B$ are such that $[A,[A,B]]=[B,[A,B]]=0$, then the following exact result holds:

[TABLE]

The algorithm we propose here constitutes in fact a direct application of the dual of the BCH theorem: the so-called Zassenhaus formula, with multiple applications in quantum mechanical systems and numerical analysis casas12eco . The problem consists essentially in finding matrices (operators) $C_{1},C_{2},\ldots$ such that ${\rm e}^{A+B}={\rm e}^{A}\,{\rm e}^{B}\,{\rm e}^{C_{1}}\,{\rm e}^{C_{2}}\cdots$, with $C_{i}$ depending only on nested commutators of $A$ and $B$.

Expressions like (1) can be useful in the perturbative treatment of quantum problems where exponentials of sums of non-commuting skew-Hermitian operators frequently appear galindo90qme .

2 Zassenhaus formula

To establish the Zassenhaus formula we consider two non commuting indeterminate variables $X$, $Y$ and the free Lie algebra generated by them, $\mathcal{L}(X,Y)$. This, roughly speaking, can be viewed as the set of linear combinations of all commutators that can be constructed with $X$ and $Y$. The elements of $\mathcal{L}(X,Y)$ are called Lie polynomials postnikov94lga . A free Lie algebra is a universal object, so that results formulated in $\mathcal{L}(X,Y)$ are valid for any (finite- or infinite-dimensional) Lie algebra munthe-kaas99cia .

Let us suppose then that $X,Y\in\mathcal{L}(X,Y)$. The Zassenhaus formula establishes that the exponential ${\rm e}^{X+Y}$ can be uniquely decomposed as

[TABLE]

where $C_{k}(X,Y)\in\mathcal{L}(X,Y)$ is a homogeneous Lie polynomial in $X$ and $Y$ of degree $k$ magnus54ote ; suzuki77otc ; weyrauch09ctb ; wilcox67eoa ; witschel75ooe . The first terms read explicitly

[TABLE]

A recursive algorithm has been proposed in casas12eco for obtaining the terms $C_{n}$ up to a prescribed value of $n$ directly in terms of the minimum number of independent commutators involving $n$ operators $X$ and $Y$. The procedure, in addition, can be easily implemented in a symbolic algebra system without any special requirement, beyond the linearity property of the commutator. It reads as follows:

[TABLE]

Here $[k/n]$ denotes the integer part of $k/n$ and the “ad” operator is defined by

[TABLE]

Whereas the factorization (8) is well defined in the free Lie algebra $\mathcal{L}(X,Y)$, it has only a finite radius of convergence when $X$ and $Y$ are $n\times n$ real or complex matrices. Specifically,

[TABLE]

only in a certain subset of the plane $(\|X\|,\|Y\|)$ bayen79otc ; suzuki77otc . As a matter of fact, by bounding appropriately the terms $f_{n,k}$ and also the $C_{n}$, i.e., by showing that

[TABLE]

and analyzing (numerically) the convergence of the series $\sum_{n=2}^{\infty}\delta_{n}$, it can be shown that the convergence domain contains the region $\|X\|+\|Y\|<1.054$, and extends to the points $(\|X\|,0)$ and $(0,\|Y\|)$ with arbitrarily large values of $\|X\|$ or $\|Y\|$ casas12eco . In practical applications, however, the infinite product (8) is truncated at some $n$ and then one takes the approximation

[TABLE]

When the Zassenhaus formula is applied to $\exp(\pm i(X+Y))$, one gets

[TABLE]

respectively, where

[TABLE]

and $C_{n}$ is determined by algorithm (10). In more detail,

[TABLE]

3 The algorithm

Expansions (2), together with (4), allow us to design a recursive procedure and obtain expressions for $\cos(X+Y)$ and $\sin(X+Y)$ in terms of the sine and cosine of $X$ and $Y$. Since

[TABLE]

all we have to do is to insert the factorizations (2) in these expressions and collect terms up to the order $n$ considered. Specifically, let us first introduce

[TABLE]

and, for $n\geq 2$,

[TABLE]

Then it is clear that

[TABLE]

Thus, up to $n=2$, one has the approximations

[TABLE]

which reproduce, with $C_{2}$ given by (9), expressions (1) (with the replacement of $X$, $Y$ by $A$ and $B$, respectively), whereas analogously

[TABLE]

The general algorithm can then be established as follows:

[TABLE]

Moreover, it is possible to establish the convergence of the procedure as follows. From (16) we have

[TABLE]

and so

[TABLE]

in the convergence domain of the Zassenhaus formula (11), in particular when $\|X\|+\|Y\|<1.054$. By applying a similar argument, it is also true that

[TABLE]

in the same domain.

The recursion (17) can be easily programmed with a symbolic algebra package in conjunction with algorithm (10) to generate the terms $C_{n}$ and thus produce approximations to $\cos(X+Y)$ and $\sin(X+Y)$ up to the desired order $n$. In particular, up to $n=4$ we have

[TABLE]

4 Examples

Next we collect two particular examples to illustrate the use of, and results obtained by, algorithm (17) to approximate $\cos(X+Y)$ and $\sin(X+Y)$.

Example 1.

Pauli matrices play an important role in many quantum mechanical problems. They are defined by

[TABLE]

and form a basis of $\mathfrak{su}(2)$, the Lie algebra of $2\times 2$ skew-Hermitian traceless matrices. They verify

[TABLE]

so that their commutators are given by

[TABLE]

where $\epsilon_{jkl}$ denotes the Levi–Civita symbol. It can be shown that

[TABLE]

where $a=\|\boldsymbol{a}\|=\sqrt{a_{1}^{2}+a_{2}^{2}+a_{3}^{2}}$ and $\boldsymbol{\sigma}=(\sigma_{1},\sigma_{2},\sigma_{3})$ galindo90qme .

Consider a parameter $\varepsilon>0$ and let us take $X=\sigma_{1}$ and $Y=\sigma_{3}$. Then, direct application of (21) shows that

[TABLE]

with $\lambda=\sqrt{1+\beta^{2}}$. On the other hand, algorithm (17) applied to this case renders

[TABLE]

with (rather involved) explicit expressions for the real functions $f_{n}^{[C]}$, $f_{n}^{[S]}$, $g_{n}^{[C]}$, $g_{n}^{[S]}$. Notice that a non-vanishing term multiplying $i\,\sigma_{2}$ appears in the expression of $\Psi_{n}^{[C]}$, contrary to the exact solution (22). It turns out, however, that $g_{n}^{[C]}(\varepsilon,\lambda)$ goes to zero when $n\rightarrow\infty$. Moreover, if a series expansion in powers of $\varepsilon$ of these functions is computed, then we reproduce the exact expressions (22) up to the order considered. Thus, in particular, up to order $\varepsilon^{8}$ we get

[TABLE]

Example 2.

For our second example we consider two $10\times 10$ matrices $A$ and $B$ whose elements are random numbers in the range $(0,1)$ and normalized so that $\|A\|_{2}=\|B\|_{2}=1$. We are therefore outside the convergence domain for the Zassenhaus formula guaranteed by casas12eco . Then we compute numerically $\cos(A+B)$ via $X={\rm e}^{i(A+B)}$, $C=\mathrm{Re}\,X$ (with Mathematica) and $\Psi_{n}^{[C]}(A,B)$ as given by algorithm (17) for several values of $n$. Finally we determine the error $\log(\|\cos(A+B)-\Psi_{n}^{[C]}(A,B)\|)$ and represent this value as a function of $n$. In this way we obtain Figure 1. We clearly observe how the error decays exponentially with $n$. In other words, algorithm (17) provides a convergent expansion for $\cos(A+B)$ well beyond the domain obtained in casas12eco . A similar conclusion is achieved if one instead considers $\log(\|\sin(A+B)-\Psi_{n}^{[S]}(A,B)\|)$.

Although algorithm (17) is used here to approximate numerically $\cos(A+B)$, it is by no means intended to be used as a practical alternative to existing numerical procedures to compute the cosine of a matrix, but rather as an analytical tool in perturbative treatments. This being said, it could also be the case that for certain matrices $A$, $B$, computing the cosine and sine is a trivial task, whereas the evaluation of $\cos(A+B)$ and $\sin(A+B)$ is much more involved from a numerical point of view. The idea is then similar to splitting methods in the integration of differential equations blanes16aci : use $\cos(A)$, $\sin(A)$, $\cos(B)$ and $\sin(B)$ to approximate $\cos(A+B)$ and $\sin(A+B)$. In this situation, our procedure could be also competitive with other methods also from the numerical point of view.

5 Generalizations

Algorithm (17) can be applied of course to get other generalized trigonometric identities involving sums and products of the cosine and sine of $X+Y$. For the sake of illustration, we next collect the expansions of $\cos(X-Y)-\cos(X+Y)$ and $\sin(X-Y)+\sin(X+Y)$ up to $n=4$ obtained with our procedure. Specifically,

[TABLE]

and

[TABLE]

Notice that if $X$ and $Y$ commute, then $C_{n}=0$ for all $n\geq 2$ and the usual expressions

[TABLE]

are recovered.

In the trigonometric expansions obtained with algorithm (17) all the successive commutators appear to the right. This of course is due to the form of the Zassenhaus formula (8). There exists, however, an alternative, “left-oriented” expression of this formula, namely

[TABLE]

with different but related exponents casas12eco :

[TABLE]

It is then clear that, by using (24) a similar algorithm can be designed to get alternative expansions for $\cos(X+Y)$ and $\sin(X+Y)$, this time with commutators appearing to the left. Also invariant expressions with respect to the interchange $X\leftrightarrow Y$ can be easily generated by just considering a symmetrized version of the previous expansions. Thus, for instance, from the first line in eq. (1) we also get

[TABLE]

and thus

[TABLE]

Acknowledgements.

This work has been partially supported by Universitat Jaume I trough the project P1-1B2015-16. The second author also acknowledges Ministerio de Economía y Competitividad (Spain) for financial support through projects MTM2013-46553-C3 and MTM2016-77660-P (AEI/FEDER, UE).