There is no need to approximate the integral by Riemann sums in order to determine the expectation and variance of the integral $X:=\int_0^T W_s^2 \, ds$. Note that this integral is a "standard" Riemann integral, i.e. for each point $\omega \in \Omega$ the integral
$$X(\omega) = \int_0^T W_s(\omega)^2 \, ds$$
is defined as a Riemann integral. Everything is well-defined since $s \mapsto W_s^2(\omega)$ is continuous, hence Riemann-integrable.
An application of Fubini's theorem shows
$$\mathbb{E}(X) = \mathbb{E} \left( \int_0^T W_s^2 \, ds \right) = \int_0^T \mathbb{E}(W_s^2) \, ds.$$
Since $W_s$ is Gaussian with mean $0$ and variance $s$, we have $\mathbb{E}(W_s^2) = s$. Consequently, $$\mathbb{E}(X)=T^2/2$$
Calculating the variance is slightly more complicated. Recall that
$$\text{var}(X) = \mathbb{E}(X^2) - (\mathbb{E}(X))^2 \tag{1}$$
and, since we already know $\mathbb{E}(X)$, it suffices to calculate $\mathbb{E}(X^2)$. To this end, note that
$$X^2 = \left( \int_0^T W_s^2 \, ds \right) \left( \int_0^T W_t^2 \, dt \right) = \int_0^T \int_0^T W_s^2 W_t^2 \, ds \, dt.$$
Applying Fubini's theorem another time, we get
$$\mathbb{E}(X^2) = \int_0^T \int_0^T \mathbb{E}(W_s^2 W_t^2) \, ds \, dt = 2 \int_0^T \int_0^t \mathbb{E}(W_s^2 W_t^2) \, ds \, dt. \tag{2}$$
Using that $W_t^2-t$ is a martingale, it is not difficult to see that
$$\mathbb{E}(W_s^2 W_t^2) = st +2s^2 \qquad \text{for all} \, \, s \leq t.$$
Plugging this into $(2)$ yields $\mathbb{E}(X^2)$, and hence, by $(1)$, the variance of $X$.