\frac{e^{n \mu \beta}+e^{- n \mu \beta }}{n^2} K_2 (y z_0)
- \frac{\alpha}{2 \pi^2 \beta}  \sum_{n=1}^{\infty} \frac{e^{n \mu \beta}+e^{- n \mu \beta }}{n} \int\limits_{z_0}^{\infty} dz
\frac{z^2}{\sqrt{z^2+\alpha^2}} K_1 (y z),\end{aligned}$$

where $K_n(x)$ is the McDonald function of order $n$, $ y = n \beta$, $z_0= m \sqrt{1-s b}$ and $\alpha=s m b/2$.

To compute the integral in the second term of Eq. (\[Grand-Potential-sst2\])

$$\begin{aligned}
I = \int\limits_{z_0}^{\infty} dz
\frac{z^2}{\sqrt{z^2+\alpha^2}} K_1 (y z),\end{aligned}$$

we follow the procedure described in Appendix A. Finally $\Omega_{st} (s)$ reads

$$\begin{aligned}
\label{Grand-Potential-sst3}
\Omega_{st}(s)= - \frac{z_0^2}{2 \pi^2 \beta^2} (1+\frac{\alpha}{\sqrt{z_0^2 + \alpha^2}}) \sum_{n=1}^{\infty} \frac{e^{n \mu \beta}+e^{- n \mu \beta }}{n^2} K_2 (y z_0)
- \frac{\alpha z_0^2}{ \pi^2 \beta^2 \sqrt{z_{0}^2 + \alpha^2}} \\ \times \sum_{n=1}^{\infty} \frac{e^{n \mu \beta}+e^{- n \mu \beta }}{n^2} \sum_{w=1}^{\infty} \frac{(-1)^w(2 w -1)!!}{(z_{0}^2 + \alpha^2)^w} \left(\frac{z_0}{y}\right)^w K_{-(w+2)}(y z_0).\nonumber\end{aligned}$$

In the low temperature limit $T \ll m$ (which for $m=2 m_e\cong1 MeV$ means $T \ll 10^{10} K$), $\Omega_{st}(-1)$ and $\Omega_{st}(0)$ vanish. Therefore, in this limit $\Omega_{st} \simeq \Omega_{st}(1)$. This means that all the particles are in the spin ground state $s=1$. The leading term of $\Omega_{st}(1)$ is the first one in Eq.(\[Grand-Potential-sst3\]). Since it admits a further simplification, for the assumed low temperatures the statistical part of the thermodynamical 