Is it possible to calculate the upper bound of neutron star mass under idealized conditions? The precise upper bound of neutron star mass is unknown because the equation of state (EOS) of nuclear matter is still poorly understood. However, there are still theoretical constraints on the stiffness of nuclear matter. For example, the bulk modulus Ks can’t exceed $c^2\rho$, otherwise the sound velocity will exceed the vacuum speed of light. Second, even for ultrarelativistic gas, the Ks is no higher than $\frac{c^2\rho}{3}$, which means the sound velocity is only $\sqrt{\frac{1}{3}}$ times the vacuum speed of light. Can we calculate the mass limit based on these idealized conditions? And how close are neutron stars to these idealized materials?
 A: Clearly the whole neutron star cannot be governed by equations of state that have these forms, because the pressure and density must decrease to zero as you get to the surface (by definition) and so the outer layers will be governed by more sensible equations of state that are reasonably well understood.
However, there is a general approach that can be taken, which is to assume an equation of state of the form
$$ P = P_o +(\rho - \rho_0)c^2$$
for densities greater than some density $\rho_0$ and to use a more realistic, and more certain, equation of state at lower densities.
This was done by Rhoades & Ruffini (1974) who used the so-called Harrison-Wheeler equation of state up to a density $\rho_0 = 4.6\times 10^{17}$ kg/m$^3$ and assumed a maximal equation of state where the speed of sound was $c$ above that. Their result was a maximum mass of $M_{\rm max}<3.2 M_{\odot}$.
An alternative treatment by Nauenberg & Chapline (1973) did something similar but used a different equation of state up to $\rho_0$ and allowed the speed of sound to be a free parameter above that. i.e.
$$ P = P_0 + (\rho - \rho_0)v^2\ .$$
Their result for $v=c$ and using the Bethe, Baym & Pethick equation of state at $\rho< 5\times 10^{17}$ kg/m$^3$ was $M_{\rm max}<3.6 M_{\odot}$. They don't give a plot to allow a straightforward estimate of the maximum mass for $v = c/\sqrt{3}$, but the results in a table suggest it is around $M_{\rm max}<2.2-2.5 M_{\odot}$.
Note that your upper limit equation of state for ultrarelativistic fermions assumes they are non-interacting ideal gases, which is a very soft equation of state. This cannot be true for real neutron stars and using such an equation of state for the entire neutron star would lead to $M_{\rm max}<0.75 M_{\odot}$, as demonstrated by Oppenheimer & Volkhoff (1939).
