Let $P=\frac 1 k \sum_i N(\theta_i, \sigma^2 I)$ and $Q=\frac 1 k \sum_i N(\mu_i, \sigma^2 I)$ be two mixtures of Gaussians in $\mathbb{R}^d$. A well known fact is that the Wasserstein-$2$ loss between $P$ and $Q$ can be written using the dual formulation as $$ W_2^2(P,Q)=\sup_{(f,g):f(x)+g(y) \leq \|x-y\|^2} \int f(x) P(x) dx+\int g(y) Q(y) dy\\ =\int f(x) P(x) dx+\int f^c(y) Q(y) dy, $$ where $f^c$ denotes the $c$-transform of the function $f$. It is also known that the optimal $f$ is related to the optimal Monge transportation map $T$ and the Kantorovich convex potential $u$ via $$ T(x) = x- \frac{\nabla f(x)}{2}, \quad T(x) = \nabla u(x) \Rightarrow \nabla^2 f(x) \preccurlyeq 2 I, $$ since $\nabla^2 u(x) \succcurlyeq 0$ from the convexity of $u$.

My question is regarding if any other gradient and Hessian information about the function $f$ is known in the literature. In particular I am looking for bounds of the type $\nabla^2 f(x) \succcurlyeq -\alpha I$ for some $\alpha \geq 0$. This is equivalent to a smoothness condition on the convex potential $u$, i.e $u$ is $\alpha$ smooth. Please point me to the relevant references if you are aware of any.