Consider a macroscopic free energy functional of the form $$\mathcal{F}_\beta(\mu):= \frac{1}{\beta}\int_{\mathbb{R}^d}\log(\mu)\mu dx + \int_{\mathbb{R}^d}V(x)\mu(x)dx + \iint_{(\mathbb{R}^d)^2}g(x-y)\mu(x)\mu(y)dxdy.$$ Above, $\beta \in (0,\infty]$, $\mu$ is a probability density on $\mathbb{R}^d$. Suppose that the confining potential $V:\mathbb{R}^d\rightarrow\mathbb{R}$ is $C^2$ and uniformly convex (i.e., $\mathrm{Hess} V\geq \lambda$ for some $\lambda>0$) and that the interaction potential $g: \mathbb{R}^d$ is globally Lipschitz (i.e., $\|\nabla g\|_{L^\infty}<\infty$).
Question. Is it possible to choose $\beta, V, g$, subject to the above constraints, such that $\mathcal{F}_{\beta}$ has a non-minimizing critical point?
Certainly, if $g$ is small (e.g., $\|g\|_{\dot{C}^2}\ll 1$), then the uniform convexity of $V$ will dominate, ensuring there is a unique critical point of $\mathcal{F}_{\beta}$ which is a minimizer. However, when $g$ is not small, convexity is lost, so it seems in principle a non-minimizing critical point could exist in this regime. But I can't readily think of a counterexample. Note that on the torus $\mathbb{T}^d$ with $V=0$, it is easy to cook up interaction potentials with non-minimizing critical points (e.g., $g(x) = \cos(2\pi x)$ for $d=1$).
This is an elaboration on Fedja's answer below.
Fix $T>0$ and apply to $1_{[-T,T]}$ the Ornstein-Uhlenbeck (OU) semigroup with stationary measure $dG(y):=(2\pi/\lambda)^{-1/2}e^{-\lambda|x|^2/2}dx$: $$f_{\epsilon}(x) := \int_{\mathbb{R}}1_{[-T,T]}(e^{-t}x+\sqrt{1-e^{-2t}}y)dG(y).$$ The OU flow preserves log-concavity, so $V_\epsilon(x):=\frac{\lambda}{2}|x|^2 -\log(f_\epsilon(x))$ is $\lambda$-convex. Also, set $V_0(x):=\frac{\lambda}{2}|x|^2-\log(1_{[-T,T]}(x))$.
The Fourier transform $\mathscr{F}(e^{-V_0})$of $e^{-V_0(x)}$ cannot be nonnegative. Otherwise, we would have $\frac{1}{2\pi}\int_{\mathbb{R}}\left|\mathscr{F}(e^{-V_0})\right|d\xi = e^{-V_0(0)}<\infty$, which would imply that $e^{-V_0}$ is continuous, a contradiction. Let $\alpha=\min \mathscr{F}(e^{-V_0}) < 0$. Choose $\epsilon$ sufficiently small so that $$\|\mathscr{F}(e^{-V_0}) - \mathscr{F}(e^{-V_\epsilon})\|_{L^\infty}\leq \|e^{-V_0}-e^{-V_\epsilon}\|_{L^1} < \frac{|\alpha|}{2}.$$ Then $\mathscr{F}(e^{-V_\epsilon})$ is strictly negative at some point. Since it is also continuous and strictly positive at the origin, by the intermediate value theorem, there exists a point $a$ at which $\mathscr{F}(e^{-V_\epsilon})(a)=0$. We fix such an $a$ for the remainder of the argument.
As explained by Fedja, for any $A\in\mathbb{R}$, one critical point of the free energy $$\mathcal{F}_A(\mu) = \int_{\mathbb{R}}\log(\mu)\mu dx + \int_{\mathbb{R}}V_\epsilon\mu dx + \frac{A}{2}\iint_{\mathbb{R}^2}\cos(a(x-y))\mu(x)\mu(y)dxdy$$ is given by $\mu_0=\frac{1}{Z_\epsilon}e^{-V_\epsilon}$, where $Z_\epsilon$ is the partition function. We argue that for a suitable range of $A$, this critical point is non-minimizing. Hereafter to, $\epsilon$ is fixed, so we write $V=V_\epsilon$ and omit the dependence in all our notation.
For $0<\delta\ll 1$, let $\mu_\delta = \frac{1}{Z_\delta}e^{\delta\cos(ax)}e^{-V}$, where $Z_\delta$ is the partition function. From the definition of $\mu_\delta$, \begin{align} \mathcal{F}_A(\mu_\delta) &= \int_{\mathbb{R}}\log\frac{e^{\delta\cos(ax)}}{Z_\delta}\mu_\delta dx+\frac12\iint_{\mathbb{R}^2}\cos(a(x-y))\mu_\delta(x)\mu_\delta(y)dxdy\\ &= \frac{\delta}{Z_\delta}\int_{\mathbb{R}}\cos(ax)e^{\delta\cos(ax)}e^{-V} dx - \log Z_\delta + \frac{A}{2Z_\delta^2}\iint_{\mathbb{R}^2}\cos(a(x-y))e^{\delta(\cos(ax)+\cos(ay))-V(x)-V(y)}dxdy. \end{align} We Taylor expand $$e^{\delta\cos(ax)}=1+\delta\cos(ax)+\frac{\delta}{2}\cos^2(ax) + O(\delta^3).$$ Substituting this expansion and using that $\mathscr{F}(e^{-V})(a)=0$, we find $$\frac{\delta}{Z_\delta}\int_{\mathbb{R}}\cos(ax)e^{\delta\cos(ax)}e^{-V} dx=\frac{\delta^2}{Z_\delta}\int_{\mathbb{R}}\cos^2(ax)e^{-V} + \frac{O(\delta^3)}{Z_\delta},$$ \begin{align} \log(Z_\delta) &= \log\left(\int_{\mathbb{R}}[1+\delta\cos(ax)+\frac{\delta^2}{2}\cos^2(ax)+O(\delta^3)]e^{-V}dx\right) \\ &=\log\left(\int_{\mathbb{R}}e^{-V}dx + \frac{\delta^2}{2}\int_{\mathbb{R}}\cos^2(ax)e^{-V}dx + O(\delta^3)\right), \end{align} and \begin{align} &\frac{A}{2Z_\delta^2}\iint_{\mathbb{R}^2}\cos(a(x-y))e^{\delta(\cos(ax)+\cos(ay))-V(x)-V(y)}dxdy \\ &=\frac{A\delta^2}{2Z_\delta^2}\iint_{\mathbb{R}^2}\cos(a(x-y))\cos(ax)\cos(ay)e^{-V(x)-V(y)}dxdy + \frac{AO(\delta^3)}{Z_\delta^2}\\ &=\frac{A\delta^2}{2Z_\delta^2}\left(\int_{\mathbb{R}}\cos^2(ax)e^{-V}dx\right)^2 + \frac{AO(\delta^3)}{Z_\delta^2}. \end{align} Choose $$A=A_\delta = -2Z_\delta\left(\int_{\mathbb{R}}\cos^2(ax)e^{-V}dx\right)^{-1}.$$ From inspection, one sees that there is a universal $\delta_0>0$, such that for all $\delta\in[0,\delta_0]$, $Z_\delta \in [\frac{2Z}{3},\frac{3Z}{2}]$. Hence, $A_\delta$ is uniformly bounded for $\delta\in [0,\delta_0]$. With our choice of $A_\delta$, a little bookkeeping shows that $$\mathcal{F}_{A_\delta}(\mu_\delta)=-\log\left(\int_{\mathbb{R}}e^{-V}dx + \frac{\delta^2}{2}\int_{\mathbb{R}}\cos^2(ax)e^{-V}dx + O(\delta^3)\right)+\frac{O(\delta^3)}{Z_\delta}+\frac{A_\delta O(\delta^3)}{Z_\delta^2}.$$ Since $\int_{\mathbb{R}}\cos^2(ax)e^{-V}dx>0$, we can now take $\delta\in (0,\delta_0]$ sufficiently small so that the preceding RHS is $<$ $$-\log\left(\int_{\mathbb{R}}e^{-V}dx\right)=\mathcal{F}_{A_\delta}(\mu_0).$$ In other words, we have found two distinct probability densities $\mu_0,\mu_\delta$ such that $$\mathcal{F}_{A_\delta}(\mu_\delta)<\mathcal{F}_{A_\delta}(\mu_0).$$ Thus, $\mu_0$ is a non-minimizing critical point.
