In the derivation of the LSZ reduction formula we write asymptotic states using ladder operators at $\pm \infty$. For example if we consider for a real scalar field $\phi$ then we have the formula
$\sqrt{2\omega_p}\big(a^{\dagger}_p(\infty)-a^{\dagger}_p(-\infty)\big)=-i\int d^4x\space e^{-ip\cdot x}(\partial^2+m^2)\phi(x)$
For a particlar scattering problem, say ${p_1,p_2 \rightarrow p_3,p_4}$ with $|i\rangle=\sqrt{2\omega_{p_1}}a^{\dagger}_{p_1}(-\infty)\sqrt{2\omega_{p_2}}a^{\dagger}_{p_2}(-\infty)|0\rangle$ and $|f\rangle=\sqrt{2\omega_{p_3}}a^{\dagger}_{p_3}(\infty)\sqrt{2\omega_{p_4}}a^{\dagger}_{p_4}(\infty)|0\rangle$, then we have scattering matrix amplitude $\langle f|S| i\rangle= \sqrt{2\omega_{p_1}}\sqrt{2\omega_{p_2}}\sqrt{2\omega_{p_3}}\sqrt{2\omega_{p_4}}\langle0|a_{p_3}(\infty) a_{p_4}(\infty)a_{p_1}^{\dagger}(-\infty)a_{p_2}^{\dagger}(-\infty)|0\rangle$.
We do the usual trick of adding the time order operator and adding in a bunch of zero terms to get $\sqrt{2\omega_{p_1}}\sqrt{2\omega_{p_2}}\sqrt{2\omega_{p_3}}\sqrt{2\omega_{p_4}}\langle0|T\{[a_{p_3}(\infty)- a_{p_3}(-\infty)][a_{p_4}(\infty)-a_{p_4}(-\infty)][a_{p_1}^{\dagger}(-\infty)-a_{p_1}^{\dagger}(\infty)][a_{p_2}^{\dagger}(-\infty)-a_{p_2}^{\dagger}(\infty)]\}|0\rangle$
Substituting in our formula at the start of this post we have $\langle f|S|i \rangle= i\int d^4x_1\space e^{-ip_1\cdot x_1}i\int d^4x_2\space e^{-ip_2\cdot x_2}i\int d^4x_3\space e^{ip_3\cdot x_3}i\int d^4x_4\space e^{ip_4\cdot x_4}\langle 0|T\{(\partial_1^2+m^2)\phi(x_1)(\partial_2^2+m^2)\phi(x_2)(\partial_3^2+m^2)\phi(x_3)(\partial_4^2+m^2)\phi(x_4)\}|0\rangle$
The next step is to 'pull' the differential operators out of the time ordering function, leaving the 4 point correlation function. The time ordering operator and differential operator don't commute so why is this ok? Why do the 'contact terms' not contribute to the S matrix?
