QFT generating functional and Green function and propagator I am confused about why does the generating functional gives the propagator by differentiation, and why that propagator is the Green function.
I understand how to take the functional derivative like this:

*

*we start with
$$
Z[J]=  N \int \mathcal D \psi  e^{i \int d^4 x L_0 + J(x)\psi(x) }\tag{1}
$$
and after some calculation we integrate out the functional determinant
$$
Z[J]=  N \det(A) \int \mathcal D \psi  e^{i \int d^4 x d^4 y J(x) G(x-y) J(y) }.\tag{2}
$$


*we take the functional derivative of the first expression, say to get two-point function
$$
\frac{\delta^2 Z}{\delta J(y_1) \delta J(y_2）} \sim N \int \mathcal D \psi  \psi(y_1) \psi(y_2) e^{i \int d^4 x d^4 y J(x) G(x-y) J(y) }\tag{3}
$$
and we take the functional derivative of the first expression, and we see the result is $<G>$.


*we assemble the whole thing at last by setting $J=0$, and somehow we have
$$\begin{align}
\left\langle 0|T[\psi(y_1)\psi(y_2)]|0\right\rangle \sim N \int \mathcal D \psi  \psi(y_1) \psi(y_2)  e^{i \int d^4 x L_0 } \sim  \frac{\delta^2 Z[J]}{\delta J(y_1) \delta J(y_2）}.
\end{align} \tag{4}$$


*Compare the result, we say $$G=<T\psi_1 \psi_2>.\tag{5}$$
What confuses me is mainly step 3. Why can we write
$$\begin{aligned}
 \int \mathcal{D} \psi(t) \psi\left(t_{1}\right) \psi\left(t_{2}\right) e^{i S[x(t)]} =\left\langle 0\left|T\left(\psi\left(t_{1}\right) \psi\left(t_{2}\right)\right)\right| 0\right\rangle~?
\end{aligned}\tag{6}$$
Seems to me that this equal sign is not obvious.
 A: You seem to be mixing definitions with computational techniques. Start by recalling statistical mechanics, where the partition function is by construction a sum over states weighted by $\exp(-\beta H)$ in the context of path integration in quantum mechanics the argument uses the evolution operator $\exp(- i H t)$ and sums over all paths, that is we have
$$\langle t_f, x_f | t_i, x_i\rangle = \int {\cal D}[x(t)] e^{-i S[x(t)]},$$
which can be obtained by breaking down to amplitude into many intermediate steps evolved for finite time intervals (you will find this in any book about path integrals). Then we immediately generalize to fields in the obvious way and consider the special case of vacuum-to-vacuum amplitudes, since we build all asymptotic particle states from vacuum states, thus
$$\langle 0 | 0 \rangle = \int {\cal D}[\phi] e^{-i S[\phi]} \equiv Z$$
Know if we understand the above in a probabilistic sense and call it partition function, you realize that the vacuum expectation value for a local generic operator ${\cal O}$ is given by
$$\langle 0 | {\cal O}(x_1,x_2,\dots) | 0 \rangle = \int {\cal D}[\phi] {\cal O}(x_1,x_2,\dots) \,e^{-i S[\phi]}.$$
All of the above is definitions and standard constructions, so now the question is simply how can you compute expectation values. You are free to choose any technique you want, but a popular one for the case of simple field insertions is to add a fictitious current $J$, as the OP writes,
$$Z[J]\equiv \int {\cal D}[\phi] e^{-i S[\phi] + i\int {\rm d}^4 x J(x)\phi(x)} ,$$
and then take derivatives and evaluate at $J=0$, so mathematically what you get is
$$\langle 0 | T\phi(x_1)\phi(x_2) | 0 \rangle = \frac{\delta^2 Z[J]}{\delta J(x_1) \delta J(x_2)}\Bigg|_{J=0} = \int {\cal D}[\phi] \phi(x_1)\phi(x_2) \,e^{-i S[\phi]}.$$
(For the origin of the time ordered operator I recommend you look into QFT books, e.g Schwartz's Quantum Field Theory and the Standard Model)
As you can see there is no need whatsoever at this point to deal with functional determinants.
To complete the story we must relate this to Green's functions. One important remark at this point is obvious, we cannot compute these objects exactly in the practice, so one must do some sort of approximation or expansion. But for the sake of the explanation, the Green's functions are solutions to the Green's equation, which in this context are the equations of motion, either for the free theory or the exact ones, or about some interesting background fields. So there are many kinds of Green's functions but for us the OP's second equation defines what the OP calls Green's functions, which seems to be the exact propagator of a possibly interacting theory. Now we can apply the above method of taking derivatives
\begin{align}
\langle 0 | T\phi(x_1)\phi(x_2) | 0 \rangle 
&= \frac{\delta^2 Z[J]}{\delta J(x_1) \delta J(x_2)}\Bigg|_{J=0} \\
&= \frac{\delta^2}{\delta J(x_1) \delta J(x_2)}\Bigg|_{J=0}{\cal N} \sqrt{\frac{1}{\det G^{-1}}} e^{\frac{1}{2}\int d^4 x\, d^4 x' J(x)G(x,x')J(x')}\\
&\propto G(x_1,x_2)
\end{align}
where the definition of $G$ was used in the second line (however corrected, the integration over paths disappears in Eq. (2)).
The last point is to connect to propagators, but their definition coincides with that of the two point time ordered correlation function unless we employ some perturbative treatment. If we think of expanding around the free theory then the propagator corresponding to the non-interacting two-point time ordered correlation function is the Feynman propagator which differs from the the exact propagator or the dressed-propagator, etc. I hope I brought some light into the matter.
A: I don't know if I understood your question, but it seems to be: "why the generating functional $Z[J]$ defined by the path integral in equation (1) has the property that it generates exactly the time-ordered correlation functions of the underlying QFT?". Here I sketch a way to justify that based on the Dyson-Schwinger (DS) functional differential equation.
The basic idea is this: the fundamental objects in QFT are the Green's functions. They are, for example, the way in which one is able to evaluate the ${\cal S}$-matrix by the LSZ prescription. Let us consider a scalar field for simplicity. In that case, the object of interest is $$G_n(x_1,\dots, x_n)=\langle \Omega |T\{\phi(x_1)\cdots \phi(x_n)\}|\Omega\rangle\tag{1}.$$
Now, suppose we know all the Green's functions (1). If that is the case, we may encode them in a generating functional:
$$Z[J]=\sum_{n=0}^\infty \dfrac{1}{n!}\int d^dx_1\cdots d^dx_n  j(x_1)\cdots j(x_n)\langle\Omega|T\{\phi(x_1)\cdots \phi(x_n)\}|\Omega\rangle\tag{2}.$$
The so-defined $Z[J]$, by definition, has the property that: $$\dfrac{\delta^nZ[J]}{\delta j(x_1)\cdots \delta j(x_n)}\bigg|_{j=0}=\langle \Omega|T\{\phi(x_1)\cdots \phi(x_n)\}|\Omega\rangle\tag{3},$$
and so, by (2) and (3) we see that knowledge of $Z[J]$ is equivalent to knowledge of the Green's functions. But how can we take advantage of that? Is there some way to obtain $Z[J]$ without knowing the Green's functions in the first place?
Well, it turns out that one may prove that $Z[J]$ must obey one functional differential equation known as the Dyson-Schwinger equation. A very clear exposition of how you may derive that $Z[J]$ obeys this equation is provided in "Path integral methods in quantum field theory" by R. J. Rivers in sections 1.1 - 1.3, and I encourage you to study the details. Here I just quote the result, given in equation (1.29) of the book, in order to address your question:
$$\left[\left(\dfrac{\delta S[\phi]}{\delta \phi(x)}\right)\bigg|_{-i\frac{\delta}{\delta J(x)}}+J(x)\right]Z[J]=0\tag{4}$$
where the first term in bracket means that we take the functional derivative of the action $S[\phi]$ and replace $\phi(x)$ by the operator $-i\frac{\delta}{\delta J(x)}$. Now we can formally solve (4) by a functional Fourier transform. In other words, let us define $\hat{Z}[\varphi]$ by $$Z[J]=\int \mathfrak{D}\varphi \hat{Z}[\varphi]e^{i\int d^dy J(y)\varphi(y)}\tag{5}.$$
Now we plug the ansatz in (4). We observe that $\frac{\delta S[\phi]}{\delta \phi(x)}$ with $\phi(x)$ replaced by $-i\frac{\delta}{\delta J(x)}$ is a polynomial in such functional derivatives. When they hit (5) each of them is effectively replaced by $\varphi(x)$. On the other hand, we can get $J(x)\hat{Z}[\varphi]$ by differentiating $\hat{Z}[\varphi]$ with respect to $\varphi(x)$. So we can write $$\left[\left(\dfrac{\delta S[\phi]}{\delta \phi(x)}\right)\bigg|_{-i\frac{\delta}{\delta J(x)}}+J(x)\right]Z[J]=\int \mathfrak{D}\varphi \bigg[\left(\dfrac{\delta S[\varphi]}{\delta \varphi(x)}\right)\hat{Z}[\varphi]-i\dfrac{\delta \hat{Z}[\varphi]}{\delta \varphi(x)}\bigg]e^{i\int d^dy J(y)\varphi(y)}.\tag{6}$$
This means we get the Fourier space equation $$\dfrac{\delta \hat{Z}[\varphi]}{\delta \varphi(x)}=i\dfrac{\delta S[\varphi]}{\delta\varphi(x)}\hat{Z}[\varphi]\tag{7}.$$
This equation immediately implies one exponential solution
$$\hat{Z}[\varphi]= N e^{iS[\varphi]}\tag{8},$$
which given (5) implies that we must have $$Z[J]=N\int \mathfrak{D}\varphi e^{iS[\varphi]+i\int d^dx J(x)\varphi(x)}\tag{9}.$$
The logic here is the following: the path integral representation of $Z[J]$ (9) gives you the Green's functions of the theory because the generating functional $Z[J]$ of such Green's functions is constrained to solve the DS equation (4) which turns out to have (9) as a solution.
One final caveat is the following. One may evaluate (9) exactly in the free theory and use perturbation theory afterwards. Still, in the free theory, solving (9) demands finding the inverse of the wave operator appearing in the equations of motion of the free theory. For example $\Box+m^2$ for scalar fields. Now, in Lorentzian signature $(\Box+m^2)$ admits more than one inverse, so a prescription is required. The key then is to recall that we want $Z[J]$ to generate time-ordered correlation functions. This is an extra, implicit, boundary condition, that one must take into account when solving (9). It translates eventually into the $i\epsilon$ prescription to the propagator.
