The HRUA algorithm is described in Stadlober's thesis "Sampling from
Poisson, binomial and hypergeometric distributions: ratio of uniforms
as a simple and fast alternative."  The pseudo-code for the algorithm
begins on page 82, where the parameter c is defined to be

    c = sqrt((N - n)*n*p*q/(N - 1) + 1/2)

`N` is the population size, `n` is the sample size, `p` = `M/N` is the
fraction of the first type (i.e. the fraction of "good" values or
successes) and `q` is `1 - p`.

The C code takes advantage of the symmetry and uses

    m = min(sample, popsize - sample);

as the sample size to be computed, where `popsize` is (not surprisingly)
the population size (i.e. `good + bad`), and `sample` is the requested
sample size.

The formula for `c` above is translated in the C code as

    d7 = sqrt((double)(popsize - m) * sample * d4 * d5 / (popsize - 1) + 0.5);

where `d4` and `d5` are `p` and `q`, resp.

The problem is that `sample` should have been replaced by `m` in that
formula. It should be

    d7 = sqrt((double)(popsize - m) * m * d4 * d5 / (popsize - 1) + 0.5);

I haven't thoroughly analyzed the implications of this mistake.  It will
only have an effect when `sample` is greater than half of `popsize`.
(When `sample` is less than or equal to half of `popsize`, `m == sample`.)
A few tests of the expected values of the distribution don't show
any obvious biases, but I haven't run it through a rigorous statistical
testing procedure.

It appears that it simply means that the bound for the rejection
method is too big, so the acceptance ratio of the method is much
lower than it should be.  This can be seen in a test such as:

    In [9]: %timeit np.random.hypergeometric(10, 190, 189, size=2500000)
    2.07 s ± 1.25 ms per loop (mean ± std. dev. of 7 runs, 1 loop each)

    In [10]: %timeit np.random.hypergeometric(10, 190, 11, size=2500000)
    1.14 s ± 1.42 ms per loop (mean ± std. dev. of 7 runs, 1 loop each)

Those two computations should take the same amount of time.
After the fixing the formula, they do:

    In [4]: %timeit np.random.hypergeometric(10, 190, 189, size=2500000)
    1.14 s ± 1.36 ms per loop (mean ± std. dev. of 7 runs, 1 loop each)

    In [5]: %timeit np.random.hypergeometric(10, 190, 11, size=2500000)
    1.14 s ± 694 µs per loop (mean ± std. dev. of 7 runs, 1 loop each)