Simplifying cubic and quartic interaction operators

Question

the probability of decay of a particle into two due to cubic in $\hat{x}$ interaction is given by $\langle f \mid \hat{x}^3\mid i\rangle$. The $\hat{x}^3$ term is written in the basis of ladder operators and takes the form: $$\sum_{k,k',q}\langle f \mid(a_k +a_k^\dagger) (a_{k'} + a_{k'}^\dagger) (a_q + a_q^\dagger)\mid i\rangle \delta (E_f-E_i)$$ where the states are number states of the form $\mid n_k \rangle \otimes \mid n_{k'} \rangle \otimes\mid n_q \rangle$. The are other factors of matrix elements and c-numbers. But let's get to the point. How do I simplify this expression and minimize the effort of considering all 8 terms? Properties like symmetry, cases where $k=k'=q$, etc., and other conditions might simplify the expression and make the calculation easier. Any words of wisdom or reference to sources like books, papers, and review articles will be a great help. Thank you.

Edit: @doublefelix pointed out that no commutator was prescribed. So let the commutator be $$[a_k,a_q^\dagger]=\delta_{k,q}$$.

Answer 1

You can reach some nice simplifications by expanding out the bra and kets as

$$|i \rangle =|n_k \rangle \otimes |n_{k'} \rangle \otimes |n_q \rangle $$

$$|f \rangle =|m_k \rangle \otimes |m_{k'} \rangle \otimes |m_q \rangle $$

Since each $a$ acts only on bras or kets for its own particle, the summand then reads

$$\langle m_k |(a_k +a_k^\dagger)| n_k \rangle \langle m_{k'} |(a_{k'} + a_{k'}^\dagger)| n_{k'} \rangle \langle m_q | (a_q + a_q^\dagger)|m_q \rangle \delta (E_f-E_i)$$

Consider just the first factor:

$$\langle m_k |(a_k +a_k^\dagger)| n_k \rangle$$

If you have $a |n_k \rangle = c_1 |n_{k-1} \rangle$ and $a^{\dagger} |n_k \rangle = c_2 |n_{k+1} \rangle$ (what $c_1$ and $c_2$ are can again depend on conventions of how $a, a^{\dagger}$ are defined, so it's best if I leave those for you to fill), then this becomes

$$ = c_1 \delta_{m_k, n_{k-1}} + c_2 \delta_{m_k, n_{k+1}}$$

and similar can be done for each of the three factors.

Lastly, notice that depending on what the dependence between $E_f$ and the particle numbers are, you can potentially express that in terms of a delta of the particle numbers as well. I would need the function for $E_f$ to do that step, though depending on your application you may or may not need to.