Problem
Show that $\hat{L}^2=\hat{r}^2\hat{p}^2 + i\hbar \hat{r} \cdot \hat{p} - (\hat{r} \cdot \hat{p})^2 $ using Einstein notation and commutation relations in QM.
Attempt
We know that $\hat{L}_i = \epsilon_{ijk} \hat{x}_j\hat{p}_k$ as the definition of orbital angular momentum. We also know that $\hat{L}^2 = \hat{L}_i\hat{L}_i$. Therefore, we have $$ \hat{L}_i\hat{L}_i = \epsilon_{ijk} \hat{x}_j\hat{p}_k\,\epsilon_{ilm} \hat{x}_l\hat{p}_m = \epsilon_{ijk}\epsilon_{ilm} \hat{x}_j\hat{p}_k\hat{x}_l\hat{p}_m $$ Using the fact that $\epsilon_{ijk}\epsilon_{ilm} = \delta_{jl}\delta_{km}-\delta_{jm}\delta_{kl}$ gives us $$ \hat{L}_i\hat{L}_i = \delta_{jl}\delta_{km} \hat{x}_j\hat{p}_k\hat{x}_l \hat{p}_m - \delta_{jm}\delta_{kl} \hat{x}_j\hat{p}_k\hat{x}_l \hat{p}_m \\ = \hat{x}_l\hat{p}_m\hat{x}_l\hat{p}_m - \hat{x}_j \hat{p}_k \hat{x}_k \hat{p}_j $$ We now make use of the Heisenberg uncertainty relationship $$ [\hat{x}_l,\hat{p}_m] = \hat{x}_l\hat{p}_m-\hat{p}_m\hat{x}_l = i\hbar \delta_{lm} \Leftrightarrow \hat{x}_l\hat{p}_m = i\hbar \delta_{lm} + \hat{p}_m\hat{x}_l $$ which therefore gives us $$ \hat{L}_i\hat{L}_i = (i\hbar \delta_{lm} + \hat{p}_m\hat{x}_l) \hat{x}_l\hat{p}_m - \hat{x}_j\hat{p}_k\hat{x}_k\hat{p}_j = i\hbar \hat{x}_l\hat{p}_l + \hat{p}_m\ hat{x}_l \hat{x}_l\hat{p}_m - \hat{x}_j \hat{p}_k\hat{x}_k\hat{p}_j $$ where $\hat{x}_l\hat{p}_l = \hat{r} \cdot \hat{p}$ as $\hat{r} = \hat{x}_i$.
All that is left to show is that $$\hat{p}_m\hat{x}_l\hat{x}_l\hat{p}_m - \hat{x}_j \hat{p}_k\hat{x}_k\hat{p}_j = \hat{r}^2\hat{p}^2 - (\hat{r} \cdot \hat{p})^2 $$ but here's where I'm struggling and don't know how to move forward. I can't really see how to manipulate those expression to get what I want. I'd be glad if anyone could propose any ideas or tips on how to continue solving this problem.
Thanks in advance!