[proofplan]
We compute how each of the three terms in the energy transforms under the scaling
\begin{align*}
u_\lambda(t,x)=\lambda^{2/(p-1)}u(\lambda t,\lambda x).
\end{align*}
The time-derivative and spatial-gradient terms acquire the same scaling exponent, and the nonlinear potential term is checked by direct algebra to acquire that same exponent. The energy is invariant exactly when this common exponent is zero, and solving that equation gives
\begin{align*}
p=\frac{n+2}{n-2}.
\end{align*}
[/proofplan]
[step:Declare the scaled derivatives and the energy at fixed time]
Fix $\lambda > 0$ and $t \in \mathbb{R}$. We carry out the computation for smooth finite-energy functions $u: \mathbb{R} \times \mathbb{R}^n \to \mathbb{R}$, so that the chain rule applies pointwise. Define the scaling exponent $a \in \mathbb{R}$ by
\begin{align*}
a := \frac{2}{p-1}.
\end{align*}
The scaled function is the map
\begin{align*}
u_\lambda: \mathbb{R} \times \mathbb{R}^n \to \mathbb{R}.
\end{align*}
Its value at $(t,x) \in \mathbb{R} \times \mathbb{R}^n$ is
\begin{align*}
u_\lambda(t,x) = \lambda^a u(\lambda t,\lambda x).
\end{align*}
For a finite-energy function $v:\mathbb R\times\mathbb R^n\to\mathbb R$, write
\begin{align*}
E[v](s):=\int_{\mathbb R^n}\left(\frac{1}{2}|\partial_t v(s,x)|^2+\frac{1}{2}|\nabla v(s,x)|^2+\frac{1}{p+1}|v(s,x)|^{p+1}\right)\,d\mathcal L^n(x).
\end{align*}
By the chain rule, for each $x \in \mathbb{R}^n$,
\begin{align*}
\partial_t u_\lambda(t,x) = \lambda^{a+1}(\partial_t u)(\lambda t,\lambda x).
\end{align*}
Similarly, the spatial gradient satisfies
\begin{align*}
\nabla u_\lambda(t,x) = \lambda^{a+1}(\nabla u)(\lambda t,\lambda x).
\end{align*}
Therefore
\begin{align*}
E[u_\lambda](t) = \int_{\mathbb{R}^n} \left( \frac{1}{2}\lambda^{2a+2}|\partial_t u(\lambda t,\lambda x)|^2 + \frac{1}{2}\lambda^{2a+2}|\nabla u(\lambda t,\lambda x)|^2 + \frac{1}{p+1}\lambda^{a(p+1)}|u(\lambda t,\lambda x)|^{p+1} \right)\, d\mathcal{L}^n(x).
\end{align*}
[/step]
[step:Change variables in the quadratic energy terms]
Use the change of variables $y = \lambda x$. Under this substitution, $x = \lambda^{-1}y$ and the [Lebesgue measure](/page/Lebesgue%20Measure) transforms as
\begin{align*}
d\mathcal{L}^n(x) = \lambda^{-n}\, d\mathcal{L}^n(y).
\end{align*}
Thus the kinetic term scales as
\begin{align*}
\int_{\mathbb{R}^n} \frac{1}{2}|\partial_t u_\lambda(t,x)|^2\, d\mathcal{L}^n(x)
= \lambda^{2a+2-n}\int_{\mathbb{R}^n} \frac{1}{2}|\partial_t u(\lambda t,y)|^2\, d\mathcal{L}^n(y).
\end{align*}
The same change of variables gives the gradient term scaling
\begin{align*}
\int_{\mathbb{R}^n} \frac{1}{2}|\nabla u_\lambda(t,x)|^2\, d\mathcal{L}^n(x)
= \lambda^{2a+2-n}\int_{\mathbb{R}^n} \frac{1}{2}|\nabla u(\lambda t,y)|^2\, d\mathcal{L}^n(y).
\end{align*}
[guided]
Define $a$ by
\begin{align*}
a:=\frac{2}{p-1},
\end{align*}
so that $u_\lambda(t,x)=\lambda^a u(\lambda t,\lambda x)$. The quadratic part of the energy contains two terms: the kinetic term involving $\partial_t u$ and the gradient term involving $\nabla u$. Both derivatives gain the same extra factor of $\lambda$ because both time and space are rescaled by the same factor $\lambda$.
For the kinetic term, the chain rule gives
\begin{align*}
\partial_t u_\lambda(t,x) = \lambda^{a+1}(\partial_t u)(\lambda t,\lambda x).
\end{align*}
Therefore
\begin{align*}
|\partial_t u_\lambda(t,x)|^2 = \lambda^{2a+2}|\partial_t u(\lambda t,\lambda x)|^2.
\end{align*}
Now set $y=\lambda x$. The map $x \mapsto y=\lambda x$ sends $\mathbb{R}^n$ bijectively onto $\mathbb{R}^n$, and Lebesgue measure transforms by
\begin{align*}
d\mathcal{L}^n(x)=\lambda^{-n}\,d\mathcal{L}^n(y).
\end{align*}
Hence
\begin{align*}
\int_{\mathbb{R}^n} \frac{1}{2}|\partial_t u_\lambda(t,x)|^2\, d\mathcal{L}^n(x)
= \lambda^{2a+2-n}\int_{\mathbb{R}^n} \frac{1}{2}|\partial_t u(\lambda t,y)|^2\, d\mathcal{L}^n(y).
\end{align*}
The gradient term is identical in structure. Since
\begin{align*}
\nabla u_\lambda(t,x)=\lambda^{a+1}(\nabla u)(\lambda t,\lambda x),
\end{align*}
we obtain
\begin{align*}
\int_{\mathbb{R}^n} \frac{1}{2}|\nabla u_\lambda(t,x)|^2\, d\mathcal{L}^n(x)
= \lambda^{2a+2-n}\int_{\mathbb{R}^n} \frac{1}{2}|\nabla u(\lambda t,y)|^2\, d\mathcal{L}^n(y).
\end{align*}
Thus both quadratic terms scale with the same exponent $2a+2-n$.
[/guided]
[/step]
[step:Change variables in the potential energy term]
For the potential term, the same substitution $y=\lambda x$ gives
\begin{align*}
\int_{\mathbb{R}^n} \frac{1}{p+1}|u_\lambda(t,x)|^{p+1}\, d\mathcal{L}^n(x)
= \lambda^{a(p+1)-n}\int_{\mathbb{R}^n} \frac{1}{p+1}|u(\lambda t,y)|^{p+1}\, d\mathcal{L}^n(y).
\end{align*}
Because $a=2/(p-1)$,
\begin{align*}
a(p+1)-n = \frac{2(p+1)}{p-1}-n.
\end{align*}
Also
\begin{align*}
2a+2-n = \frac{4}{p-1}+2-n.
\end{align*}
These exponents agree, since
\begin{align*}
\frac{2(p+1)}{p-1}-n = \frac{2p+2}{p-1}-n = 2+\frac{4}{p-1}-n.
\end{align*}
Thus every term in the energy scales by the common factor
\begin{align*}
\lambda^\gamma,
\end{align*}
where
\begin{align*}
\gamma := \frac{4}{p-1}+2-n.
\end{align*}
Consequently
\begin{align*}
E[u_\lambda](t)=\lambda^\gamma E[u](\lambda t).
\end{align*}
[/step]
[step:Solve the zero-exponent condition]
The energy is invariant under the scaling for all $\lambda>0$ and all finite-energy functions $u$ precisely when
\begin{align*}
\lambda^\gamma = 1
\end{align*}
for every $\lambda>0$. For necessity, choose any nonzero compactly supported smooth function $u: \mathbb{R} \times \mathbb{R}^n \to \mathbb{R}$ independent of $t$ with nonzero spatial gradient or nonzero value; then $E[u](\lambda t)>0$ for some $t \in \mathbb{R}$, so invariance forces $\lambda^\gamma=1$. This holds for every $\lambda>0$ if and only if $\gamma=0$. Therefore energy invariance is equivalent to
\begin{align*}
\frac{4}{p-1}+2-n=0.
\end{align*}
Solving,
\begin{align*}
\frac{4}{p-1}=n-2.
\end{align*}
Since $n\geq 3$, we have $n-2>0$, and hence
\begin{align*}
p-1=\frac{4}{n-2}.
\end{align*}
Thus
\begin{align*}
p=1+\frac{4}{n-2}=\frac{n+2}{n-2}.
\end{align*}
Therefore the conserved energy is invariant under the scaling exactly for the exponent
\begin{align*}
p=\frac{n+2}{n-2}.
\end{align*}
This proves the scaling criterion.
[/step]