I'll take a look, but I don't know if I'll be able to reproduce the problem. The call np.random.hypergeometric([1,2,2], [2], [2]) ultimately runs the legacy version, not the new version implemented in #13761. To see if the new code also causes a problem, try:

>>> rng = np.random.generator.Generator(np.random.PCG64(123))
>>> rng.hypergeometric([1, 2, 2], [2], [2])

The platform being ppc64(be) and those huge values of sample visible in the stack trace suggests there might be an endianess issue somewhere.

The new code runs no problem.

It does seem to be some kind of cython/multiiter endianness issue. Here's some simple cython code that seems to be failing (I don't write much cython so I haven't figured out the exact problem yet).

This is debug code I added in discrete_broadcast_iii in numpy/random/common.pyx:

    it = np.PyArray_MultiIterNew4(randoms, a_arr, b_arr, c_arr)
    with lock:
        for i in range(n):
            a_val = (<int64_t*>np.PyArray_MultiIter_DATA(it, 1))[0]
            b_val = (<int64_t*>np.PyArray_MultiIter_DATA(it, 2))[0]
            c_val = (<int64_t*>np.PyArray_MultiIter_DATA(it, 3))[0]
            print(c_val, c_arr)

The last line prints 8866477368 [2], so something seems to go wrong with PyArray_MultiIter_DATA here.

All right, I think I see the bug. The hypergeometric code assumes that npy_long is int64, but on ppc64(be) it is int32. I think that's because (confusingly) npy_long maps to python-int, which maps to c-int, which is int32 on ppc64, even though c-long is an int64. (edit: wrong, see comment below)

So the following code in hypergeometric is converting the inputs to int32, but later code assumes it is int64.

        # This cast to long is required to ensure that the values are inbounds
        ongood = <np.ndarray>np.PyArray_FROM_OTF(ngood, np.NPY_LONG, np.NPY_ALIGNED)
        onbad = <np.ndarray>np.PyArray_FROM_OTF(nbad, np.NPY_LONG, np.NPY_ALIGNED)
        onsample = <np.ndarray>np.PyArray_FROM_OTF(nsample, np.NPY_LONG, np.NPY_ALIGNED)

The fix might be as simple as replacing NPY_LONG by NPY_INT64, but the comment about being "inbounds" is throwing me off since I'm not sure what it means. Would that fix make sense @WarrenWeckesser ?

Hmm, npy_long is used all over the legacy random code. Is this perhaps a cython bug?

I think my description of the size of npy_long above wasn't quite right.

For reasons I don't really remember, my gcc on my ppc64 system is in 32-bit ELF mode for executables/libraries, which is why long is 4 bytes. I think if it was in 64-bit mode long should be 8 bytes. In any case we should support 32-bit executables on ppc64.

I don't think it's a cython bug since long is well documented to be 4 bytes on 32-bit systems.

npy_long maps to python-int, which maps to c-int, which is int32 on ppc64, even though c-long is an int64.

This sounds surprising to me. I was under the impression that python2's int mapped to C long, meaning that npy_long really is just long. C long is unfortunately unspecified as either int32 or int64, but I thought the former was only true in practice on windows.

C long is unfortunately unspecified as either int32 or int64, but I thought the former was only true in practice on windows.

I think it is also 4 bytes in 32-bit executables on linux 64-bit systems, and is probably controlled by compiler flags too.

I think the reason my ppc64(be) system is using 32-bit mode has to do with ubuntu support. I remember spending a while investigating ubuntu versions available on this system (a VM provided by OSUOSL, who did the ubuntu install for us), and ubuntu only provided 32-bit binaries for big-endian ppc when I asked for it in 2016.

I think ubuntu used to support "powerpc" (which was 32-bit big-endian) but dropped support 2 years ago. It currently supports ppc64el (little endian, 64-bit) used by new IBM computers.