Is this function concave? A question concernig a paper on the Kuramoto model In the paper "DYNAMICAL ASPECTS OF MEAN FIELD PLANE ROTATORS AND THE KURAMOTO MODEL" by L. Bertini, G. Giacomin, AND K. Pakdaman we read:

$$r:=\Psi(2Kr),\quad\text{with}\quad\Psi(x):=\dfrac{\int_{\mathbb{S}}\cos(\theta)\exp(x\cos(\theta))\,\mathrm{d}\theta}{\int_{\mathbb{S}}\exp(x\cos(\theta))\,\mathrm{d}\theta}.\tag{1.14}$$ In general, there is more than one solution to $(1.14)$: in fact, there can be at most two, more precisely there is only the trivial solution $r=0$ for $K\leq1$ and there is also a second solution $r>0$ if $K>1$. This is because $\Psi'(0)=1$ and because $\Psi(\cdot):[0,\infty)\to[0,1)$ is strictly concave $[16]$. In terms of stationary solutions, this means that for $K\leq 1$ only

Define 
$$
\Psi(x) = \frac{\int_{0}^{2\pi}\exp\{x\cos\theta\}\cos \theta\,d \theta}{\int_{0}^{2\pi}\exp\{x\cos\theta\}\, d \theta}.
$$
The first derivative of $\psi$ is 
$$\partial_{x}\Psi(x) = \frac{\int_{0}^{2\pi}\exp\{x\cos\theta\}\, d \theta\int_{0}^{2\pi}\exp\{x\cos\theta\}\cos^{2}\theta\, d \theta - \big(\int_{0}^{2\pi}\exp\{x\cos\theta\}\cos\theta\, d \theta\big)^{2}}{(\int_{0}^{2\pi}\exp\{x\cos\theta\}\, d \theta)^{2}}
$$
Rewrite this expression and use Holder  inequality (with p = q = 2) to get:
$$
\int_{0}^{2\pi}\exp\{x\cos\theta\}\cos\theta\, d \theta 
= \int_{0}^{2\pi}\exp\bigg\{\frac{1}{2}x\cos\theta\bigg\}\exp\bigg\{\frac{1}{2}x\cos\theta\bigg\}\cos\theta\, d \theta
\leq \Big(\int_{0}^{2\pi}\exp\{x\cos\theta\}\, d \theta\Big)^{1/2}\Big(\int_{0}^{2\pi}\exp\{x\cos\theta\}\cos^{2}\theta\, d \theta\Big)^{1/2}
$$
Taking the square on both sides yields
$$\Big(\int_{0}^{2\pi}\exp\{x\cos\theta\}\cos\theta\, d \theta\Big)^{2} \leq \int_{0}^{2\pi}\exp\{x\cos\theta\}\, d \theta\int_{0}^{2\pi}\exp\{x\cos\theta\}\cos^{2}\theta\, d \theta
$$
which proves that $\partial \Psi$ is non-negative.
But to prove concavity, we need to look at the sign of the second derivative. let's rewrite:
$$
\partial_{x}\Psi(x) = \text{I} - \text{II}
= \frac{\int f''(x, \theta)\,d \theta}{\int f(x, \theta)\,d \theta} - \frac{(\int f'(x, \theta)\,d \theta)^{2}}{(\int f(x, \theta)\,d \theta)^{2}}
$$where $f(x, \theta) = \exp\{x\cos\theta\}$ and the prime refers to a derivative taken with respect to $x$. The integrals are always from $0$ to $2\pi$. Taking the derivative of the first term gives
$$
\text{I}' = \frac{\int f(x, \theta)\,d \theta\int f'''(x, \theta)\,d \theta - \int f''(x, \theta)\,d \theta\int f'(x, \theta)\,d \theta}{(\int f(x, \theta)\,d \theta)^{2}}
$$
while for the second
$$\text{II}' = 2\frac{\int f'(x, \theta)\,d \theta}{\int f(x, \theta)\,d \theta}\times\frac{\int f(x, \theta)\,d \theta\int f''(x, \theta)\,d \theta - (\int f'(x, \theta)\,d \theta)^{2}}{(\int f(x, \theta)\,d \theta)^{2}}.
$$
Using a common denominator we can write the difference as 
$$
\partial_{x}^{2}\Psi(x) = \frac{1}{\big(\int f(x, \theta)\,d \theta\big)^{3}}\Big[\Big(\int f(x, \theta)\,d \theta\Big)^{2}\int f'''(x, \theta)\,d \theta \\
- 3\int f(x, \theta)\,d \theta\int f'(x, \theta)\,d \theta\int f''(x, \theta)\,d \theta + 2 \Big(\int f'(x, \theta)\,d \theta\Big)^{3} \Big]
$$
But how do I check that $\partial_x^2 \Psi(x) <0$?
 A: Here is an attempt:
First, let
$$2c =  \int e^{x \cos \theta} \,\mathrm d\theta,$$
then the desired expression becomes
$$\psi''(x) = \int f'''(x, \theta)\,\frac{\mathrm d \theta}{2c} 
- 3\int f'(x, \theta)\,\frac{\mathrm d \theta}{2c}\int f''(x, \theta)\,\frac{\mathrm d \theta}{2c} + 2 \left(\int f'(x, \theta)\,\frac{\mathrm d \theta}{2c}\right)^3.$$
Now writting the change of variables using 
$$\arccos_1: (-1,1) \to (0,\pi),\\
\arccos_2: (-1,1) \to (\pi,2\pi),$$
since
\begin{align*}
\int_0^{2\pi} (\cos(\theta) )^k  e^{x \cos \theta} \frac{\mathrm d \theta}{2c} &=\int_0^{\pi} (\cos(\theta) )^k  e^{x \cos \theta} \frac{\mathrm d \theta}{2c} +\int_\pi^{2\pi} (\cos(\theta) )^k  e^{x \cos \theta} \frac{\mathrm d \theta}{2c}\\
&=\int_{-1}^1 u^k e^{x u} \frac{\mathrm du}{c\sqrt{1 - u^2}},
\end{align*}
in particular
$$1 =\int_{-1}^1  e^{x u} \frac{\mathrm du}{c\sqrt{1 - u^2}},$$
Let 
$$\mathrm d\nu(u) = \frac{1}{c\sqrt{1-u^2}}\mathrm du,$$
we end up with
$$\psi''(x) = \int u^3e^{xu}\,\mathrm d \nu(u) 
- 3\int ue^{xu}\,\mathrm d \nu(u)\int u^2e^{xu}\,\mathrm d \nu(u) + 2 \left(\int ue^{xu}\,\mathrm d \nu(u)\right)^3.$$
Now let $m = m(x) = \int ue^{xu}\,\mathrm d \nu(u)$, we obtain now
$$\psi''(x) = \int u^3e^{xu}\,\mathrm d \nu(u) - 3m\int u^2e^{xu}\,\mathrm d \nu(u) + 2 m^3. \tag{*}$$
Now note that (*) equals to 
$$\int (u-m)^3e^{xu}\,\mathrm d \nu(u).$$
Indeed, since $\int 3m^2 u e^{xu}\,\mathrm d \nu(u) = 3m^3$,
\begin{align*}
\int (u-m)^3e^{xu}\,\mathrm d \nu(u) &= \int (u^3-3mu^2 + 3m^2 u - m^3)e^{xu}\,\mathrm d \nu(u) \\
&= \int u^3e^{xu}\,\mathrm d \nu(u) - 3m\int u^2e^{xu}\,\mathrm d \nu(u) + 2 m^3.
\end{align*}
Now the last step look up the paper from  P. A. Pearce, "Mean-field bounds on the magnetization for ferromagnetic spin models."

2. ONE-COMPONENT BOUNDED SPINS
We show that the mean-field magnetization is an upper bound whenever the single-spin measure $\nu$ lies in a certain class $\mathscr{P}$. To defined $\mathscr{P}$, let $m=m(k)$ depend on $\nu$ as prescribed by $(11)$. The class $\mathscr{P}$ is then the set of all even probability measures on $\mathbb{R}$ with compact support, such that $$\int d\nu(\sigma)e^{k\sigma}(m-\sigma)^p\geqslant0\tag{18}$$

where (11) is 

$$m=\left.\int \mathrm d\nu(\sigma)\sigma e^{k\sigma}\,\middle/\,\int \mathrm d\nu(\sigma)e^{k\sigma}\right..\tag{11}$$

And more important:

Lemma 4. Let $\nu$ be an even probability measure with support on $[-1,1]$, and suppose $\nu$ is absolutely continuous, i.e., $\mathrm d\nu(\sigma)=f(\sigma)\,\mathrm d\sigma$, with $f$ nondecreasing on $[0,1]$. Then $\nu\in\mathscr{P}$.

The concludsion is that 
$$\int (u-m)^3e^{xu}\,\mathrm d \nu(u) \leq 0.$$
To get stric inequality it suffices to note that the bounds used in the paper from Pearce are not equalities for our $\nu$.
