The gauge covariant derivative and its substitution I was wondering wether it would make a difference (in general) if I were to were introduce the gauge covariant derivative
$$D_\mu=\partial_\mu+ieA_\mu$$
In the Lagrangian density and then derive the equations of motion. OR if I were to derive the equations of motion and then introduce the gauge covariant derivative.
I checked this for the Klein-Gordon and Dirac equations where it made no difference at all.
My guess is that it would make no difference at all since we introduce this one to make our system (action, equations of motion, ...) gauge-invariant. So it doesn't really matter at which point I'd introduce the gauge covariant derivative.
But I'm unsure if this reasoning is correct ...
To make the question really basic I guess it comes to the question: "Do we always need to introduce the local symmetries in the action, or can we wait with that untill the equations of motion ?"

EDIT: I'm also a big believer of the fact that you just substitute the covariant derivative into the action to introduce your interacitons (this is the way to go since Noether's theorem works with the action and this theorem guarantees the conserved quantities).
Now I'm studying "Quarks & Leptons: An Introductory Course in Modern Particle Physics" by Halzen & Martin, and when they study QED of spin-0 particles (chapter 4 in my edition which is the 1984 edition apperently) they introduce the gauge covariant derivative in the EQUATIONS OF MOTION in stead of the action.
Somewhere this has some logic to me since if your system is invariant under some kind of symmetry, the the equations of motion should be invariant for this kind of transformation as well (otherwise two solutions related by a symmetry-transformation wouldn't yield the same result). I think this is what they're proving here (arxiv.org/abs/0907.2301), but I'm not sure to which extent this goes ?
Of course (as said in a comment to one of the answers), if you want to give your vectorfield $A_\mu$ a general dynamics you need Yang-Mills theory ans some kind of kinetic term for this field, in this case you need a Lagrangian to "unify" both dynamics.
EDIT2: @QMechanic's answer also gives a very nice scematic version of the question (just for reference).

 A: I am not sure of what you mean by "make a difference".
The covariant derivative is introduced in the Lagrangian density to add Local Gauge Symmetry to the Action. A given field theory is described by its Action. Hence if the Action does not have a particular symmetry, you cannot introduce the symmetry later. 
By "making a difference", if you mean a difference to the physical consequences, then it will make a difference if you do not introduce the covariant derivative in the Lagrangian density itself. Taking your example of the complex scalar field theory, the complete Lagrangian density is,
$(1/2)(D_{\mu}\phi)^{*}D^{\mu}\phi-(1/4)F_{\mu\nu}F^{\mu\nu}-V(\phi,\phi^{*})$
The first term gives the kinetic term for your scalar field and an interaction term with $\phi, \phi^{*}, A^{\mu}$ coupled. The second term is the kinetic term for $A^{\mu}$, while the term is an interaction term with $\phi$ and $\phi^{*}$.
If you do not have the covariant derivative in your Lagrangian density, you do not have the interaction between the complex scalar field and the gauge field, and hence it is a lot of difference to the physics. In QED that would mean that Compton scattering and a lot that happens in nature will not happen.  
A: OP is pondering if extremization of the action commutes with the minimal coupling (MC) recipe, i.e., if the following diagram would commute:  
$$ \begin{array}{ccc}
\text{Lagrangian density} && \text{Lagrangian density} \cr
{\cal L}(\phi(x),\partial\phi(x),A(x),F(x),x) &\stackrel{\text{MC}}{\longrightarrow} 
& {\cal L}(\phi(x),D\phi(x),A(x),F(x),x) \cr\cr
\text{extremization}\downarrow &&\downarrow \text{extremization}\cr\cr
{\cal E}(\phi(x),\partial\phi(x),\partial^2\phi(x),A(x), &\stackrel{\text{MC}}{\longrightarrow} 
& {\cal E}(\phi(x),D\phi(x),D^2\phi(x),A(x),\cr
F(x),\partial F(x),x)=0 &  
& F(x),D F(x),x)=0\cr
\text{EOMs} && \text{EOMs} 
\end{array} $$
Here $\phi(x)$ is shorthand for various matter fields. 
No. There are several problems with this proposal. E.g. the extremization wrt. $A_{\mu}$ of (scalar) QED before minimal coupling would produce Maxwell's equations without matter source terms. A later minimal coupling procedure would not repair this. This fact is related to Sourav's answer.
