4 i d d t δa ( t ) = [ − Δ$_{C}$ + N 〈 U 〉 − iκ $^{]}$δa ( t ) + Nα$_{ss}$ ∫ ϕ ( x ) U ( x ) $^{[}$δ ˜ Ψ( x, t ) + δ ˜ Ψ † ( x, t ) $^{]}$dx + iξ ( t ) , (7a) i ℏ ∂ ∂t δ ˜ Ψ( x, t ) = [ ℏ $^{2}$Δ 2 m + V$_{ext}$ ( x ) + ℏ | α$_{ss}$ | $^{2}$U ( x ) − µ + gNϕ $^{2}$( x )$^{]}$ δ ˜ Ψ( x, t ) + ℏ U ( x ) ϕ ( x ) $^{[}$α ∗ $_{ss}$δa ( t ) + α$_{ss}$δa † ( t ) ] + gNϕ $^{2}$( x ) $^{[}$δ ˜ Ψ † ( x, t ) + δ ˜ Ψ( x, t ) $^{]}$. (7b) i ∂ ∂t R ( t ) = M R ( t ) + iZ ( t ) , (8a) joints δa † and δ ˜ Ψ † ; this approach is familiar from the Bogoliubov-de Gennes theory of superfluidity. Adopting the Bogoliubov-de Gennes way we gather the fluctuations into the following column vector R = ( δa, δa † , δ ˜ Ψ , δ ˜ Ψ † ) $^{T}$, where the superscript T stands for transposition, and noises to the other column vector Z = ( ξ, ξ † , 0 , 0) $^{T}$. Eqs. (7) can now be cast into the closed form where the zeroth order terms cancel, since they fulfill Eq. (4b) and (5) with d α$_{ss}$/ d t = 0. We have introduced δ ˜ Ψ( x, t ) = N − 1 / $^{2}$δ Ψ( x, t ) e iµt/ ℏ . A closer look on Eqs. (7) reveals that the time evolu- tion of the annihilation operators couple to those of the creation operators. It is a consequence of the complex nature of the photonic and particle fields. In order to diagonalize Eqs. (7) one can choose two equivalent ways. Either one can separate the complex quantities into real and imaginary parts and study their time evolution; this is the approach mainly used in optomechanical studies. The other way is to diagonalize the set of equations not just containing δa and δ ˜ Ψ, but also their hermitian ad- with eigenvectors of M , i.e., with ω$_{k}$ being the corresponding eigenvalue of M . The fluctuation operator R can be expanded with the help of the eigenvectors (if they form a complete set) Since δa is not independent of δa † and similarly δ ˜ Ψ is also not independent of δ ˜ Ψ † , the matrix M has an important symmetry property. It is a consequence that the effect of hermitian conjugation of R can be obtained with a linear transformation C that swaps the first row with the second one and simultaneously the third one with the fourth, so R † = C R . It directly follows from Eq. (8a) and from this symmetry property that M =    A 0 Nα$_{ss}$X Nα$_{ss}$X 0 − A ∗ − Nα ∗ $_{ss}$X − Nα ∗ $_{ss}$X α ∗ $_{ss}$Y ( x ) α$_{ss}$Y ( x ) ℏ − $^{1}$[ H$_{0}$ + gNϕ $^{2}$( x )] ℏ $^{1}$gNϕ$^{2}$( x ) − α $_{ss}$Y ( x ) − α$_{ss}$Y ( x ) − ℏ − $^{1}$gNϕ$^{2}$( x ) − ℏ − $^{1}$[ H$_{0}$ + gNϕ $^{2}$( x )]   $_{}$, (8b) H$_{0}$ = − ℏ $^{2}$Δ 2 m + V$_{ext}$ ( x ) + ℏ | α$_{ss}$ ( x ) | $^{2}$U ( x ) − µ + gNϕ $^{2}$( x ) , (8c) X · f ( x ) = ∫ ϕ ( x ) U ( x ) f ( x ) dx, (8d) Y ( x ) = U ( x ) ϕ ( x ) , (8e) A = − Δ$_{C}$ + N 〈 U 〉 − iκ. (8f) M r ( k ) = ω$_{k}$ r ( k $^{)}$, (10) R ( t ) = ∑ k ρ$_{k}$ ( t ) r ( k $^{)}$, (11) C · M · C = − M ∗ . (9) ρ$_{k}$ ( t ) = $^{(}$l ( k $^{)}$, R ( t ) $^{)}$, (12) M † l ( k ) = ω ∗ k l ( k $^{)}$. (13) with ρ$_{k}$ being the operator, or expansion coefficient, of the quasi-normal mode k . The operator ρ$_{k}$ is given by: with ( · , · ) is the Euclidean scalar product, and l ( k ) is the left eigenvector of M , defined as In order to study the correlations of the fluctuations one has to determine the time evolution of the fluctua- tion operators first. One has to introduce quasi-normal modes that diagonalize Eq. (8a) and therefore have a simple time evolution. Let us denote by r ( k ) the right