Distribution of Variance of draws from Multivariate Normal Distribution Let the vector $\boldsymbol x$ be a draw of $n$ values from a multivariate normal distribution with zero mean.
$$
\boldsymbol x \sim \mathcal{N}(\boldsymbol 0, \Sigma)
$$
It may be assumed that $\Sigma$ is standardised.
Let the scalar $y$ be the variance of $\boldsymbol x$:
$$
y = var(\boldsymbol x)
$$
Then, what is the distribution of $y$, given $n$ and $\Sigma$?
$$
y \sim ?(n, \Sigma)
$$

Context: an autoregressive process (AR1) generates time series where values close in time are more correlated than values further apart. Generating $n$ such values is equivalent to sampling from a multivariate normal distribution with a particular covariance structure, namely $\Sigma_{ij} = \rho^{|t_i - t_j|}$. What is the variance $y$ of the values of these time series, and how is that distributed when multiple time series are generated?

You mentioned simulation in one of your comments.  This is easily done and the R code below does so for the particular case where the vector has $10$ elements, $\mu=0$ and $\Sigma$ is the identity matrix, but this can easily be changed for other $\mu$ and $\Sigma$.  Note that the variance formula used is the $\frac1{n-1}$ version.
As whuber commented, in this particular case you get a chi-square distribution with $10-1=9$ degrees of freedom, with values multiplied by $\frac19$, so this distribution of variances has mean $1$ and variance $\frac29$.
library(MASS)
library(matrixStats)
set.seed(2022)
mu <- rep(0,10)
Sigma <- diag(10)  
sims <- rowVars(mvrnorm(10^5, mu, Sigma)) 
mean(sims)
# 0.9987011
var(sims)  
# 0.2225065

To illustrate that this particular case does have that distribution, here are the cdf and density from the simulation (black) compared to the theoretical curves (red).  The overlap is so good that it is difficult to distinguish the curves
plot.ecdf(sims)
curve(pchisq(9*x, df=9), add=TRUE, col="red")


plot(density(sims))
curve(9*dchisq(9*x, df=9), add=TRUE, col="red")


