Tis, a good question. Two related questions arise from it. The first one is, will the hydrogen all be used up in finite time? The second related one is, will star formation completely stop in finite time? They sound related, but the first result doesn't necessarily imply the same result for the second, or vice versa. I.E. a low but nonzero gas density might possibly not allow further star formation, and maybe we could have no hydrogen, but have other types of gas (or even solid objects) still collect into stellar mass objects.
I don't know for sure the answers. The rate of star formation (and hydrogen consumption) could decline slowly enough as to never formally reach zero. Or not.
We do know that a lot of gas gets blown out of galaxies by massive stars, supernova, and black hole activity, and becomes intergalactic gas -usually staying within the galaxy cluster. On a long time scale this should eventually fall back into the cluster's galaxies. So I would think the star formation rate would have a very long tail.
The WKB approximation states that in one dimension, the tunneling probability $P$ can be approximated as
$\ln P=-\frac{2\sqrt{2m}}{\hbar}\int_a^b \sqrt{V-E} dx$ ,
where the limits of integration $a$ and $b$ are the classical turning points, $m$ is the reduced mass, the electrical potential $V$ is a function of $x$, and $E$ is the total energy. Setting $V=kq_1q_2/x$, we have for the integral
$I=\int_a^b \sqrt{V-E} dx$
$=\frac{kq_1k_2}{\sqrt{E}}\int_{A}^1\sqrt{u^{-1}-1} du$ ,
where $A=a/b$. The indefinite integral equals $-u\sqrt{u^{-1}-1}+\tan^{-1}\sqrt{u^{-1}-1}$, and for $A\ll 1$ the definite integral is then $\pi/2$. The result is
$\ln P=-\frac{\pi kq_1q_2}{\hbar}\sqrt{\frac{2m}{E}} $ .
This result was obtained in Gamow 1938, and $G=-(1/2)\ln P$ is referred to as the Gamow factor or Gamow-Sommerfeld factor.
The fact that the integral $\int_A^1\ldots$ can be approximated as $\int_0^1\ldots$ tells us that the right-hand tail of the barrier dominates, i.e., it is hard for the nuclei to travel through the very long stretch of $\sim 1$ nm over which the motion is only mildly classically forbidden, but if they can do that, it's relatively easy for them to penetrate the highly classically forbidden region at $x\sim1$ fm. Surprisingly, the result can be written in a form that depends only on $m$ and $E$, but not on $a$, i.e., we don't even have to know the range of the strong nuclear force in order to calculate the result.
The generic WKB expression depends on $E$ through an expression of the form $V-E$, which might have led us to believe that with a 1 MeV barrier, it would make little difference whether $E$ was 1 eV or 1 keV, and fusion would be just as likely in trees and houses as in the sun. But because the tunneling probability is dominated by the tail of the barrier, not its peak, the final result ends up depending on $1/\sqrt{E}$.
Because $P$ increases extremely rapidly as a function of $E$, fusion is dominated by nuclei whose energies lie in the tails of the Maxwellian distribution. There is a narrow range of energies, known as the Gamow window, in which the product of $P$ and the Maxwell distribution is large enough to contibute significantly to the rate of fusion.
Gamow and Teller, Phys. Rev. 53 (1938) 608
Best Answer
We can only give an answer on the basis of what we currently know about cosmological parameters. If indeed these have been correctly estimated, and that the cosmological constant is constant, then the universe will continue to expand at an accelerating rate.
Given that about half the baryons in the universe currently exist outside galaxies in the "warm/hot intergalactic medium", then it seems to me quite likely that these will never form part of stars and that the increasingly rarefied universe will cease to form (many) stars in the future.
In that scenario, isolated protons would always be the most common nuclei.
One complication is that the answer may be different for that part of the universe that we can see.
Gravitationally bound systems like the Milky Way, the local group and cluster will likely not take part in the accelerating expansion. Thus the baryons within this region could continue to be processed in stars, even if they do not constitute the majority of baryons in the universe. How to estimate the timescale for this processing? A lower limit is probably the free-fall timescale of the local supercluster, which has a mass of about $10^{15}$ M$_{\odot}$ in a diameter of 33 Mpc. The free-fall timescale is $\sim (G \rho)^{-1/2} \simeq 100$ billion years. However, this time must be increased to allow that most stars formed are low-mass red dwarfs with timescales to turn hydrogen into helium of a further $10^{12}$ years. So within our local group of galaxies, I would estimate $10^{11}$ years to incorporate all the gas in our local cluster into stars and a further $10^{12}$ years for them to fuse most of the hydrogen.