Concerning the bit about:
Why do we have the two-particle interaction "sandwiched" between the field operators, but the annihilation/creation operators do not follow the same pattern?
They do follow the same pattern, only when using field operators you are dealing with a continuous variable (i.e. $\textbf{x}$) which makes easier to write the relations in another way. To see this consider for example the kinetic energy operator $T$.
Discrete basis case:
Let $c_k^\dagger, c_k$ are the creation/annihilation operators of momentum eigenstates, and let $a_i^\dagger,a_i$ creation/annihilation operators for a generic (complete) set of states $\{ | i \rangle \}_i$. We have
$$ \tag{D} T =
\sum_{\textbf{k}} \frac{| \textbf{k} |^2}{2m} c_\textbf{k}^\dagger c_\textbf{k} = \sum_{i,j} t_{ij} a_i^\dagger a_j
$$
In both cases you are kind of expanding $T$ in terms of its matrix elements between Fock states. The "kind of" here is important because these are not really expansions in terms of an orthonormal basis for a number of reasons, for example because $a_i^\dagger a_j$ are not projection operators. However, it is still true (as you can readily verify using the (anti)commutation relations) that
$$ \langle \textbf{k} | T | \textbf{k}' \rangle
\equiv \langle 0 | c_\textbf{k} T c_{\textbf{k}'}^\dagger | 0 \rangle
= \delta_{\textbf{k}\textbf{k}'} \frac{ | \textbf{k} |^2}{2m}$$
$$ \langle i | T | j \rangle
\equiv \langle 0 | a_i T a_j^\dagger | 0 \rangle
= t_{ij}$$
Note that the "expanding $T$" interpretation starts to break down when you consider many-particle states, and that in these cases you start appreciating the difference between fermions and bosons. For example for two particle states you have for bosons:
$$ \langle \textbf{k}\textbf{q} | T | \textbf{k}\textbf{q} \rangle
\equiv \langle 0 | c_\textbf{q}c_\textbf{k} T c_\textbf{k}^\dagger c_\textbf{q}^\dagger | 0 \rangle
= \left( \frac{| \textbf{k} |^2}{2m} + \frac{| \textbf{q} |^2}{2m} \right) (1+\delta_{\textbf{k}\textbf{q}})$$
while for fermions:
$$ \langle \textbf{k}\textbf{q} | T | \textbf{k}\textbf{q} \rangle
\equiv \langle 0 | c_\textbf{q}c_\textbf{k} T c_\textbf{k}^\dagger c_\textbf{q}^\dagger | 0 \rangle
= \left( \frac{| \textbf{k} |^2}{2m} + \frac{| \textbf{q} |^2}{2m} \right) (1-\delta_{\textbf{k}\textbf{q}})$$
Continuous basis case
Using the field operators $\psi(\textbf{x})$ and $\psi^\dagger(\textbf{x})$ we write $T$ as:
$$ \tag{C} T = \int d^3 x \psi^\dagger(\textbf{x}) \left (\frac{-1}{2m} \nabla^2 \right) \psi(\textbf{x}) $$
so how is this similar to something like (D)?
To better see this lets assume we are in a 1D space (so $\textbf{x} \approx x_n$, $\psi(\textbf{x}) \approx \psi_n $) and switch to discrete varibles. We then have
$$ \nabla^2 \psi(\textbf{x}) \approx ( \psi_{n+2} - 2 \psi_{n+1} + \psi_n ) $$
and using this $T$ becomes
$$ \tag{R} T = \frac{-1}{2m} \sum_n \psi_n^\dagger ( \psi_{n+2} - 2 \psi_{n+1} + \psi_n ) = \sum_{n,m} \psi^\dagger_n T_{nm} \psi_m = \sum_{n,m} T_{nm} \psi^\dagger_n \psi_m $$
where we have defined
$$ T_{nm} \equiv \frac{-1}{2m} ( \delta_{n,m+2} - 2 \delta_{n,m+1} + \delta_{n,m} ) $$
as you can see (R) is again a form like (D), and the kinetic operator is no longer "sandwiched" between field operators (but it has acquired a more cumbersome form).
