Question about source terms in scalar quantum field theory I'm having a bit of a mental block when trying to interpret the inhomogeneous Klein-Gordon equation $$(\Box +m^{2})\phi(x,t)=j(x,t)$$
In particular, how does one interpret the term on the right-hand side of the equation, $j(x,t)$ as a source term for the scalar field $\phi(x,t)$? 
Is it analogous to the situation in classical electromagnetism? For example, the flux of an electric field through an closed surface, $S$ $\bigl(\oint_{S}(\mathbf{E}\cdot\mathbf{n})\,dS$, where $\mathbf{n}$ is the unit normal vector to the surface $S\bigr)$ is proportional to the total charge density $\rho$ in the enclosed volume $V$. This can be written in differential form as $$\nabla\cdot\mathbf{E}=\frac{\rho}{\varepsilon_{0}}$$ and so the charge density $\rho$ is interpreted as a source for a spatially varying electric field. 
 A: The source terms have different interpretations depending on context.
Within classical field theory, they are called sources because they affect the solutions of the equation; in particular, they force the field to be nonzero. Take the equations of electrostatics, for example:
$$\nabla \cdot \mathbf{E} = \rho \qquad \nabla \times \mathbf{E} = 0$$
If $\rho = 0$, then $\mathbf{E} = 0$ is a possible solution. Of course there are others, but if, for example, we fix as boundary conditions that the field go to zero at infinity, then that's the only one. On the other hand, if $\rho$ is nonzero, $\mathbf{E} = 0$ is not a solution. In this case the interpretation of the "sources" is pretty easy, because we know what charges and currents are, and it is intuitive to say that they are the source of electromagnetic fields.
The Klein-Gordon equation doesn't correspond to a familiar classical situation, but we deal with sources in the same way. The free equation $(\partial^2 + m^2)\phi = 0$ is a kind of wave equation; if we solve by Fourier transform, the equation says that the wave vectors satisfy the dispersion relation $\omega^2 = \mathbf{k}^2 + m^2$. But again, if $j$ is nonzero the $\phi = 0$ solution disappears. The physical intuition is trickier because I can't give one unless I know what $\phi$ represents, but we can take a special case: if $m=0$ we recover the usual wave equation, and we know that the electromagnetic potentials satisfy that one. So we could interpret $j$ as being the "charge" for a potential $\phi$ that satisfies a slightly weird dispersion relation. You can try plugging in simple functions in place of $j$ (for example, one that varies sinusoidally in time) and see what you get.
In quantum field theory, there are two ways of thinking about the sources that I know of. One is to interpret them as background fields. This is an approach taken, for example, by Peskin and Schroeder when discussing electron scattering of a proton. If we make the approximation that the proton is stationary during scattering, we can treat its field as classical, and only quantize the electron field. This is done by adding to the free electron Lagrangian a term $\mathcal{L}_{\text{int}} = e \bar{\psi} \gamma_\mu \psi A^\mu_\text{cl}$, where $A^\mu_\text{cl}$ is a fixed classical field. This adds certain kinds of vertices to the Feynman diagrams, similar to the ones in fully quantum QED. One can derive Rutherford's formula with this formalism. This is covered, for example, in exercise 4.4 of Peskin and Schroeder, or section 6.4 of Weinberg (warning: if you read Weinberg you might end up more confused than when you started!).
The other use of sources is inside a path integral. In this formalism, the time-ordered correlation function of a number of fields is given by Feynman's path integral:
$$\langle 0 | T \phi(x_1) \cdots \phi(x_n) | 0\rangle = \int \mathcal{D}\phi\  e^{i\int d^4x \mathcal{L}} \phi(x_1) \cdots \phi(x_n)$$
The way to calculate this is to add to the Lagrangian a $J\phi$ term, with $J$ some fixed function (we don't actually care what it is) called a classical source or just a source. If we define the generating functional
$$Z[J] = \int \mathcal{D}\phi\  e^{i\int d^4x (\mathcal{L} + J\phi)}$$
we can insert factors of $\phi$ in the integral by taking functional derivatives with respect to $J$ and then setting $J=0$ (this is why we don't care what $J$ is). For example, take the 1-point function:
$$\langle 0 | \phi(x_1) | 0 \rangle = \int \mathcal{D}\phi\ e^{i\int d^4x \mathcal{L}} \phi(x_1) = \frac{1}{Z[0]} \frac{1}{i} \frac{\delta}{\delta J(x_1)} \int \mathcal{D}\phi\ e^{i\int d^4x (\mathcal{L} + J\phi)}$$
If you know $Z[J]$, you can calculate any correlation function. This approach is used extensively in any QFT book that uses path integrals. See for example Srednicki, in particular sections 6 to 9.
A: Generally, source terms are what stands in the rhs of the differential equation $$\hat{\Theta}_x f(x) = j(x),$$ where $\hat{\Theta}_x$ is some linear differential operator (in your case $\hat{\Theta}_x = \Box_x + m^2$).
A generic solution of this equation is a sum of any solution of the homogeneous equation $\hat{\Theta}_x f(x) = 0$ (which usually is a superposition of plane waves) and any particular solution of the nonhomogeneous equation. Moreover, this particular solution is given by the Green's function:
$$ f(x) = f_0 (x) + \int d^4 y \, G(x, y) j(y) $$
where $G(x, y)$ is the Green's function of the differential operator $\hat{\Theta}_x$ (which means that it satisfies $\hat{\Theta}_x G(x, y) = \delta^{(4)}(x - y)$ as a distribution), and $f_0 (x)$ is any solution of the homogeneous equation (a superposition of plane waves).
