Non-Abelian gauge theory action with external currents: gauge invariance?

I am currently trying to work out the gauge invariance of the Yang-Mills action coupled to an external source and I noticed something weird.

The Yang-Mills action in presence of an external source $J^a_\mu$ reads

$S = \intop_x - \frac{1}{4} F^{a\mu\nu}F^a_{\mu\nu} + J^{a\mu} A^a_\mu = S_{YM} + S_{J}$.

The pure Yang-Mills part $S_{YM}$ is invariant under gauge transformations $V(x) \in \text{SU(N)}$ via

$\mathcal{A}_\mu = A^a_\mu t^a \rightarrow V \mathcal{A}_\mu V^\dagger - i V \partial_\mu V^\dagger$.

However, $S_J$, the part where we coupling to the external source happens, is not gauge-invariant. I know that the source should transform gauge covariantly via

$\mathcal{J_\mu} = J^a_\mu t^a \rightarrow V \mathcal{J}_\mu V^\dagger$,

but this does not render $S_J$ gauge-invariant. In fact I obtain

$S_J \rightarrow S_J - 2i\intop_x tr(V \mathcal{J_\mu} \partial^\mu V^\dagger)$,

which does not vanish. What's exactly happening here?

Note: if we look at the Abelian case and repeat the same steps we could actually make the action invariant by requiring charge conservation $\partial_\mu J^\mu=0$.

I know that the QCD action is gauge-invariant and in that case I obtain an extra term from the matter fields which cancels against terms from $S_J$. Does this mean that Yang-Mills with external sources is meaningless? Does Yang-Mills only make sense if the current $J^a_\mu$ comes from "dynamical" fields like quarks?

It is not possible to couple a conventional $c$-number source to a quantized non-abelian gauge field and maintain gauge invariance. For a current $J^\mu=\lambda^a J^\mu_a$, gauge invariance requires covariant current conservation $$0= \nabla_\mu J^\mu\equiv \partial_\mu J^\mu + [A_\mu,J^\mu]=0$$ and this requires a $c$-number(the first term) to equal a term containing an the operator $A_\mu$. It is because of this problem that external sources are introduced into a gauge theory as Wilson Loops $$W= P\exp\{ i\int A_\mu dx^\mu\}$$ where $P$ denoters a path-ordered integral of the matrix valued ($A_\mu = A_\mu^a \lambda^a$) gauge field.
Though mike stone have given a excellent answer, I want to collect some basic facts here. For Yang-Mills with external current $$S_{\rm YM}= \int d^4x\Big[-\frac{1}{4}F^a_{\mu\nu}F^{\mu\nu a}-J^{\mu a}A_\mu^a\Big]$$ The gauge transformation is $$A_\mu(x)\,\,\, \longrightarrow \,\,\,\Omega A_\mu(x) \Omega^{-1}+\frac{i}{g} \Omega \partial_\mu\Omega^{-1}$$ and its infinitesimal form is $$A_\mu(x) \,\,\,\longrightarrow \,\,\, A_\mu(x)-\partial_\mu\omega +ig[A_\mu(x),\,\omega(x)]$$ The classical equation of motion is \begin{align} D_\mu F^{\mu\nu }=J^{\nu} \,\,\,\,\,\hbox{or}\,\,\,\,\,\, D_\mu F^{\mu\nu a }=J^{\nu a} \end{align} which implies the gauge transformation of the external current: \begin{align} J^{\mu}\,\,\,\longrightarrow\,\,\,\Omega J^{\mu} \Omega^{-1} \end{align} since \begin{align} F^{\mu\nu}\,\,\,\longrightarrow\,\,\,\Omega F^{\mu\nu} \Omega^{-1}, \,\,\,\,\,\,\,\,\, D_\alpha F^{\mu\nu}\,\,\,\longrightarrow\,\,\,\Omega D_\alpha F^{\mu\nu} \Omega^{-1} \end{align} For the external current, we have \begin{align} D_\nu J^{\nu}=D_\nu D_\mu F^{\mu\nu }=-\frac{1}{2}[D_\mu,\,D_\nu]F^{\mu\nu} =\frac{ig}{2}[F_{\mu\nu},\,F^{\mu\nu}]=0 \end{align}
Gauge invariance of $${\rm tr} (j^\mu A_\mu)$$ requires \begin{align} \partial_\mu J^\mu=0 \end{align} so we need independent conditions \begin{align} \partial_\mu J^\mu=0,\,\,\,\,\,\,\,\,\, [A_\mu,\,J^\mu]=0 \end{align} Thus at some extend we make the external currrent term gauge invariance.