How to calculate the limiting magnitude of Hubble? I've found various sources on the web stating that the limiting apparent magnitude of the Hubble Space Telescope is about 31 or 31.5. Wikipedia lists that here and it has a cite to a pdf about it.

What I want to know is how is that calculated?

I'm aware of two methods found at the bottom of a table here, but they seems to be general calculations for ground-based telescopes, because if I plug in Hubble's numbers, I get something around 21. The difference between an apparent magnitude of 21 and 31 is actually something like ten thousands times (in terms of brightness).

Also, one key factor missing in their formulas is exposure time. I'm pretty sure the limiting magnitude has something to do with exposure time. Hubble Ultra Deep Field, for example, had an exposure time of about 1 million seconds.

So to sum up, how would we calculate the limiting magnitude of a space-based telescope, in visible light?
The table you linked to gives limiting magnitudes for direct observations through a telescope with the human eye, so it's definitely not what you want to use.

The quoted number for HST is an empirical one, determined from the actual "Extreme Deep Field" data (total exposure time ~ 2 million seconds) after the fact; the Illingworth et al. PDF you linked to explains how it was done. Briefly, they decided the limit would be defined by a signal-to-noise ratio of 5, and then figured out what the noise level in the final combined images was by integrating the total flux in many circular apertures of diameter 0.35 arc seconds placed in regions of blank sky background, and then computing the standard deviation ($\sigma$) of these measurements. Five times $\sigma$ was then the limiting flux. Knowing the photometric calibration (the conversion between observed counts and apparent magnitude; see the paper for specific details), they converted this $S/N = 5$ limit into an apparent magnitude.

If, instead, you want try to estimate in advance what the limiting magnitude should be, then you have to consider several factors. Simplifying a bit, the equation for long-exposure, background-limited $S/N$ is
$$
S/N \approx \frac{s t}{ \sqrt{ s t + n_{\rm pix}(b t)} }
$$
where $s$ = detected flux from the source (star, galaxy, etc.) in electrons (i.e., detected photons) per second, $t$ is the integration time in seconds , $b$ is the flux from the background in electrons per pixel, and $n_{\rm pix}$ is the number of pixels you integrate over for the measurement. $s$ can be computed from the combination of the object's apparent magnitude, the diameter of the telescope's mirror and the overall efficiency of the telescope + detector system. $b$ can be determined in a similar fashion, using the estimated background brightness (which, for a telescope like HST, is mostly sunlight reflected from dust in the inner solar system). $n_{\rm pix}$ is a number chosen as a compromise between maximizing the amount of light from the source (higher $n_{\rm pix}$) and minimizing the noise from the background (lower $n_{\rm pix}$). This depends on the telescope's resolution (better resolution = smaller images of stars and other point-like sources = fewer pixels needed); the empirical calculations by Illingworth et al. used an aperture size amounting to about 100 pixels.