Deriving effective potential energy from the Lagrangian of a two-body system I'm having some issues understanding how the effective potential energy of a two-body system is derived from the Lagrangian of the system. Specifically my issue is with one step... 
Suppose we are analyzing the system in the centre of mass frame with a reduced mass $\mu$, radial separation $r$, and angular speed $\dot{\phi}$. Then the Lagrangian can be expressed as: 
$$\mathscr{L}= \frac{1}{2}\mu(\dot{r}^2+r^2\dot{\phi}^2)-U(r)$$
I understand that $\phi$ is a cyclic coordinate, since $\frac{\partial \mathscr{L}}{\partial \phi}=0$
, and thus the angular momentum $l$ is conserved. Since 
$l=\mu r^2 \dot{\phi}$
, we can use this to substitute in for 
$\dot{\phi}$
 and make the Lagrangian one-dimensional. 
My problem is this: if we substitute in the 
$\dot{\phi}$
 before using the Euler-Lagrange equation to derive the equation of motion, we get the wrong result. There is a minus sign where there shouldn't be! $$\begin{align}
\mathscr{L}
&= \frac{1}{2}\mu\dot{r}^2+\frac{1}{2}\mu r^2 \biggl(\frac{l}{\mu r^2}\biggr)^2-U(r) \\
&= \frac{1}{2}\mu\dot{r}^2+\frac{l^2}{2\mu r^2}-U(r) \ \Rightarrow \ U_{eff}(r) = U(r) - \frac{l^2}{2\mu r^2}
\end{align}\\ 
\frac{d}{dt} \frac{\partial \mathscr{L}}{\partial \dot{r}}=\frac{\partial \mathscr{L}}{\partial r} \ \Rightarrow\ \ \mu\ddot{r}=\color{red}{-}\frac{l^2}{\mu r^3}-U'(r)$$
 Whereas if we reverse the order of these steps, we get the correct result. $$\mathscr{L}= \frac{1}{2}\mu\dot{r}^2+\frac{1}{2}\mu r^2\dot{\phi}^2-U(r) \\ 
\begin{align}\frac{d}{dt} \frac{\partial \mathscr{L}}{\partial \dot{r}}=\frac{\partial \mathscr{L}}{\partial r} \Rightarrow\mu\ddot{r}
&= \mu r \dot{\phi}^2 -U'(r)\\
&= \mu r \biggl(\frac{l}{\mu r^2}\biggr)^2 -U'(r)\\
&= \frac{l^2}{\mu r^3} -U'(r)\\
\end{align}$$
What's going on here? Why does the order of these operations matter and what does this mean physically?
POST-EDIT: 
In response to being marked as a duplicate of: Lagrangian of an effective potential. I did not fully understand the answer to this question as it is beyond the scope of what I've learned (as an undergraduate sophomore).
As for: How can you solve this "paradox"? Central potential, the title isn't very descriptive of what's being asked, so I didn't find this post. However, the accepted answer is very succinct and given with less advanced theory. Thank you for bringing it to my attention!
 A: This is because of how the partial derivatives are defined. When you consider the partial derivative with respect to r you don’t have to consider the chain rule in $ \phi_{(r)} $. This is because you consider r and $ \phi $ as independent coordinates of the system. You dont apply the chain rule to independent coordinates.
Try applying the chain rule to $ \frac{\partial}{\partial r} \left( \frac{1}{2} \mu r^2 \dot{\phi}_{(r)}^2 \right) $ with $ \dot{\phi}=\frac{l}{\mu r^2} $ and see that you will arrive to the same wrong answer.
A: Lagrange’s equations are an independent set of equations, which is obtained by deriving the Lagrangian in a particular way;
For the different equations to be consistent to one another, the function you start from in deriving must be the same for every equation. 
What you did is basically plugging the result of the first equation of motion (the one in $\phi$) in the Lagrangian before you calculated the second equation of motion (the one in $r$), de facto modifying it adding further components in $r$, that otherwise would not have taken part in the derivation, and that obviously altered the final result;
