+1 to indicate spin up and -1 to indicate spin down). Each energy level has a degeneracy given by $D=2\mathcal{A}B/\phi$, where $\mathcal{A}$ is the sheet area of graphene, $\phi=h/e$ is the magnetic unit flux and the factor of 2 takes into account the valley degeneracy [@Goerbig2011].

For a Fermi energy $\mu>0$ at zero temperature, the valence band is full while the conduction band is partially filled. We will write the conduction energy levels as $\varepsilon_{m}=\alpha\sqrt{n_{m}B}-s_{m}\mu_{\mathrm{B}}B,$ where we have introduced the decreasing energy sorting index $m=0,\:1,\,2,\ldots$, so $n_{m}$ gives the LL and $s_{m}$ the spin for the $m$ position. In general, the mixing of the LL depends on the magnetic field [@Chang1997], but for usual values there is no spin mixing so $n_{m}=\sqrt{m/2-\left[1-\left(-1\right)^{m}\right]/4}$ and $s_{m}=\left(-1\right)^{m}$. At a given $\mu>0$, all conduction energy levels $m=0,\:1,\,2,\ldots,\,f$ are filled, where $f$ is such that $\varepsilon_{f}<\mu\leq\varepsilon_{f+1}$. Then the grand potential at zero temperature is

$$\Omega=\Omega_{V}+\sum_{m=0}^{f}D\left(\varepsilon_{m}-\mu\right),\label{GP}$$

where $\Omega_{V}$ is the grand potential due to the filled valence band. From equation (\[GP\]) we shall expect the oscillatory contribution coming only from the last term, associated with the change in the last energy level as $B$ is modified. On the other hand, the first term $\Omega_{V}$ makes a 