I've read several proofs of Central Limit Theorem and they all seemed inaccurate to me, because they drop last members of Taylor series, whereas those members are not infinitesimal.
The classical proof of Central Limit Theorem via characteristic functions is based on the fact that $e^{-\frac{x^2}{2}}$ is eigenfunction of Fourier transform.
If $\xi$ is a random variable with $E\xi=0$ and $Var(\xi)=1$, we substitute it with another random variable $\nu=\xi/\sqrt{n}$. We look into Fourier spectrum of its probability density function:
$\varphi_{\nu}(t)=\int\limits _{x=-\infty}^{\infty}e^{itx}f_{\nu}(x)dx=\int\limits _{x=-\infty}^{\infty}(e^{i\cdot0\cdot x}+itxe^{i\cdot0\cdot x}+\frac{i^{2}t^{2}x^{2}}{2!}e^{i\cdot0\cdot x}+\frac{i^{3}t^{3}x^{3}}{3!}e^{i\cdot0\cdot x}+...)f_{\nu}(x)dx=$
$=\int\limits _{x=-\infty}^{\infty}(e^{i\cdot0\cdot x}+itxe^{i\cdot0\cdot x}+\frac{i^{2}t^{2}x^{2}}{2!}e^{i\cdot0\cdot x}+\frac{i^{3}t^{3}x^{3}}{3!}e^{i\cdot0\cdot x}+...)f_{\xi}(\sqrt{n}x)\sqrt{n}dx=\int\limits _{y=-\infty}^{\infty}(e^{i\cdot0\cdot \frac{y}{\sqrt{n}}}+it\frac{y}{\sqrt{n}}e^{i\cdot0\cdot \frac{y}{\sqrt{n}}}+\frac{i^{2}t^{2}y^2}{2!\cdot n}e^{i\cdot0\cdot \frac{y}{\sqrt{n}}}+\frac{i^{3}t^{3}y^{3}}{3! \cdot n\sqrt{n}}e^{i\cdot0\cdot \frac{y}{\sqrt{n}}}+...)f_{\xi}(y)dy$
And then into Fourier spectrum $\hat{f}_{{S}_{\nu,n}}(t)$ of p.d.s. of sum $f_{S_{\nu,n}}(x)$ of $n$ instances of $\nu$, which is n-th power of $\varphi_{\nu}(t)$. Each $t$-th harmonic of $\hat{f}_{{S}_{\nu,n}}(t)$ in the interval $t=[-C,C]$ has amplitude, close to $e^{-t^2/2}$ (with allowed error $\epsilon$) and that interval $C(\epsilon, n)$ grows infinitely with growing $n$:
To complete the proof, we need to show that starting from some $n$ amplitudes of harmonics outside that interval $t=[-C,C]$ can be considered negligible and dropped without significantly affecting the calculation of Fourier synthesis $f_{S_{\nu,n}}(x)$. In other words, this picture should not take place:
I don't know hot to do that.
I tried Parseval's identity for p.d.f. of $\nu$ and its characteristic function: $\int\limits_{x=-\infty}^{\infty}f^2_\nu(x)dx = \int\limits_{t=-\infty}^{\infty}\hat{f}^2_\nu(t)dt$.
As $\sqrt{n}$ grows larger, $\nu$ becomes more of a Dirac delta-function, thus its square integral gets infinitely large:
At the same time, central part of the spectrum of $\nu$ is $e^{-t^2/2n}$, and its square integral on $[-\infty, \infty]$ converges to $\int\limits_{-\infty}^{\infty}e^{-\frac{t^2}{n}}dt=\sqrt{\frac{n}{2}}$. Ok, this goes to infinity as n goes to infinity and allows remaining harmonics to have infinitesimal amplitudes, but doesn't guarantee that. Dead end.
EDIT: As Yemon Choi suggested, Levy theorem proves that higher harmonics have infinitesimal amplitudes for large enough $t$ as $n$ grows. That's it.
Convolution of the Scaled Probability Measure
Suppose that $f$ is a probability density with mean $0$ and variance $1$.
Consider the Fourier Transform of $f$ $$ \widehat{f}(\xi)=\int_{\mathbb{R}}f(x)\,e^{-2\pi ix\xi}\,\mathrm{d}x\tag{1} $$ The assumptions we've made allow us to say that $$ \widehat{f}(0)=1\qquad\widehat{f}\vphantom{f}'(0)=0\qquad\widehat{f}\vphantom{f}''(0)=-4\pi^2\tag{2} $$ Equations $(2)$, say that near $\xi=0$ $$ \widehat{f}(\xi)=1-4\pi^2\xi^2+o(\xi^2)\tag{3} $$ Define $$ f_n(x)=\overbrace{f\ast f\ast f\ast\cdots\ast f}^{\text{$n$ copies of $f$}}(x\sqrt{n})\sqrt{n}\tag{4} $$ Then $$ \begin{align} \widehat{f_n}(\xi) &=\widehat{f}(\xi/\sqrt{n})^n\\[3pt] &=\left(1-2\pi^2\xi^2/n+o(\xi^2/n)\right)^n\\[3pt] &=e^{-2\pi^2\xi^2}e^{\psi_n(\xi)}\tag{5} \end{align} $$ where $\psi_n(\xi)=no(\xi^2/n)$. That is, on any compact set $\psi_n(\xi)\to0$ uniformly.
Weak Convergence of the Convolution
Let $\phi(x)$ be a function so that $\widehat{\phi}\in L^1$ (e.g. Schwartz class). Choose an $\epsilon\gt0$. There is an $M$ so that $$ \int_{|\xi|\gt M}\left|\widehat{\phi}(\xi)\right|\,\mathrm{d}\xi\le\frac\epsilon4\tag{6} $$ Furthermore, since $\psi_n(\xi)\to0$ uniformly on $|\xi|\le M$, there is an $N$ so that for $n\ge N$ $$ \begin{align} \int_{|\xi|\le M}\left|\widehat{f_n}(\xi)-e^{-2\pi^2\xi^2}\right|\left|\widehat{\phi}(\xi)\right|\,\mathrm{d}\xi &=\int_{|\xi|\le M}e^{-2\pi^2\xi^2}\left|1-e^{\psi_n(\xi)}\right|\left|\widehat{\phi}(\xi)\right|\,\mathrm{d}\xi\\ &\le\frac\epsilon2\tag{7} \end{align} $$ Therefore, when testing against a function $\phi$ so that $\widehat{\phi}\in L^1$, $$ \begin{align} \left|\,\int_{\mathbb{R}}\left[f_n(x)-\frac1{\sqrt{2\pi}}e^{-x^2/2}\right]\overline{\phi(x)}\,\mathrm{d}x\,\right| &=\left|\,\int_{\mathbb{R}}\left[\widehat{f_n}(\xi)-e^{-2\pi^2\xi^2}\right]\overline{\widehat{\phi}(\xi)}\,\mathrm{d}\xi\,\right|\\ &\le\int_{|\xi|\le M}e^{-2\pi^2\xi^2}\left|1-e^{\psi_n(\xi)}\right|\left|\widehat{\phi}(\xi)\right|\,\mathrm{d}\xi\\ &+\int_{|\xi|\gt M}2\left|\widehat{\phi}(\xi)\right|\,\mathrm{d}\xi\\[6pt] &\le\epsilon\tag{8} \end{align} $$ Since $\epsilon\gt0$ was arbitrary, $(8)$ says that, at least weakly, $$ \lim_{n\to\infty}f_n(x)=\frac1{\sqrt{2\pi}}e^{-x^2/2}\tag{9} $$