$\newcommand\R{\mathbb R}\newcommand{\Ga}{\Gamma}$The integral in question is
\begin{equation*}
\begin{aligned}
\int_0^1 dr\,\int_{rS_{n-1}}dx\,f(x)
&=\int_0^1 dr\,r^{2m+n-1}\int_{S_{n-1}}du\,f(u) \\
&=\frac1{2m+n}\,|S_{n-1}|Ef(U_n),
\end{aligned}
\end{equation*}
where
\begin{equation*}
f(x):=\Big(\sum_{j=1}^n a_j x_j^2\Big)^m
\end{equation*}
(so that $f$ is $2m$-homogeneous), $|S_{n-1}|$ is the surface area of the unit sphere $S_{n-1}$ in $\R^n$, and $U_n$ is a random vector uniformly distributed on $S_{n-1}$.
So, the problem reduces to finding the asymptotics of $Ef(U_n)$. Without loss of generality,
\begin{equation*}
U_n=\frac G{\|G\|},
\end{equation*}
where $G=(G_1,\dots,G_n)$ is a standard Gaussian vector in $\R^n$ and $\|\cdot\|$ is the Euclidean norm.
Let us normalize the $a_j$'s by assuming that
\begin{equation*}
\sum_{j=1}^n a_j=1 \tag{0}
\end{equation*}
and consider the following two extreme cases:
Case 1: $a_1=1$, $a_2=\cdots=a_n=0$ and
Case 2: $a_1=\cdots=a_n=1/n$.
In Case 1,
\begin{equation*}
f(U_n)=\Big(\frac{G_1^2}{\|G\|^2}\Big)^m,
\end{equation*}
and $\dfrac{G_1^2}{\|G\|^2}$ has the beta distribution with parameters $1/2,(n-1)/2$. So, here
\begin{equation*}
Ef(U_n)=\frac{\Ga(n/2)\Ga(m+1/2)}{\Ga(1/2)\Ga(m+n/2)}.
\end{equation*}
So, if e.g. $m=n\to\infty$, then the integral in question is
\begin{equation*}
\sim\frac1{3n}\,\sqrt6\,\Big(\frac2{3\sqrt3}\Big)^n\,|S_{n-1}|. \tag{1}
\end{equation*}
So, here the coefficient of $|S_{n-1}|$ decreases exponentially.
In the trivial Case 2, $f(U_n)=1/n^m$, and hence for $m=n\to\infty$ the integral in question is
\begin{equation*}
=\frac1{3n^{n+1}}\,|S_{n-1}|,
\end{equation*}
which is very different from the asymptotics for Case 1 given in (1).
Between these two extreme cases, any intermediate asymptotics (or absense of any asymptotics) should be possible. As was noted in comments, the asymptotics in general (if any) will very much depend on how much the $a_j$'s differ from one another.
Addendum 1:
To be more specific, consider the following setting, intermediate between Cases 1 and 2:
$a_1=\cdots=a_k=1/k$, $a_{k+1}=\cdots=a_n=0$ for natural $k\in[2,n-1]$, where $k$ is allowed to vary with $n$ and $m$.
Then
\begin{equation*}
f(U_n)=\frac1{k^m}\,\Big(\frac{G_1^2+\cdots+G_k^2}{\|G\|^2}\Big)^m,
\end{equation*}
and $\dfrac{G_1^2+\cdots+G_k^2}{\|G\|^2}$ has the beta distribution with parameters $k/2,(n-k)/2$. So, here
\begin{equation*}
Ef(U_n)=\frac1{k^m}\,\prod_{r=0}^{m-1}\frac{k/2+r}{n/2+r}. \tag{2}
\end{equation*}
So, if $k/m^2\to\infty$, then
\begin{equation*}
Ef(U_n)\sim\frac1{n^m}, \tag{3}
\end{equation*}
which is the same asymptotics as for $m$ fixed.
However, if $n/m^2\to\infty$ and $k\sim cm^2$ for some $c\in(0,\infty)$, then
\begin{equation*}
Ef(U_n)\sim\frac{e^{1/c}}{n^m},
\end{equation*}
now with the additional factor $e^{1/c}$ as compared to (3).
If we now let $n/m^2\to\infty$ and $k\sim cm$ for some $c\in(0,\infty)$, then an additional exponentially growing (with $m$) factor (as compared to (3)) will appear.
Addendum 2: Letting $g(a_1,\dots,a_n):=Ef(U_n)$, we see that the function $g$ is symmetric and convex, and hence Schur convex -- see e.g. Theorem A on p. 258. So, in view of condition (0),
\begin{equation*}
Ef(U_n)\ge g(1/n,\dots,1/n)=\frac1{n^m}. \tag{$\clubsuit$}
\end{equation*}
That is, the smallest value of $Ef(U_n)$ is attained in Case 2 ($a_1=\cdots=a_n=1/n$), considered above.
Addendum 3: For completeness, consider also the case when $m$ is fixed (even though $m\to\infty$ in the OP).
Then, by (2),
\begin{equation*}
Ef(U_n)\sim\frac1{n^m}\,\prod_{r=0}^{m-1}\Big(1+\frac{2r}k\Big)
\end{equation*}
if $a_1=\cdots=a_k=1/k$ and $a_{k+1}=\cdots=a_n=0$ for a fixed natural $k\in[2,n-1]$ and $n\to\infty$. We see that, even when $m$ is fixed, the asymptotics depends on $k$ and, more generally, on how much the $a_j$'s differ from one another.
Addendum 4: Here we complement the lower bound on $Ef(U_n)$ given by ($\clubsuit$) in Addendum 2 by providing a matching upper bound on $Ef(U_n)$ that implies the following:
If the $a_i$'s are uniformly small in the sense that
\begin{equation*}
a:=\max_{i=1}^n a_i\to0
\end{equation*}
and, moreover, $m$ is at most moderately large in the sense that
\begin{equation*}
m^3 a\to0, \tag{4}
\end{equation*}
then
\begin{equation*}
Ef(U_n)\sim\frac1{n^m} \tag{$\heartsuit$}
\end{equation*}
(as $n\to\infty$).
Indeed, the random point $(x_1^2,\dots,x_n^2)$ has the Dirichlet distribution with parameters $1/n,\dots,1/n$, and the Dirichlet distribution has the negative association (NA) property. So, by (say) Theorem 2,
\begin{equation*}
Ef(U_n)\le E\Big(\sum_{j=1}^n a_j Y_j\Big)^m, \tag{5}
\end{equation*}
where the $Y_j$'s are iid random variables each with the beta distribution with parameters $1/2,(n-1)/2$. Denoting now the $L^m$ norm by $\|\cdot\|_m$ and using Minkowski's inequality, we get
\begin{equation*}
(Ef(U_n))^{1/m}\le\Big\|\sum_{j=1}^n a_j Y_j\Big\|_m
\le \Big\|\sum_{j=1}^n a_j EY_j\Big\|_m +\Big\|\sum_{j=1}^n a_j Z_j\Big\|_m,
\end{equation*}
where $Z_j:=Y_j-EY_j$. Since $EY_j=1/n$, we have
\begin{equation*}
\Big\|\sum_{j=1}^n a_j EY_j\Big\|_m
= \sum_{j=1}^n a_j EY_j=\frac1n. \tag{6}
\end{equation*}
Note also that $Var\,Y_j\sim2/n^2$ and $\|Z_j\|_m\le2\|Y_j\|_m\ll m/n$; here and in what follows, $A\ll B$ means $A\le CB$ for some universal real constant $C>0$.
Using now an appropriate version of Rosenthal's inequality (see e.g. Theorem 6.1), we get
\begin{equation*}
\begin{aligned}
\Big\|\sum_{j=1}^n a_j Z_j\Big\|_m
&\ll\frac1n\,(m^2 a^{1-1/m}+m^{1/2} a^{1/2}) \\
&=\frac1n\,\frac1m\,((m^3 a)^{1-1/m}m^{3/m}+(m^3 a)^{1/2})
=o\Big(\frac1n\,\frac1m\Big),
\end{aligned}
\end{equation*}
by (4).
So, by (5) and (6),
\begin{equation*}
(Ef(U_n))\le\frac1{n^m}\Big(1+\frac{o(1)}m\Big)^m
=\frac{1+o(1)}{n^m}.
\end{equation*}
Now ($\heartsuit$) follows, in view of ($\clubsuit$).
