# How to find a common eigenstate of commuting operators?

I have multiple different operators in matrix form and I need to find their common eigenstates. The challenge is that the common eigenstate is in a superposition of multiple states and isn't just a single eigenvector. Can someone give me an idea of how to do this? How would I figure out what a common eigenstate in the form of a linear combination of vectors is for multiple matrices? For context, I am using Python and NumPy.

Superposition is a linear combination of vectors which is a single vector. It may not be equal to some basis vector (from computational basis, for example), but this is not a big deal in general.

As for question, I suppose operators are normal (so we can apply the spectral theorem).

In the simplest case if operator $$A$$ has eigenvalue $$\lambda$$ with multiplicity 1, then the corresponding eigenvector $$v$$ will also be eigenvector for any commuting operator $$B$$. Since $$Av=\lambda v$$ and $$AB=BA$$ then $$ABv = BAv= B\lambda v \implies A(Bv) = \lambda (Bv) \implies Bv = \mu v$$ Here we obtained that $$Bv$$ is also eigenvector of $$A$$ with eigenvalue $$\lambda$$ and since it has multiplicity 1 it must be collinear to $$v$$.
So such eigenvector $$v$$ will automatically be common eigenvector for a set of commuting with $$A$$ operators.

In general case $$A$$ may not have eigenvalue with multiplicity 1. In this case $$A$$ has invariant subspaces (eigenspaces) $$H_{\lambda_i}, \text{dim}H_{\lambda_i}>1$$ such that $$Av = \lambda_i v, \forall v \in H_{\lambda_i}$$.They split space $$H = H_{\lambda_1} \oplus ... \oplus H_{\lambda_k}$$.
Similarly to the simple case, we can show that if $$AB=BA$$ then $$H_{\lambda_i}$$ are also invariant for $$B: B(H_{\lambda_i}) \subset H_{\lambda_i}$$. Now we can consider restrictions of operator $$B$$ on the subspaces $$H_{\lambda_i}$$. Every restriction has its own eigendecomposition. But every eigenspace of every such restriction of $$B$$ will be common eigenspace of $$A$$ and $$B$$ since $$H_{\lambda_i}$$ is eigenspace of $$A$$.
Actually, if we take the intersection of eigendecompositions $$\{H_{\lambda_i}\}$$ and $$\{H_{\mu_i}\}$$ correspondingly for $$A$$ and $$B$$ then we will obtain the decomposition of $$H$$ into common eigenspaces of $$A$$ and $$B$$.

Though it may not be that efficient to calculate intersection of all eigendecompositions for operators $$A_i$$.

I guess we can just pick some random numbers $$r_i$$ and find eigenvector of $$\sum_i r_i A_i$$. With high probability it will be the common eigenvector for commuting operators $$A_i$$.