I would like to propagate a two-mode squeezed vacuum state through a collection of beamsplitters. I'd like to represent this in the continuous variable quantum information formalism, so that I start with a covariance matrix, $\gamma$ and apply the many beamsplitter symplectic matrices, $S_1, S_2, S_3, ...$.
If I apply $S_1$ first, $\gamma_1=S_1^T\gamma S_1$, this of course makes $\gamma_1$ a diagonal matrix, so that I think it would then be improper to apply $S_2$. What is the appropriate way to do this in order to get a physically correct final covariance matrix?