Why is the contribution of a path in Feynmans path integral formalism $\sim e^{(i/\hbar)S[x(t)]}$?

There are already several good answers. Here I will only answer the very last question, i.e., if the Boltzmann factor in the path integral is $f(S(t_f,t_i))$, with action $$S(t_f,t_i)~=~\int_{t_i}^{t_f} dt \ L(t)\tag{1},$$ why is the function $f:\mathbb{R}\to\mathbb{C}$ an exponential function, and not something else?

Well, since the Feynman "sum over histories" propagator should have the group property

$$\begin{align} K(x_3,t_3;x_1,t_1) ~=~&\cr \int_{-\infty}^{\infty}\mathrm{d}x_2 \ K(x_3,t_3;x_2,t_2)& K(x_2,t_2;x_1,t_1),\end{align}\tag{2}$$

one must demand that

$$\begin{align}f(S(t_3,t_2)&f(S(t_2,t_1)) \cr~=~& f(S(t_3,t_1)) \cr~=~& f(S(t_3,t_2)+S(t_2,t_1)).\end{align}\tag{3}$$ In the last equality of eq. (3) we used the additivity of the action (1). Eq. (3) implies that

$$f(0)~=~f(S(t_1,t_1)) ~=~ 1.\tag{4}$$ (The other possibility $f\equiv 0$ is physically un-acceptable.)

So the question boils down to: How many continuous functions $f:\mathbb{R}\to\mathbb{C}$ satisfy $$f(s)f(s^{\prime}) ~=~f(s+s^{\prime})\quad\text{and}\quad f(0) ~=~1~?\tag{5}$$

Answer: The exponential function!

Proof (ignoring some mathematical technicalities): If $s$ is infinitesimally small, then one may Taylor expand

$$\begin{align}f(s) ~=~& f(0) + f^{\prime}(0)s +{\cal O}(s^{2}) \cr ~=~& 1+cs+{\cal O}(s^{2}) \end{align}\tag{6}$$

with some constant $c:=f^{\prime}(0)$. Then one calculates

$$\begin{align} f(s) ~=~&\lim_{n\to\infty}f(\frac{s}{n})^n\cr ~=~&\lim_{n\to\infty}\left(1+\frac{cs}{n}+o(\frac{1}{n})\right)^n\cr ~=~&e^{cs},\end{align} \tag{7}$$

i.e., the exponential function! $\Box$

You start by writing down the probability to find a particle at $y$ at time $t$ when it was at $x$ at time $0$, denoted as $K(y,t;x,0)$. You get this by solving the Schrödinger equation with the initial condition $\psi(y,0) = \delta(y-x)$. Then, $K(y,t;x,0) = \psi(y,t)$. Thus, to solve this, we need to know the time development of the initial condition $\psi(y,0)$.

Let us start with the simple example of a free particle. This is easiest solved in momentum-representation, obtained by Fourier-transforming $\psi(y,t)$:

$$\psi(y,t) = \frac{1}{\sqrt{2\pi\hbar}} \int dp \exp(ipy/\hbar) \tilde \psi(p,t)$$ For $\tilde \psi$, the Schrödinger equation gives $$\tilde \psi(p,t) = \frac{1}{\sqrt{2\pi\hbar}} \exp \left(-\frac{i}{\hbar} \left[\frac{p^2 t}{2m} - px\right]\right)$$ This can be inserted back into the equation for $\psi(y,t)$. The integral over $p$ can be solved exactly. The final result is $$K_\text{free}(y,t;x,0) = \sqrt{\frac{m}{2\pi i\hbar t}} \exp\left(\frac{im(x-y)^2}{2\hbar t}\right)$$

Next step: The solution of the Schrödinger equation can generally be written as $$|\psi, t\rangle = \exp\left(-\frac{iHt}{\hbar}\right) |\psi,0\rangle$$ with $H$ being the Hamiltonian of your system. Writing $H = T+V$, the general formula for $K$ becomes $$K(y,t;x,0) = \langle y \mid \exp(-\frac{i(T+V)t}{\hbar}) \mid x \rangle$$ We use the Trotter-Kato Formula (which holds under certain conditions which I won't go into detail at this point. It allows us to write $$K(y,t;x,0) = \lim_{N\rightarrow \infty} \langle y \mid \left[ \exp(-\frac{iTt}{N\hbar}) \exp(-\frac{iVt}{N\hbar})\right]^N \mid x\rangle$$ We insert the unity operator, decomposed as $1 = \int dx | x \rangle \langle x |$ $N-1$ times, which gives us $$K(y,t;x,0) = \int dx_1 dx_2 \dots dx_{N-1} \prod_{j=0}{N-1} \langle x_{j+1} \mid \exp(-iTt/N\hbar) \exp(-iVt/N\hbar) \mid x_j \rangle$$ Note that $V$ as an operator acting on $|x\rangle$ gives just $V(x) |x\rangle$. And $\langle x_{j+1} | \exp(-iTt/N\hbar) | x_j \rangle$ gives us just the contribution of a free particle, i.e. $$\sqrt{\frac{mN}{2\pi i\hbar t}} \exp\left(\frac{imN}{2\hbar t}(x_{j+1} - x_j)\right)^2$$. If we abbreviate $\tau = t/N$, we can write: $$K(y,t;x,0) = \lim_{N\rightarrow \infty} \int dx_1 dx_2 \dots dx_{N-1} \left( \frac{m}{2\pi i\hbar \tau}\right)^{N/2} \times$$ $$\exp \left(\frac{i\tau}{\hbar} \sum_{j=0}^{N-1} \left[ \frac{m}{2}\left(\frac{x_{j+1}-x_j}{\tau}\right)^2 - V(x_j)\right]\right)$$

The next step is to see the values $x_j$ as points of a certain path $x(t')$ evaluated at points $t' = t_j = j\tau = jt/N$. If $\tau$ is small, we write $$\sum_{j=0}^{N-1} \tau f(t_j) \rightarrow \int f(t') dt'$$ $$\frac{x_{j+1} - x_j}{\tau} \rightarrow \dot x(t')$$ where the dot denotes the time-derivative.

The argument of the exponential then becomes $$\frac{i}{\hbar} \int_0^t dt' \left( \frac{m\dot x(t')^2}{2} - V(x(t'))\right)$$ You will have no trouble identifying the integrand as the Lagrangian $L = T-V$. The integral itself, therefore, is the classical action.

Thus, the formula we have for $K$ can be interpreted as the sum over all possible paths from $(x,0)$ to $(y,t)$ of the function $\exp\left(\frac{i}{\hbar} S(t,0)\right)$ of the classical action.

The interpretation of this was given in other answers: The classical path is that which minimizes the action, i.e. the action is stationary for the classical path. In your path-integral formula, this path will have a large contribution, as all paths that vary only slightly from the classical path will still have pretty much the same phase factor as the classical one, leading to constructive interference of those paths. For paths far from the classical path, the action will vary greater among the paths, so that there all possible phases occur, which will ultimately cancel each outer out.

Reference A lecture on advanced quantum mechanics given by Prof. Crispin Gardiner. Lecture notes are, unfortunately, not freely available. It was a good lecture :)

If you accept that Quantum Mechanics is built upon the fact that you sum complex amplitudes of processes (see this previous Question/Answers about this fact) you would expect that a sum over multiple paths behaves like a sum of different complex phases: $$M \sim \sum e^{i*\text{phase}}$$

Applying the variational principle to the phase, you see that the paths which vary their phase the least will contribute the most to the sum (because the others will average each other). Add the fact that you want the classical path to be the main contribution (because we want to match classical physics, this is the correspondence princple), and that the classical path is the path where the action $S$ varies the least, you can identify the phase with the action and get $\text{phase} \sim S[x(t)]$. Then you get $1 / \hbar$ as an experimental constant.

I'm not sure if this is a satisfactory answer, but most of the "strangeness" here comes from the QM superposition principle in the first place anyway. Note that the variational principle in classical mechanics was known and used before QM was invented and had the teleological property of "sniffing out" paths of least action. In the QM path-integral method this is at least explained from a more local point of view.