Diffeomorphisms preserving harmonic functions

Solution 1:

These are harmonic morphisms: maps $\varphi$ that "preserve Laplace's equation $\Delta u=0$" upon precomposition by $\varphi$. There is a vast literature on them (they are more generally defined between arbitrary Riemannian manifolds using the Laplace-Beltrami operator).

In particular with Euclidean spaces, we have this characterization: a map is a harmonic morphism if and only if it is itself harmonic and "horizontally weakly conformal," meaning its components' gradients are orthogonal and equal length. Finding all harmonic morphisms between arbitrary open subsets of Euclidean spaces is an open problem. See Proposition $1.10$ in Wood's reference article Harmonic Morphisms Between Riemannian Manifolds. Wood & Helein have a broader-purpose article Harmonic Maps in the text Handbook of Global Analysis.