I would like to provide some details on the answers by ofer zeitouni and Noam D. Elkies.
Selberg'sBy Selberg's formula seems a bit different:
$$\int_0^1 \cdots \int_0^1 V(a)^{2\beta} \prod_{i=0}^{n-1} da_i=
\prod_{j=0}^{n-1} \frac{(\Gamma(1+j\beta))^2 \cdot \Gamma(1+(j+1)\beta)}
{\Gamma(2+(n+j-1)\beta)\cdot \Gamma(1+\beta)}=:A(n,\beta);$$
see e.g. Theorem 3.2 in http://www.sciencedirect.com/science/article/pii/S0196885801907985 . However, this does not affect the easy calculation of the limitas presented in the following formula:answer by ofer zeitouni, one easily finds
$$(1)\qquad M_n=\lim_{\beta\to\infty} A(n,\beta)^{1/(2\beta)}=
\prod _{j=0}^{n-1} \frac{j^j (j+1)^{(j+1)/2}}{(j+n-1)^{(j+n-1)/2}}
$$
(with $0^0:=1$), which follows immediately from the observation that for any real $a$ and any real $c>0$
$$\Gamma(a+c\beta)^{1/(2\beta)}\sim(c\beta/e)^{c/2}
$$
as $\beta\to\infty$.
I asked for an upper bound on $M_n$, which would be asymptotic to $M_n$, thinking that an explicit expression for $M_n$ would not be possible. However, as is now clear from the answers by ofer zeitouni and Noam D. Elkies, such an expression is not so hard to obtain, and formula (1) presents this expression.
On the other hand, the relation between the logarithmic energy
$E_{n,a}:=2\sum_{0\le i<j\le n-1}\ln(a_j-a_i)$ and its ostensibly more general version $n^2E_\mu$ with $E_\mu:=\int\ln|x-y|\mu(dx) \mu(dy)$
for probability measures $\mu$ supported on $[0,1]$ seems unclear, I guess because of the singularities on the diagonal. Namely, one would expect that
$$(2)\qquad 2\ln M_n=\max_{0=a_0<\dots<a_{n-1}=1}E_{n,a}
\le n^2 \max_\mu E_\mu.
$$
However, it is not clear if $2E_{n,a}=2\sum_{0\le i<j\le n-1}\ln(a_j-a_i)$ could be written as
$n^2E_\mu=n^2\int\ln|x-y|\mu(dx) \mu(dy)$
for some probability measure $\mu$ on $[0,1]$.
In fact, quite surprisingly to me, the inequality in (2) is false, at least for large enough $n$. Indeed, it is not hard to show based on formula (1) that for some real $c\in(0,\infty)$ $$(3)\qquad M_n=(c+o(1))m_n,\quad\text{where}\quad m_n:=2^{-n^2} \sqrt{(n-1)!}\,(8e)^{n/2} n^{3/8} $$ $$(4)\qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad >>2^{-n^2}=\text{(?)}=\exp\Big\{\frac{n^2}2\,\max_\mu E_\mu\Big\}. $$ The asymptotics $M_n=(c+o(1))m_n$ follows because $$\ln\frac{M_{n+1}}{M_n} =-n \ln2+\frac{n-1}{2}\,\ln (n-1)+\frac{n+1}{2}\,\ln (n+1)-\frac{2n-1}{2}\,\ln(2 n-1)$$ $$=d_n+O(1/n^2)$$ for large $n$, where $d_n:=\frac{3}{8 n}-2n\ln2+\frac{1}{2}\ln n+\frac12(1+\ln2)$.
A curious corollary to (3)--(4) is that for each large enough natural $n$ there is some $a=(a_0,\dots,a_{n-1})$ with $0=a_0<\dots<a_{n-1}=1$ such that $2E_{n,a}$ cannot be approximated by (let alone written as) $n^2E_{\mu_k}$ for any sequence $(\mu_k)$ of probability measures on $[0,1]$.
Yet, it also follows from (3)--(4) that the logarithmic asymptotics $\ln M_n\sim\frac{n^2}2 \max_\mu E_\mu$ holds.
One can also see that $\ln(M_{n+1}/M_n)<d_n$ for $n\ge4$, and hence for each $k\ge4$ and all $n>k$ one has the upper bound $M_k\exp\sum_{j=k}^{n-1}d_j$ on $M_n$, which is asymptotic to $M_n$ as $n>k\to\infty$. In particular, for $n\ge4$ one has $M_n<\tilde c_4 m_n$, where $m_n$ is as before and $\tilde c_4:=\frac{512}{25} \sqrt{\frac{2}{15}} e^{(6 \gamma -43)/16}=0.631\dots$ and $\gamma$ is the Euler constant.