[proofplan] We compute the first variation of each term in $\mathcal E$ at a smooth positive density $\rho$ along a smooth mass-preserving perturbation. The internal and potential terms give $U'(\rho)$ and $V$, while the evenness of $W$ allows the interaction variation to symmetrise to $W*\rho$. Thus the formal variational derivative is $U'(\rho)+V+W*\rho$. The formal Wasserstein gradient-flow rule then sets the velocity equal to the negative spatial gradient of this variational derivative, and substituting this velocity into the continuity equation gives the displayed PDE. [/proofplan]

[step:Compute the first variation along a mass-preserving perturbation]Fix a time $t\in(0,T)$ and abbreviate $\rho=\rho_t$. Let \begin{align*} \sigma:\mathbb R^n\to\mathbb R \end{align*} be a smooth admissible perturbation with sufficient decay and zero total mass: \begin{align*} \int_{\mathbb R^n}\sigma(x)\,d\mathcal L^n(x)=0. \end{align*} For $\varepsilon$ in an interval around $0$, define \begin{align*} \rho_\varepsilon:\mathbb R^n&\to(0,\infty) \end{align*} by $\rho_\varepsilon(x)=\rho(x)+\varepsilon\sigma(x)$, assuming $\rho_\varepsilon$ remains positive and admissible for the formal computation. The zero-mass condition gives \begin{align*} \int_{\mathbb R^n}\rho_\varepsilon(x)\,d\mathcal L^n(x)=1, \end{align*} so $\rho_\varepsilon$ is a mass-preserving variation. By the assumed justification for differentiating under the integral sign and by the chain rule applied to $U\in C^2((0,\infty);\mathbb R)$, \begin{align*} \left.\frac{d}{d\varepsilon}\right|_{\varepsilon=0}\int_{\mathbb R^n}U(\rho_\varepsilon(x))\,d\mathcal L^n(x)=\int_{\mathbb R^n}U'(\rho(x))\sigma(x)\,d\mathcal L^n(x). \end{align*} Similarly, differentiating the potential term gives \begin{align*} \left.\frac{d}{d\varepsilon}\right|_{\varepsilon=0}\int_{\mathbb R^n}V(x)\rho_\varepsilon(x)\,d\mathcal L^n(x)=\int_{\mathbb R^n}V(x)\sigma(x)\,d\mathcal L^n(x). \end{align*}[/step]

[guided]Fix $t\in(0,T)$ and write $\rho=\rho_t$ to focus on the variational calculation at one density. A Wasserstein tangent vector preserves total mass, so we test the first variation against a smooth admissible perturbation \begin{align*} \sigma:\mathbb R^n\to\mathbb R \end{align*} satisfying \begin{align*} \int_{\mathbb R^n}\sigma(x)\,d\mathcal L^n(x)=0. \end{align*} For $\varepsilon$ sufficiently close to $0$, define \begin{align*} \rho_\varepsilon:\mathbb R^n&\to(0,\infty) \end{align*} by $\rho_\varepsilon(x)=\rho(x)+\varepsilon\sigma(x)$, with admissibility and positivity included in the formal hypotheses. Then \begin{align*} \int_{\mathbb R^n}\rho_\varepsilon(x)\,d\mathcal L^n(x)=\int_{\mathbb R^n}\rho(x)\,d\mathcal L^n(x)+\varepsilon\int_{\mathbb R^n}\sigma(x)\,d\mathcal L^n(x)=1. \end{align*} Thus the variation stays in the affine hyperplane of probability densities. For the internal energy, the map $\varepsilon\mapsto U(\rho(x)+\varepsilon\sigma(x))$ is differentiable for each $x$ because $U\in C^2((0,\infty);\mathbb R)$ and $\rho_\varepsilon(x)>0$. The theorem assumes enough decay and domination to move the $\varepsilon$-derivative through the integral. Therefore the chain rule gives \begin{align*} \left.\frac{d}{d\varepsilon}\right|_{\varepsilon=0}\int_{\mathbb R^n}U(\rho_\varepsilon(x))\,d\mathcal L^n(x)=\int_{\mathbb R^n}U'(\rho(x))\sigma(x)\,d\mathcal L^n(x). \end{align*} For the potential term, the integrand depends linearly on $\rho_\varepsilon$. Differentiation under the integral sign gives \begin{align*} \left.\frac{d}{d\varepsilon}\right|_{\varepsilon=0}\int_{\mathbb R^n}V(x)\rho_\varepsilon(x)\,d\mathcal L^n(x)=\int_{\mathbb R^n}V(x)\sigma(x)\,d\mathcal L^n(x). \end{align*} These two computations identify the internal and confinement contributions to the variational derivative as $U'(\rho)$ and $V$.[/guided]

[step:Use evenness of $W$ to identify the interaction variation] Define the interaction functional $\mathcal W$ on admissible densities by \begin{align*} \mathcal W[\eta]=\frac12\int_{\mathbb R^n}\int_{\mathbb R^n}W(x-y)\eta(x)\eta(y)\,d\mathcal L^n(y)\,d\mathcal L^n(x). \end{align*} Using the assumed differentiability under the double integral, \begin{align*} \left.\frac{d}{d\varepsilon}\right|_{\varepsilon=0}\mathcal W[\rho_\varepsilon]=\frac12\int_{\mathbb R^n}\int_{\mathbb R^n}W(x-y)\sigma(x)\rho(y)\,d\mathcal L^n(y)\,d\mathcal L^n(x)+\frac12\int_{\mathbb R^n}\int_{\mathbb R^n}W(x-y)\rho(x)\sigma(y)\,d\mathcal L^n(y)\,d\mathcal L^n(x). \end{align*} In the second double integral, interchange the names of the variables $x$ and $y$. Since $\mathcal L^n\otimes\mathcal L^n$ is invariant under this coordinate swap and $W$ is even, \begin{align*} W(y-x)=W(-(x-y))=W(x-y). \end{align*} Hence the second double integral equals the first. Therefore \begin{align*} \left.\frac{d}{d\varepsilon}\right|_{\varepsilon=0}\mathcal W[\rho_\varepsilon]=\int_{\mathbb R^n}\int_{\mathbb R^n}W(x-y)\sigma(x)\rho(y)\,d\mathcal L^n(y)\,d\mathcal L^n(x). \end{align*} By the definition of convolution, \begin{align*} \left.\frac{d}{d\varepsilon}\right|_{\varepsilon=0}\mathcal W[\rho_\varepsilon]=\int_{\mathbb R^n}(W*\rho)(x)\sigma(x)\,d\mathcal L^n(x). \end{align*} [/step]

[step:Identify the formal variational derivative of $\mathcal E$] Combining the three first-variation formulas, we obtain \begin{align*} \left.\frac{d}{d\varepsilon}\right|_{\varepsilon=0}\mathcal E[\rho_\varepsilon]=\int_{\mathbb R^n}\left(U'(\rho(x))+V(x)+(W*\rho)(x)\right)\sigma(x)\,d\mathcal L^n(x). \end{align*} Thus the formal variational derivative of $\mathcal E$ at $\rho$ is \begin{align*} \frac{\delta\mathcal E}{\delta\rho}=U'(\rho)+V+W*\rho. \end{align*} Since $t\in(0,T)$ was arbitrary, this gives for the curve $(\rho_t)_{t\in(0,T)}$ the time-dependent potential \begin{align*} \Phi_t:\mathbb R^n&\to\mathbb R \end{align*} defined by \begin{align*} \Phi_t(x)=U'(\rho_t(x))+V(x)+(W*\rho_t)(x). \end{align*} By hypothesis, $\Phi_t$ is $C^1$ in the spatial variable. [/step]

[step:Substitute the Wasserstein gradient velocity into the continuity equation] In the formal $2$-Wasserstein Riemannian calculus, the gradient-flow velocity associated to the first variation $\Phi_t=\delta\mathcal E/\delta\rho_t$ is \begin{align*} v_t:\mathbb R^n&\to\mathbb R^n \end{align*} defined by \begin{align*} v_t(x)=-\nabla\Phi_t(x)=-\nabla\left(U'(\rho_t(x))+V(x)+(W*\rho_t)(x)\right). \end{align*} The formal gradient-flow curve is therefore required to satisfy the continuity equation \begin{align*} \partial_t\rho_t+\nabla\cdot(\rho_t v_t)=0. \end{align*} Substituting the displayed formula for $v_t$ gives \begin{align*} \partial_t\rho_t+\nabla\cdot\left(-\rho_t\nabla\left(U'(\rho_t)+V+W*\rho_t\right)\right)=0. \end{align*} Equivalently, \begin{align*} \partial_t\rho_t=\nabla\cdot\left(\rho_t\nabla\left(U'(\rho_t)+V+W*\rho_t\right)\right). \end{align*} This is exactly the asserted formal aggregation-diffusion Wasserstein gradient-flow equation. [/step]