[proofplan]
Run the harmonic map heat flow starting from the given smooth map $u_0$. The Eells-Sampson heat-flow theorem gives global smooth existence, preservation of the initial homotopy class along finite time intervals, and compactness of large-time slices when the target is compact with nonpositive sectional curvature. A smooth subsequential limit is harmonic because the tension field tends to zero along a large-time sequence, and it remains homotopic to $u_0$ by combining the heat-flow homotopy with a short geodesic homotopy to the limit.
[/proofplan]
[step:Run harmonic map heat flow from the initial map]
Let $\tau(u)$ denote the tension field of a smooth map $u:M\to N$, viewed as a smooth section of $u^*TN$. Apply the [Eells-Sampson Heat Flow Theorem](/page/Eells-Sampson%20Heat%20Flow%20Theorem) to the closed Riemannian manifold $(M,g)$, the compact Riemannian manifold $(N,h)$ with nonpositive sectional curvature, and the smooth initial map $u_0:M\to N$ from the theorem statement. It gives a smooth map $u:M\times [0,\infty)\to N$, where $u(x,t)$ denotes the value at $(x,t)$, satisfying the harmonic map [heat equation](/page/Heat%20Equation)
\begin{align*}
\partial_t u(x,t)=\tau(u)(x,t), \qquad u(x,0)=u_0(x),
\end{align*}
for every $x\in M$ and $t\geq 0$.
[guided]
The natural way to find a harmonic representative is to deform $u_0$ by the negative gradient flow of the energy. For a smooth map $v:M\to N$, let $v^*TN\to M$ denote the pullback tangent bundle, and let $\Gamma(v^*TN)$ denote the [vector space](/page/Vector%20Space) of smooth sections of this bundle. The Euler-Lagrange operator for the energy is the tension field $\tau(v)\in \Gamma(v^*TN)$, and harmonicity is the equation $\tau(v)=0$.
We apply the [Eells-Sampson Heat Flow Theorem](/page/Eells-Sampson%20Heat%20Flow%20Theorem). Its hypotheses are exactly the geometric hypotheses in the theorem statement: $(M,g)$ is closed, $(N,h)$ is compact, and every sectional curvature of $h$ is nonpositive. Therefore the theorem supplies a global smooth solution $u:M\times [0,\infty)\to N$, where $u(x,t)$ denotes the value at $(x,t)$, to
\begin{align*}
\partial_t u(x,t)=\tau(u)(x,t), \qquad u(x,0)=u_0(x).
\end{align*}
This is the point where closedness of $M$ and compactness plus nonpositive curvature of $N$ are used: they are the assumptions that prevent finite-time breakdown and provide the estimates needed at large time.
[/guided]
[/step]
[step:Choose large-time slices with vanishing tension]
For $t\geq 0$, define the energy function $E:[0,\infty)\to [0,\infty)$ by
\begin{align*}
E(t)=\frac{1}{2}\int_M |\nabla u(x,t)|_{g,h}^2\,d\operatorname{vol}_g(x),
\end{align*}
where $d\operatorname{vol}_g$ is the Riemannian volume measure on $M$ and $|\nabla u|_{g,h}$ is the Hilbert-Schmidt norm induced by $g$ and $h$. Let $\mathcal{L}^1$ denote one-dimensional Lebesgue measure on $[0,\infty)$. The Eells-Sampson energy identity gives
\begin{align*}
E(T)-E(0)=-\int_0^\top\int_M |\tau(u)(x,t)|_h^2\,d\operatorname{vol}_g(x)\,d\mathcal{L}^1(t)
\end{align*}
for every $T>0$. Since $E(T)\geq 0$, the finite-time identity gives
\begin{align*}
\int_0^\top\int_M |\tau(u)(x,t)|_h^2\,d\operatorname{vol}_g(x)\,d\mathcal{L}^1(t)\leq E(0)
\end{align*}
for every $T>0$. Applying the [Monotone Convergence Theorem](/page/Monotone%20Convergence%20Theorem) to the nonnegative functions $\mathbb{1}_{[0,T]}(t)|\tau(u)(x,t)|_h^2$ on $M\times [0,\infty)$ with respect to $d\operatorname{vol}_g(x)\otimes d\mathcal{L}^1(t)$ yields
\begin{align*}
\int_0^\infty\int_M |\tau(u)(x,t)|_h^2\,d\operatorname{vol}_g(x)\,d\mathcal{L}^1(t)\leq E(0)<\infty.
\end{align*}
Hence there exists a sequence $(t_j)_{j=1}^\infty$ in $[0,\infty)$ with $t_j\to\infty$ and
\begin{align*}
\int_M |\tau(u)(x,t_j)|_h^2\,d\operatorname{vol}_g(x)\to 0.
\end{align*}
[/step]
[step:Extract a smooth large-time limit]
For each $j\in\mathbb{N}$, define the large-time slice $u_j:M\to N$ by $u_j(x)=u(x,t_j)$ for $x\in M$. We now use the large-time compactness conclusion of the [Eells-Sampson Heat Flow Theorem](/page/Eells-Sampson%20Heat%20Flow%20Theorem) as a substantial external input: for heat-flow solutions with closed domain and compact nonpositively curved target, the large-time slices have a smoothly convergent subsequence. Therefore there is a subsequence, still denoted $(u_j)_{j=1}^\infty$, and a smooth map $u_\infty:M\to N$ such that $u_j\to u_\infty$ in $C^k(M,N)$ for every integer $k\geq 0$, where $C^k(M,N)$ denotes convergence of the maps and their coordinate derivatives through order $k$ uniformly on $M$ in any finite smooth atlas on $N$. Since the tension operator is a second-order smooth differential operator in local coordinates, $C^2$ convergence implies
\begin{align*}
\tau(u_j)\to \tau(u_\infty)
\end{align*}
uniformly on $M$. Here $L^2$ convergence of $\tau(u_j)$ to zero means convergence to zero in the norm
\begin{align*}
\|\tau(u_j)\|_{L^2(M)}^2=\int_M |\tau(u_j)(x)|_h^2\,d\operatorname{vol}_g(x).
\end{align*}
Combining the [uniform convergence](/page/Uniform%20Convergence) of $\tau(u_j)$ to $\tau(u_\infty)$ with this $L^2$ convergence gives $\tau(u_\infty)=0$ on $M$. Thus $u_\infty$ is harmonic.
[/step]
[step:Keep the limiting map in the original homotopy class]
For each $T>0$, define the finite-time heat-flow homotopy $H_T:M\times [0,1]\to N$ by $H_T(x,s)=u(x,sT)$ for $(x,s)\in M\times [0,1]$. This map is smooth, satisfies $H_T(x,0)=u_0(x)$ and $H_T(x,1)=u(x,T)$, and therefore shows that $u(\cdot,T)$ is homotopic to $u_0$.
It remains to compare $u_j$ with $u_\infty$. Since $N$ is compact, its [injectivity radius](/page/Injectivity%20Radius) is a positive number $\rho>0$. Since $u_j\to u_\infty$ uniformly on $M$, choose $j_0\in\mathbb{N}$ such that
\begin{align*}
d_h(u_{j_0}(x),u_\infty(x))<\rho
\end{align*}
for every $x\in M$, where $d_h$ is the Riemannian distance induced by $h$. Define the short geodesic homotopy $G:M\times [0,1]\to N$ by
\begin{align*}
G(x,s)=\exp_{u_{j_0}(x)}\left(s\,\exp_{u_{j_0}(x)}^{-1}(u_\infty(x))\right)
\end{align*}
for $(x,s)\in M\times [0,1]$, where $\exp_p:T_pN\to N$ is the Riemannian [exponential map](/page/Exponential%20Map). The injectivity-radius bound makes $\exp_{u_{j_0}(x)}^{-1}(u_\infty(x))$ well-defined and smoothly dependent on $x$. Hence $G$ is a homotopy from $u_{j_0}$ to $u_\infty$. Since $u_{j_0}=u(\cdot,t_{j_0})$ is homotopic to $u_0$, transitivity of homotopy gives that $u_\infty$ is homotopic to $u_0$.
[/step]
[step:Conclude existence of the harmonic representative]
The map $u_\infty:M\to N$ is smooth, satisfies $\tau(u_\infty)=0$, and is homotopic to $u_0$. Therefore $u_\infty$ is a smooth harmonic representative in the homotopy class of $u_0$, which proves the theorem.
[/step]