I want to calculate the eigenvalues of a matrix $H=A+B$, where $A$ and $B$ are Kronecker/tensor products of Pauli matrices, that is, $$A=\bigotimes_{i=1}^{2N^2}P_i\qquad B=\bigotimes_{i=1}^{2N^2}Q_i$$ where $P_i\in\{I,X\}$ and $Q_i\in\{I,Z\}$, $\forall i,j=1,...,2N^2$ and where $I=[[1,0],[0,1]]$, $X=[[0,1],[1,0]]$ and $Z=[[1,0],[0,-1]]$. Also, the number of $X's$ and $Z's$ in $A$ and $B$ come in pairs, so that $A$ and $B$ commute (example: $A=I\otimes I\otimes X\otimes X\otimes I\otimes I\otimes X\otimes X$, and $B=Z\otimes Z\otimes I\otimes I\otimes I\otimes I\otimes Z\otimes Z$). My argument to calculate the eigenvalues is the following: If we first consider $C$, $D$ random matrices of the same dimension $n$, and thanks to the properties of the Kronecker product, we know that the eigenvalues of $M=C\otimes D$, are $eig(M)=\{c_1d_1,...,c_1d_n,c_2d_1,...,c_2d_n,...,c_nd_1,...,c_nd_n\}$, where $c_i$ and $d_j$ are the eigenvalues of $C$ and $D$, respectively. We can extrapolate this result to calculate the eigenvalues of $A$ and $B$ as all the possible combinations of products of the eigenvalues of $A_i$ and $B_i$, respectively. Also, $[A,B]=[A,H]=[B,H]=0$, so all the terms commute. Because $H,A$ and $B$ are also normal, they can be simultaneously diagonalised. Meaning that the eigenvalues of $H$, $h_i$, $\forall i=1,...,2^{2N^2}$, can be expressed as \begin{equation} h_i=a_i+b_i,\quad\forall i=1,...,2^{2N^2}, \end{equation} where $a_i$ and $b_i$ are eigenvalues of $A$ and $B$, respectively. With the two above results, one should be able to calculate $h_i$ as $$h_i=a_i+b_i=p_{1,j}p_{2,k}...p_{2N^2,l}+q_{1,j}q_{2,k}...q_{2N^2,l}$$ where $j+k+...+l=i$, $\forall i=1,...,2^{2N^2}$, where $p_{m,i}$ and $q_{m,i}$ are eigenvalues of $P_m$ and $Q_m$, respectively, $\forall m=1,...,2N^2$. However, I'm a bit concerned about the order of the eigenvalues of $A$ and $B$. As there are many ways to diagonalize every matrix (in the sense that they can be ordered in many ways), I am not sure that the sum made to obtain $h_i$ will give me an appropriate result. So, my question is, how can I make sure that the eigenvalues obtained using this method are the proper ones? Also, what if we consider the generalisation: $$H=A+B=\sum_{i=1}^{N^2}A_i+\sum_{i=1}^{N^2}B_i=\sum_{i=1}^{N^2}\bigotimes_{j=1}^{2N^2}P_{i,j}+\sum_{i=1}^{N^2}\bigotimes_{j=1}^{2N^2}Q_{i,j}$$, where $P_{i,j}\in\{I,X\}$, $\forall i,j$ and $Q_{i,j}\in\{I,Z\}$, $\forall i,j$. That is, $A$ is a sum of matrices which are Kronecker products of Pauli matrices, and same for $B$. My question, again, is how to calculate the eigenvalues and their degenerations. Thanks in advance