Matrix exponentiation of a Kronecker product of Pauli matrices I need to numerically compute the matrix-exponential of a Kronecker product of Pauli matrices (including the identity).
For example,
$$
\exp( X \otimes Y \otimes I \otimes Z \;\otimes \;... )
$$
or generally
$$
\exp \bigotimes\limits_j \sigma_j, \; \; \text{where} \;\; \sigma_j \in \{I,X,Y,Z\}.
$$
I can construct the Pauli product easily enough. However, I want to avoid implementing a numerical routine for exponentiating general (or just square) complex matrices, since the Hermitian & unitary matrix resulting from Pauli products is very particular. For example, it's clear that the resulting matrix will only contain the elements $\{ \pm 1, \pm i, 0 \}$.
Surely this begs an analytic form, or at least a significantly simplified numerical routine, for computing the matrix exponential!
 A: Since Pauli products and their Kronecker products are involutory, this is actually straight forward.
Since
$$
\sigma^n = \begin{cases} I, &n \text{ even} \\ \sigma, &n \text{ odd} \end{cases}
$$
and
$$
(A \otimes B) \; (C \otimes D) = (A \, C) \otimes (B \, D),
$$
then
$$ 
\implies
\left(\bigotimes_j \sigma_j\right)^n = \bigotimes_j \begin{cases} I, &n \text{ even} \\ \sigma_j, &n \text{ odd} \end{cases}.
$$
By the Taylor expansion of the exponential, we have
$$
\exp \bigotimes_j \sigma_j = \sum\limits_{n=0}^\infty \frac{1}{n!} \left( \bigotimes_j \sigma_j \right)^n
= \sum\limits_{n=0}^\infty \frac{1}{n!} \bigotimes_j\begin{cases} I, &n \text{ even} \\ \sigma_j, &n \text{ odd} \end{cases}
$$
and separating the terms by parity,
$$
= \sum\limits_{n=0}^\infty \frac{1}{(2n)!} \bigotimes_j I  \;\;+ \;\;\sum\limits_{n=0}^\infty \frac{1}{(2n+1)!} \bigotimes_j \sigma_j
$$
and "recognising" expansions, we conclude
$$
\exp \bigotimes_j \sigma_j = \cosh(1) \bigotimes_j I + \sinh(1) \bigotimes_j \sigma_j.
$$
