/ area(A) = 2 * umax * (vmax - vmin) / area(pdf)``,
where `area(pdf)` is the integral of `pdf` (which is equal to one if the
probability density function is used but can take on other values if a
function proportional to the density is used). The equality holds since
the area of ``A`` is equal to ``0.5 * area(pdf)`` (Theorem 7.1 in [1]_).
If the sampling fails to generate a single random variate after 50000
iterations (i.e. not a single draw is in ``A``), an exception is raised.

If the bounding rectangle is not correctly specified (i.e. if it does not
contain ``A``), the algorithm samples from a distribution different from
the one given by `pdf`. It is therefore recommended to perform a
test such as `~scipy.stats.kstest` as a check.

References
----------
.. [1] L. Devroye, "Non-Uniform Random Variate Generation",
   Springer-Verlag, 1986.

.. [2] W. Hoermann and J. Leydold, "Generating generalized inverse Gaussian
   random variates", Statistics and Computing, 24(4), p. 547--557, 2014.

.. [3] A.J. Kinderman and J.F. Monahan, "Computer Generation of Random
   Variables Using the Ratio of Uniform Deviates",
   ACM Transactions on Mathematical Software, 3(3), p. 257--260, 1977.

Examples
--------
>>> import numpy as np
>>> from scipy import stats

>>> from scipy.stats.sampling import RatioUniforms
>>> rng = np.random.default_rng()

Simulate normally distributed random variables. It is easy to compute the
bounding rectangle explicitly in that case. For simplicity, we drop the
normalization factor of the density.

>>> f = lambda x: np.exp(-x**2 / 2)
>>> v = np.sqrt(f(np.sqrt(2))) * np.sqrt(2)
>>> umax = np.sqrt(f(0))
>>> gen = RatioUniforms(f, umax=umax, vmin=-v, vmax=v, random_state=rng)
>>> r = gen.rvs(size=2500)

The K-S test confirms that the random variates are indeed normally
distributed (normality is not rejected at 5% significance level):

>>> stats.kstest(r, 'norm')[1]
0.250634764150542

The exponential distribution provides another example where the bounding
rectangle can be determined explicitly.

>>> gen = RatioUniforms(lambda x: np.exp(-x), umax=1, vmin=0,
...                     vmax=2*np.exp(-1), random_state=rng)
>>> r = gen.rvs(1000)
>>> stats.kstest(r, 'expon')[1]
0.21121052054580314

r