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? Say, that I don't sum two eigenvalues a_i and b_i that don't "correspond to the same row"? I don't know if my question makes sense or if I'm properly explaining myself.

Any indication will be much appreciated, also if you know about any easier method to obtain the eigenvalues. Thanks and sorry for the notation, I hope it is not too confusing.