The Maximum of $X_1,\dots,X_n. \sim$ i.i.d. Standardnormals converges to the Standard Gumbel Distribution according to Extreme Value Theory.
How can we show that?
We have
$$P(\max X_i \leq x) = P(X_1 \leq x, \dots, X_n \leq x) = P(X_1 \leq x) \cdots P(X_n \leq x) = F(x)^n $$
We need to find/choose $a_n>0,b_n\in\mathbb{R}$ sequences of constants such that: $$\Phi\left(a_n x+b_n\right)^n\rightarrow^{n\rightarrow\infty} G(x) = e^{-\exp(-x)}$$
Can you solve it or find it in literature?
There are some examples pg.6/71, but not for the Normal case.
$$\Phi\left(a_n x+b_n\right)^n=\left(\frac{1}{\sqrt{2\pi}}\int_{-\infty}^{a_n x+b_n} e^{-\frac{y^2}{2}}dy\right)^n\rightarrow e^{-\exp(-x)}$$
Use the asymptotic expansion
$$\Phi(x) = 1-\frac{e^{-x^2/2}}{\sqrt{2\pi}}\left(\frac{1}{x}+ \ldots\right)$$
If you can select $a_n$ and $b_n$ such that as $n \rightarrow \infty,$
$$\frac{e^{-(a_n x + b_n)^2/2}}{(a_n x +b_n)\sqrt{2\pi}}\sim \frac{e^{-x}}{n}, $$
then, as desired,
$$\lim_{n \rightarrow \infty}\Phi(a_nx+b_n)^n =\lim_{n \rightarrow \infty}\left(1-\frac{e^{-x}}{n}\right)^n= \exp(-e^{-x}).$$
First select $a_n = 1/b_n$, leaving $b_n$ unspecified for now , but with the requirement $b_n \rightarrow \infty.$
Then
$$\frac{e^{-(a_nx+b_n)^2/2}}{(a_nx+b_n)\sqrt{2\pi}} = \frac{e^{-x^2/2b_n^2}e^{-x}}{x/b_n^2+1}\frac{e^{-b_n^2/2}}{b_n\sqrt{2\pi}}. $$
This a good start as
$$\lim_{n \rightarrow \infty} \frac{e^{-x^2/2b_n^2}e^{-x}}{x/b_n^2+1}=e^{-x}.$$
Now using the inverse normal distribution, choose $b_n = \Phi^{-1}(1 - 1/n)$, which satisfies the requirement, $b_n \rightarrow \infty$, and
$$\frac{e^{-b_n^2/2}}{b_n\sqrt{2\pi}} \sim 1-\Phi[\Phi^{-1}(1-1/n)]= \frac1{n}.$$
Putting everything together, we have as $n \rightarrow \infty$,
$$\frac{e^{-(a_n x + b_n)^2/2}}{(a_n x +b_n)\sqrt{2\pi}}\sim \frac{e^{-x}}{n}. $$