diff --git a/.github/workflows/build-with-clang.yml b/.github/workflows/build-with-clang.yml index be95089..99b5bc3 100644 --- a/.github/workflows/build-with-clang.yml +++ b/.github/workflows/build-with-clang.yml @@ -5,8 +5,8 @@ on: branches: [master] jobs: - generate-coverage: - name: Generate coverage and push to Coveralls.io + build-with-clang: + name: Build project with IntelLLVM clang compiler runs-on: ubuntu-latest env: @@ -14,7 +14,7 @@ jobs: steps: - name: Cancel Previous Runs - uses: styfle/cancel-workflow-action@0.11.0 + uses: styfle/cancel-workflow-action@0.12.1 with: access_token: ${{ github.token }} @@ -33,13 +33,13 @@ jobs: sudo apt-get install intel-oneapi-mkl-devel - name: Setup Python - uses: actions/setup-python@v4 + uses: actions/setup-python@v5 with: - python-version: '3.11' + python-version: '3.12' architecture: x64 - name: Checkout repo - uses: actions/checkout@v3 + uses: actions/checkout@v4 with: fetch-depth: 0 @@ -57,10 +57,9 @@ jobs: run: | source /opt/intel/oneapi/setvars.sh echo $CMPLR_ROOT - export CC=$CMPLR_ROOT/../latest/bin-llvm/clang - export CXX=$CMPLR_ROOT/../latest/bin-llvm/clang++ - echo "CC = ${CC} CXX=${CXX}" - ls -l ${CC} ${CXX} + export CC=$CMPLR_ROOT/bin/icx + export CXX=$CMPLR_ROOT/bin/icpx + export CFLAGS="${CFLAGS} -fno-fast-math -O2" python setup.py develop - name: Run mkl_random tests diff --git a/.github/workflows/conda-package.yml b/.github/workflows/conda-package.yml index acbb8d7..444abc9 100644 --- a/.github/workflows/conda-package.yml +++ b/.github/workflows/conda-package.yml @@ -6,16 +6,16 @@ env: PACKAGE_NAME: mkl_random MODULE_NAME: mkl_random VER_SCRIPT1: "import json; f = open('ver.json', 'r'); j = json.load(f); f.close(); " - VER_SCRIPT2: "d = j['dpctl'][0]; print('='.join((d[s] for s in ('version', 'build'))))" + VER_SCRIPT2: "d = j['mkl_random'][0]; print('='.join((d[s] for s in ('version', 'build'))))" jobs: build_linux: runs-on: ubuntu-latest strategy: matrix: - python: [3.9] + python: ["3.9", "3.10", "3.11", "3.12"] steps: - - uses: actions/checkout@v2 + - uses: actions/checkout@v4 with: fetch-depth: 0 @@ -23,7 +23,7 @@ jobs: run: | echo "pkgs_dirs: [~/.conda/pkgs]" >> ~/.condarc - name: Cache conda packages - uses: actions/cache@v2 + uses: actions/cache@v4 env: CACHE_NUMBER: 0 # Increase to reset cache with: @@ -40,7 +40,7 @@ jobs: run: conda install conda-build - name: Build conda package run: | - CHANNELS="-c intel -c defaults --override-channels" + CHANNELS="-c conda-forge -c intel --override-channels" VERSIONS="--python ${{ matrix.python }}" TEST="--no-test" @@ -50,7 +50,7 @@ jobs: $CHANNELS \ conda-recipe - name: Upload artifact - uses: actions/upload-artifact@v2 + uses: actions/upload-artifact@v4 with: name: ${{ env.PACKAGE_NAME }} ${{ runner.os }} Python ${{ matrix.python }} path: /usr/share/miniconda/conda-bld/linux-64/${{ env.PACKAGE_NAME }}-*.tar.bz2 @@ -60,20 +60,20 @@ jobs: strategy: matrix: - python: ['3.9'] + python: ['3.9', '3.10', '3.11', '3.12'] env: conda-bld: C:\Miniconda\conda-bld\win-64\ steps: - - uses: actions/checkout@v2 + - uses: actions/checkout@v4 with: fetch-depth: 0 - - uses: conda-incubator/setup-miniconda@v2 + - uses: conda-incubator/setup-miniconda@v3 with: auto-activate-base: true activate-environment: "" - name: Cache conda packages - uses: actions/cache@v3 + uses: actions/cache@v4 env: CACHE_NUMBER: 3 # Increase to reset cache with: @@ -86,9 +86,9 @@ jobs: - name: Install conda-build run: conda install conda-build - name: Build conda package - run: conda build --no-test --python ${{ matrix.python }} -c intel -c defaults --override-channels conda-recipe + run: conda build --no-test --python ${{ matrix.python }} -c intel -c conda-forge --override-channels conda-recipe - name: Upload artifact - uses: actions/upload-artifact@v2 + uses: actions/upload-artifact@v4 with: name: ${{ env.PACKAGE_NAME }} ${{ runner.os }} Python ${{ matrix.python }} path: ${{ env.conda-bld }}${{ env.PACKAGE_NAME }}-*.tar.bz2 @@ -99,16 +99,16 @@ jobs: strategy: matrix: - python: [3.9] + python: ['3.9', '3.10', '3.11', '3.12'] experimental: [false] runner: [ubuntu-latest] continue-on-error: ${{ matrix.experimental }} env: - CHANNELS: -c intel -c defaults --override-channels + CHANNELS: -c conda-forge -c intel --override-channels steps: - name: Download artifact - uses: actions/download-artifact@v2 + uses: actions/download-artifact@v4 with: name: ${{ env.PACKAGE_NAME }} ${{ runner.os }} Python ${{ matrix.python }} - name: Add conda to system path @@ -135,7 +135,7 @@ jobs: run: | echo "pkgs_dirs: [~/.conda/pkgs]" >> ~/.condarc - name: Cache conda packages - uses: actions/cache@v2 + uses: actions/cache@v4 env: CACHE_NUMBER: 0 # Increase to reset cache with: @@ -151,14 +151,14 @@ jobs: . $CONDA/etc/profile.d/conda.sh CHANNELS="-c $GITHUB_WORKSPACE/channel ${{ env.CHANNELS }}" export PACKAGE_VERSION=$(python -c "${VER_SCRIPT1} ${VER_SCRIPT2}") - conda create -n test_mkl_random $PACKAGE_NAME=${PACKAGE_VERSION} nose python=${{ matrix.python }} $CHANNELS + conda create -n test_mkl_random $PACKAGE_NAME=${PACKAGE_VERSION} pytest python=${{ matrix.python }} $CHANNELS # Test installed packages conda list - name: Run tests run: | . $CONDA/etc/profile.d/conda.sh conda activate test_mkl_random - nosetests -v mkl_random + pytest -vv --pyargs mkl_random test_windows: needs: build_windows @@ -166,19 +166,19 @@ jobs: strategy: matrix: - python: ['3.9'] + python: ['3.9', '3.10', '3.11', '3.12'] experimental: [false] runner: [windows-latest] continue-on-error: ${{ matrix.experimental }} env: - CHANNELS: -c intel -c defaults --override-channels + CHANNELS: -c conda-forge -c intel --override-channels steps: - name: Download artifact - uses: actions/download-artifact@v2 + uses: actions/download-artifact@v4 with: name: ${{ env.PACKAGE_NAME }} ${{ runner.os }} Python ${{ matrix.python }} - - uses: conda-incubator/setup-miniconda@v2 + - uses: conda-incubator/setup-miniconda@v3 with: auto-activate-base: true activate-environment: "" @@ -205,7 +205,7 @@ jobs: conda create -n test_mkl_random ${{ env.PACKAGE_NAME }}=%PACKAGE_VERSION% python=${{ matrix.python }} -c ${{ env.GITHUB_WORKSPACE }}/channel ${{ env.CHANNELS }} --only-deps --dry-run > lockfile more lockfile - name: Cache conda packages - uses: actions/cache@v3 + uses: actions/cache@v4 env: CACHE_NUMBER: 3 # Increase to reset cache with: @@ -224,10 +224,10 @@ jobs: FOR /F "tokens=* USEBACKQ" %%F IN (`python -c "%SCRIPT%"`) DO ( SET PACKAGE_VERSION=%%F ) - conda create -n test_mkl_random ${{ env.PACKAGE_NAME }}=%PACKAGE_VERSION% nose python=${{ matrix.python }} -c ${{ env.GITHUB_WORKSPACE }}/channel ${{ env.CHANNELS }} + conda create -n test_mkl_random ${{ env.PACKAGE_NAME }}=%PACKAGE_VERSION% pytest python=${{ matrix.python }} -c ${{ env.GITHUB_WORKSPACE }}/channel ${{ env.CHANNELS }} # Test installed packages conda list - name: Run tests run: | conda activate -n test_mkl_random - nosetests -v ${{ env.MODULE_NAME }} + pytest -v --pyargs ${{ env.MODULE_NAME }} diff --git a/conda-recipe/bld.bat b/conda-recipe/bld.bat index c639b91..1b04207 100644 --- a/conda-recipe/bld.bat +++ b/conda-recipe/bld.bat @@ -1,5 +1,5 @@ @rem Remember to source the compiler set MKLROOT=%CONDA_PREFIX% -%PYTHON% setup.py install --old-and-unmanageable +%PYTHON% setup.py install if errorlevel 1 exit 1 diff --git a/conda-recipe/build.sh b/conda-recipe/build.sh index 0429491..34241cc 100644 --- a/conda-recipe/build.sh +++ b/conda-recipe/build.sh @@ -1,4 +1,5 @@ #!/bin/bash -x export CFLAGS="-I$PREFIX/include $CFLAGS" -MKLROOT=$CONDA_PREFIX $PYTHON setup.py install --old-and-unmanageable +export MKLROOT=$CONDA_PREFIX +$PYTHON setup.py install diff --git a/conda-recipe/meta.yaml b/conda-recipe/meta.yaml index 9cffa23..9883238 100644 --- a/conda-recipe/meta.yaml +++ b/conda-recipe/meta.yaml @@ -13,7 +13,6 @@ build: ignore_run_exports: - blas - requirements: build: - {{ compiler('c') }} @@ -24,16 +23,17 @@ requirements: - mkl-devel - cython - numpy-base + - pip run: - python - - {{ pin_compatible('mkl', min_pin="x.x", max_pin="x.x") }} - - {{ pin_compatible('numpy') }} + - {{ pin_compatible('mkl', min_pin="x.x", max_pin="x") }} + - {{ pin_compatible('numpy', min_pin="x.x", max_pin="x") }} test: commands: - - nosetests -v mkl_random + - pytest --pyargs mkl_random requires: - - nose + - pytest - mkl-service - numpy imports: diff --git a/mkl_random/mklrand.pyx b/mkl_random/mklrand.pyx index 809e296..93c7e26 100644 --- a/mkl_random/mklrand.pyx +++ b/mkl_random/mklrand.pyx @@ -74,113 +74,113 @@ cdef extern from "randomkit.h": NONDETERM = 9 ARS5 = 10 - void irk_fill(void *buffer, size_t size, irk_state *state) nogil + void irk_fill(void *buffer, size_t size, irk_state *state) noexcept nogil void irk_dealloc_stream(irk_state *state) void irk_seed_mkl(irk_state * state, unsigned int seed, irk_brng_t brng, unsigned int stream_id) void irk_seed_mkl_array(irk_state * state, unsigned int * seed_vec, int seed_len, irk_brng_t brng, unsigned int stream_id) irk_error irk_randomseed_mkl(irk_state * state, irk_brng_t brng, unsigned int stream_id) - int irk_get_stream_size(irk_state * state) nogil + int irk_get_stream_size(irk_state * state) noexcept nogil void irk_get_state_mkl(irk_state * state, char * buf) int irk_set_state_mkl(irk_state * state, char * buf) - int irk_get_brng_mkl(irk_state *state) nogil - int irk_get_brng_and_stream_mkl(irk_state *state, unsigned int * stream_id) nogil - int irk_leapfrog_stream_mkl(irk_state *state, int k, int nstreams) nogil - int irk_skipahead_stream_mkl(irk_state *state, long long int nskips) nogil + int irk_get_brng_mkl(irk_state *state) noexcept nogil + int irk_get_brng_and_stream_mkl(irk_state *state, unsigned int * stream_id) noexcept nogil + int irk_leapfrog_stream_mkl(irk_state *state, int k, int nstreams) noexcept nogil + int irk_skipahead_stream_mkl(irk_state *state, long long int nskips) noexcept nogil cdef extern from "mkl_distributions.h": - void irk_double_vec(irk_state *state, npy_intp len, double *res) nogil - void irk_uniform_vec(irk_state *state, npy_intp len, double *res, double dlow, double dhigh) nogil - - void irk_normal_vec_BM1(irk_state *state, npy_intp len, double *res, double mean, double sigma) nogil - void irk_normal_vec_BM2(irk_state *state, npy_intp len, double *res, double mean, double sigma) nogil - void irk_normal_vec_ICDF(irk_state *state, npy_intp len, double *res, double mean, double sigma) nogil - - void irk_standard_normal_vec_BM1(irk_state *state, npy_intp len, double *res) nogil - void irk_standard_normal_vec_BM2(irk_state *state, npy_intp len, double *res) nogil - void irk_standard_normal_vec_ICDF(irk_state *state, npy_intp len, double *res) nogil - - void irk_standard_exponential_vec(irk_state *state, npy_intp len, double *res) nogil - void irk_exponential_vec(irk_state *state, npy_intp len, double *res, double scale) nogil - - void irk_standard_cauchy_vec(irk_state *state, npy_intp len, double *res) nogil - void irk_standard_gamma_vec(irk_state *state, npy_intp len, double *res, double shape) nogil - void irk_gamma_vec(irk_state *state, npy_intp len, double *res, double shape, double scale) nogil - - void irk_beta_vec(irk_state *state, npy_intp len, double *res, double p, double q) nogil - - void irk_chisquare_vec(irk_state *state, npy_intp len, double *res, double df) nogil - void irk_standard_t_vec(irk_state *state, npy_intp len, double *res, double df) nogil - - void irk_rayleigh_vec(irk_state *state, npy_intp len, double *res, double sigma) nogil - void irk_pareto_vec(irk_state *state, npy_intp len, double *res, double alp) nogil - void irk_power_vec(irk_state *state, npy_intp len, double *res, double alp) nogil - void irk_weibull_vec(irk_state *state, npy_intp len, double *res, double alp) nogil - void irk_f_vec(irk_state *state, npy_intp len, double *res, double df_num, double df_den) nogil - void irk_noncentral_chisquare_vec(irk_state *state, npy_intp len, double *res, double df, double nonc) nogil - void irk_laplace_vec(irk_state *state, npy_intp len, double *res, double loc, double scale) nogil - void irk_gumbel_vec(irk_state *state, npy_intp len, double *res, double loc, double scale) nogil - void irk_logistic_vec(irk_state *state, npy_intp len, double *res, double loc, double scale) nogil - void irk_wald_vec(irk_state *state, npy_intp len, double *res, double mean, double scale) nogil - void irk_lognormal_vec_ICDF(irk_state *state, npy_intp len, double *res, double mean, double scale) nogil - void irk_lognormal_vec_BM(irk_state *state, npy_intp len, double *res, double mean, double scale) nogil - void irk_vonmises_vec(irk_state *state, npy_intp len, double *res, double mu, double kappa) nogil - - void irk_noncentral_f_vec(irk_state *state, npy_intp len, double *res, double df_num, double df_den, double nonc) nogil - void irk_triangular_vec(irk_state *state, npy_intp len, double *res, double left, double mode, double right) nogil - - void irk_geometric_vec(irk_state *state, npy_intp len, int *res, double p) nogil - void irk_negbinomial_vec(irk_state *state, npy_intp len, int *res, double a, double p) nogil - void irk_binomial_vec(irk_state *state, npy_intp len, int *res, int n, double p) nogil - void irk_multinomial_vec(irk_state *state, npy_intp len, int *res, int n, int d, double *pvec) nogil - void irk_hypergeometric_vec(irk_state *state, npy_intp len, int *res, int ls, int ss, int ms) nogil - - void irk_poisson_vec_PTPE(irk_state *state, npy_intp len, int *res, double lam) nogil - void irk_poisson_vec_POISNORM(irk_state *state, npy_intp len, int *res, double lam) nogil - void irk_poisson_vec_V(irk_state *state, npy_intp len, int *res, double *lam_vec) nogil - - void irk_zipf_long_vec(irk_state *state, npy_intp len, long *res, double alpha) nogil - void irk_logseries_vec(irk_state *state, npy_intp len, int *res, double theta) nogil + void irk_double_vec(irk_state *state, npy_intp len, double *res) noexcept nogil + void irk_uniform_vec(irk_state *state, npy_intp len, double *res, double dlow, double dhigh) noexcept nogil + + void irk_normal_vec_BM1(irk_state *state, npy_intp len, double *res, double mean, double sigma) noexcept nogil + void irk_normal_vec_BM2(irk_state *state, npy_intp len, double *res, double mean, double sigma) noexcept nogil + void irk_normal_vec_ICDF(irk_state *state, npy_intp len, double *res, double mean, double sigma) noexcept nogil + + void irk_standard_normal_vec_BM1(irk_state *state, npy_intp len, double *res) noexcept nogil + void irk_standard_normal_vec_BM2(irk_state *state, npy_intp len, double *res) noexcept nogil + void irk_standard_normal_vec_ICDF(irk_state *state, npy_intp len, double *res) noexcept nogil + + void irk_standard_exponential_vec(irk_state *state, npy_intp len, double *res) noexcept nogil + void irk_exponential_vec(irk_state *state, npy_intp len, double *res, double scale) noexcept nogil + + void irk_standard_cauchy_vec(irk_state *state, npy_intp len, double *res) noexcept nogil + void irk_standard_gamma_vec(irk_state *state, npy_intp len, double *res, double shape) noexcept nogil + void irk_gamma_vec(irk_state *state, npy_intp len, double *res, double shape, double scale) noexcept nogil + + void irk_beta_vec(irk_state *state, npy_intp len, double *res, double p, double q) noexcept nogil + + void irk_chisquare_vec(irk_state *state, npy_intp len, double *res, double df) noexcept nogil + void irk_standard_t_vec(irk_state *state, npy_intp len, double *res, double df) noexcept nogil + + void irk_rayleigh_vec(irk_state *state, npy_intp len, double *res, double sigma) noexcept nogil + void irk_pareto_vec(irk_state *state, npy_intp len, double *res, double alp) noexcept nogil + void irk_power_vec(irk_state *state, npy_intp len, double *res, double alp) noexcept nogil + void irk_weibull_vec(irk_state *state, npy_intp len, double *res, double alp) noexcept nogil + void irk_f_vec(irk_state *state, npy_intp len, double *res, double df_num, double df_den) noexcept nogil + void irk_noncentral_chisquare_vec(irk_state *state, npy_intp len, double *res, double df, double nonc) noexcept nogil + void irk_laplace_vec(irk_state *state, npy_intp len, double *res, double loc, double scale) noexcept nogil + void irk_gumbel_vec(irk_state *state, npy_intp len, double *res, double loc, double scale) noexcept nogil + void irk_logistic_vec(irk_state *state, npy_intp len, double *res, double loc, double scale) noexcept nogil + void irk_wald_vec(irk_state *state, npy_intp len, double *res, double mean, double scale) noexcept nogil + void irk_lognormal_vec_ICDF(irk_state *state, npy_intp len, double *res, double mean, double scale) noexcept nogil + void irk_lognormal_vec_BM(irk_state *state, npy_intp len, double *res, double mean, double scale) noexcept nogil + void irk_vonmises_vec(irk_state *state, npy_intp len, double *res, double mu, double kappa) noexcept nogil + + void irk_noncentral_f_vec(irk_state *state, npy_intp len, double *res, double df_num, double df_den, double nonc) noexcept nogil + void irk_triangular_vec(irk_state *state, npy_intp len, double *res, double left, double mode, double right) noexcept nogil + + void irk_geometric_vec(irk_state *state, npy_intp len, int *res, double p) noexcept nogil + void irk_negbinomial_vec(irk_state *state, npy_intp len, int *res, double a, double p) noexcept nogil + void irk_binomial_vec(irk_state *state, npy_intp len, int *res, int n, double p) noexcept nogil + void irk_multinomial_vec(irk_state *state, npy_intp len, int *res, int n, int d, double *pvec) noexcept nogil + void irk_hypergeometric_vec(irk_state *state, npy_intp len, int *res, int ls, int ss, int ms) noexcept nogil + + void irk_poisson_vec_PTPE(irk_state *state, npy_intp len, int *res, double lam) noexcept nogil + void irk_poisson_vec_POISNORM(irk_state *state, npy_intp len, int *res, double lam) noexcept nogil + void irk_poisson_vec_V(irk_state *state, npy_intp len, int *res, double *lam_vec) noexcept nogil + + void irk_zipf_long_vec(irk_state *state, npy_intp len, long *res, double alpha) noexcept nogil + void irk_logseries_vec(irk_state *state, npy_intp len, int *res, double theta) noexcept nogil # random integers madness - void irk_discrete_uniform_vec(irk_state *state, npy_intp len, int *res, int low, int high) nogil - void irk_discrete_uniform_long_vec(irk_state *state, npy_intp len, long *res, long low, long high) nogil - void irk_rand_bool_vec(irk_state *state, npy_intp len, npy_bool *res, npy_bool low, npy_bool high) nogil - void irk_rand_uint8_vec(irk_state *state, npy_intp len, npy_uint8 *res, npy_uint8 low, npy_uint8 high) nogil - void irk_rand_int8_vec(irk_state *state, npy_intp len, npy_int8 *res, npy_int8 low, npy_int8 high) nogil - void irk_rand_uint16_vec(irk_state *state, npy_intp len, npy_uint16 *res, npy_uint16 low, npy_uint16 high) nogil - void irk_rand_int16_vec(irk_state *state, npy_intp len, npy_int16 *res, npy_int16 low, npy_int16 high) nogil - void irk_rand_uint32_vec(irk_state *state, npy_intp len, npy_uint32 *res, npy_uint32 low, npy_uint32 high) nogil - void irk_rand_int32_vec(irk_state *state, npy_intp len, npy_int32 *res, npy_int32 low, npy_int32 high) nogil - void irk_rand_uint64_vec(irk_state *state, npy_intp len, npy_uint64 *res, npy_uint64 low, npy_uint64 high) nogil - void irk_rand_int64_vec(irk_state *state, npy_intp len, npy_int64 *res, npy_int64 low, npy_int64 high) nogil - - void irk_long_vec(irk_state *state, npy_intp len, long *res) nogil + void irk_discrete_uniform_vec(irk_state *state, npy_intp len, int *res, int low, int high) noexcept nogil + void irk_discrete_uniform_long_vec(irk_state *state, npy_intp len, long *res, long low, long high) noexcept nogil + void irk_rand_bool_vec(irk_state *state, npy_intp len, npy_bool *res, npy_bool low, npy_bool high) noexcept nogil + void irk_rand_uint8_vec(irk_state *state, npy_intp len, npy_uint8 *res, npy_uint8 low, npy_uint8 high) noexcept nogil + void irk_rand_int8_vec(irk_state *state, npy_intp len, npy_int8 *res, npy_int8 low, npy_int8 high) noexcept nogil + void irk_rand_uint16_vec(irk_state *state, npy_intp len, npy_uint16 *res, npy_uint16 low, npy_uint16 high) noexcept nogil + void irk_rand_int16_vec(irk_state *state, npy_intp len, npy_int16 *res, npy_int16 low, npy_int16 high) noexcept nogil + void irk_rand_uint32_vec(irk_state *state, npy_intp len, npy_uint32 *res, npy_uint32 low, npy_uint32 high) noexcept nogil + void irk_rand_int32_vec(irk_state *state, npy_intp len, npy_int32 *res, npy_int32 low, npy_int32 high) noexcept nogil + void irk_rand_uint64_vec(irk_state *state, npy_intp len, npy_uint64 *res, npy_uint64 low, npy_uint64 high) noexcept nogil + void irk_rand_int64_vec(irk_state *state, npy_intp len, npy_int64 *res, npy_int64 low, npy_int64 high) noexcept nogil + + void irk_long_vec(irk_state *state, npy_intp len, long *res) noexcept nogil ctypedef enum ch_st_enum: MATRIX = 0 PACKED = 1 DIAGONAL = 2 - void irk_multinormal_vec_ICDF(irk_state *state, npy_intp len, double *res, int dim, double *mean_vec, double *ch, ch_st_enum storage_mode) nogil - void irk_multinormal_vec_BM1(irk_state *state, npy_intp len, double *res, int dim, double *mean_vec, double *ch, ch_st_enum storage_mode) nogil - void irk_multinormal_vec_BM2(irk_state *state, npy_intp len, double *res, int dim, double *mean_vec, double *ch, ch_st_enum storage_mode) nogil + void irk_multinormal_vec_ICDF(irk_state *state, npy_intp len, double *res, int dim, double *mean_vec, double *ch, ch_st_enum storage_mode) noexcept nogil + void irk_multinormal_vec_BM1(irk_state *state, npy_intp len, double *res, int dim, double *mean_vec, double *ch, ch_st_enum storage_mode) noexcept nogil + void irk_multinormal_vec_BM2(irk_state *state, npy_intp len, double *res, int dim, double *mean_vec, double *ch, ch_st_enum storage_mode) noexcept nogil -ctypedef void (* irk_cont0_vec)(irk_state *state, npy_intp len, double *res) nogil -ctypedef void (* irk_cont1_vec)(irk_state *state, npy_intp len, double *res, double a) nogil -ctypedef void (* irk_cont2_vec)(irk_state *state, npy_intp len, double *res, double a, double b) nogil -ctypedef void (* irk_cont3_vec)(irk_state *state, npy_intp len, double *res, double a, double b, double c) nogil +ctypedef void (* irk_cont0_vec)(irk_state *state, npy_intp len, double *res) noexcept nogil +ctypedef void (* irk_cont1_vec)(irk_state *state, npy_intp len, double *res, double a) noexcept nogil +ctypedef void (* irk_cont2_vec)(irk_state *state, npy_intp len, double *res, double a, double b) noexcept nogil +ctypedef void (* irk_cont3_vec)(irk_state *state, npy_intp len, double *res, double a, double b, double c) noexcept nogil -ctypedef void (* irk_disc0_vec)(irk_state *state, npy_intp len, int *res) nogil -ctypedef void (* irk_disc0_vec_long)(irk_state *state, npy_intp len, long *res) nogil -ctypedef void (* irk_discnp_vec)(irk_state *state, npy_intp len, int *res, int n, double a) nogil -ctypedef void (* irk_discdd_vec)(irk_state *state, npy_intp len, int *res, double n, double p) nogil -ctypedef void (* irk_discnmN_vec)(irk_state *state, npy_intp len, int *res, int n, int m, int N) nogil -ctypedef void (* irk_discd_vec)(irk_state *state, npy_intp len, int *res, double a) nogil -ctypedef void (* irk_discd_long_vec)(irk_state *state, npy_intp len, long *res, double a) nogil -ctypedef void (* irk_discdptr_vec)(irk_state *state, npy_intp len, int *res, double *a) nogil +ctypedef void (* irk_disc0_vec)(irk_state *state, npy_intp len, int *res) noexcept nogil +ctypedef void (* irk_disc0_vec_long)(irk_state *state, npy_intp len, long *res) noexcept nogil +ctypedef void (* irk_discnp_vec)(irk_state *state, npy_intp len, int *res, int n, double a) noexcept nogil +ctypedef void (* irk_discdd_vec)(irk_state *state, npy_intp len, int *res, double n, double p) noexcept nogil +ctypedef void (* irk_discnmN_vec)(irk_state *state, npy_intp len, int *res, int n, int m, int N) noexcept nogil +ctypedef void (* irk_discd_vec)(irk_state *state, npy_intp len, int *res, double a) noexcept nogil +ctypedef void (* irk_discd_long_vec)(irk_state *state, npy_intp len, long *res, double a) noexcept nogil +ctypedef void (* irk_discdptr_vec)(irk_state *state, npy_intp len, int *res, double *a) noexcept nogil cdef int r = _import_array() diff --git a/mkl_random/src/mkl_distributions.cpp b/mkl_random/src/mkl_distributions.cpp index c942207..38be263 100644 --- a/mkl_random/src/mkl_distributions.cpp +++ b/mkl_random/src/mkl_distributions.cpp @@ -25,7 +25,7 @@ OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. */ -#include /* for NULL */ +#include /* for nullptr */ #include /* for ULONG_MAX */ #include #include /* fmod, fabs */ @@ -60,7 +60,7 @@ void irk_double_vec(irk_state *state, npy_intp len, double *res) { - int err; + int err = 0; const double d_zero = 0.0, d_one = 1.0; if (len < 1) @@ -81,7 +81,7 @@ void irk_double_vec(irk_state *state, npy_intp len, double *res) void irk_uniform_vec(irk_state *state, npy_intp len, double *res, const double low, const double high) { - int err; + int err = 0; if (len < 1) return; @@ -101,7 +101,7 @@ void irk_uniform_vec(irk_state *state, npy_intp len, double *res, const double l void irk_standard_normal_vec_ICDF(irk_state *state, npy_intp len, double *res) { - int err; + int err = 0; const double d_zero = 0.0, d_one = 1.0; if (len < 1) @@ -122,7 +122,7 @@ void irk_standard_normal_vec_ICDF(irk_state *state, npy_intp len, double *res) void irk_normal_vec_ICDF(irk_state *state, npy_intp len, double *res, const double loc, const double scale) { - int err; + int err = 0; if (len < 1) return; @@ -142,7 +142,7 @@ void irk_normal_vec_ICDF(irk_state *state, npy_intp len, double *res, const doub void irk_standard_normal_vec_BM1(irk_state *state, npy_intp len, double *res) { - int err; + int err = 0; const double d_zero = 0.0, d_one = 1.0; if (len < 1) @@ -163,7 +163,7 @@ void irk_standard_normal_vec_BM1(irk_state *state, npy_intp len, double *res) void irk_normal_vec_BM1(irk_state *state, npy_intp len, double *res, const double loc, const double scale) { - int err; + int err = 0; if (len < 1) return; @@ -183,7 +183,7 @@ void irk_normal_vec_BM1(irk_state *state, npy_intp len, double *res, const doubl void irk_standard_normal_vec_BM2(irk_state *state, npy_intp len, double *res) { - int err; + int err = 0; const double d_zero = 0.0, d_one = 1.0; if (len < 1) @@ -204,7 +204,7 @@ void irk_standard_normal_vec_BM2(irk_state *state, npy_intp len, double *res) void irk_normal_vec_BM2(irk_state *state, npy_intp len, double *res, const double loc, const double scale) { - int err; + int err = 0; if (len < 1) return; @@ -224,7 +224,7 @@ void irk_normal_vec_BM2(irk_state *state, npy_intp len, double *res, const doubl void irk_standard_exponential_vec(irk_state *state, npy_intp len, double *res) { - int err; + int err = 0; const double d_zero = 0.0, d_one = 1.0; if (len < 1) @@ -245,7 +245,7 @@ void irk_standard_exponential_vec(irk_state *state, npy_intp len, double *res) void irk_exponential_vec(irk_state *state, npy_intp len, double *res, const double scale) { - int err; + int err = 0; const double d_zero = 0.0; if (len < 1) @@ -266,7 +266,7 @@ void irk_exponential_vec(irk_state *state, npy_intp len, double *res, const doub void irk_standard_cauchy_vec(irk_state *state, npy_intp len, double *res) { - int err; + int err = 0; const double d_zero = 0.0, d_one = 1.0; if (len < 1) @@ -287,7 +287,7 @@ void irk_standard_cauchy_vec(irk_state *state, npy_intp len, double *res) void irk_standard_gamma_vec(irk_state *state, npy_intp len, double *res, const double shape) { - int err; + int err = 0; const double d_zero = 0.0, d_one = 1.0; if (len < 1) @@ -308,7 +308,7 @@ void irk_standard_gamma_vec(irk_state *state, npy_intp len, double *res, const d void irk_gamma_vec(irk_state *state, npy_intp len, double *res, const double shape, const double scale) { - int err; + int err = 0; const double d_zero = 0.0; if (len < 1) @@ -330,10 +330,10 @@ void irk_gamma_vec(irk_state *state, npy_intp len, double *res, const double sha /* X ~ Z * (G*(2/df))**-0.5 */ void irk_standard_t_vec(irk_state *state, npy_intp len, double *res, const double df) { - int err; + int err = 0; const double d_zero = 0.0, d_one = 1.0; double shape = df / 2; - double *sn = NULL; + double *sn = nullptr; if (len < 1) return; @@ -352,7 +352,7 @@ void irk_standard_t_vec(irk_state *state, npy_intp len, double *res, const doubl vmdInvSqrt(len, res, res, VML_HA); sn = (double *)mkl_malloc(len * sizeof(double), 64); - assert(sn != NULL); + assert(sn != nullptr); err = vdRngGaussian(VSL_RNG_METHOD_GAUSSIAN_ICDF, state->stream, len, sn, d_zero, d_one); assert(err == VSL_STATUS_OK); @@ -364,7 +364,7 @@ void irk_standard_t_vec(irk_state *state, npy_intp len, double *res, const doubl /* chisquare(df) ~ G(df/2, 2) */ void irk_chisquare_vec(irk_state *state, npy_intp len, double *res, const double df) { - int err; + int err = 0; const double d_zero = 0.0, d_two = 2.0; double shape = 0.5 * df; @@ -386,7 +386,8 @@ void irk_chisquare_vec(irk_state *state, npy_intp len, double *res, const double /* P ~ U^(-1/a) - 1 = */ void irk_pareto_vec(irk_state *state, npy_intp len, double *res, const double alp) { - int i, err; + int err = 0; + npy_intp i = 0; const double d_zero = 0.0, d_one = 1.0; double neg_rec_alp = -1.0 / alp; @@ -408,14 +409,14 @@ void irk_pareto_vec(irk_state *state, npy_intp len, double *res, const double al vmdPowx(len, res, neg_rec_alp, res, VML_HA); DIST_PRAGMA_VECTOR - for (i = 0; i < len; i++) + for (i = 0; i < len; ++i) res[i] -= 1.0; } /* W ~ E^(1/alp) */ void irk_weibull_vec(irk_state *state, npy_intp len, double *res, const double alp) { - int err; + int err = 0; const double d_zero = 0.0, d_one = 1.0; double rec_alp = 1.0 / alp; @@ -439,7 +440,7 @@ void irk_weibull_vec(irk_state *state, npy_intp len, double *res, const double a /* pow(1 - exp(-E(1))), 1./a) == pow(U, 1./a) */ void irk_power_vec(irk_state *state, npy_intp len, double *res, const double alp) { - int err; + int err = 0; const double d_zero = 0.0, d_one = 1.0; double rec_alp = 1.0 / alp; @@ -464,7 +465,8 @@ void irk_power_vec(irk_state *state, npy_intp len, double *res, const double alp /* scale * sqrt(2.0 * E(1)) */ void irk_rayleigh_vec(irk_state *state, npy_intp len, double *res, const double scale) { - int i, err; + int err = 0; + npy_intp i = 0; const double d_zero = 0.0, d_two = 2.0; if (len < 1) @@ -484,13 +486,13 @@ void irk_rayleigh_vec(irk_state *state, npy_intp len, double *res, const double vmdSqrt(len, res, res, VML_HA); DIST_PRAGMA_VECTOR - for (i = 0; i < len; i++) + for (i = 0; i < len; ++i) res[i] *= scale; } void irk_beta_vec(irk_state *state, npy_intp len, double *res, const double a, const double b) { - int err; + int err = 0; const double d_zero = 0.0, d_one = 1.0; if (len < 1) @@ -512,10 +514,10 @@ void irk_beta_vec(irk_state *state, npy_intp len, double *res, const double a, c /* F(df_num, df_den) ~ G( df_num/2, 2/df_num) / G(df_den/2, 2/df_den)) */ void irk_f_vec(irk_state *state, npy_intp len, double *res, const double df_num, const double df_den) { - int err; + int err = 0; const double d_zero = 0.0; double shape = 0.5 * df_num, scale = 2.0 / df_num; - double *den = NULL; + double *den = nullptr; if (len < 1) return; @@ -532,7 +534,7 @@ void irk_f_vec(irk_state *state, npy_intp len, double *res, const double df_num, assert(err == VSL_STATUS_OK); den = (double *)mkl_malloc(len * sizeof(double), 64); - assert(den != NULL); + assert(den != nullptr); shape = 0.5 * df_den; scale = 2.0 / df_den; @@ -549,7 +551,8 @@ void irk_f_vec(irk_state *state, npy_intp len, double *res, const double df_num, */ void irk_noncentral_chisquare_vec(irk_state *state, npy_intp len, double *res, const double df, const double nonc) { - int i, err; + int err = 0; + npy_intp i = 0; const double d_zero = 0.0, d_one = 1.0, d_two = 2.0; double shape, loc; @@ -573,7 +576,7 @@ void irk_noncentral_chisquare_vec(irk_state *state, npy_intp len, double *res, c err = vdRngGamma(VSL_RNG_METHOD_GAMMA_GNORM_ACCURATE, state->stream, len, res, shape, d_zero, d_two); nvec = (double *)mkl_malloc(len * sizeof(double), 64); - assert(nvec != NULL); + assert(nvec != nullptr); loc = sqrt(nonc); err = vdRngGaussian(VSL_RNG_METHOD_GAUSSIAN_ICDF, state->stream, len, nvec, loc, d_one); @@ -597,7 +600,7 @@ void irk_noncentral_chisquare_vec(irk_state *state, npy_intp len, double *res, c double lambda; int *pvec = (int *)mkl_malloc(len * sizeof(int), 64); - assert(pvec != NULL); + assert(pvec != nullptr); lambda = 0.5 * nonc; err = viRngPoisson(VSL_RNG_METHOD_POISSON_PTPE, state->stream, len, pvec, lambda); @@ -607,15 +610,15 @@ void irk_noncentral_chisquare_vec(irk_state *state, npy_intp len, double *res, c if (0.125 * len > sqrt(lambda)) { - int *idx = NULL; - double *tmp = NULL; + int *idx = nullptr; + double *tmp = nullptr; idx = (int *)mkl_malloc(len * sizeof(int), 64); - assert(idx != NULL); + assert(idx != nullptr); DIST_PRAGMA_VECTOR - for (i = 0; i < len; i++) - idx[i] = i; + for (i = 0; i < len; ++i) + idx[i] = (int)i; std::sort(idx, idx + len, [pvec](int i1, int i2) { return pvec[i1] < pvec[i2]; }); @@ -623,13 +626,14 @@ void irk_noncentral_chisquare_vec(irk_state *state, npy_intp len, double *res, c /* allocate workspace to store samples of gamma, enough to hold entire output */ tmp = (double *)mkl_malloc(len * sizeof(double), 64); - assert(tmp != NULL); + assert(tmp != nullptr); for (i = 0; i < len;) { - int k, j, cv = pvec[idx[i]]; + int cv = pvec[idx[i]]; + npy_intp k = 0, j = 0; - for (j = i + 1; (j < len) && (pvec[idx[j]] == cv); j++) + for (j = i + 1; (j < len) && (pvec[idx[j]] == cv); ++j) { } @@ -639,8 +643,8 @@ void irk_noncentral_chisquare_vec(irk_state *state, npy_intp len, double *res, c assert(err == VSL_STATUS_OK); DIST_PRAGMA_VECTOR - for (k = i; k < j; k++) - res[idx[k]] = tmp[k - i]; + for (k = 0; k < j - i; ++k) + res[idx[k + i]] = tmp[k]; i = j; } @@ -651,7 +655,7 @@ void irk_noncentral_chisquare_vec(irk_state *state, npy_intp len, double *res, c else { - for (i = 0; i < len; i++) + for (i = 0; i < len; ++i) { err = vdRngGamma(VSL_RNG_METHOD_GAMMA_GNORM_ACCURATE, state->stream, 1, res + i, shape + pvec[i], d_zero, d_two); @@ -675,7 +679,7 @@ void irk_noncentral_chisquare_vec(irk_state *state, npy_intp len, double *res, c void irk_laplace_vec(irk_state *state, npy_intp len, double *res, const double loc, const double scale) { - int err; + int err = 0; if (len < 1) return; @@ -695,7 +699,7 @@ void irk_laplace_vec(irk_state *state, npy_intp len, double *res, const double l void irk_gumbel_vec(irk_state *state, npy_intp len, double *res, const double loc, const double scale) { - int err; + int err = 0; if (len < 1) return; @@ -716,7 +720,8 @@ void irk_gumbel_vec(irk_state *state, npy_intp len, double *res, const double lo /* Logistic(loc, scale) ~ loc + scale * log(u/(1.0 - u)) */ void irk_logistic_vec(irk_state *state, npy_intp len, double *res, const double loc, const double scale) { - int i, err; + int err = 0; + npy_intp i = 0; const double d_one = 1.0, d_zero = 0.0; if (len < 1) @@ -735,17 +740,17 @@ void irk_logistic_vec(irk_state *state, npy_intp len, double *res, const double /* can MKL optimize computation of the logit function p \mapsto \ln(p/(1-p)) */ DIST_PRAGMA_VECTOR - for (i = 0; i < len; i++) + for (i = 0; i < len; ++i) res[i] = log(res[i] / (1.0 - res[i])); DIST_PRAGMA_VECTOR - for (i = 0; i < len; i++) + for (i = 0; i < len; ++i) res[i] = loc + scale * res[i]; } void irk_lognormal_vec_ICDF(irk_state *state, npy_intp len, double *res, const double mean, const double sigma) { - int err; + int err = 0; const double d_zero = 0.0, d_one = 1.0; if (len < 1) @@ -766,7 +771,7 @@ void irk_lognormal_vec_ICDF(irk_state *state, npy_intp len, double *res, const d void irk_lognormal_vec_BM(irk_state *state, npy_intp len, double *res, const double mean, const double sigma) { - int err; + int err = 0; const double d_zero = 0.0, d_one = 1.0; if (len < 1) @@ -788,9 +793,10 @@ void irk_lognormal_vec_BM(irk_state *state, npy_intp len, double *res, const dou /* direct transformation method */ void irk_wald_vec(irk_state *state, npy_intp len, double *res, const double mean, const double scale) { - int i, err; + int err = 0; + npy_intp i = 0; const double d_zero = 0., d_one = 1.0; - double *uvec = NULL; + double *uvec = nullptr; double gsc = sqrt(0.5 * mean / scale); if (len < 1) @@ -811,7 +817,7 @@ void irk_wald_vec(irk_state *state, npy_intp len, double *res, const double mean vmdSqr(len, res, res, VML_HA); DIST_PRAGMA_VECTOR - for (i = 0; i < len; i++) + for (i = 0; i < len; ++i) { if (res[i] <= 2.0) { @@ -824,13 +830,13 @@ void irk_wald_vec(irk_state *state, npy_intp len, double *res, const double mean } uvec = (double *)mkl_malloc(len * sizeof(double), 64); - assert(uvec != NULL); + assert(uvec != nullptr); err = vdRngUniform(VSL_RNG_METHOD_UNIFORM_STD_ACCURATE, state->stream, len, uvec, d_zero, d_one); assert(err == VSL_STATUS_OK); DIST_PRAGMA_VECTOR - for (i = 0; i < len; i++) + for (i = 0; i < len; ++i) { if (uvec[i] * (1.0 + res[i]) <= res[i]) res[i] = mean / res[i]; @@ -855,10 +861,11 @@ void irk_wald_vec(irk_state *state, npy_intp len, double *res, const double mean static void irk_vonmises_vec_small_kappa(irk_state *state, npy_intp len, double *res, const double mu, const double kappa) { - int i, err, n, size; + int err = 0; + npy_intp n = 0, i = 0, size = 0; double rho_over_kappa, rho, r, s_kappa, Z, W, Y, V; - double *Uvec = NULL, *Vvec = NULL; - float *VFvec = NULL; + double *Uvec = nullptr, *Vvec = nullptr; + float *VFvec = nullptr; const double d_zero = 0.0, d_one = 1.0; assert(0. < kappa <= 1.0); @@ -871,9 +878,9 @@ irk_vonmises_vec_small_kappa(irk_state *state, npy_intp len, double *res, const s_kappa = (1 + rho * rho) / (2 * rho_over_kappa); Uvec = (double *)mkl_malloc(len * sizeof(double), 64); - assert(Uvec != NULL); + assert(Uvec != nullptr); Vvec = (double *)mkl_malloc(len * sizeof(double), 64); - assert(Vvec != NULL); + assert(Vvec != nullptr); for (n = 0; n < len;) { @@ -883,7 +890,7 @@ irk_vonmises_vec_small_kappa(irk_state *state, npy_intp len, double *res, const err = vdRngUniform(VSL_RNG_METHOD_UNIFORM_STD_ACCURATE, state->stream, size, Vvec, d_zero, d_one); assert(err == VSL_STATUS_OK); - for (i = 0; i < size; i++) + for (i = 0; i < size; ++i) { Z = cos(Uvec[i]); V = Vvec[i]; @@ -903,7 +910,7 @@ irk_vonmises_vec_small_kappa(irk_state *state, npy_intp len, double *res, const assert(err == VSL_STATUS_OK); DIST_PRAGMA_VECTOR - for (i = 0; i < len; i++) + for (i = 0; i < len; ++i) { double mod, resi; @@ -919,11 +926,12 @@ irk_vonmises_vec_small_kappa(irk_state *state, npy_intp len, double *res, const static void irk_vonmises_vec_large_kappa(irk_state *state, npy_intp len, double *res, const double mu, const double kappa) { - int i, err, n, size; + int err = 0; + npy_intp i = 0, n = 0, size = 0; double r_over_two_kappa, recip_two_kappa; double s_minus_one, hpt, r_over_two_kappa_minus_one, rho_minus_one; - double *Uvec = NULL, *Vvec = NULL; - float *VFvec = NULL; + double *Uvec = nullptr, *Vvec = nullptr; + float *VFvec = nullptr; const double d_zero = 0.0, d_one = 1.0; assert(kappa > 1.0); @@ -938,9 +946,9 @@ irk_vonmises_vec_large_kappa(irk_state *state, npy_intp len, double *res, const s_minus_one = rho_minus_one * (0.5 * rho_minus_one / (1 + rho_minus_one)); Uvec = (double *)mkl_malloc(len * sizeof(double), 64); - assert(Uvec != NULL); + assert(Uvec != nullptr); Vvec = (double *)mkl_malloc(len * sizeof(double), 64); - assert(Vvec != NULL); + assert(Vvec != nullptr); for (n = 0; n < len;) { @@ -950,8 +958,7 @@ irk_vonmises_vec_large_kappa(irk_state *state, npy_intp len, double *res, const err = vdRngUniform(VSL_RNG_METHOD_UNIFORM_STD_ACCURATE, state->stream, size, Vvec, d_zero, d_one); assert(err == VSL_STATUS_OK); - DIST_PRAGMA_VECTOR - for (i = 0; i < size; i++) + for (i = 0; i < size; ++i) { double sn, cn, sn2, cn2; double neg_W_minus_one, V, Y; @@ -985,7 +992,7 @@ irk_vonmises_vec_large_kappa(irk_state *state, npy_intp len, double *res, const assert(err == VSL_STATUS_OK); DIST_PRAGMA_VECTOR - for (i = 0; i < len; i++) + for (i = 0; i < len; ++i) { double mod, resi; @@ -1019,8 +1026,8 @@ void irk_vonmises_vec(irk_state *state, npy_intp len, double *res, const double void irk_noncentral_f_vec(irk_state *state, npy_intp len, double *res, const double df_num, const double df_den, const double nonc) { - int i; - double *den = NULL, fctr; + npy_intp i; + double *den = nullptr, fctr; if (len < 1) return; @@ -1040,7 +1047,7 @@ void irk_noncentral_f_vec(irk_state *state, npy_intp len, double *res, const dou den = (double *)mkl_malloc(len * sizeof(double), 64); - if (den == NULL) + if (den == nullptr) return; irk_noncentral_chisquare_vec(state, len, den, df_den, nonc); @@ -1051,13 +1058,14 @@ void irk_noncentral_f_vec(irk_state *state, npy_intp len, double *res, const dou fctr = df_den / df_num; DIST_PRAGMA_VECTOR - for (i = 0; i < len; i++) + for (i = 0; i < len; ++i) res[i] *= fctr; } void irk_triangular_vec(irk_state *state, npy_intp len, double *res, const double x_min, const double x_mode, const double x_max) { - int i, err; + int err = 0; + npy_intp i = 0; const double d_zero = 0.0, d_one = 1.0; double ratio, lpr, rpr; @@ -1092,7 +1100,7 @@ void irk_triangular_vec(irk_state *state, npy_intp len, double *res, const doubl if (ratio <= 0) { DIST_PRAGMA_VECTOR - for (i = 0; i < len; i++) + for (i = 0; i < len; ++i) { /* U and 1 - U are equal in distribution */ res[i] = x_max - sqrt(res[i] * rpr); @@ -1101,7 +1109,7 @@ void irk_triangular_vec(irk_state *state, npy_intp len, double *res, const doubl else if (ratio >= 1) { DIST_PRAGMA_VECTOR - for (i = 0; i < len; i++) + for (i = 0; i < len; ++i) { res[i] = x_min + sqrt(res[i] * lpr); } @@ -1109,7 +1117,7 @@ void irk_triangular_vec(irk_state *state, npy_intp len, double *res, const doubl else { DIST_PRAGMA_VECTOR - for (i = 0; i < len; i++) + for (i = 0; i < len; ++i) { double ui = res[i]; res[i] = (ui > ratio) ? x_max - sqrt((1.0 - ui) * rpr) : x_min + sqrt(ui * lpr); @@ -1119,7 +1127,7 @@ void irk_triangular_vec(irk_state *state, npy_intp len, double *res, const doubl void irk_binomial_vec(irk_state *state, npy_intp len, int *res, const int n, const double p) { - int err; + int err = 0; if (len < 1) return; @@ -1145,7 +1153,7 @@ void irk_binomial_vec(irk_state *state, npy_intp len, int *res, const int n, con void irk_multinomial_vec(irk_state *state, npy_intp len, int *res, const int n, const int k, const double *pvec) { - int err; + int err = 0; if (len < 1) return; @@ -1171,7 +1179,7 @@ void irk_multinomial_vec(irk_state *state, npy_intp len, int *res, const int n, void irk_geometric_vec(irk_state *state, npy_intp len, int *res, const double p) { - int err; + int err = 0; if (len < 1) return; @@ -1208,7 +1216,7 @@ void irk_geometric_vec(irk_state *state, npy_intp len, int *res, const double p) void irk_negbinomial_vec(irk_state *state, npy_intp len, int *res, const double a, const double p) { - int err; + int err = 0; if (len < 1) return; @@ -1229,7 +1237,7 @@ void irk_negbinomial_vec(irk_state *state, npy_intp len, int *res, const double void irk_hypergeometric_vec(irk_state *state, npy_intp len, int *res, const int lot_s, const int sampling_s, const int marked_s) { - int err; + int err = 0; if (len < 1) return; @@ -1251,7 +1259,7 @@ void irk_hypergeometric_vec(irk_state *state, npy_intp len, int *res, const int void irk_poisson_vec_PTPE(irk_state *state, npy_intp len, int *res, const double lambda) { - int err; + int err = 0; if (len < 1) return; @@ -1271,7 +1279,7 @@ void irk_poisson_vec_PTPE(irk_state *state, npy_intp len, int *res, const double void irk_poisson_vec_POISNORM(irk_state *state, npy_intp len, int *res, const double lambda) { - int err; + int err = 0; if (len < 1) return; @@ -1291,7 +1299,7 @@ void irk_poisson_vec_POISNORM(irk_state *state, npy_intp len, int *res, const do void irk_poisson_vec_V(irk_state *state, npy_intp len, int *res, double *lambdas) { - int err; + int err = 0; if (len < 1) return; @@ -1312,9 +1320,10 @@ void irk_poisson_vec_V(irk_state *state, npy_intp len, int *res, double *lambdas void irk_zipf_long_vec(irk_state *state, npy_intp len, long *res, const double a) { - int i, err, n_accepted, batch_size; + int err = 0; + npy_intp i = 0, n_accepted = 0, batch_size = 0; double T, U, V, am1, b; - double *Uvec = NULL, *Vvec = NULL; + double *Uvec = nullptr, *Vvec = nullptr; long X; const double d_zero = 0.0, d_one = 1.0; @@ -1333,9 +1342,9 @@ void irk_zipf_long_vec(irk_state *state, npy_intp len, long *res, const double a b = pow(2.0, am1); Uvec = (double *)mkl_malloc(len * sizeof(double), 64); - assert(Uvec != NULL); + assert(Uvec != nullptr); Vvec = (double *)mkl_malloc(len * sizeof(double), 64); - assert(Vvec != NULL); + assert(Vvec != nullptr); for (n_accepted = 0; n_accepted < len;) { @@ -1346,7 +1355,7 @@ void irk_zipf_long_vec(irk_state *state, npy_intp len, long *res, const double a assert(err == VSL_STATUS_OK); DIST_PRAGMA_VECTOR - for (i = 0; i < batch_size; i++) + for (i = 0; i < batch_size; ++i) { U = d_one - Uvec[i]; V = Vvec[i]; @@ -1371,9 +1380,10 @@ void irk_zipf_long_vec(irk_state *state, npy_intp len, long *res, const double a void irk_logseries_vec(irk_state *state, npy_intp len, int *res, const double theta) { - int i, err, n_accepted, batch_size; + int err = 0; + npy_intp i = 0, n_accepted = 0, batch_size = 0; double q, r, V; - double *Uvec = NULL, *Vvec = NULL; + double *Uvec = nullptr, *Vvec = nullptr; int result; const double d_zero = 0.0, d_one = 1.0; @@ -1391,9 +1401,9 @@ void irk_logseries_vec(irk_state *state, npy_intp len, int *res, const double th r = log(d_one - theta); Uvec = (double *)mkl_malloc(len * sizeof(double), 64); - assert(Uvec != NULL); + assert(Uvec != nullptr); Vvec = (double *)mkl_malloc(len * sizeof(double), 64); - assert(Vvec != NULL); + assert(Vvec != nullptr); for (n_accepted = 0; n_accepted < len;) { @@ -1404,7 +1414,7 @@ void irk_logseries_vec(irk_state *state, npy_intp len, int *res, const double th assert(err == VSL_STATUS_OK); DIST_PRAGMA_VECTOR - for (i = 0; i < batch_size; i++) + for (i = 0; i < batch_size; ++i) { V = Vvec[i]; if (V >= theta) @@ -1450,7 +1460,7 @@ void irk_logseries_vec(irk_state *state, npy_intp len, int *res, const double th /* samples discrete uniforms from [low, high) */ void irk_discrete_uniform_vec(irk_state *state, npy_intp len, int *res, const int low, const int high) { - int err; + int err = 0; if (len < 1) return; @@ -1470,9 +1480,9 @@ void irk_discrete_uniform_vec(irk_state *state, npy_intp len, int *res, const in void irk_discrete_uniform_long_vec(irk_state *state, npy_intp len, long *res, const long low, const long high) { - int err; + int err = 0; unsigned long max; - int i; + npy_intp i = 0; if (len < 1) return; @@ -1489,7 +1499,7 @@ void irk_discrete_uniform_long_vec(irk_state *state, npy_intp len, long *res, co if (max == 0) { DIST_PRAGMA_VECTOR - for (i = 0; i < len; i++) + for (i = 0; i < len; ++i) res[i] = low; return; @@ -1498,13 +1508,13 @@ void irk_discrete_uniform_long_vec(irk_state *state, npy_intp len, long *res, co if (max <= (unsigned long)INT_MAX) { int *buf = (int *)mkl_malloc(len * sizeof(int), 64); - assert(buf != NULL); + assert(buf != nullptr); err = viRngUniform(VSL_RNG_METHOD_UNIFORM_STD, state->stream, len, buf, -1, (const int)max); assert(err == VSL_STATUS_OK); DIST_PRAGMA_VECTOR - for (i = 0; i < len; i++) + for (i = 0; i < len; ++i) res[i] = low + ((long)buf[i]) + 1L; mkl_free(buf); @@ -1512,7 +1522,7 @@ void irk_discrete_uniform_long_vec(irk_state *state, npy_intp len, long *res, co else { unsigned long mask = max; - unsigned long *buf = NULL; + unsigned long *buf = nullptr; int n_accepted; /* Smallest bit mask >= max */ @@ -1526,7 +1536,7 @@ void irk_discrete_uniform_long_vec(irk_state *state, npy_intp len, long *res, co #endif buf = (unsigned long *)mkl_malloc(len * sizeof(long), 64); - assert(buf != NULL); + assert(buf != nullptr); n_accepted = 0; while (n_accepted < len) @@ -1536,7 +1546,7 @@ void irk_discrete_uniform_long_vec(irk_state *state, npy_intp len, long *res, co err = viRngUniformBits64(VSL_RNG_METHOD_UNIFORM_STD, state->stream, batchSize, (unsigned MKL_INT64 *)buf); assert(err == VSL_STATUS_OK); - for (k = 0; k < batchSize; k++) + for (k = 0; k < batchSize; ++k) { unsigned long value = buf[k] & mask; if (value <= max) @@ -1552,7 +1562,7 @@ void irk_discrete_uniform_long_vec(irk_state *state, npy_intp len, long *res, co void irk_ulong_vec(irk_state *state, npy_intp len, unsigned long *res) { - int err; + int err = 0; if (len < 1) return; @@ -1576,20 +1586,21 @@ void irk_ulong_vec(irk_state *state, npy_intp len, unsigned long *res) void irk_long_vec(irk_state *state, npy_intp len, long *res) { - npy_intp i; + npy_intp i = 0; unsigned long *ulptr = (unsigned long *)res; irk_ulong_vec(state, len, ulptr); DIST_PRAGMA_VECTOR - for (i = 0; i < len; i++) + for (i = 0; i < len; ++i) res[i] = (long)(ulptr[i] >> 1); } void irk_rand_bool_vec(irk_state *state, npy_intp len, npy_bool *res, const npy_bool lo, const npy_bool hi) { - int err, i; - int *buf = NULL; + int err = 0; + npy_intp i = 0; + int *buf = nullptr; if (len < 1) return; @@ -1605,7 +1616,7 @@ void irk_rand_bool_vec(irk_state *state, npy_intp len, npy_bool *res, const npy_ if (lo == hi) { DIST_PRAGMA_VECTOR - for (i = 0; i < len; i++) + for (i = 0; i < len; ++i) res[i] = lo; return; @@ -1613,13 +1624,13 @@ void irk_rand_bool_vec(irk_state *state, npy_intp len, npy_bool *res, const npy_ assert((lo == 0) && (hi == 1)); buf = (int *)mkl_malloc(len * sizeof(int), 64); - assert(buf != NULL); + assert(buf != nullptr); err = viRngUniform(VSL_RNG_METHOD_UNIFORM_STD, state->stream, len, buf, (const int)lo, (const int)hi + 1); assert(err == VSL_STATUS_OK); DIST_PRAGMA_VECTOR - for (i = 0; i < len; i++) + for (i = 0; i < len; ++i) res[i] = (npy_bool)buf[i]; mkl_free(buf); @@ -1627,8 +1638,9 @@ void irk_rand_bool_vec(irk_state *state, npy_intp len, npy_bool *res, const npy_ void irk_rand_uint8_vec(irk_state *state, npy_intp len, npy_uint8 *res, const npy_uint8 lo, const npy_uint8 hi) { - int err, i; - int *buf = NULL; + int err = 0; + npy_intp i = 0; + int *buf = nullptr; if (len < 1) return; @@ -1644,7 +1656,7 @@ void irk_rand_uint8_vec(irk_state *state, npy_intp len, npy_uint8 *res, const np if (lo == hi) { DIST_PRAGMA_VECTOR - for (i = 0; i < len; i++) + for (i = 0; i < len; ++i) res[i] = lo; return; @@ -1652,13 +1664,13 @@ void irk_rand_uint8_vec(irk_state *state, npy_intp len, npy_uint8 *res, const np assert(lo < hi); buf = (int *)mkl_malloc(len * sizeof(int), 64); - assert(buf != NULL); + assert(buf != nullptr); err = viRngUniform(VSL_RNG_METHOD_UNIFORM_STD, state->stream, len, buf, (const int)lo, (const int)hi + 1); assert(err == VSL_STATUS_OK); DIST_PRAGMA_VECTOR - for (i = 0; i < len; i++) + for (i = 0; i < len; ++i) res[i] = (npy_uint8)buf[i]; mkl_free(buf); @@ -1666,8 +1678,9 @@ void irk_rand_uint8_vec(irk_state *state, npy_intp len, npy_uint8 *res, const np void irk_rand_int8_vec(irk_state *state, npy_intp len, npy_int8 *res, const npy_int8 lo, const npy_int8 hi) { - int err, i; - int *buf = NULL; + int err = 0; + npy_intp i = 0; + int *buf = nullptr; if (len < 1) return; @@ -1683,7 +1696,7 @@ void irk_rand_int8_vec(irk_state *state, npy_intp len, npy_int8 *res, const npy_ if (lo == hi) { DIST_PRAGMA_VECTOR - for (i = 0; i < len; i++) + for (i = 0; i < len; ++i) res[i] = lo; return; @@ -1691,13 +1704,13 @@ void irk_rand_int8_vec(irk_state *state, npy_intp len, npy_int8 *res, const npy_ assert(lo < hi); buf = (int *)mkl_malloc(len * sizeof(int), 64); - assert(buf != NULL); + assert(buf != nullptr); err = viRngUniform(VSL_RNG_METHOD_UNIFORM_STD, state->stream, len, buf, (const int)lo, (const int)hi + 1); assert(err == VSL_STATUS_OK); DIST_PRAGMA_VECTOR - for (i = 0; i < len; i++) + for (i = 0; i < len; ++i) res[i] = (npy_int8)buf[i]; mkl_free(buf); @@ -1705,8 +1718,9 @@ void irk_rand_int8_vec(irk_state *state, npy_intp len, npy_int8 *res, const npy_ void irk_rand_uint16_vec(irk_state *state, npy_intp len, npy_uint16 *res, const npy_uint16 lo, const npy_uint16 hi) { - int err, i; - int *buf = NULL; + int err = 0; + npy_intp i = 0; + int *buf = nullptr; if (len < 1) return; @@ -1722,7 +1736,7 @@ void irk_rand_uint16_vec(irk_state *state, npy_intp len, npy_uint16 *res, const if (lo == hi) { DIST_PRAGMA_VECTOR - for (i = 0; i < len; i++) + for (i = 0; i < len; ++i) res[i] = lo; return; @@ -1730,13 +1744,13 @@ void irk_rand_uint16_vec(irk_state *state, npy_intp len, npy_uint16 *res, const assert(lo < hi); buf = (int *)mkl_malloc(len * sizeof(int), 64); - assert(buf != NULL); + assert(buf != nullptr); err = viRngUniform(VSL_RNG_METHOD_UNIFORM_STD, state->stream, len, buf, (const int)lo, (const int)hi + 1); assert(err == VSL_STATUS_OK); DIST_PRAGMA_VECTOR - for (i = 0; i < len; i++) + for (i = 0; i < len; ++i) res[i] = (npy_uint16)buf[i]; mkl_free(buf); @@ -1744,8 +1758,9 @@ void irk_rand_uint16_vec(irk_state *state, npy_intp len, npy_uint16 *res, const void irk_rand_int16_vec(irk_state *state, npy_intp len, npy_int16 *res, const npy_int16 lo, const npy_int16 hi) { - int err, i; - int *buf = NULL; + int err = 0; + npy_intp i = 0; + int *buf = nullptr; if (len < 1) return; @@ -1761,7 +1776,7 @@ void irk_rand_int16_vec(irk_state *state, npy_intp len, npy_int16 *res, const np if (lo == hi) { DIST_PRAGMA_VECTOR - for (i = 0; i < len; i++) + for (i = 0; i < len; ++i) res[i] = lo; return; @@ -1769,13 +1784,13 @@ void irk_rand_int16_vec(irk_state *state, npy_intp len, npy_int16 *res, const np assert(lo < hi); buf = (int *)mkl_malloc(len * sizeof(int), 64); - assert(buf != NULL); + assert(buf != nullptr); err = viRngUniform(VSL_RNG_METHOD_UNIFORM_STD, state->stream, len, buf, (const int)lo, (const int)hi + 1); assert(err == VSL_STATUS_OK); DIST_PRAGMA_VECTOR - for (i = 0; i < len; i++) + for (i = 0; i < len; ++i) res[i] = (npy_int16)buf[i]; mkl_free(buf); @@ -1783,7 +1798,7 @@ void irk_rand_int16_vec(irk_state *state, npy_intp len, npy_int16 *res, const np void irk_rand_uint32_vec(irk_state *state, npy_intp len, npy_uint32 *res, const npy_uint32 lo, const npy_uint32 hi) { - int err; + int err = 0; unsigned int intm = INT_MAX; if (len < 1) @@ -1820,7 +1835,7 @@ void irk_rand_uint32_vec(irk_state *state, npy_intp len, npy_uint32 *res, const assert(err == VSL_STATUS_OK); DIST_PRAGMA_VECTOR - for (i = 0; i < len; i++) + for (i = 0; i < len; ++i) res[i] += shft; } else @@ -1832,7 +1847,7 @@ void irk_rand_uint32_vec(irk_state *state, npy_intp len, npy_uint32 *res, const void irk_rand_int32_vec(irk_state *state, npy_intp len, npy_int32 *res, const npy_int32 lo, const npy_int32 hi) { - int err; + int err = 0; int intm = INT_MAX; if (len < 1) @@ -1853,7 +1868,7 @@ void irk_rand_int32_vec(irk_state *state, npy_intp len, npy_int32 *res, const np irk_rand_uint32_vec(state, len, (npy_uint32 *)res, 0U, (npy_uint32)(hi - lo)); DIST_PRAGMA_VECTOR - for (i = 0; i < len; i++) + for (i = 0; i < len; ++i) res[i] += lo; } else @@ -1866,7 +1881,8 @@ void irk_rand_int32_vec(irk_state *state, npy_intp len, npy_int32 *res, const np void irk_rand_uint64_vec(irk_state *state, npy_intp len, npy_uint64 *res, const npy_uint64 lo, const npy_uint64 hi) { npy_uint64 rng; - int i, err; + int err = 0; + npy_intp i = 0; if (len < 1) return; @@ -1892,7 +1908,7 @@ void irk_rand_uint64_vec(irk_state *state, npy_intp len, npy_uint64 *res, const if (!rng) { DIST_PRAGMA_VECTOR - for (i = 0; i < len; i++) + for (i = 0; i < len; ++i) res[i] = lo; return; @@ -1903,13 +1919,13 @@ void irk_rand_uint64_vec(irk_state *state, npy_intp len, npy_uint64 *res, const if (rng <= (npy_uint64)INT_MAX) { int *buf = (int *)mkl_malloc(len * sizeof(int), 64); - assert(buf != NULL); + assert(buf != nullptr); err = viRngUniform(VSL_RNG_METHOD_UNIFORM_STD, state->stream, len, buf, 0, (const int)rng); assert(err == VSL_STATUS_OK); DIST_PRAGMA_VECTOR - for (i = 0; i < len; i++) + for (i = 0; i < len; ++i) res[i] = lo + ((npy_uint64)buf[i]); mkl_free(buf); @@ -1917,8 +1933,8 @@ void irk_rand_uint64_vec(irk_state *state, npy_intp len, npy_uint64 *res, const else { npy_uint64 mask = rng; - npy_uint64 *buf = NULL; - int n_accepted = 0; + npy_uint64 *buf = nullptr; + npy_intp n_accepted = 0; mask |= mask >> 1; mask |= mask >> 2; @@ -1928,16 +1944,17 @@ void irk_rand_uint64_vec(irk_state *state, npy_intp len, npy_uint64 *res, const mask |= mask >> 32; buf = (npy_uint64 *)mkl_malloc(len * sizeof(npy_uint64), 64); - assert(buf != NULL); + assert(buf != nullptr); while (n_accepted < len) { - int k, batchSize = len - n_accepted; + npy_intp k = 0; + npy_intp batchSize = len - n_accepted; err = viRngUniformBits64(VSL_RNG_METHOD_UNIFORM_STD, state->stream, batchSize, (unsigned MKL_INT64 *)buf); assert(err == VSL_STATUS_OK); - for (k = 0; k < batchSize; k++) + for (k = 0; k < batchSize; ++k) { npy_uint64 value = buf[k] & mask; if (value <= rng) @@ -1953,8 +1970,8 @@ void irk_rand_uint64_vec(irk_state *state, npy_intp len, npy_uint64 *res, const void irk_rand_int64_vec(irk_state *state, npy_intp len, npy_int64 *res, const npy_int64 lo, const npy_int64 hi) { - npy_uint64 rng; - npy_intp i; + npy_uint64 rng = 0; + npy_intp i = 0; if (len < 1) return; @@ -1964,7 +1981,7 @@ void irk_rand_int64_vec(irk_state *state, npy_intp len, npy_int64 *res, const np irk_rand_uint64_vec(state, len, (npy_uint64 *)res, 0, rng); DIST_PRAGMA_VECTOR - for (i = 0; i < len; i++) + for (i = 0; i < len; ++i) res[i] = res[i] + lo; } @@ -1976,7 +1993,7 @@ const MKL_INT cholesky_storage_flags[3] = { void irk_multinormal_vec_ICDF(irk_state *state, npy_intp len, double *res, const int dim, double *mean_vec, double *ch, const ch_st_enum storage_flag) { - int err; + int err = 0; const MKL_INT storage_mode = cholesky_storage_flags[storage_flag]; err = vdRngGaussianMV(VSL_RNG_METHOD_GAUSSIANMV_ICDF, state->stream, len, res, dim, storage_mode, mean_vec, ch); @@ -1986,7 +2003,7 @@ void irk_multinormal_vec_ICDF(irk_state *state, npy_intp len, double *res, const void irk_multinormal_vec_BM1(irk_state *state, npy_intp len, double *res, const int dim, double *mean_vec, double *ch, const ch_st_enum storage_flag) { - int err; + int err = 0; const MKL_INT storage_mode = cholesky_storage_flags[storage_flag]; if (len < 1) @@ -2008,7 +2025,7 @@ void irk_multinormal_vec_BM1(irk_state *state, npy_intp len, double *res, const void irk_multinormal_vec_BM2(irk_state *state, npy_intp len, double *res, const int dim, double *mean_vec, double *ch, const ch_st_enum storage_flag) { - int err; + int err = 0; const MKL_INT storage_mode = cholesky_storage_flags[storage_flag]; if (len < 1) diff --git a/mkl_random/tests/test_random.py b/mkl_random/tests/test_random.py index 7315643..3f9babd 100644 --- a/mkl_random/tests/test_random.py +++ b/mkl_random/tests/test_random.py @@ -24,931 +24,1022 @@ # OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE # OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. -from __future__ import division, absolute_import, print_function +from typing import NamedTuple import numpy as np -from unittest import TestCase +import mkl_random as rnd from numpy.testing import ( assert_, assert_raises, assert_equal, assert_warns, suppress_warnings) -import mkl_random as rnd from numpy.compat import asbytes import sys import warnings +import pytest -class TestSeed_Intel(TestCase): - def test_scalar(self): - evs_zero_seed = { - 'MT19937' : 844, 'SFMT19937' : 857, - 'WH' : 0, 'MT2203' : 890, - 'MCG31' : 0, 'R250' : 229, - 'MRG32K3A' : 0, 'MCG59' : 0 } - for brng_algo in evs_zero_seed: - s = rnd.RandomState(0, brng = brng_algo) - assert_equal(s.get_state()[0], brng_algo) - assert_equal(s.randint(1000), evs_zero_seed[brng_algo]) - evs_max_seed = { - 'MT19937' : 635, 'SFMT19937' : 25, - 'WH' : 100, 'MT2203' : 527, - 'MCG31' : 0, 'R250' : 229, - 'MRG32K3A' : 961, 'MCG59' : 0 } - for brng_algo in evs_max_seed: - s = rnd.RandomState(4294967295, brng = brng_algo) - assert_equal(s.get_state()[0], brng_algo) - assert_equal(s.randint(1000), evs_max_seed[brng_algo]) - - def test_array(self): - s = rnd.RandomState(range(10), brng='MT19937') - assert_equal(s.randint(1000), 410) - s = rnd.RandomState(np.arange(10), brng='MT19937') - assert_equal(s.randint(1000), 410) - s = rnd.RandomState([0], brng='MT19937') - assert_equal(s.randint(1000), 844) - s = rnd.RandomState([4294967295], brng='MT19937') - assert_equal(s.randint(1000), 635) - - def test_invalid_scalar(self): - # seed must be an unsigned 32 bit integers - assert_raises(TypeError, rnd.RandomState, -0.5) - assert_raises(ValueError, rnd.RandomState, -1) - - def test_invalid_array(self): - # seed must be an unsigned 32 bit integers - assert_raises(TypeError, rnd.RandomState, [-0.5]) - assert_raises(ValueError, rnd.RandomState, [-1]) - assert_raises(ValueError, rnd.RandomState, [4294967296]) - assert_raises(ValueError, rnd.RandomState, [1, 2, 4294967296]) - assert_raises(ValueError, rnd.RandomState, [1, -2, 4294967296]) - - def test_non_deterministic(self): - rs = rnd.RandomState(brng='nondeterministic') - rs.rand(10) - rs.randint(0, 10) - -class TestBinomial_Intel(TestCase): - def test_n_zero(self): - # Tests the corner case of n == 0 for the binomial distribution. - # binomial(0, p) should be zero for any p in [0, 1]. - # This test addresses issue #3480. - zeros = np.zeros(2, dtype='int') - for p in [0, .5, 1]: - assert_(rnd.binomial(0, p) == 0) - np.testing.assert_array_equal(rnd.binomial(zeros, p), zeros) - - def test_p_is_nan(self): - # Issue #4571. - assert_raises(ValueError, rnd.binomial, 1, np.nan) - - -class TestMultinomial_Intel(TestCase): - def test_basic(self): - rnd.multinomial(100, [0.2, 0.8]) - - def test_zero_probability(self): - rnd.multinomial(100, [0.2, 0.8, 0.0, 0.0, 0.0]) - - def test_int_negative_interval(self): - assert_(-5 <= rnd.randint(-5, -1) < -1) - x = rnd.randint(-5, -1, 5) - assert_(np.all(-5 <= x)) - assert_(np.all(x < -1)) - - def test_size(self): - # gh-3173 - p = [0.5, 0.5] - assert_equal(rnd.multinomial(1, p, np.uint32(1)).shape, (1, 2)) - assert_equal(rnd.multinomial(1, p, np.uint32(1)).shape, (1, 2)) - assert_equal(rnd.multinomial(1, p, np.uint32(1)).shape, (1, 2)) - assert_equal(rnd.multinomial(1, p, [2, 2]).shape, (2, 2, 2)) - assert_equal(rnd.multinomial(1, p, (2, 2)).shape, (2, 2, 2)) - assert_equal(rnd.multinomial(1, p, np.array((2, 2))).shape, - (2, 2, 2)) - - assert_raises(TypeError, rnd.multinomial, 1, p, - np.float64(1)) - - -class TestSetState_Intel(TestCase): - def setUp(self): - self.seed = 1234567890 - self.prng = rnd.RandomState(self.seed) - self.state = self.prng.get_state() - - def test_basic(self): - old = self.prng.tomaxint(16) - self.prng.set_state(self.state) - new = self.prng.tomaxint(16) - assert_(np.all(old == new)) - - def test_gaussian_reset(self): - # Make sure the cached every-other-Gaussian is reset. - old = self.prng.standard_normal(size=3) - self.prng.set_state(self.state) - new = self.prng.standard_normal(size=3) - assert_(np.all(old == new)) - - def test_gaussian_reset_in_media_res(self): - # When the state is saved with a cached Gaussian, make sure the - # cached Gaussian is restored. - self.prng.standard_normal() - state = self.prng.get_state() - old = self.prng.standard_normal(size=3) - self.prng.set_state(state) - new = self.prng.standard_normal(size=3) - assert_(np.all(old == new)) - - def test_backwards_compatibility(self): - # Make sure we can accept old state tuples that do not have the - # cached Gaussian value. - if len(self.state) == 5: - old_state = self.state[:-2] - x1 = self.prng.standard_normal(size=16) - self.prng.set_state(old_state) - x2 = self.prng.standard_normal(size=16) - self.prng.set_state(self.state) - x3 = self.prng.standard_normal(size=16) - assert_(np.all(x1 == x2)) - assert_(np.all(x1 == x3)) - - def test_negative_binomial(self): - # Ensure that the negative binomial results take floating point - # arguments without truncation. - self.prng.negative_binomial(0.5, 0.5) - -class TestRandint_Intel(TestCase): - - rfunc = rnd.randint - - # valid integer/boolean types - itype = [np.bool_, np.int8, np.uint8, np.int16, np.uint16, - np.int32, np.uint32, np.int64, np.uint64] - - def test_unsupported_type(self): - assert_raises(TypeError, self.rfunc, 1, dtype=np.float64) - - def test_bounds_checking(self): - for dt in self.itype: - lbnd = 0 if dt is np.bool_ else np.iinfo(dt).min - ubnd = 2 if dt is np.bool_ else np.iinfo(dt).max + 1 - assert_raises(ValueError, self.rfunc, lbnd - 1, ubnd, dtype=dt) - assert_raises(ValueError, self.rfunc, lbnd, ubnd + 1, dtype=dt) - assert_raises(ValueError, self.rfunc, ubnd, lbnd, dtype=dt) - assert_raises(ValueError, self.rfunc, 1, 0, dtype=dt) - - def test_rng_zero_and_extremes(self): - for dt in self.itype: - lbnd = 0 if dt is np.bool_ else np.iinfo(dt).min - ubnd = 2 if dt is np.bool_ else np.iinfo(dt).max + 1 - tgt = ubnd - 1 - assert_equal(self.rfunc(tgt, tgt + 1, size=1000, dtype=dt), tgt) - tgt = lbnd - assert_equal(self.rfunc(tgt, tgt + 1, size=1000, dtype=dt), tgt) - tgt = (lbnd + ubnd)//2 - assert_equal(self.rfunc(tgt, tgt + 1, size=1000, dtype=dt), tgt) - - def test_in_bounds_fuzz(self): - # Don't use fixed seed - rnd.seed() - for dt in self.itype[1:]: - for ubnd in [4, 8, 16]: - vals = self.rfunc(2, ubnd, size=2**16, dtype=dt) - assert_(vals.max() < ubnd) - assert_(vals.min() >= 2) - vals = self.rfunc(0, 2, size=2**16, dtype='bool') - assert_(vals.max() < 2) - assert_(vals.min() >= 0) - - def test_repeatability(self): - import hashlib - # We use a md5 hash of generated sequences of 1000 samples - # in the range [0, 6) for all but np.bool, where the range - # is [0, 2). Hashes are for little endian numbers. - tgt = {'bool': '4fee98a6885457da67c39331a9ec336f', - 'int16': '80a5ff69c315ab6f80b03da1d570b656', - 'int32': '15a3c379b6c7b0f296b162194eab68bc', - 'int64': 'ea9875f9334c2775b00d4976b85a1458', - 'int8': '0f56333af47de94930c799806158a274', - 'uint16': '80a5ff69c315ab6f80b03da1d570b656', - 'uint32': '15a3c379b6c7b0f296b162194eab68bc', - 'uint64': 'ea9875f9334c2775b00d4976b85a1458', - 'uint8': '0f56333af47de94930c799806158a274'} - - for dt in self.itype[1:]: - rnd.seed(1234, brng='MT19937') - - # view as little endian for hash - if sys.byteorder == 'little': - val = self.rfunc(0, 6, size=1000, dtype=dt) - else: - val = self.rfunc(0, 6, size=1000, dtype=dt).byteswap() - - res = hashlib.md5(val.view(np.int8)).hexdigest() - assert_(tgt[np.dtype(dt).name] == res) - - # bools do not depend on endianess - rnd.seed(1234, brng='MT19937') - val = self.rfunc(0, 2, size=1000, dtype='bool').view(np.int8) - res = hashlib.md5(val).hexdigest() - assert_(tgt[np.dtype('bool').name] == res) - - def test_respect_dtype_singleton(self): - # See gh-7203 - for dt in self.itype: - lbnd = 0 if dt is np.bool_ else np.iinfo(dt).min - ubnd = 2 if dt is np.bool_ else np.iinfo(dt).max + 1 - - sample = self.rfunc(lbnd, ubnd, dtype=dt) - self.assertEqual(sample.dtype, np.dtype(dt)) - - for dt in (bool, int): - lbnd = 0 if dt is bool else np.iinfo(np.dtype(dt)).min - ubnd = 2 if dt is bool else np.iinfo(np.dtype(dt)).max + 1 - - # gh-7284: Ensure that we get Python data types - sample = self.rfunc(lbnd, ubnd, dtype=dt) - self.assertFalse(hasattr(sample, 'dtype')) - self.assertEqual(type(sample), dt) - - -class TestRandomDist_Intel(TestCase): - # Make sure the random distribution returns the correct value for a - # given seed. Low value of decimal argument is intended, since functional - # transformations's implementation or approximations thereof used to produce non-uniform - # random variates can vary across platforms, yet be statistically indistinguishable to the end user, - # that is no computationally feasible statistical experiment can detect the difference. - - def setUp(self): - self.seed = 1234567890 - self.brng = 'SFMT19937' - - def test_rand(self): - rnd.seed(self.seed, brng=self.brng) - actual = rnd.rand(3, 2) - desired = np.array([[0.9838694715872407, 0.019142669625580311], - [0.1767608025111258, 0.70966427633538842], - [0.518550637178123, 0.98780936631374061]]) - np.testing.assert_array_almost_equal(actual, desired, decimal=10) - - def test_randn(self): - rnd.seed(self.seed, brng=self.brng) - actual = rnd.randn(3, 2) - desired = np.array([[2.1411609928913298, -2.0717866791744819], - [-0.92778018318550248, 0.55240420724917727], - [0.04651632135517459, 2.2510674226058036]]) - np.testing.assert_array_almost_equal(actual, desired, decimal=10) - - def test_randint(self): - rnd.seed(self.seed, brng=self.brng) - actual = rnd.randint(-99, 99, size=(3, 2)) - desired = np.array([[95, -96], [-65, 41], [3, 96]]) - np.testing.assert_array_equal(actual, desired) - def test_random_integers(self): - rnd.seed(self.seed, brng=self.brng) - with suppress_warnings() as sup: - w = sup.record(DeprecationWarning) - actual = rnd.random_integers(-99, 99, size=(3, 2)) - assert_(len(w) == 1) +def test_zero_scalar_seed(): + evs_zero_seed = { + 'MT19937' : 844, 'SFMT19937' : 857, + 'WH' : 0, 'MT2203' : 890, + 'MCG31' : 0, 'R250' : 229, + 'MRG32K3A' : 0, 'MCG59' : 0 } + for brng_algo in evs_zero_seed: + s = rnd.RandomState(0, brng = brng_algo) + assert_equal(s.get_state()[0], brng_algo) + assert_equal(s.randint(1000), evs_zero_seed[brng_algo]) - desired = np.array([[96, -96], [-64, 42], [4, 97]]) - np.testing.assert_array_equal(actual, desired) +def test_max_scalar_seed(): + evs_max_seed = { + 'MT19937' : 635, 'SFMT19937' : 25, + 'WH' : 100, 'MT2203' : 527, + 'MCG31' : 0, 'R250' : 229, + 'MRG32K3A' : 961, 'MCG59' : 0 } + for brng_algo in evs_max_seed: + s = rnd.RandomState(4294967295, brng = brng_algo) + assert_equal(s.get_state()[0], brng_algo) + assert_equal(s.randint(1000), evs_max_seed[brng_algo]) - def test_random_integers_max_int(self): - # Tests whether random_integers can generate the - # maximum allowed Python int that can be converted - # into a C long. Previous implementations of this - # method have thrown an OverflowError when attempting - # to generate this integer. - with suppress_warnings() as sup: - w = sup.record(DeprecationWarning) - actual = rnd.random_integers(np.iinfo('l').max, - np.iinfo('l').max) - assert_(len(w) == 1) - desired = np.iinfo('l').max - np.testing.assert_equal(actual, desired) - - def test_random_integers_deprecated(self): - with warnings.catch_warnings(): - warnings.simplefilter("error", DeprecationWarning) - - # DeprecationWarning raised with high == None - assert_raises(DeprecationWarning, - rnd.random_integers, - np.iinfo('l').max) - - # DeprecationWarning raised with high != None - assert_raises(DeprecationWarning, - rnd.random_integers, - np.iinfo('l').max, np.iinfo('l').max) - - def test_random_sample(self): - rnd.seed(self.seed, brng=self.brng) - actual = rnd.random_sample((3, 2)) - desired = np.array([[0.9838694715872407, 0.01914266962558031], - [0.1767608025111258, 0.7096642763353884], - [0.518550637178123, 0.9878093663137406]]) - np.testing.assert_array_almost_equal(actual, desired, decimal=10) - - def test_choice_uniform_replace(self): - rnd.seed(self.seed, brng=self.brng) - actual = rnd.choice(4, 4) - desired = np.array([3, 0, 0, 2]) - np.testing.assert_array_equal(actual, desired) - def test_choice_nonuniform_replace(self): - rnd.seed(self.seed, brng=self.brng) - actual = rnd.choice(4, 4, p=[0.4, 0.4, 0.1, 0.1]) - desired = np.array([3, 0, 0, 1]) - np.testing.assert_array_equal(actual, desired) +def test_array_seed(): + s = rnd.RandomState(range(10), brng='MT19937') + assert_equal(s.randint(1000), 410) + s = rnd.RandomState(np.arange(10), brng='MT19937') + assert_equal(s.randint(1000), 410) + s = rnd.RandomState([0], brng='MT19937') + assert_equal(s.randint(1000), 844) + s = rnd.RandomState([4294967295], brng='MT19937') + assert_equal(s.randint(1000), 635) - def test_choice_uniform_noreplace(self): - rnd.seed(self.seed, brng=self.brng) - actual = rnd.choice(4, 3, replace=False) - desired = np.array([2, 1, 3]) - np.testing.assert_array_equal(actual, desired) - def test_choice_nonuniform_noreplace(self): - rnd.seed(self.seed, brng=self.brng) - actual = rnd.choice(4, 3, replace=False, - p=[0.1, 0.3, 0.5, 0.1]) - desired = np.array([3, 0, 1]) - np.testing.assert_array_equal(actual, desired) +def test_invalid_scalar_seed(): + # seed must be an unsigned 32 bit integers + pytest.raises(TypeError, rnd.RandomState, -0.5) + pytest.raises(ValueError, rnd.RandomState, -1) + + +def test_invalid_array_seed(): + # seed must be an unsigned 32 bit integers + pytest.raises(TypeError, rnd.RandomState, [-0.5]) + pytest.raises(ValueError, rnd.RandomState, [-1]) + pytest.raises(ValueError, rnd.RandomState, [4294967296]) + pytest.raises(ValueError, rnd.RandomState, [1, 2, 4294967296]) + pytest.raises(ValueError, rnd.RandomState, [1, -2, 4294967296]) + + +def test_non_deterministic_brng(): + rs = rnd.RandomState(brng='nondeterministic') + v = rs.rand(10) + assert isinstance(v, np.ndarray) + v = rs.randint(0, 10) + assert isinstance(v, int) + + +def test_binomial_n_zero(): + zeros = np.zeros(2, dtype='int') + for p in [0, .5, 1]: + assert rnd.binomial(0, p) == 0 + actual = rnd.binomial(zeros, p) + np.testing.assert_allclose(actual, zeros) + + +def test_binomial_p_is_nan(): + # Issue #4571. + pytest.raises(ValueError, rnd.binomial, 1, np.nan) + + +def test_multinomial_basic(): + rnd.multinomial(100, [0.2, 0.8]) + + +def test_multinomial_zero_probability(): + rnd.multinomial(100, [0.2, 0.8, 0.0, 0.0, 0.0]) + + +def test_multinomial_int_negative_interval(): + assert -5 <= rnd.randint(-5, -1) < -1 + x = rnd.randint(-5, -1, 5) + assert np.all(-5 <= x) + assert np.all(x < -1) + + +def test_size(): + # gh-3173 + p = [0.5, 0.5] + assert_equal(rnd.multinomial(1, p, np.uint32(1)).shape, (1, 2)) + assert_equal(rnd.multinomial(1, p, np.uint32(1)).shape, (1, 2)) + assert_equal(rnd.multinomial(1, p, np.uint32(1)).shape, (1, 2)) + assert_equal(rnd.multinomial(1, p, [2, 2]).shape, (2, 2, 2)) + assert_equal(rnd.multinomial(1, p, (2, 2)).shape, (2, 2, 2)) + assert_equal(rnd.multinomial(1, p, np.array((2, 2))).shape, + (2, 2, 2)) + + pytest.raises(TypeError, rnd.multinomial, 1, p, + np.float64(1)) + +class RngState(NamedTuple): + seed: int + prng: object + state: object + + +@pytest.fixture +def rng_state(): + seed = 1234567890 + prng = rnd.RandomState(seed) + state = prng.get_state() + return RngState(seed, prng, state) + + +def test_set_state_basic(rng_state): + sample_ref = rng_state.prng.tomaxint(16) + new_rng = rnd.RandomState() + new_rng.set_state(rng_state.state) + sample_from_new = new_rng.tomaxint(16) + assert_equal(sample_ref, sample_from_new) + + +def test_set_state_gaussian_reset(rng_state): + # Make sure the cached every-other-Gaussian is reset. + sample_ref = rng_state.prng.standard_normal(size=3) + new_rng = rnd.RandomState() + new_rng.set_state(rng_state.state) + sample_from_new = new_rng.standard_normal(size=3) + assert_equal(sample_ref, sample_from_new) + + +def test_set_state_gaussian_reset_in_media_res(rng_state): + # When the state is saved with a cached Gaussian, make sure the + # cached Gaussian is restored. + prng = rng_state.prng + _ = prng.standard_normal() + state_after_draw = prng.get_state() + sample_ref = prng.standard_normal(size=3) + new_rng = rnd.RandomState() + new_rng.set_state(state_after_draw) + sample_from_new = new_rng.standard_normal(size=3) + assert_equal(sample_ref, sample_from_new) + + +def test_set_state_backward_compatibility(rng_state): + # Make sure we can accept old state tuples that do not have the + # cached Gaussian value. + if len(rng_state.state) == 5: + state_old_format = rng_state.state[:-2] + x1 = rng_state.prng.standard_normal(size=16) + new_rng = rnd.RandomState() + new_rng.set_state(state_old_format) + x2 = new_rng.standard_normal(size=16) + new_rng.set_state(rng_state.state) + x3 = new_rng.standard_normal(size=16) + assert_equal(x1, x2) + assert_equal(x1, x3) + + +def test_set_state_negative_binomial(rng_state): + # Ensure that the negative binomial results take floating point + # arguments without truncation. + v = rng_state.prng.negative_binomial(0.5, 0.5) + assert isinstance(v, int) + + +class RandIntData(NamedTuple): + rfunc : object + itype : list + + +@pytest.fixture +def randint(): + rfunc_method = rnd.randint + integral_dtypes = [ + np.bool_, np.int8, np.uint8, np.int16, np.uint16, + np.int32, np.uint32, np.int64, np.uint64 + ] + return RandIntData(rfunc_method, integral_dtypes) + + +def test_randint_unsupported_type(randint): + pytest.raises(TypeError, randint.rfunc, 1, dtype=np.float64) + + +def test_randint_bounds_checking(randint): + for dt in randint.itype: + lbnd = 0 if dt is np.bool_ else np.iinfo(dt).min + ubnd = 2 if dt is np.bool_ else np.iinfo(dt).max + 1 + pytest.raises(ValueError, randint.rfunc, lbnd - 1, ubnd, dtype=dt) + pytest.raises(ValueError, randint.rfunc, lbnd, ubnd + 1, dtype=dt) + pytest.raises(ValueError, randint.rfunc, ubnd, lbnd, dtype=dt) + pytest.raises(ValueError, randint.rfunc, 1, 0, dtype=dt) + + +def test_randint_rng_zero_and_extremes(randint): + for dt in randint.itype: + lbnd = 0 if dt is np.bool_ else np.iinfo(dt).min + ubnd = 2 if dt is np.bool_ else np.iinfo(dt).max + 1 + tgt = ubnd - 1 + assert_equal(randint.rfunc(tgt, tgt + 1, size=1000, dtype=dt), tgt) + tgt = lbnd + assert_equal(randint.rfunc(tgt, tgt + 1, size=1000, dtype=dt), tgt) + tgt = lbnd + ((ubnd - lbnd)//2) + assert_equal(randint.rfunc(tgt, tgt + 1, size=1000, dtype=dt), tgt) + + +def test_randint_in_bounds_fuzz(randint): + # Don't use fixed seed + rnd.seed() + for dt in randint.itype[1:]: + for ubnd in [4, 8, 16]: + vals = randint.rfunc(2, ubnd, size=2**16, dtype=dt) + assert_(vals.max() < ubnd) + assert_(vals.min() >= 2) + vals = randint.rfunc(0, 2, size=2**16, dtype='bool') + assert (vals.max() < 2) + assert (vals.min() >= 0) + + +def test_randint_repeatability(randint): + import hashlib + # We use a md5 hash of generated sequences of 1000 samples + # in the range [0, 6) for all but np.bool, where the range + # is [0, 2). Hashes are for little endian numbers. + tgt = {'bool': '4fee98a6885457da67c39331a9ec336f', + 'int16': '80a5ff69c315ab6f80b03da1d570b656', + 'int32': '15a3c379b6c7b0f296b162194eab68bc', + 'int64': 'ea9875f9334c2775b00d4976b85a1458', + 'int8': '0f56333af47de94930c799806158a274', + 'uint16': '80a5ff69c315ab6f80b03da1d570b656', + 'uint32': '15a3c379b6c7b0f296b162194eab68bc', + 'uint64': 'ea9875f9334c2775b00d4976b85a1458', + 'uint8': '0f56333af47de94930c799806158a274'} + + for dt in randint.itype[1:]: + rnd.seed(1234, brng='MT19937') - def test_choice_noninteger(self): - rnd.seed(self.seed, brng=self.brng) - actual = rnd.choice(['a', 'b', 'c', 'd'], 4) - desired = np.array(['d', 'a', 'a', 'c']) - np.testing.assert_array_equal(actual, desired) + # view as little endian for hash + if sys.byteorder == 'little': + val = randint.rfunc(0, 6, size=1000, dtype=dt) + else: + val = randint.rfunc(0, 6, size=1000, dtype=dt).byteswap() - def test_choice_exceptions(self): - sample = rnd.choice - assert_raises(ValueError, sample, -1, 3) - assert_raises(ValueError, sample, 3., 3) - assert_raises(ValueError, sample, [[1, 2], [3, 4]], 3) - assert_raises(ValueError, sample, [], 3) - assert_raises(ValueError, sample, [1, 2, 3, 4], 3, - p=[[0.25, 0.25], [0.25, 0.25]]) - assert_raises(ValueError, sample, [1, 2], 3, p=[0.4, 0.4, 0.2]) - assert_raises(ValueError, sample, [1, 2], 3, p=[1.1, -0.1]) - assert_raises(ValueError, sample, [1, 2], 3, p=[0.4, 0.4]) - assert_raises(ValueError, sample, [1, 2, 3], 4, replace=False) - assert_raises(ValueError, sample, [1, 2, 3], 2, replace=False, - p=[1, 0, 0]) - - def test_choice_return_shape(self): - p = [0.1, 0.9] - # Check scalar - assert_(np.isscalar(rnd.choice(2, replace=True))) - assert_(np.isscalar(rnd.choice(2, replace=False))) - assert_(np.isscalar(rnd.choice(2, replace=True, p=p))) - assert_(np.isscalar(rnd.choice(2, replace=False, p=p))) - assert_(np.isscalar(rnd.choice([1, 2], replace=True))) - assert_(rnd.choice([None], replace=True) is None) - a = np.array([1, 2]) - arr = np.empty(1, dtype=object) - arr[0] = a - assert_(rnd.choice(arr, replace=True) is a) - - # Check 0-d array - s = tuple() - assert_(not np.isscalar(rnd.choice(2, s, replace=True))) - assert_(not np.isscalar(rnd.choice(2, s, replace=False))) - assert_(not np.isscalar(rnd.choice(2, s, replace=True, p=p))) - assert_(not np.isscalar(rnd.choice(2, s, replace=False, p=p))) - assert_(not np.isscalar(rnd.choice([1, 2], s, replace=True))) - assert_(rnd.choice([None], s, replace=True).ndim == 0) - a = np.array([1, 2]) - arr = np.empty(1, dtype=object) - arr[0] = a - assert_(rnd.choice(arr, s, replace=True).item() is a) - - # Check multi dimensional array - s = (2, 3) - p = [0.1, 0.1, 0.1, 0.1, 0.4, 0.2] - assert_(rnd.choice(6, s, replace=True).shape, s) - assert_(rnd.choice(6, s, replace=False).shape, s) - assert_(rnd.choice(6, s, replace=True, p=p).shape, s) - assert_(rnd.choice(6, s, replace=False, p=p).shape, s) - assert_(rnd.choice(np.arange(6), s, replace=True).shape, s) - - def test_bytes(self): - rnd.seed(self.seed, brng=self.brng) - actual = rnd.bytes(10) - desired = asbytes('\xa4\xde\xde{\xb4\x88\xe6\x84*2') - np.testing.assert_equal(actual, desired) - - def test_shuffle(self): - # Test lists, arrays (of various dtypes), and multidimensional versions - # of both, c-contiguous or not: - for conv in [lambda x: np.array([]), - lambda x: x, - lambda x: np.asarray(x).astype(np.int8), - lambda x: np.asarray(x).astype(np.float32), - lambda x: np.asarray(x).astype(np.complex64), - lambda x: np.asarray(x).astype(object), - lambda x: [(i, i) for i in x], - lambda x: np.asarray([[i, i] for i in x]), - lambda x: np.vstack([x, x]).T, - # gh-4270 - lambda x: np.asarray([(i, i) for i in x], - [("a", object, (1,)), - ("b", np.int32, (1,))])]: - rnd.seed(self.seed, brng=self.brng) - alist = conv([1, 2, 3, 4, 5, 6, 7, 8, 9, 0]) - rnd.shuffle(alist) - actual = alist - desired = conv([9, 8, 5, 1, 6, 4, 7, 2, 3, 0]) - np.testing.assert_array_equal(actual, desired) - - - def test_shuffle_masked(self): - # gh-3263 - a = np.ma.masked_values(np.reshape(range(20), (5,4)) % 3 - 1, -1) - b = np.ma.masked_values(np.arange(20) % 3 - 1, -1) - a_orig = a.copy() - b_orig = b.copy() - for i in range(50): - rnd.shuffle(a) - assert_equal( - sorted(a.data[~a.mask]), sorted(a_orig.data[~a_orig.mask])) - rnd.shuffle(b) - assert_equal( - sorted(b.data[~b.mask]), sorted(b_orig.data[~b_orig.mask])) - - - def test_beta(self): - rnd.seed(self.seed, brng=self.brng) - actual = rnd.beta(.1, .9, size=(3, 2)) - desired = np.array( - [[0.9856952034381025, 4.35869375658114e-08], - [0.0014230232791189966, 1.4981856288121975e-06], - [1.426135763875603e-06, 4.5801786040477326e-07]]) - np.testing.assert_array_almost_equal(actual, desired, decimal=10) - - - def test_binomial(self): - rnd.seed(self.seed, brng=self.brng) - actual = rnd.binomial(100.123, .456, size=(3, 2)) - desired = np.array([[43, 48], [55, 48], [46, 53]]) - np.testing.assert_array_equal(actual, desired) + res = hashlib.md5(val.view(np.int8)).hexdigest() + assert tgt[np.dtype(dt).name] == res + # bools do not depend on endianess + rnd.seed(1234, brng='MT19937') + val = randint.rfunc(0, 2, size=1000, dtype='bool').view(np.int8) + res = hashlib.md5(val).hexdigest() + assert (tgt[np.dtype('bool').name] == res) - def test_chisquare(self): - rnd.seed(self.seed, brng=self.brng) - actual = rnd.chisquare(50, size=(3, 2)) - desired = np.array([[50.955833609920589, 50.133178918244099], - [61.513615847062013, 50.757127871422448], - [52.79816819717081, 49.973023331993552]]) - np.testing.assert_array_almost_equal(actual, desired, decimal=7) - - - def test_dirichlet(self): - rnd.seed(self.seed, brng=self.brng) - alpha = np.array([51.72840233779265162, 39.74494232180943953]) - actual = rnd.dirichlet(alpha, size=(3, 2)) - desired = np.array([[[0.6332947001908874, 0.36670529980911254], - [0.5376828907571894, 0.4623171092428107]], - [[0.6835615930093024, 0.3164384069906976], - [0.5452378139016114, 0.45476218609838875]], - [[0.6498494402738553, 0.3501505597261446], - [0.5622024400324822, 0.43779755996751785]]]) - np.testing.assert_array_almost_equal(actual, desired, decimal=10) - - - def test_dirichlet_size(self): - # gh-3173 - p = np.array([51.72840233779265162, 39.74494232180943953]) - assert_equal(rnd.dirichlet(p, np.uint32(1)).shape, (1, 2)) - assert_equal(rnd.dirichlet(p, np.uint32(1)).shape, (1, 2)) - assert_equal(rnd.dirichlet(p, np.uint32(1)).shape, (1, 2)) - assert_equal(rnd.dirichlet(p, [2, 2]).shape, (2, 2, 2)) - assert_equal(rnd.dirichlet(p, (2, 2)).shape, (2, 2, 2)) - assert_equal(rnd.dirichlet(p, np.array((2, 2))).shape, (2, 2, 2)) - - assert_raises(TypeError, rnd.dirichlet, p, np.float64(1)) - - - def test_exponential(self): - rnd.seed(self.seed, brng=self.brng) - actual = rnd.exponential(1.1234, size=(3, 2)) - desired = np.array([[0.01826877748252199, 4.4439855151117005], - [1.9468048583654507, 0.38528493864979607], - [0.7377565464231758, 0.013779117663987912]]) - np.testing.assert_array_almost_equal(actual, desired, decimal=10) - - - def test_f(self): - rnd.seed(self.seed, brng=self.brng) - actual = rnd.f(12, 77, size=(3, 2)) - desired = np.array([[1.325076177478387, 0.8670927327120197], - [2.1190792007836827, 0.9095296301824258], - [1.4953697422236187, 0.9547125618834837]]) - np.testing.assert_array_almost_equal(actual, desired, decimal=9) - - - def test_gamma(self): - rnd.seed(self.seed, brng=self.brng) - actual = rnd.gamma(5, 3, size=(3, 2)) - desired = np.array([[15.073510060334929, 14.525495858042685], - [22.73897210140115, 14.94044782480266], - [16.327929995271095, 14.419692564592896]]) - np.testing.assert_array_almost_equal(actual, desired, decimal=7) - - - def test_geometric(self): - rnd.seed(self.seed, brng=self.brng) - actual = rnd.geometric(.123456789, size=(3, 2)) - desired = np.array([[0, 30], [13, 2], [4, 0]]) - np.testing.assert_array_equal(actual, desired) +def test_randint_respect_dtype_singleton(randint): + # See gh-7203 + for dt in randint.itype: + lbnd = 0 if dt is np.bool_ else np.iinfo(dt).min + ubnd = 2 if dt is np.bool_ else np.iinfo(dt).max + 1 - def test_gumbel(self): - rnd.seed(self.seed, self.brng) - actual = rnd.gumbel(loc=.123456789, scale=2.0, size=(3, 2)) - desired = np.array([[-8.114386462751979, 2.873840411460178], - [1.2231161758452016, -2.0168070493213532], - [-0.7175455966332102, -8.678464904504784]]) - np.testing.assert_array_almost_equal(actual, desired, decimal=7) + sample = randint.rfunc(lbnd, ubnd, dtype=dt) + assert_equal(sample.dtype, np.dtype(dt)) + for dt in (bool, int): + lbnd = 0 if dt is bool else np.iinfo(np.dtype(dt)).min + ubnd = 2 if dt is bool else np.iinfo(np.dtype(dt)).max + 1 - def test_hypergeometric(self): - rnd.seed(self.seed, brng=self.brng) - actual = rnd.hypergeometric(10.1, 5.5, 14, size=(3, 2)) - desired = np.array([[10, 9], [9, 10], [9, 10]]) - np.testing.assert_array_equal(actual, desired) + # gh-7284: Ensure that we get Python data types + sample = randint.rfunc(lbnd, ubnd, dtype=dt) + assert not hasattr(sample, 'dtype') + assert (type(sample) == dt) - # Test nbad = 0 - actual = rnd.hypergeometric(5, 0, 3, size=4) - desired = np.array([3, 3, 3, 3]) - np.testing.assert_array_equal(actual, desired) - actual = rnd.hypergeometric(15, 0, 12, size=4) - desired = np.array([12, 12, 12, 12]) - np.testing.assert_array_equal(actual, desired) +class RandomDistData(NamedTuple): + seed : int + brng : str - # Test ngood = 0 - actual = rnd.hypergeometric(0, 5, 3, size=4) - desired = np.array([0, 0, 0, 0]) - np.testing.assert_array_equal(actual, desired) - actual = rnd.hypergeometric(0, 15, 12, size=4) - desired = np.array([0, 0, 0, 0]) - np.testing.assert_array_equal(actual, desired) +@pytest.fixture +def randomdist(): + return RandomDistData(seed=1234567890, brng='SFMT19937') - def test_laplace(self): - rnd.seed(self.seed, brng=self.brng) - actual = rnd.laplace(loc=.123456789, scale=2.0, size=(3, 2)) - desired = np.array([[0.15598087210935016, -3.3424589282252994], - [-1.189978401356375, 3.0607925598732253], - [0.0030946589024587745, 3.14795824463997]]) - np.testing.assert_array_almost_equal(actual, desired, decimal=10) - - - def test_logistic(self): - rnd.seed(self.seed, brng=self.brng) - actual = rnd.logistic(loc=.123456789, scale=2.0, size=(3, 2)) - desired = np.array([[8.345015961402696, -7.749557532940552], - [-2.9534419690278444, 1.910964962531448], - [0.2719300361499433, 8.913100396613983]]) - np.testing.assert_array_almost_equal(actual, desired, decimal=10) - - - def test_lognormal(self): - rnd.seed(self.seed, brng=self.brng) - actual = rnd.lognormal(mean=.123456789, sigma=2.0, size=(3, 2)) - desired = np.array([[81.92291750917155, 0.01795087229603931], - [0.1769118704670423, 3.415299544410577], - [1.2417099625339398, 102.0631392685238]]) - np.testing.assert_array_almost_equal(actual, desired, decimal=6) - actual = rnd.lognormal(mean=.123456789, sigma=2.0, size=(3,2), - method='Box-Muller2') - desired = np.array([[0.2585388231094821, 0.43734953048924663], - [26.050836228611697, 26.76266237820882], - [0.24216420175675096, 0.2481945765083541]]) - np.testing.assert_array_almost_equal(actual, desired, decimal=7) - - - def test_logseries(self): - rnd.seed(self.seed, brng=self.brng) - actual = rnd.logseries(p=.923456789, size=(3, 2)) - desired = np.array([[18, 1], [1, 1], [5, 19]]) - np.testing.assert_array_equal(actual, desired) +# Make sure the random distribution returns the correct value for a +# given seed. Low value of decimal argument is intended, since functional +# transformations's implementation or approximations thereof used to produce non-uniform +# random variates can vary across platforms, yet be statistically indistinguishable to the end user, +# that is no computationally feasible statistical experiment can detect the difference. +def test_randomdist_rand(randomdist): + rnd.seed(randomdist.seed, brng=randomdist.brng) + actual = rnd.rand(3, 2) + desired = np.array([[0.9838694715872407, 0.019142669625580311], + [0.1767608025111258, 0.70966427633538842], + [0.518550637178123, 0.98780936631374061]]) + np.testing.assert_allclose(actual, desired, atol=1e-10, rtol=1e-10) - def test_multinomial(self): - rs = rnd.RandomState(self.seed, brng=self.brng) - actual = rs.multinomial(20, [1/6.]*6, size=(3, 2)) - desired = np.full((3, 2), 20, dtype=actual.dtype) - np.testing.assert_array_equal(actual.sum(axis=-1), desired) - expected = np.array([ - [[6, 2, 1, 3, 2, 6], [7, 5, 1, 2, 3, 2]], - [[5, 1, 8, 3, 2, 1], [4, 6, 0, 4, 4, 2]], - [[6, 3, 1, 4, 4, 2], [3, 2, 4, 2, 1, 8]]], actual.dtype) - np.testing.assert_array_equal(actual, expected) - - def test_multivariate_normal(self): - rnd.seed(self.seed, brng=self.brng) - mean = (.123456789, 10) - # Hmm... not even symmetric. - cov = [[1, 0], [1, 0]] - size = (3, 2) - actual = rnd.multivariate_normal(mean, cov, size) - desired = np.array([[[-2.42282709811266, 10.0], - [1.2267795840027274, 10.0]], - [[0.06813924868067336, 10.0], - [1.001190462507746, 10.0]], - [[-1.74157261455869, 10.0], - [1.0400952859037553, 10.0]]]) - np.testing.assert_array_almost_equal(actual, desired, decimal=10) - - # Check for default size, was raising deprecation warning - actual = rnd.multivariate_normal(mean, cov) - desired = np.array([1.0579899448949994, 10.0]) - np.testing.assert_array_almost_equal(actual, desired, decimal=10) - - # Check that non positive-semidefinite covariance raises warning - mean = [0, 0] - cov = [[1, 1 + 1e-10], [1 + 1e-10, 1]] - assert_warns(RuntimeWarning, rnd.multivariate_normal, mean, cov) - - - def test_multinormal_cholesky(self): - rnd.seed(self.seed, brng=self.brng) - mean = (.123456789, 10) - # lower-triangular cholesky matrix - chol_mat = [[1, 0], [-0.5, 1]] - size = (3, 2) - actual = rnd.multinormal_cholesky(mean, chol_mat, size, method='ICDF') - desired = np.array([[[2.26461778189133, 6.857632824379853], - [-0.8043233941855025, 11.01629429884193]], - [[0.1699731103551746, 12.227809261928217], - [-0.6146263106001378, 9.893801873973892]], - [[1.691753328795276, 10.797627196240155], - [-0.647341237129921, 9.626899489691816]]]) - np.testing.assert_array_almost_equal(actual, desired, decimal=10) - - - def test_negative_binomial(self): - rnd.seed(self.seed, brng=self.brng) - actual = rnd.negative_binomial(n=100, p=.12345, size=(3, 2)) - desired = np.array([[667, 679], [677, 676], [779, 648]]) - np.testing.assert_array_equal(actual, desired) +def test_randomdist_randn(randomdist): + rnd.seed(randomdist.seed, brng=randomdist.brng) + actual = rnd.randn(3, 2) + desired = np.array([[2.1411609928913298, -2.0717866791744819], + [-0.92778018318550248, 0.55240420724917727], + [0.04651632135517459, 2.2510674226058036]]) + np.testing.assert_allclose(actual, desired, atol=1e-10) - def test_noncentral_chisquare(self): - rnd.seed(self.seed, brng=self.brng) - actual = rnd.noncentral_chisquare(df=5, nonc=5, size=(3, 2)) - desired = np.array([[5.871334619375055, 8.756238913383225], - [17.29576535176833, 3.9028417087862177], - [5.1315133729432505, 9.942717979531027]]) - np.testing.assert_array_almost_equal(actual, desired, decimal=7) - - actual = rnd.noncentral_chisquare(df=.5, nonc=.2, size=(3, 2)) - desired = np.array([[0.0008971007339949436, 0.08948578998156566], - [0.6721835871997511, 2.8892645287699352], - [5.0858149962761007e-05, 1.7315797643658821]]) - np.testing.assert_array_almost_equal(actual, desired, decimal=7) - - - def test_noncentral_f(self): - rnd.seed(self.seed, brng=self.brng) - actual = rnd.noncentral_f(dfnum=5, dfden=2, nonc=1, - size=(3, 2)) - desired = np.array([[0.2216297348371284, 0.7632696724492449], - [98.67664232828238, 0.9500319825372799], - [0.3489618249246971, 1.5035633972571092]]) - np.testing.assert_array_almost_equal(actual, desired, decimal=7) - - def test_normal(self): - rnd.seed(self.seed, brng=self.brng) - actual = rnd.normal(loc=.123456789, scale=2.0, size=(3, 2)) - desired = np.array([[4.405778774782659, -4.020116569348963], - [-1.732103577371005, 1.2282652034983546], - [0.21648943171034918, 4.625591634211608]]) - np.testing.assert_array_almost_equal(actual, desired, decimal=7) - - rnd.seed(self.seed, brng=self.brng) - actual = rnd.normal(loc=.123456789, scale=2.0, size=(3, 2), method="BoxMuller") - desired = np.array([[0.16673479781277187, -3.4809986872165952], - [-0.05193761082535492, 3.249201213154922], - [-0.11915582299214138, 3.555636100927892]]) - np.testing.assert_array_almost_equal(actual, desired, decimal=8) - - rnd.seed(self.seed, brng=self.brng) - actual = rnd.normal(loc=.123456789, scale=2.0, size=(3, 2), method="BoxMuller2") - desired = np.array([[0.16673479781277187, 0.48153966449249175], - [-3.4809986872165952, -0.8101190082826486], - [-0.051937610825354905, 2.4088402362484342]]) - np.testing.assert_array_almost_equal(actual, desired, decimal=7) - - - def test_pareto(self): - rnd.seed(self.seed, brng=self.brng) - actual = rnd.pareto(a=.123456789, size=(3, 2)) - desired = np.array( - [[0.14079174875385214, 82372044085468.92], - [1247881.6368437486, 15.086855668610944], - [203.2638558933401, 0.10445383654349749]]) - # For some reason on 32-bit x86 Ubuntu 12.10 the [1, 0] entry in this - # matrix differs by 24 nulps. Discussion: - # http://mail.scipy.org/pipermail/numpy-discussion/2012-September/063801.html - # Consensus is that this is probably some gcc quirk that affects - # rounding but not in any important way, so we just use a looser - # tolerance on this test: - np.testing.assert_array_almost_equal_nulp(actual, desired, nulp=30) - - def test_poisson(self): - rnd.seed(self.seed, brng=self.brng) - actual = rnd.poisson(lam=.123456789, size=(3, 2)) - desired = np.array([[1, 0], [0, 0], [0, 1]]) - np.testing.assert_array_equal(actual, desired) - rnd.seed(self.seed, brng=self.brng) - actual = rnd.poisson(lam=1234.56789, size=(3, 2)) - desired = np.array([[1310, 1162], [1202, 1254], [1236, 1314]]) - np.testing.assert_array_equal(actual, desired) +def test_randomdist_randint(randomdist): + rnd.seed(randomdist.seed, brng=randomdist.brng) + actual = rnd.randint(-99, 99, size=(3, 2)) + desired = np.array([[95, -96], [-65, 41], [3, 96]]) + np.testing.assert_array_equal(actual, desired) - def test_poisson_exceptions(self): - lambig = np.iinfo('l').max - lamneg = -1 - assert_raises(ValueError, rnd.poisson, lamneg) - assert_raises(ValueError, rnd.poisson, [lamneg]*10) - assert_raises(ValueError, rnd.poisson, lambig) - assert_raises(ValueError, rnd.poisson, [lambig]*10) - - def test_power(self): - rnd.seed(self.seed, brng=self.brng) - actual = rnd.power(a=.123456789, size=(3, 2)) - desired = np.array([[0.8765841803224415, 1.2140041091640163e-14], - [8.013574117268635e-07, 0.06216255187464781], - [0.004895628723087296, 0.9054248959192386]]) - np.testing.assert_array_almost_equal(actual, desired, decimal=10) - - def test_rayleigh(self): - rnd.seed(self.seed, brng=self.brng) - actual = rnd.rayleigh(scale=10, size=(3, 2)) - desired = np.array([[1.80344345931194, 28.127692489122378], - [18.6169699930609, 8.282068232120208], - [11.460520015934597, 1.5662406536967712]]) - np.testing.assert_array_almost_equal(actual, desired, decimal=7) - - def test_standard_cauchy(self): - rnd.seed(self.seed, brng=self.brng) - actual = rnd.standard_cauchy(size=(3, 2)) - desired = np.array([[19.716487700629912, -16.608240276131227], - [-1.6117703817332278, 0.7739915895826882], - [0.058344614106131, 26.09825325697747]]) - np.testing.assert_array_almost_equal(actual, desired, decimal=9) - - def test_standard_exponential(self): - rnd.seed(self.seed, brng=self.brng) - actual = rnd.standard_exponential(size=(3, 2)) - desired = np.array([[0.016262041554675085, 3.955835423813157], - [1.7329578586126497, 0.3429632710074738], - [0.6567175951781875, 0.012265548926462446]]) - np.testing.assert_array_almost_equal(actual, desired, decimal=10) - - def test_standard_gamma(self): - rnd.seed(self.seed, brng=self.brng) - actual = rnd.standard_gamma(shape=3, size=(3, 2)) - desired = np.array([[2.939330965027084, 2.799606052259993], - [4.988193705918075, 2.905305108691164], - [3.2630929395548147, 2.772756340265377]]) - np.testing.assert_array_almost_equal(actual, desired, decimal=7) - - def test_standard_normal(self): - rnd.seed(self.seed, brng=self.brng) - actual = rnd.standard_normal(size=(3, 2)) - desired = np.array([[2.1411609928913298, -2.071786679174482], - [-0.9277801831855025, 0.5524042072491773], - [0.04651632135517459, 2.2510674226058036]]) - np.testing.assert_array_almost_equal(actual, desired, decimal=7) - - rnd.seed(self.seed, brng=self.brng) - actual = rnd.standard_normal(size=(3, 2), method='BoxMuller2') - desired = np.array([[0.021639004406385935, 0.17904143774624587], - [-1.8022277381082976, -0.4667878986413243], - [-0.08769719991267745, 1.1426917236242171]]) - np.testing.assert_array_almost_equal(actual, desired, decimal=7) - - def test_standard_t(self): - rnd.seed(self.seed, brng=self.brng) - actual = rnd.standard_t(df=10, size=(3, 2)) - desired = np.array([[-0.783927044239963, 0.04762883516531178], - [0.7624597987725193, -1.8045540288955506], - [-1.2657694296239195, 0.307870906117017]]) - np.testing.assert_array_almost_equal(actual, desired, decimal=10) - - def test_triangular(self): - rnd.seed(self.seed, brng=self.brng) - actual = rnd.triangular(left=5.12, mode=10.23, right=20.34, - size=(3, 2)) - desired = np.array([[18.764540652669638, 6.340166306695037], - [8.827752689522429, 13.65605077739865], - [11.732872979633328, 18.970392754850423]]) - np.testing.assert_array_almost_equal(actual, desired, decimal=10) - - def test_uniform(self): - rnd.seed(self.seed, brng=self.brng) - actual = rnd.uniform(low=1.23, high=10.54, size=(3, 2)) - desired = np.array([[10.38982478047721, 1.408218254214153], - [2.8756430713785814, 7.836974412682466], - [6.057706432128325, 10.426505200380925]]) - np.testing.assert_array_almost_equal(actual, desired, decimal=10) - - def test_uniform_range_bounds(self): - fmin = np.finfo('float').min - fmax = np.finfo('float').max - - func = rnd.uniform - np.testing.assert_raises(OverflowError, func, -np.inf, 0) - np.testing.assert_raises(OverflowError, func, 0, np.inf) - # this should not throw any error, since rng can be sampled as fmin*u + fmax*(1-u) - # for 0 -np.pi) and np.all(r <= np.pi)) - - def test_hypergeometric_range(self): - # Test for ticket #921 - assert_(np.all(rnd.hypergeometric(3, 18, 11, size=10) < 4)) - assert_(np.all(rnd.hypergeometric(18, 3, 11, size=10) > 0)) - - # Test for ticket #5623 - args = [ - (2**20 - 2, 2**20 - 2, 2**20 - 2), # Check for 32-bit systems - (2 ** 30 - 1, 2 ** 30 - 2, 2 ** 30 - 1) - ] - for arg in args: - assert_(rnd.hypergeometric(*arg) > 0) - - def test_logseries_convergence(self): - # Test for ticket #923 - N = 1000 - rnd.seed(0, brng='MT19937') - rvsn = rnd.logseries(0.8, size=N) - # these two frequency counts should be close to theoretical - # numbers with this large sample - # theoretical large N result is 0.49706795 - freq = np.sum(rvsn == 1) / float(N) - msg = "Frequency was %f, should be > 0.45" % freq - assert_(freq > 0.45, msg) - # theoretical large N result is 0.19882718 - freq = np.sum(rvsn == 2) / float(N) - msg = "Frequency was %f, should be < 0.23" % freq - assert_(freq < 0.23, msg) - - def test_permutation_longs(self): - rnd.seed(1234, brng='MT19937') - a = rnd.permutation(12) - rnd.seed(1234, brng='MT19937') - b = rnd.permutation(long(12)) - assert_array_equal(a, b) - - def test_randint_range(self): - # Test for ticket #1690 - lmax = np.iinfo('l').max - lmin = np.iinfo('l').min - try: - rnd.randint(lmin, lmax) - except: - raise AssertionError - - def test_shuffle_mixed_dimension(self): - # Test for trac ticket #2074 - for t in [[1, 2, 3, None], - [(1, 1), (2, 2), (3, 3), None], - [1, (2, 2), (3, 3), None], - [(1, 1), 2, 3, None]]: - rnd.seed(12345, brng='MT2203') - shuffled = np.array(list(t), dtype=object) - rnd.shuffle(shuffled) - expected = np.array([t[0], t[2], t[1], t[3]], dtype=object) - assert_array_equal(shuffled, expected) - - def test_call_within_randomstate(self): - # Check that custom RandomState does not call into global state - m = rnd.RandomState() - res = np.array([5, 7, 5, 4, 5, 5, 6, 9, 6, 1]) - for i in range(3): - rnd.seed(i) - m.seed(4321, brng='SFMT19937') - # If m.state is not honored, the result will change - assert_array_equal(m.choice(10, size=10, p=np.ones(10)/10.), res) - - def test_multivariate_normal_size_types(self): - # Test for multivariate_normal issue with 'size' argument. - # Check that the multivariate_normal size argument can be a - # numpy integer. - rnd.multivariate_normal([0], [[0]], size=1) - rnd.multivariate_normal([0], [[0]], size=np.int_(1)) - rnd.multivariate_normal([0], [[0]], size=np.int64(1)) - - def test_beta_small_parameters(self): - # Test that beta with small a and b parameters does not produce - # NaNs due to roundoff errors causing 0 / 0, gh-5851 - rnd.seed(1234567890, brng='MT19937') - x = rnd.beta(0.0001, 0.0001, size=100) - assert_(not np.any(np.isnan(x)), 'Nans in rnd.beta') - - def test_choice_sum_of_probs_tolerance(self): - # The sum of probs should be 1.0 with some tolerance. - # For low precision dtypes the tolerance was too tight. - # See numpy github issue 6123. - rnd.seed(1234, brng='MT19937') - a = [1, 2, 3] - counts = [4, 4, 2] - for dt in np.float16, np.float32, np.float64: - probs = np.array(counts, dtype=dt) / sum(counts) - c = rnd.choice(a, p=probs) - assert_(c in a) - assert_raises(ValueError, rnd.choice, a, p=probs*0.9) - - def test_shuffle_of_array_of_different_length_strings(self): - # Test that permuting an array of different length strings - # will not cause a segfault on garbage collection - # Tests gh-7710 - rnd.seed(1234, brng='MT19937') - - a = np.array(['a', 'a' * 1000]) - - for _ in range(100): - rnd.shuffle(a) - - # Force Garbage Collection - should not segfault. - import gc - gc.collect() - - - def test_shuffle_of_array_of_objects(self): - # Test that permuting an array of objects will not cause - # a segfault on garbage collection. - # See gh-7719 - rnd.seed(1234, brng='MT19937') - a = np.array([np.arange(4), np.arange(4)]) - - for _ in range(1000): - rnd.shuffle(a) - - # Force Garbage Collection - should not segfault. - import gc - gc.collect() - - - def test_non_central_chi_squared_df_one(self): - a = rnd.noncentral_chisquare(df = 1.0, nonc=2.3, size=10**4) - assert(a.min() > 0.0) +import pytest +import gc + + +def test_VonMises_range(): + # Make sure generated random variables are in [-pi, pi]. + # Regression test for ticket #986. + for mu in np.linspace(-7., 7., 5): + r = rnd.vonmises(mu, 1, 50) + assert_(np.all(r > -np.pi) and np.all(r <= np.pi)) + + +def test_hypergeometric_range(): + # Test for ticket #921 + assert_(np.all(rnd.hypergeometric(3, 18, 11, size=10) < 4)) + assert_(np.all(rnd.hypergeometric(18, 3, 11, size=10) > 0)) + + # Test for ticket #5623 + args = [ + (2**20 - 2, 2**20 - 2, 2**20 - 2), # Check for 32-bit systems + (2 ** 30 - 1, 2 ** 30 - 2, 2 ** 30 - 1) + ] + for arg in args: + assert_(rnd.hypergeometric(*arg) > 0) + + +def test_logseries_convergence(): + # Test for ticket #923 + N = 1000 + rnd.seed(0, brng='MT19937') + rvsn = rnd.logseries(0.8, size=N) + # these two frequency counts should be close to theoretical + # numbers with this large sample + # theoretical large N result is 0.49706795 + freq = np.sum(rvsn == 1) / float(N) + msg = "Frequency was %f, should be > 0.45" % freq + assert_(freq > 0.45, msg) + # theoretical large N result is 0.19882718 + freq = np.sum(rvsn == 2) / float(N) + msg = "Frequency was %f, should be < 0.23" % freq + assert_(freq < 0.23, msg) + + +def test_permutation_longs(): + rnd.seed(1234, brng='MT19937') + a = rnd.permutation(12) + rnd.seed(1234, brng='MT19937') + b = rnd.permutation(long(12)) + assert_array_equal(a, b) + + +def test_randint_range(): + # Test for ticket #1690 + lmax = np.iinfo('l').max + lmin = np.iinfo('l').min + try: + rnd.randint(lmin, lmax) + except: + raise AssertionError + + +def test_shuffle_mixed_dimension(): + # Test for trac ticket #2074 + for t in [[1, 2, 3, None], + [(1, 1), (2, 2), (3, 3), None], + [1, (2, 2), (3, 3), None], + [(1, 1), 2, 3, None]]: + rnd.seed(12345, brng='MT2203') + shuffled = np.array(list(t), dtype=object) + rnd.shuffle(shuffled) + expected = np.array([t[0], t[2], t[1], t[3]], dtype=object) + assert_array_equal(shuffled, expected) + + +def test_call_within_randomstate(): + # Check that custom RandomState does not call into global state + m = rnd.RandomState() + res = np.array([5, 7, 5, 4, 5, 5, 6, 9, 6, 1]) + for i in range(3): + rnd.seed(i) + m.seed(4321, brng='SFMT19937') + # If m.state is not honored, the result will change + assert_array_equal(m.choice(10, size=10, p=np.ones(10)/10.), res) + + +def test_multivariate_normal_size_types(): + # Test for multivariate_normal issue with 'size' argument. + # Check that the multivariate_normal size argument can be a + # numpy integer. + rnd.multivariate_normal([0], [[0]], size=1) + rnd.multivariate_normal([0], [[0]], size=np.int_(1)) + rnd.multivariate_normal([0], [[0]], size=np.int64(1)) + + +def test_beta_small_parameters(): + # Test that beta with small a and b parameters does not produce + # NaNs due to roundoff errors causing 0 / 0, gh-5851 + rnd.seed(1234567890, brng='MT19937') + x = rnd.beta(0.0001, 0.0001, size=100) + assert_(not np.any(np.isnan(x)), 'Nans in rnd.beta') + + +def test_choice_sum_of_probs_tolerance(): + # The sum of probs should be 1.0 with some tolerance. + # For low precision dtypes the tolerance was too tight. + # See numpy github issue 6123. + rnd.seed(1234, brng='MT19937') + a = [1, 2, 3] + counts = [4, 4, 2] + for dt in np.float16, np.float32, np.float64: + probs = np.array(counts, dtype=dt) / sum(counts) + c = rnd.choice(a, p=probs) + assert_(c in a) + assert_raises(ValueError, rnd.choice, a, p=probs*0.9) + + +def test_shuffle_of_array_of_different_length_strings(): + # Test that permuting an array of different length strings + # will not cause a segfault on garbage collection + # Tests gh-7710 + rnd.seed(1234, brng='MT19937') + + a = np.array(['a', 'a' * 1000]) + + for _ in range(100): + rnd.shuffle(a) + + # Force Garbage Collection - should not segfault. + gc.collect() + + +def test_shuffle_of_array_of_objects(): + # Test that permuting an array of objects will not cause + # a segfault on garbage collection. + # See gh-7719 + rnd.seed(1234, brng='MT19937') + a = np.array([np.arange(4), np.arange(4)]) + + for _ in range(1000): + rnd.shuffle(a) + + # Force Garbage Collection - should not segfault. + gc.collect() + + +def test_non_central_chi_squared_df_one(): + a = rnd.noncentral_chisquare(df = 1.0, nonc=2.3, size=10**4) + assert(a.min() > 0.0) diff --git a/setup.py b/setup.py index e01a444..083b211 100644 --- a/setup.py +++ b/setup.py @@ -146,6 +146,14 @@ def extensions(): python_requires = '>=3.7', setup_requires=["Cython",], install_requires = ["numpy >=1.16"], + packages=[ + "mkl_random", + ], + package_data={ + "mkl_random" : [ + "tests/*.*", + ] + }, keywords=["MKL", "VSL", "true randomness", "pseudorandomness", "Philox", "MT-19937", "SFMT-19937", "MT-2203", "ARS-5", "R-250", "MCG-31",],