Skip to content

Commit

Permalink
BUG: Sync upstream big fixes
Browse files Browse the repository at this point in the history 
Sync with 1.19 bug fixes
  • Loading branch information
@bashtage
bashtage authored and Kevin Sheppard committed Jun 3, 2020
1 parent 437bdf5 commit 4c1239d
Show file tree
Hide file tree
Showing 5 changed files with 205 additions and 39 deletions.
31 changes: 17 additions & 14 deletions randomgen/generator.pyx
View file
Original file line number Diff line number Diff line change
Expand Up @@ -3341,7 +3341,7 @@ cdef class Generator:
"""
return disc(&random_negative_binomial, &self._bitgen, size, self.lock, 2, 0,
n, "n", CONS_POSITIVE_NOT_NAN,
p, "p", CONS_BOUNDED_0_1,
p, "p", CONS_BOUNDED_GT_0_1,
0.0, "", CONS_NONE)

def poisson(self, lam=1.0, size=None):
Expand Down Expand Up @@ -3918,8 +3918,7 @@ cdef class Generator:
else:
raise ValueError("covariance is not positive-semidefinite.")

x = np.dot(x, np.sqrt(s)[:, None] * v)
x += mean
x = mean + x @ (u * np.sqrt(s)).T
x.shape = tuple(final_shape)
return x

Expand Down Expand Up @@ -4019,6 +4018,9 @@ cdef class Generator:
cdef int64_t ni
cdef np.broadcast it

if np.ndim(pvals) != 1:
raise ValueError("pvals must be 1d array")

d = len(pvals)
on = <np.ndarray>np.PyArray_FROM_OTF(n, np.NPY_INT64, api.NPY_ARRAY_ALIGNED)
parr = <np.ndarray>np.PyArray_FROM_OTF(pvals, np.NPY_DOUBLE, api.NPY_ARRAY_ALIGNED | api.NPY_ARRAY_C_CONTIGUOUS)
Expand Down Expand Up @@ -4085,23 +4087,23 @@ cdef class Generator:
Parameters
----------
alpha : array
Parameter of the distribution (k dimension for sample of
dimension k).
alpha : sequence of floats, length k
Parameter of the distribution (length ``k`` for sample of
length ``k``).
size : int or tuple of ints, optional
Output shape. If the given shape is, e.g., ``(m, n, k)``, then
Output shape. If the given shape is, e.g., ``(m, n)``, then
``m * n * k`` samples are drawn. Default is None, in which case a
single value is returned.
vector of length ``k`` is returned.
Returns
-------
samples : ndarray,
The drawn samples, of shape (size, alpha.ndim).
The drawn samples, of shape ``(size, k)``.
Raises
-------
ValueError
If any value in alpha is less than or equal to zero
If any value in ``alpha`` is less than or equal to zero.
Notes
-----
Expand Down Expand Up @@ -4175,7 +4177,7 @@ cdef class Generator:
cdef double acc, invacc

k = len(alpha)
alpha_arr = <np.ndarray>np.PyArray_FROM_OTF(alpha, np.NPY_DOUBLE, api.NPY_ARRAY_ALIGNED | api.NPY_ARRAY_C_CONTIGUOUS)
alpha_arr = <np.ndarray>np.PyArray_FROMANY(alpha, np.NPY_DOUBLE, 1, 1, np.NPY_ARRAY_ALIGNED | np.NPY_ARRAY_C_CONTIGUOUS)
if np.any(np.less_equal(alpha_arr, 0)):
raise ValueError("alpha <= 0")
alpha_data = <double*>np.PyArray_DATA(alpha_arr)
Expand Down Expand Up @@ -4254,15 +4256,15 @@ cdef class Generator:
# Fast, statically typed path: shuffle the underlying buffer.
# Only for non-empty, 1d objects of class ndarray (subclasses such
# as MaskedArrays may not support this approach).
x_ptr = <char*><size_t>x.ctypes.data
x_ptr = <char*><size_t>np.PyArray_DATA(x)
stride = x.strides[0]
itemsize = x.dtype.itemsize
# As the array x could contain python objects we use a buffer
# of bytes for the swaps to avoid leaving one of the objects
# within the buffer and erroneously decrementing it's refcount
# when the function exits.
buf = np.empty(itemsize, dtype=np.int8) # GC'd at function exit
buf_ptr = <char*><size_t>buf.ctypes.data
buf_ptr = <char*><size_t>np.PyArray_DATA(buf)
with self.lock:
# We trick gcc into providing a specialized implementation for
# the most common case, yielding a ~33% performance improvement.
Expand Down Expand Up @@ -4383,7 +4385,8 @@ cdef class Generator:
return arr

arr = np.asarray(x)

if arr.ndim < 1:
raise IndexError("x must be an integer or at least 1-dimensional")
# shuffle has fast-path for 1-d
if arr.ndim == 1:
# Return a copy if same memory
Expand Down
31 changes: 18 additions & 13 deletions randomgen/mtrand.pyx
View file
Original file line number Diff line number Diff line change
Expand Up @@ -3095,14 +3095,14 @@ cdef class RandomState:
Samples are drawn from a negative binomial distribution with specified
parameters, `n` successes and `p` probability of success where `n`
is > 0 and `p` is in the interval [0, 1].
is > 0 and `p` is in the interval (0, 1].
Parameters
----------
n : float or array_like of floats
Parameter of the distribution, > 0.
p : float or array_like of floats
Parameter of the distribution, >= 0 and <=1.
Parameter of the distribution. Must satisfy 0 < p <= 1.
size : int or tuple of ints, optional
Output shape. If the given shape is, e.g., ``(m, n, k)``, then
``m * n * k`` samples are drawn. If size is ``None`` (default),
Expand Down Expand Up @@ -3160,7 +3160,7 @@ cdef class RandomState:
"""
out = disc(&legacy_negative_binomial, &self._aug_state, size, self.lock, 2, 0,
n, "n", CONS_POSITIVE,
p, "p", CONS_BOUNDED_0_1,
p, "p", CONS_BOUNDED_GT_0_1,
0.0, "", CONS_NONE)
# Match historical output type
return int64_to_long(out)
Expand Down Expand Up @@ -3830,6 +3830,9 @@ cdef class RandomState:
cdef long *mnix
cdef long ni

if np.ndim(pvals) != 1:
raise ValueError("pvals must be 1d array")

d = len(pvals)
parr = <np.ndarray>np.PyArray_FROM_OTF(pvals, np.NPY_DOUBLE, api.NPY_ARRAY_ALIGNED | api.NPY_ARRAY_C_CONTIGUOUS)
check_array_constraint(parr, "pvals", CONS_BOUNDED_0_1)
Expand Down Expand Up @@ -3873,23 +3876,23 @@ cdef class RandomState:
Parameters
----------
alpha : array
Parameter of the distribution (k dimension for sample of
dimension k).
alpha : sequence of floats, length k
Parameter of the distribution (length ``k`` for sample of
length ``k``).
size : int or tuple of ints, optional
Output shape. If the given shape is, e.g., ``(m, n, k)``, then
Output shape. If the given shape is, e.g., ``(m, n)``, then
``m * n * k`` samples are drawn. Default is None, in which case a
single value is returned.
vector of length ``k`` is returned.
Returns
-------
samples : ndarray,
The drawn samples, of shape (size, alpha.ndim).
The drawn samples, of shape ``(size, k)``.
Raises
-------
ValueError
If any value in alpha is less than or equal to zero
If any value in ``alpha`` is less than or equal to zero.
Notes
-----
Expand Down Expand Up @@ -3963,7 +3966,7 @@ cdef class RandomState:
cdef double acc, invacc

k = len(alpha)
alpha_arr = <np.ndarray>np.PyArray_FROM_OTF(alpha, np.NPY_DOUBLE, api.NPY_ARRAY_ALIGNED | api.NPY_ARRAY_C_CONTIGUOUS)
alpha_arr = <np.ndarray>np.PyArray_FROMANY(alpha, np.NPY_DOUBLE, 1, 1, np.NPY_ARRAY_ALIGNED | np.NPY_ARRAY_C_CONTIGUOUS)
if np.any(np.less_equal(alpha_arr, 0)):
raise ValueError("alpha <= 0")
alpha_data = <double*>np.PyArray_DATA(alpha_arr)
Expand Down Expand Up @@ -4042,15 +4045,15 @@ cdef class RandomState:
# Fast, statically typed path: shuffle the underlying buffer.
# Only for non-empty, 1d objects of class ndarray (subclasses such
# as MaskedArrays may not support this approach).
x_ptr = <char*><size_t>x.ctypes.data
x_ptr = <char*><size_t>np.PyArray_DATA(x)
stride = x.strides[0]
itemsize = x.dtype.itemsize
# As the array x could contain python objects we use a buffer
# of bytes for the swaps to avoid leaving one of the objects
# within the buffer and erroneously decrementing it's refcount
# when the function exits.
buf = np.empty(itemsize, dtype=np.int8) # GC'd at function exit
buf_ptr = <char*><size_t>buf.ctypes.data
buf_ptr = <char*><size_t>np.PyArray_DATA(buf)
with self.lock:
# We trick gcc into providing a specialized implementation for
# the most common case, yielding a ~33% performance improvement.
Expand Down Expand Up @@ -4127,6 +4130,8 @@ cdef class RandomState:
return arr

arr = np.asarray(x)
if arr.ndim < 1:
raise IndexError("x must be an integer or at least 1-dimensional")

# shuffle has fast-path for 1-d
if arr.ndim == 1:
Expand Down
30 changes: 24 additions & 6 deletions randomgen/src/distributions/distributions.c
View file
Original file line number Diff line number Diff line change
Expand Up @@ -1138,9 +1138,7 @@ static NPY_INLINE uint64_t bounded_lemire_uint64(bitgen_t *bitgen_state,
if (leftover < rng_excl) {
/* `rng_excl` is a simple upper bound for `threshold`. */

const uint64_t threshold = -rng_excl % rng_excl;
/* Same as: threshold=((uint64_t)(0x10000000000000000ULLL - rng_excl)) %
* rng_excl; */
const uint64_t threshold = (UINT64_MAX - rng) % rng_excl;

while (leftover < threshold) {
m = ((__uint128_t)next_uint64(bitgen_state)) * rng_excl;
Expand All @@ -1164,7 +1162,7 @@ static NPY_INLINE uint64_t bounded_lemire_uint64(bitgen_t *bitgen_state,
if (leftover < rng_excl) {
/* `rng_excl` is a simple upper bound for `threshold`. */

const uint64_t threshold = -rng_excl % rng_excl;
const uint64_t threshold = (UINT64_MAX - rng) % rng_excl;
/* Same as:threshold=((uint64_t)(0x10000000000000000ULLL - rng_excl)) %
* rng_excl; */

Expand Down Expand Up @@ -1320,6 +1318,14 @@ uint64_t random_bounded_uint64(bitgen_t *bitgen_state, uint64_t off,
return off;
} else if (rng <= 0xFFFFFFFFUL) {
/* Call 32-bit generator if range in 32-bit. */
if (rng == 0xFFFFFFFFUL) {
/*
* The 32-bit Lemire method does not handle rng=0xFFFFFFFF, so we'll
* call next_uint32 directly. This also works when use_masked is True,
* so we handle both cases here.
*/
return off + (uint64_t) next_uint32(bitgen_state);
}
if (use_masked) {
return off + buffered_bounded_masked_uint32(bitgen_state, (uint32_t)rng,
(uint32_t)mask, NULL, NULL);
Expand Down Expand Up @@ -1433,10 +1439,21 @@ void random_bounded_uint64_fill(bitgen_t *bitgen_state, uint64_t off,
out[i] = off;
}
} else if (rng <= 0xFFFFFFFFUL) {
uint32_t buf = 0;
/*
* The 32-bit Lemire method does not handle rng=0xFFFFFFFF, so we'll
* call next_uint32 directly. This also works when use_masked is True,
* so we handle both cases here.
*/
if (rng == 0xFFFFFFFFUL) {
for (i = 0; i < cnt; i++) {
out[i] = off + (uint64_t) next_uint32(bitgen_state);
}
}
/* Call 32-bit generator if range in 32-bit. */
else {
uint32_t buf = 0;
int bcnt = 0;

/* Call 32-bit generator if range in 32-bit. */
if (use_masked) {
/* Smallest bit mask >= max */
uint64_t mask = gen_mask(rng);
Expand All @@ -1451,6 +1468,7 @@ void random_bounded_uint64_fill(bitgen_t *bitgen_state, uint64_t off,
buffered_bounded_lemire_uint32(bitgen_state, (uint32_t)rng, &bcnt, &buf);
}
}
}
} else if (rng == 0xFFFFFFFFFFFFFFFFULL) {
/* Lemire64 doesn't support rng = 0xFFFFFFFFFFFFFFFF. */
for (i = 0; i < cnt; i++) {
Expand Down
13 changes: 7 additions & 6 deletions randomgen/tests/test_generator_117.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -25,9 +25,7 @@
NP_LT_1174 = LooseVersion(np.__version__) < v1174
NP_LT_1174_OR_GT_118 = NP_LT_1174 or LooseVersion(np.__version__) > v118

pytestmark = pytest.mark.skipif(
NP_LT_1174_OR_GT_118, reason="Only test 1.17.4 to 1.18.x"
)
pytestmark = pytest.mark.skipif(NP_LT_1174, reason="Only test 1.17.4 to 1.18.x")


def positive_param():
Expand Down Expand Up @@ -88,15 +86,18 @@ def above_1():

def triangular():
out = product(*[positive_param() for _ in range(3)])
out = [(l, l + m, l + m + r) for l, m, r in out]
out = [(lft, lft + mid, lft + mid + rgt) for lft, mid, rgt in out]
return out


def uniform():
lo = positive_param()
low = positive_param()
high = positive_param()
scale = positive_param()
return [(l / l + h + s, (l + h) / (l + h + s)) for l, h, s in zip(lo, high, scale)]
return [
(lo / lo + hi + sc, (lo + hi) / (lo + hi + sc))
for lo, hi, sc in zip(low, high, scale)
]


def integers():
Expand Down

0 comments on commit 4c1239d

Please sign in to comment.