Deriving a Useful Solution of the free Schrödinger equation 
How does one derive the fact that
  $$\psi(t,x) = (\tfrac{2 \pi \hbar t}{m})^{-d/2}\int_{\mathbb{R}^d} e^{im\tfrac{(x-y)^2}{2\hbar t}}\psi_0(y)dy$$
  is a solution of the time-dependent free Schrödinger equation
  $$i\hbar \tfrac{\partial }{\partial t}\psi(t,x) = -\tfrac{\hbar^2}{2m} \Delta \psi(t,x), \\ \psi(0,x) = \psi_0(x)~?$$

It may follow from the operator solution of Schrödinger
$$\psi(t,x) = e^{\tfrac{i}{\hbar}H_ot}\psi_0(x)$$
i.e. start with
$$ e^{\tfrac{i}{\hbar}H_ot}\psi_0(x) $$
and insert a delta function or something and end up with
$$\psi(t,x) = e^{\tfrac{i}{\hbar}H_ot}\psi_0(x) = (\tfrac{2 \pi \hbar t}{m})^{-d/2}\int_{\mathbb{R}^d} e^{im\tfrac{(x-y)^2}{2\hbar t}}\psi_0(y)dy$$
as page 3 of S. Mazzucchi, Mathematical Feynman Path Integrals and Their Applications, seems to imply (readable on amazon's search feature). 
 A: This question is about Guassian wave-packet propagation and the corresponding Green's function in ordinary quantum mechanics.
Assuming $\hbar=m=1$ for simplicity, consider the solution (with its initial condition) to the following Schrodinger equation:
$$i\partial_t G=-\partial^2_x G \\
G(t=0,x)=\delta(x)
$$
Now assume the Fourier ansatz for $G$:
$$G(t,x)=\int dk A(t,k)e^{ikx}$$
Putting it back to the equation above determines the form of $A$:
$$i\partial_t A=k^2A$$
therefore
$$A(t,k)=C(k)e^{-ik^2t}$$
which gives
$$
G(t,x)=\int dk C(k)e^{-ik^2t+ikx}
$$
Imposing the initial condition now yields:
$$
G(t=0,x)=\int dk C(k)e^{ikx}=\delta(x) \implies C(k)=1/2\pi \\
G(t,x)=\int \frac{dk}{2\pi} e^{-ik^2t+ikx}
$$
The last integral is the usual Gaussian integral, which can be evaluated by making a complete square out of the integrand:
$$
-ik^2t+ikx = -it\left[k^2-\frac{kx}{t}\right] = -it\left[\left(k-\frac{x}{2t}\right)^2-\frac{x^2}{4t^2}\right]
$$
therefore the integral will evaluate to 
$$
G(t,x) \propto \frac{1}{\sqrt{t}}\exp\left(\frac{ix^2}{4t}\right)
$$
Now if you imagine $\psi_0$ is made up of weighted $\delta$ contributions, i.e.:
$$
\psi_0(x)=\int dy \psi_0(y) \delta(y-x)
$$
then you arrive at your original $\psi$ function in your question. In plain language, the $\psi$ function integrates the result of propagating the delta functions which make up the initial condition $\psi_0$. Mathematically:
$$
\psi(t,x)=\int dy G(t,x-y)\psi_0(y)
$$
