  We study nonnegative optimizers of a Gagliardo-Nirenberg type inequality
$$\iint_{\mathbb{R}^N \times \mathbb{R}^N} \frac{|u(x)|^p\,|u(y)|^p}{|x -
y|^{N-\alpha}} dx\, dy\le C\Big(\int_{{\mathbb R}^N}|u|^2 dx\Big)^{p\theta}
\Big(\int_{{\mathbb R}^N}|u|^q dx\Big)^{2p(1-\theta)/q},$$ that involves the
nonlocal Riesz energy with $0<\alpha<N$, $p>\frac{N+\alpha}{N}$,
$q>\frac{2Np}{N+\alpha}$ and $\theta=\frac{(N+\alpha)q-2Np}{Np(q-2)}$. For
$p=2$, the equivalent problem has been studied in connection with the
Keller-Segel diffusion-aggregation models in the past few decades. The general
case $p\neq 2$ considered here appears in the study of Thomas-Fermi limit
regime for the Choquard equations with local repulsion. We establish optimal
ranges of parameters for the validity of the above interpolation inequality,
discuss the existence and qualitative properties of the nonnegative maximizers,
and in some special cases estimate the optimal constant. For $p=2$ it is known
that the maximizers are H\"older continuous and compactly supported on a ball.
We show that for $p<2$ the maximizers are smooth functions supported on
$\mathbb{R}^N$, while for $p>2$ the maximizers consist of a characteristic
function of a ball and a nonconstant nonincreasing H\"older continuous function
supported on the same ball. We use these qualitative properties of the
maximizers to establish the validity of the Thomas-Fermi approximations for the
Choquard equations with local repulsion. The results are verified numerically
with extensive examples.
