The implementation of the Pauli matrices and their algebra were reviewed for Maple 2019, including the algebraic manipulation of nested commutators, resulting in faster computations using simpler and more flexible input. As it frequently happens, improvements of this type suddenly transform research problems presented in the literature as untractable in practice, into tractable.

As an illustration, we tackle below the derivation of the coefficients entering the Zassenhaus formula shown in section 4 of [1] for the Pauli matrices up to order 10 (results in the literature go up to order 5). The computation presented can be reused to compute these coefficients up to any desired higher-order (hardware limitations may apply). A number of examples which exploit this formula and its dual, the Baker-Campbell-Hausdorff formula, occur in connection with the Weyl prescription for converting a classical function to a quantum operator (see sec. 5 of [1]), as well as when solving the eigenvalue problem for classes of mathematical-physics partial differential equations [2].

References

[1] R.M. Wilcox. "Exponential Operators and Parameter Differentiation in Quantum Physics", Journal of Mathematical Physics, V.8, 4, (1967.    [2] S. Steinberg. "Applications of the lie algebraic formulas of Baker, Campbell, Hausdorff, and Zassenhaus to the calculation of explicit solutions of partial differential equations". Journal of Differential Equations, V.26, 3, 1977.    [3] K. Huang. "Statistical Mechanics". John Wiley & Sons, Inc. 1963, p217, Eq.(10.60).

Formulation of the problem

The Zassenhaus formula expresses ${\ⅇ}^{\mathrm{\λ}\left(A+B\right)}$ as an infinite product of exponential operators involving nested commutators of increasing complexity

$$\begin{gathered} {\ⅇ}^{\mathrm{\lambda }\left(A+B\right)} \= {\ⅇ}^{\mathrm{\lambda }A}\cdot {\ⅇ}^{\mathrm{\lambda }B}\cdot {\ⅇ}^{{\mathrm{\lambda }}^2{C}_2}\cdot {\ⅇ}^{{\mathrm{\lambda }}^3{C}_3}\cdot ..\. \\ \= {\ⅇ}^{\mathrm{\lambda }A}\cdot {\ⅇ}^{\mathrm{\lambda }B}\cdot {\ⅇ}^{-\frac{1}{2}{\mathrm{\lambda }}^2{\left[A,B\right]}_{-}}\cdot {\ⅇ}^{\frac{1}{6}{\mathrm{\lambda }}^3\left({\left[A,{\left[A,B\right]}_{-}\right]}_{-}+2\cdot {\left[B,{\left[A,B\right]}_{-}\right]}_{-}\right)}\cdot ..\. \end{gathered}$$

Given $A$, $B$ and their commutator $E={\left[A,B\right]}_{-}$, if $A$ and $B$ commute with $E$, $C_n\=0$ for $n\ge 3$ and the Zassenhaus formula reduces to the product of the first three exponentials above. The interest here is in the general case, when ${\left[A,E\right]}_{-}\ne 0$ and ${\left[B,E\right]}_{-}\ne 0$, and the goal is to compute the Zassenhaus coefficients $C_n$in terms of $A$, $B$ for arbitrary n. Following [1], in that general case, differentiating the Zassenhaus formula with respect to $\mathrm{\λ}$ and multiplying from the right by ${\ⅇ}^{-\mathrm{\λ} \left(A+B\right)}$ one obtains

$A+B=A+{\ⅇ}^{\mathrm{\lambda }A}B{\ⅇ}^{-\mathrm{\lambda }A}+{\ⅇ}^{\mathrm{\lambda }A}+{\ⅇ}^{\mathrm{\lambda }B}2\mathrm{\lambda }{C}_2{\ⅇ}^{-\mathrm{\λ}B}{\ⅇ}^{-\mathrm{\lambda }A}+ ..\.$

This is an intricate formula, which however (see eq.(4.20) of [1]) can be represented in abstract form as

$0=\left(\sum _{n=1}^{\mathrm{\infty }}\frac{{\mathrm{\lambda }}^n}{n!}\left\{A^n\,B\right\}\right)+2\mathrm{\lambda }\left(\sum _{m=0}^{\mathrm{\infty }}\sum _{n=0}^{\mathrm{\infty }}\frac{{\mathrm{\lambda }}^{n+m}}{n!m!}\left\{A^m\,B^n\,C_2\right\}\right)+3{\mathrm{\lambda }}^2\left(\sum _{k=0}^{\mathrm{\infty }}\sum _{m=0}^{\mathrm{\infty }}\sum _{n=0}^{\mathrm{\infty }}\frac{{\mathrm{\lambda }}^{n+m+k}}{n!m!k!}\left\{A^k\,B^m\,{C_2}^n\,C_3\right\}\right)\+ ..\.$

from where an equation to be solved for each $C_n$ is obtained by equating to 0 the coefficient of ${\mathrm{\lambda }}^{n-1}$. In this formula, the repeated commutator bracket is defined inductively in terms of the standard commutator ${\left[A,B\right]}_{-}$by

$\left\{A^0\,B^{}\right\}=B, \left\{A^{n+1}\,B\right\}={\left[A,\left\{A^n\,B^{}\right\}\right]}_{-}$

$\left\{A^0\,B^n\,{C_j}^{}\right\}=\left\{B^n\,{C_j}^{}\right\}, \left\{A^m\,B^n\,{C_j}^{}\right\}={\left[A,\left\{A^{m-1}\,B^n\,{C_j}^{}\right\}\right]}_{-}$

and higher-order repeated-commutator brackets are similarly defined. For example, taking the coefficient of $\mathrm{\λ}$ and ${\mathrm{\λ}}^2$ and respectively solving each of them for $C_2$ and $C_3$ one obtains

$C_2=-\frac{1}{2}{\left[A,B\right]}_{-}$

$C_3=\frac{1}{6}{\left[A,{\left[A,B\right]}_{-}\right]}_{-}+\frac{1}{3}{\left[B,{\left[A,B\right]}_{-}\right]}_{-}$

This method is used in [3] to treat quantum deviations from the classical limit of the partition function for both a Bose-Einstein and Fermi-Dirac gas. The complexity of the computation of $C_n$ grows rapidly and in the literature only the coefficients up to $C_5$ have been published. Taking advantage of developments in the Physics package for Maple 2019, below we show the computation up to $C_{10}$ and provide a compact approach to compute them up to arbitrary order.

Computing up to ${\mathit{C}}_{\mathbf{10}}$

Set the signature of spacetime such that its space part is equal to +++ and use lowercaselatin letters to represent space indices. Set also $A$, $B$ and $C_n$ to represent quantum operators

To illustrate the computation up to $C_{10}$, a convenient example, where the commutator algebra is closed, consists of taking $A$ and $B$ as Pauli Matrices which, multiplied by the imaginary unit, form a basis for the $\mathrm{\𝔰\𝔲}\left(2\right)$group, which in turn exponentiate to the relevant Special Unitary Group $\mathrm{SU}\left(2\right)$. The algebra for the Pauli matrices involves a commutator and an anticommutator

Assign now $A$ and $B$ to two Pauli matrices, for instance

Next, to extract the coefficient of ${\mathrm{\lambda }}^n$ from

$0=\left(\sum _{n=1}^{\mathrm{\infty }}\frac{{\mathrm{\lambda }}^n}{n!}\left\{A^n\,B\right\}\right)+2\mathrm{\lambda }\left(\sum _{m=0}^{\mathrm{\infty }}\sum _{n=0}^{\mathrm{\infty }}\frac{{\mathrm{\lambda }}^{n+m}}{n!m!}\left\{A^m\,B^n\,C_2\right\}\right)+3{\mathrm{\lambda }}^2\left(\sum _{k=0}^{\mathrm{\infty }}\sum _{m=0}^{\mathrm{\infty }}\sum _{n=0}^{\mathrm{\infty }}\frac{{\mathrm{\lambda }}^{n+m+k}}{n!m!k!}\left\{A^k\,B^m\,{C_2}^n\,C_3\right\}\right)\+..\.$

to solve it for $C_{n+1}$ we note that each term has a factor ${\mathrm{\lambda }}^m$ multiplying a sum, so we only need to take into account the first $n\+1$ terms (sums) and in each sum replace $\mathrm{\∞}$ by the corresponding $n-m$. For example, given $C_2=-\frac{1}{2}{\left[A,B\right]}_{-}\,$to compute $C_3$ we only need to compute these first three terms:

$0=\left(\sum _{n=1}^2\frac{{\mathrm{\lambda }}^n}{n!}\left\{A^n\,B\right\}\right)+2\mathrm{\lambda }\left(\sum _{m=0}^1\sum _{n=0}^1\frac{{\mathrm{\lambda }}^{n+m}}{n!m!}\left\{A^m\,B^n\,C_2\right\}\right)+3{\mathrm{\lambda }}^2\left(\sum _{k=0}^0\sum _{m=0}^0\sum _{n=0}^0\frac{{\mathrm{\lambda }}^{n+m+k}}{n!m!k!}\left\{A^k\,B^m\,{C_2}^n\,C_3\right\}\right)$

then solving for $C_3$ one gets $C_3=\frac{1}{3}{\left[B,{\left[A,B\right]}_{-}\right]}_{-}+\frac{1}{6}{\left[A,{\left[A,B\right]}_{-}\right]}_{-}$.

Also, since to compute $C_n$ we only need the coefficient of ${\mathrm{\lambda }}^{n-1}$, it is not necessary to compute all the terms of each multiple-sum. One way of restricting the multiple-sums to only one power of $\mathrm{\lambda }$ consists of using multi-index summation, available in the Physics package (see Physics:-Library:-Add). For that purpose, redefine sum to extend its functionality with multi-index summation

Now we can represent the same computation of $C_3$ without multiple sums and without computing unnecessary terms as

$0=\left(\sum _{n=1}\frac{{\mathrm{\lambda }}^n}{n!}\left\{A^n\,B\right\}\right)+2\mathrm{\lambda }\left(\sum _{n\+m=1}\frac{{\mathrm{\lambda }}^{n+m}}{n!m!}\left\{A^m\,B^n\,C_2\right\}\right)+3{\mathrm{\lambda }}^2\left(\sum _{n\+m\+k=0}\frac{{\mathrm{\lambda }}^{n+m+k}}{n!m!k!}\left\{A^k\,B^m\,{C_2}^n\,C_3\right\}\right)$

Finally, we need a computational representation for the repeated commutator bracket 

$\left\{A^0\,B\right\}=B, \left\{A^{n+1}\,B\right\}={\left[A,\left\{A^n\,B^{}\right\}\right]}_{-}$

One way of representing this commutator bracket operation is defining a procedure, say F, with a cache to avoid recomputing lower order nested commutators, as follows

For example,

We can set now the value of $C_2$

and enter the formula that involves only multi-index summation

from where we compute $C_3$ by solving for it the coefficient of ${\mathrm{\lambda }}^2$, and since due to the multi-index summation this expression already contains ${\mathrm{\λ}}^2$ as a factor,

In order to generalize the formula for H for higher powers of $\mathrm{\λ}$, the right-hand side of the multi-index summation limit can be expressed in terms of an abstract N, and H transformed into a mapping:

Now we have

The following is already equal to (11)

In this way, we can reproduce the results published in the literature for the coefficients of Zassenhaus formula up to $C_5$ by adding two more multi-index sums to (13). Unassign $C$ first

We compute now up to $C_5$ in one go

The nested-commutator expression solved in the last step for $C_5$ is

With everything understood, we want now to extend these results generalizing them into an approach to compute an arbitrarily large coefficient $C_n$, then use that generalization to compute all the Zassenhaus coefficients up to $C_{10}$. To type the formula for H for higher powers of $\mathrm{\λ}$ is however prone to typographical mistakes. The following is a program, using the Maple programming language, that produces these formulas for an arbitrary integer power of $\mathrm{\λ}$:

This Formula program uses a sequence of summation indices with as much indices as the order of the coefficient $C_n$ we want to compute, in this case we need 10 of them

To avoid interference of the results computed in the loop (17), unassign $C$ again

Now the formulas typed by hand, used lines above to compute each of $C_2$, $C_3$ and $C_5$, are respectively constructed by the computer

Construct then the formula for $C_{10}$ and make it be a mapping with respect to N, as done for $C_5$ after (16)

Compute now the coefficients of the Zassenhaus formula up to $C_{10}$ all in one go