Here is a good asymptotics for this:
Let $f(x)=\frac{e^{-x^2/2}}{\sqrt{2\pi}}$ be the density and $\Phi(x)=\int_{-\infty}^x f(u)du$ be the CDF of the standard normal random variable $X$. Let $Z=Z(n)=\max(X_1,\dots,X_n)$ where $X_i$ are i.i.d.\ normal ${\cal N}(0,1)$. We are interested in $\mathbb{E}Z(n)$.
Since the CDF of $Z$ is given by $F_Z(z)=\mathbb{P}(Z\le z)=\mathbb{P}(X_1\le z,\dots,X_n\le z)=\mathbb{P}(X\le z)^n=\Phi(z)^n$, hence its density is $\frac{d}{dz}\left[\Phi(z)^n\right]=n f(z)\Phi(z)^{n-1}$ and
\begin{align*}
\mathbb{E}Z(n)&=\int_{-\infty}^\infty z\cdot n f(z)\Phi(z)^{n-1} dz
=\int_{-\infty}^0 z\cdot \frac{d}{dz}\left[\Phi(z)^n\right] dz
-\int_0^{\infty} z\cdot \frac{d}{dz}\left[1-\Phi(z)^n\right] dz
\\
&=-\int_{-\infty}^0 \Phi(z)^n dz
+\int_0^{\infty} 1-\Phi(z)^n dz
\end{align*}
using the integration parts and noting that $\lim_{z\to-\infty}z\cdot \Phi(z)^n=0$ and $ \lim_{z\to\infty}z\cdot \left[1-\Phi(z)^n\right]=0$ by e.g.\ L'Hôpital's rule. Next,
\begin{align*}
0\le \int_{-\infty}^0 \Phi(z)^n dz
\le
\int_{-\infty}^{-1} e^{zn} dz
+\int_{-1}^{0} \Phi(0)^n dz
= \frac{1}{n}e^{-n} + 2^{-n}\to 0
\end{align*}
since $f(z)<e^z$ (and hence $\Phi(z)<e^{z}$ too) for $x\le -1$, and $\Phi(z)\le \Phi(0)=1/2$ for $z\le 0$.
Now it only remains to estimate $\int_0^{\infty} 1-\Phi(z)^n dz$.
Observe that
$$
\frac{1-\Phi(z)}{f(z)/z}\to 1\text{ as }z\to+\infty.
$$
since for $z> 0$
\begin{align*}
1-\Phi(z)&=\int_z^\infty \frac{e^{-x^2/2}}{\sqrt{2\pi}}dx \le\int_z^\infty \frac{x}{z} \cdot \frac{ e^{-x^2/2}}{\sqrt{2\pi}}dx =\frac 1z\cdot \frac{ e^{-z^2/2}}{\sqrt{2\pi}};
\\
1-\Phi(z) &\ge \int_z^{z+1} \frac{e^{-x^2/2}}{\sqrt{2\pi}}dx \ge \int_z^{z+1} \frac{x}{z+1} \cdot \frac{ e^{-x^2/2}}{\sqrt{2\pi}}dx =
\frac {1-e^{-1/2-z}}{1+1/z}\cdot \frac{ e^{-z^2/2}}{z\sqrt{2\pi}},
\end{align*}
so we can write
$$
1-\Phi(z)=c_z\frac{ e^{-z^2/2}}{z\sqrt{2\pi}}\quad\text{for }z>0
$$
where $c_z\to 1$ as $z\to +\infty$.
Let $y=y(n)=\sqrt{2\ln n}$.
For $z\ge y$ we have
%$1-\Phi(z)\ll \frac 1n$,
\begin{align*}
\int_{y}^{\infty} 1-\Phi(z)^n dz&=
\int_{y}^{\infty} 1-\left[1- c_z\frac{ e^{-z^2/2}}{z\sqrt{2\pi}}\right]^n dz
\le
\int_{y}^{\infty} n c_z\frac{ e^{-z^2/2}}{z\sqrt{2\pi}} dz
\\ &\le
\int_{y}^{\infty} n c_z \frac{z}{y^2} \frac{e^{-z^2/2}}{\sqrt{2\pi}} dz
\le
\frac{n\max_{z\ge y} c_z}{y^2\sqrt{2\pi}} \int_{y}^{\infty} z e^{-z^2/2} dz
=O\left(\frac1{\ln n}\right).
\end{align*}
Let us estimate now $\int_{0}^{y} \Phi(z)^n dz$.
Fix $\varepsilon>0$ small, and note that if $1\le z\le y-\varepsilon$ then
$$
1-\Phi(z)\ge \left(\min_{1\le z\le y} c_z\right) \frac{n^{-1}\, e^{\varepsilon y-\frac{\varepsilon^2}2}}{(y-\varepsilon)\sqrt{2\pi}}
=
\frac 1n\, \frac{c_1}{1-o(1)} \frac{e^{\varepsilon y}}{y}\gg \frac 1n
$$
for some $c_1(\varepsilon)>0$ and where $o(1)\to 0$ as $n\to\infty$; note also that $e^{\varepsilon y}/y\to\infty $ as $y\to\infty$.
This yields $\Phi(z)^n\le \exp\left(-(c_1+o(1))e^{\varepsilon y}/y\right)$. Hence
\begin{align*}
\int_0^{y-\varepsilon}\Phi(z)^n dz&\le
\int_0^1 \Phi(1)^n dz+
\int_1^{y-\varepsilon} \Phi(z)^n dz
\le 0.85^n+\int_1^{y-\varepsilon}
\exp\left(-(c_1+o(1))e^{\varepsilon y}/y\right) dz
\\ &
\le 0.85^n+ y\exp\left(-(c_1+o(1))e^{\varepsilon y}/y\right) \to 0\qquad\text{ as }n,y\to\infty.
\end{align*}
Consequently, since $\int_{y-\varepsilon}^y \Phi(z)^n dz\le \varepsilon$, we conclude
$$
\limsup_{n\to\infty} |\mathbb{E}Z(n)-y(n)|\le \varepsilon
$$
Since $\varepsilon>0$ is arbitrary, this means
$$
\mathbb{E}Z(n)=\sqrt{2\ln n}+o(1).
$$
