The random variable $Z = \max_{i=1}^n(X_i)$ is known as order statistics, and is sometimes denoted as $X_{n:n}$.
The cumulative density function of $Z$ is easy to find:
$$
F_Z(z) = \mathbb{P}\left(Z \leqslant z\right) = \mathbb{P}\left(\max_{i=1}^n(X_i) \leqslant z\right) = \mathbb{P}\left( X_1 \leqslant z, X_2 \leqslant z, \ldots, X_n \leqslant z\right)
$$
using independence:
$$
F_Z(z) = \left(F_X(z)\right)^n
$$
Thus the density function is
$$
f_Z(z) = n f_X(z) F_X^{n-1}(z)
$$
In particular, it follows that $Z$ is not normal.
Expected values of $Z$ are known in closed form for $n=1,2,3,4,5$ (asking Mathematica):
In[31]:= Table[
Mean[OrderDistribution[{NormalDistribution[m, s], n}, n]], {n, 1, 5}]
Out[31]= {m, m + s/Sqrt[Pi], m + (3 s)/(2 Sqrt[Pi]),
m + (6 s ArcTan[Sqrt[2]])/Pi^(3/2),
m - (5 s)/(2 Sqrt[Pi]) + (15 s ArcTan[Sqrt[2]])/Pi^(3/2)}
Large $n$ asymptotics is discussed here.
The distribution of the mean and variance of a normal rv is very well known:
$$\sqrt n \left( \begin{array}{c}
\overline X - \mu \\
{S^2} - {\sigma ^2}
\end{array} \right) \sim \ \left(\begin{array}{c}
\mathcal{N}(0,1) \\
\sigma^2\left(\frac{\sqrt{n}\chi^2_{n-1}}{n-1}-1\right)
\end{array} \right)$$
Best Answer
This is the Marsaglia polar method, a variant of the Box–Muller transform.
The very largest (and smallest) values of the result occur when the pair you're calling $(u,v)$ is as close to $0$ as possible. It looks like there is an option for a $32$-bit and $64$-bit version. With $b$ bits available, the smallest positive nonzero value is $1/2^{b-1}$ (this is then cast to a floating point value, which shouldn't cause problems). When $u = 1/2^{b-1}$ and $v = 0$, we get $s = 1/2^{2b-2}$; the final generated value (assuming a standard Gaussian) is $2\sqrt{(b-1)\ln 2}$. The most negative value should similarly be $-2\sqrt{(b-1)\ln 2}$.
This is about $9.27094$ for the $32$-bit version and about $13.2164$ for the $64$-bit version.