It seems to be "folklore" knowledge that all the (source) symmetries of the $d$-Dirichlet energy are conformal maps.
Specifically, I have found this nice proof for the following claim:
Proposition: Let $M,N$ be oriented $n$-dimensional Riemannian manifolds, let $f:M \to N$ be a (local) diffeomorphism. Suppose that for any $h \in C^{\infty}(N)$, $$ E_N^p(h)=E_M^p(h \circ f)$$ where $E^p$ is the $p$-energy: $E_N^p(h)=\int_N \|dh\|^p \operatorname{Vol}_N$. (we assume $p\ge 1$).
Then,
If $p \neq n$, $\,f$ is an isometry. If $p = n, \,f$ is conformal.
Are there any other proofs around? While the proof I mentioned above is nice, there are some "magical choices" involved, and I wonder if there is a simpler proof (with less computations).
Edit: (The pointwise problem is easy-the challenge is to pass from global to local)
Using michelle's comment, let's try to build an elementary proof:
Suppose $f:M \to N$ is an orientation-preserving diffeomorphism:
We want $$\int_M \|dh \circ df\|^p \operatorname{Vol}_M=\int_N \|dh\|^p \operatorname{Vol}_N=\int_M f^*(\|dh \|^p \operatorname{Vol}_N)=\int_M (\|dh \|^p \circ f ) f^*\operatorname{Vol}_N=\int_M (\|dh \|^p \circ f ) \det df \operatorname{Vol}_M. \tag{1}$$
So, $$ \int_M \|dh \circ df\|^p \operatorname{Vol}_M=\int_M (\|dh \|^p \circ f ) \det df \operatorname{Vol}_M \, \, \text{for every} \, h \in C^{\infty}(N).$$
Since $h$ can be arbitrary one can hope to deduce equality of the integrands, i.e
$$ \|dh \circ df\|^p = (\|dh \|^p \circ f ) \det df \, \, \text{for every} \, h \in C^{\infty}(N). \tag{2}$$
The main challenge: Passing from equation $(1)$ to equation $(2)$; Something like selecting $h$ to be very "localized" and using a Lebesgue differentiation-type theorem might work. We can try selecting $h$ which are quickly decaying, and are zero outside a small set. But then at the "transition phase" their differential will have a non-negligible integral.
(Note that if we assume that $f$ restricted to arbitrary small open subsets preserve the energy, then the proof is trivial.)
Writing equation $(2)$ explicitly, in a given point $x \in M$, we get
$$ \|dh_{f(x)} \circ df_x\|^p = (\|dh_{f(x)} \|^p ) \det df_x $$
Since $dh_{f(x)}$ can be chosen at will, this means (thinking now in matrix terms, writing $df_x=A \in M_n,dh_{f(x)}=B \in \mathbb{R}^n$) that
$$ \|BA\|^p=\|B\|^p\det A \tag{3} \, \, \text{for every} \, B \in \mathbb{R}^n.$$
Taking $B=e_i$, we get $$\|A_{i\rightarrow}\|^p=\det A, \tag{4}$$ i.e the norms of all rows are identical.
Taking $B=e_i+e_j$, we get $\|A_{i\rightarrow}+A_{j\rightarrow}\|^p=(\sqrt{2})^p\det A=(\sqrt{2}\|A_{i\rightarrow}\|)^p$. This implies $A_{i\rightarrow},A_{j\rightarrow}$ are orthogonal; Indeed
$$ 2\|A_{i\rightarrow}\|^2=\|A_{i\rightarrow}+A_{j\rightarrow}\|^2=2\|A_{i\rightarrow}\|^2+2\langle A_{i\rightarrow}, A_{j\rightarrow} \rangle.$$
So the rows of $A$ are orthogonal and have the same norm, thus $A$ is conformal. Suppose $A \in \lambda \text{SO}(n)$, and plug this into equation $(3)$:
$$ \|B\|^p\lambda^p=\|B\|^p\lambda^n \Rightarrow \lambda^p=\lambda^n.$$
So, if $p \neq n$, $\lambda=1$ and $A$ is an isometry. If $p=n$, we can only conclude $A$ is conformal.
