Does Noether's theorem apply to constrained system? The Lagrangian of a constrained system will be
$$L-\lambda_1f_1-\lambda_2f_2-...\lambda_kf_k.$$
If a transformation will not affect the constrained Lagrangian, the there is some corresponding conservative quantities. E.g. If a the transformation is in position, then momentum conservative.
But consider the simple case that two balls $m_1,m_2$ are connected by a rod with length $l$ placed in the empty space. Then the Lagrangian is $$L=\frac{1}{2}(m_1\dot{\textbf{r}}_1^2+m_2\dot{\textbf{r}}_2^2)-\lambda(l^2-(\textbf{r}_1-\textbf{r}_2)^2).$$
With some Newtonian intuition, this system has apparently conservative momenta. But now if we transform a bit in $\textbf{r}_1\to \textbf{r}_1+\delta$, but $\textbf{r}_2$ is related. Are sure that $L$ is invariant under symmetric transformation of positions?
So I am think may be Noether's theorem does not hold for constrained system?
 A: When you are dealing with constrained systems, you need to make sure that the symmetries you propose respect this constraint. In the example that you have given, the transformation $r_1 \to r_1 + \delta r_1$ is not a symmetry of the system because it does not respect the constraint. That is, you cannot move only one of the balls by $\delta r_1$ without also moving the other ball because of the rod between them. Mathematically, let's decompose the Lagrangian as $L = L_{dynamics} + L_{constraint}$ then we have: 
$$ \delta L_{constraint} = 2 \lambda\, (\delta \vec r_1) \cdot (\vec r_1 -\vec r_2 )\neq 0$$
For example, if we varied both of the coordinates, you see that $r_i \to r_i + \delta r_i$ is only a symmetry of the system iff
$$ \delta L_{constraint} = 2 \lambda\, (\delta \vec r_1) \cdot (\vec r_1 -\vec r_2 ) - 2 \lambda\, (\delta \vec r_2) \cdot (\vec r_1 -\vec r_2 ) = 0 \iff \delta \vec r_1 = \delta \vec r_2$$
telling you that the only transnational symmetry permitted by the constraint is the one, where you display both of the balls the same amount. Intuitively this also makes sense. Note that this doesn't yet mean that it is a symmetry of the system but that it is the only kind of permitted symmetry; you need to additionally check that $\delta L_{dynmaics}$ is also zero or a total derivative.
A: *

*Noether's theorem doesn't discriminate between dynamical and auxiliary variables (meaning: with and without the presence of their corresponding time derivatives, respectively). As long as the extended action functional $S$ has a (quasi)symmetry, there is a conservation law. 

*In OP's concrete case, the infinitesimal symmetry is given by
$$ \delta {\bf r}_1~=~\varepsilon~=~\delta {\bf r}_2, \qquad  \delta\lambda~=~0~~=~\delta t, $$
and it leads to total momentum conservation.
A: To calculate the equations of motion and the constraint forces, we use the Euler Lagrange function:
\begin{align*}
 &\boxed{\frac{\partial L}{\partial \vec{r}}-\frac{d}{dt}\left(\frac{\partial L}{\partial \dot{\vec{r}}}\right)=\left(\frac{\partial \vec{R}}{\partial \vec{r}}\right)^T\vec{f}_a+\left(\frac{\partial \vec{F}_c}{\partial \vec{r}}\right)^T\vec{\lambda}}&(1)
\end{align*}
With:
$\vec{r}$ Unconstraint coordinates
$\vec{F}_c$ Constraint equations
$\vec{\lambda}$ Generalized constraint forces
$\vec{f}_a$ External forces
To get the EOM's for generalized coordinates, we transfer the velocity vector $\vec{\dot{r}}$ to:
$\vec{\dot{r}}\mapsto J\,\vec{\dot{q}}$
with $J$ the jacobi matrix and $\vec{q}$ the generalized coordinates.
to eliminate the constraint froces from equation (1) , we multiply equation (1) from the left with $J^T$  .
D’Alembert approach:
$J^T\,\left(\frac{\partial \vec{F}_c}{\partial \vec{r}}\right)^T=0\qquad(2)$
Conclusion:
Each transformation of the coordinates $\vec{r}\mapsto R\,\vec{r}$ must fulfill
equation (2)  
