diff --git a/include/mbedtls/bignum.h b/include/mbedtls/bignum.h
index 073b4a4..e938232 100644
--- a/include/mbedtls/bignum.h
+++ b/include/mbedtls/bignum.h
@@ -871,6 +871,37 @@
                      int (*f_rng)(void *, unsigned char *, size_t),
                      void *p_rng );
 
+/** Generate a random number uniformly in a range.
+ *
+ * This function generates a random number between \p min inclusive and
+ * \p N exclusive.
+ *
+ * The procedure complies with RFC 6979 §3.3 (deterministic ECDSA)
+ * when the RNG is a suitably parametrized instance of HMAC_DRBG
+ * and \p min is \c 1.
+ *
+ * \note           There are `N - min` possible outputs. The lower bound
+ *                 \p min can be reached, but the upper bound \p N cannot.
+ *
+ * \param X        The destination MPI. This must point to an initialized MPI.
+ * \param min      The minimum value to return.
+ *                 It must be nonnegative.
+ * \param N        The upper bound of the range, exclusive.
+ *                 In other words, this is one plus the maximum value to return.
+ *                 \p N must be strictly larger than \p min.
+ * \param f_rng    The RNG function to use. This must not be \c NULL.
+ * \param p_rng    The RNG parameter to be passed to \p f_rng.
+ *
+ * \return         \c 0 if successful.
+ * \return         #MBEDTLS_ERR_MPI_ALLOC_FAILED if a memory allocation failed.
+ * \return         Another negative error code on failure.
+ */
+int mbedtls_mpi_random( mbedtls_mpi *X,
+                        mbedtls_mpi_sint min,
+                        const mbedtls_mpi *N,
+                        int (*f_rng)(void *, unsigned char *, size_t),
+                        void *p_rng );
+
 /**
  * \brief          Compute the greatest common divisor: G = gcd(A, B)
  *
diff --git a/library/bignum.c b/library/bignum.c
index c20c6b7..aa4b1ce 100644
--- a/library/bignum.c
+++ b/library/bignum.c
@@ -2432,6 +2432,59 @@
     return( ret );
 }
 
+int mbedtls_mpi_random( mbedtls_mpi *X,
+                        mbedtls_mpi_sint min,
+                        const mbedtls_mpi *N,
+                        int (*f_rng)(void *, unsigned char *, size_t),
+                        void *p_rng )
+{
+    /* SEC1 3.2.1: Generate X such that 1 <= n < N */
+    int ret = MBEDTLS_ERR_MPI_BAD_INPUT_DATA;
+    int count = 0;
+    unsigned cmp = 0;
+    size_t n_bits = mbedtls_mpi_bitlen( N );
+    size_t n_bytes = ( n_bits + 7 ) / 8;
+
+    /*
+     * Match the procedure given in RFC 6979 §3.3 (deterministic ECDSA)
+     * when f_rng is a suitably parametrized instance of HMAC_DRBG:
+     * - use the same byte ordering;
+     * - keep the leftmost n_bits bits of the generated octet string;
+     * - try until result is in the desired range.
+     * This also avoids any bias, which is especially important for ECDSA.
+     */
+    do
+    {
+        MBEDTLS_MPI_CHK( mbedtls_mpi_fill_random( X, n_bytes, f_rng, p_rng ) );
+        MBEDTLS_MPI_CHK( mbedtls_mpi_shift_r( X, 8 * n_bytes - n_bits ) );
+
+        /*
+         * Each try has at worst a probability 1/2 of failing (the msb has
+         * a probability 1/2 of being 0, and then the result will be < N),
+         * so after 30 tries failure probability is a most 2**(-30).
+         *
+         * For most curves, 1 try is enough with overwhelming probability,
+         * since N starts with a lot of 1s in binary, but some curves
+         * such as secp224k1 are actually very close to the worst case.
+         */
+        if( ++count > 30 )
+        {
+            ret = MBEDTLS_ERR_MPI_NOT_ACCEPTABLE;
+            goto cleanup;
+        }
+
+        ret = mbedtls_mpi_lt_mpi_ct( X, N, &cmp );
+        if( ret != 0 )
+        {
+            goto cleanup;
+        }
+    }
+    while( mbedtls_mpi_cmp_int( X, min ) < 0 || cmp != 1 );
+
+cleanup:
+    return( ret );
+}
+
 /*
  * Modular inverse: X = A^-1 mod N  (HAC 14.61 / 14.64)
  */
diff --git a/library/ecp.c b/library/ecp.c
index 98ada82..3e5cebc 100644
--- a/library/ecp.c
+++ b/library/ecp.c
@@ -3091,62 +3091,6 @@
 }
 #endif /* MBEDTLS_ECP_MONTGOMERY_ENABLED */
 
-#if defined(MBEDTLS_ECP_SHORT_WEIERSTRASS_ENABLED)
-MBEDTLS_STATIC_TESTABLE
-int mbedtls_mpi_random( mbedtls_mpi *X,
-                        mbedtls_mpi_sint min,
-                        const mbedtls_mpi *N,
-                        int (*f_rng)(void *, unsigned char *, size_t),
-                        void *p_rng )
-{
-    /* SEC1 3.2.1: Generate X such that 1 <= n < N */
-    int ret = MBEDTLS_ERR_MPI_BAD_INPUT_DATA;
-    int count = 0;
-    unsigned cmp = 0;
-    size_t n_bits = mbedtls_mpi_bitlen( N );
-    size_t n_bytes = ( n_bits + 7 ) / 8;
-
-    /*
-     * Match the procedure given in RFC 6979 §3.3 (deterministic ECDSA)
-     * when f_rng is a suitably parametrized instance of HMAC_DRBG:
-     * - use the same byte ordering;
-     * - keep the leftmost n_bits bits of the generated octet string;
-     * - try until result is in the desired range.
-     * This also avoids any bias, which is especially important for ECDSA.
-     */
-    do
-    {
-        MBEDTLS_MPI_CHK( mbedtls_mpi_fill_random( X, n_bytes, f_rng, p_rng ) );
-        MBEDTLS_MPI_CHK( mbedtls_mpi_shift_r( X, 8 * n_bytes - n_bits ) );
-
-        /*
-         * Each try has at worst a probability 1/2 of failing (the msb has
-         * a probability 1/2 of being 0, and then the result will be < N),
-         * so after 30 tries failure probability is a most 2**(-30).
-         *
-         * For most curves, 1 try is enough with overwhelming probability,
-         * since N starts with a lot of 1s in binary, but some curves
-         * such as secp224k1 are actually very close to the worst case.
-         */
-        if( ++count > 30 )
-        {
-            ret = MBEDTLS_ERR_MPI_NOT_ACCEPTABLE;
-            goto cleanup;
-        }
-
-        ret = mbedtls_mpi_lt_mpi_ct( X, N, &cmp );
-        if( ret != 0 )
-        {
-            goto cleanup;
-        }
-    }
-    while( mbedtls_mpi_cmp_int( X, min ) < 0 || cmp != 1 );
-
-cleanup:
-    return( ret );
-}
-#endif /* MBEDTLS_ECP_SHORT_WEIERSTRASS_ENABLED */
-
 /*
  * Generate a private key
  */
diff --git a/library/ecp_invasive.h b/library/ecp_invasive.h
index 1ca65fd..71c7702 100644
--- a/library/ecp_invasive.h
+++ b/library/ecp_invasive.h
@@ -76,39 +76,6 @@
 
 #endif /* MBEDTLS_ECP_MONTGOMERY_ENABLED */
 
-#if defined(MBEDTLS_ECP_SHORT_WEIERSTRASS_ENABLED)
-/** Generate a random number uniformly in a range.
- *
- * This function generates a random number between \p min inclusive and
- * \p N exclusive.
- *
- * The procedure complies with RFC 6979 §3.3 (deterministic ECDSA)
- * when the RNG is a suitably parametrized instance of HMAC_DRBG
- * and \p min is \c 1.
- *
- * \note           There are `N - min` possible outputs. The lower bound
- *                 \p min can be reached, but the upper bound \p N cannot.
- *
- * \param X        The destination MPI. This must point to an initialized MPI.
- * \param min      The minimum value to return.
- *                 It must be nonnegative.
- * \param N        The upper bound of the range, exclusive.
- *                 In other words, this is one plus the maximum value to return.
- *                 \p N must be strictly larger than \p min.
- * \param f_rng    The RNG function to use. This must not be \c NULL.
- * \param p_rng    The RNG parameter to be passed to \p f_rng.
- *
- * \return         \c 0 if successful.
- * \return         #MBEDTLS_ERR_MPI_ALLOC_FAILED if a memory allocation failed.
- * \return         Another negative error code on failure.
- */
-int mbedtls_mpi_random( mbedtls_mpi *X,
-                        mbedtls_mpi_sint min,
-                        const mbedtls_mpi *N,
-                        int (*f_rng)(void *, unsigned char *, size_t),
-                        void *p_rng );
-#endif /* MBEDTLS_ECP_SHORT_WEIERSTRASS_ENABLED */
-
 #endif /* MBEDTLS_TEST_HOOKS && MBEDTLS_ECP_C */
 
 #endif /* MBEDTLS_ECP_INVASIVE_H */
diff --git a/tests/suites/test_suite_ecp.data b/tests/suites/test_suite_ecp.data
index c007cad..5f92ca4 100644
--- a/tests/suites/test_suite_ecp.data
+++ b/tests/suites/test_suite_ecp.data
@@ -312,105 +312,6 @@
 ECP generate Montgomery key: Curve448, not enough entropy
 genkey_mx_known_answer:447:"4f0102030405060708090a0b0c0d0e0f101112131415161718191a1b1c1d1e1f202122232425262728292a2b2c2d2e2f30313233343536":""
 
-MPI random in range: 1..4
-mpi_random_many:1:"04":1000
-
-MPI random in range: 1..5
-mpi_random_many:1:"05":1000
-
-MPI random in range: 1..6
-mpi_random_many:1:"06":1000
-
-MPI random in range: 1..7
-mpi_random_many:1:"07":1000
-
-MPI random in range: 1..8
-mpi_random_many:1:"08":1000
-
-MPI random in range: 1..9
-mpi_random_many:1:"09":1000
-
-MPI random in range: 1..10
-mpi_random_many:1:"0a":1000
-
-MPI random in range: 1..11
-mpi_random_many:1:"0b":1000
-
-MPI random in range: 1..12
-mpi_random_many:1:"0c":1000
-
-MPI random in range: 1..255
-mpi_random_many:1:"ff":100
-
-MPI random in range: 1..256
-mpi_random_many:1:"0100":100
-
-MPI random in range: 1..257
-mpi_random_many:1:"0101":100
-
-MPI random in range: 1..272
-mpi_random_many:1:"0110":100
-
-MPI random in range: 1..2^64-1
-mpi_random_many:1:"ffffffffffffffff":100
-
-MPI random in range: 1..2^64
-mpi_random_many:1:"010000000000000000":100
-
-MPI random in range: 1..2^64+1
-mpi_random_many:1:"010000000000000001":100
-
-MPI random in range: 1..2^64+2^63
-mpi_random_many:1:"018000000000000000":100
-
-MPI random in range: 1..2^65-1
-mpi_random_many:1:"01ffffffffffffffff":100
-
-MPI random in range: 1..2^65
-mpi_random_many:1:"020000000000000000":100
-
-MPI random in range: 1..2^65+1
-mpi_random_many:1:"020000000000000001":100
-
-MPI random in range: 1..2^65+2^64
-mpi_random_many:1:"030000000000000000":100
-
-MPI random in range: 1..2^66+2^65
-mpi_random_many:1:"060000000000000000":100
-
-MPI random in range: 1..2^71-1
-mpi_random_many:1:"7fffffffffffffffff":100
-
-MPI random in range: 1..2^71
-mpi_random_many:1:"800000000000000000":100
-
-MPI random in range: 1..2^71+1
-mpi_random_many:1:"800000000000000001":100
-
-MPI random in range: 1..2^71+2^63
-mpi_random_many:1:"c00000000000000000":100
-
-MPI random in range: 1..2^72-1
-mpi_random_many:1:"ffffffffffffffffff":100
-
-MPI random in range: 1..2^72
-mpi_random_many:1:"01000000000000000000":100
-
-MPI random in range: 1..2^72+1
-mpi_random_many:1:"01000000000000000001":100
-
-MPI random in range: 1..2^72+2^63
-mpi_random_many:1:"01800000000000000000":100
-
-MPI random in range: 0..4
-mpi_random_many:0:"04":1000
-
-MPI random in range: 2..4
-mpi_random_many:1:"04":1000
-
-MPI random in range: 3..4
-mpi_random_many:1:"04":1000
-
 ECP read key #1 (short weierstrass, too small)
 depends_on:MBEDTLS_ECP_DP_SECP192R1_ENABLED
 mbedtls_ecp_read_key:MBEDTLS_ECP_DP_SECP192R1:"00":MBEDTLS_ERR_ECP_INVALID_KEY:0
diff --git a/tests/suites/test_suite_ecp.function b/tests/suites/test_suite_ecp.function
index fadc300..934598d 100644
--- a/tests/suites/test_suite_ecp.function
+++ b/tests/suites/test_suite_ecp.function
@@ -16,42 +16,6 @@
     mbedtls_ecp_point_free( x );    \
     mbedtls_ecp_point_init( x );
 
-#if defined(MBEDTLS_TEST_HOOKS) && defined(MBEDTLS_ECP_SHORT_WEIERSTRASS_ENABLED)
-/* Test whether bytes represents (in big-endian base 256) a number B that
- * is "significantly" above a power of 2, which is defined as follows.
- * Let n be the integer such that 2^n <= B < 2^{n+1}. B is significantly
- * above a power of 2 if (B - 2^n) / 2^n is not negligible. "Negligible"
- * is defined as having a negligible chance that if you draw an integer
- * in the range [1, B-1] K times, the number will always be less than 2^n,
- * where K is the iteration count passed to genkey_sw_many.
- */
-static int is_significantly_above_a_power_of_2( data_t *bytes )
-{
-    const uint8_t *p = bytes->x;
-    size_t len = bytes->len;
-    unsigned x;
-    while( len > 0 && p[0] == 0 )
-    {
-        ++p;
-        --len;
-    }
-    if( len == 0 )
-        return( 0 );
-    else if( len == 1 )
-        x = p[0];
-    else
-        x = ( p[0] << 8 ) | p[1];
-
-    if( x <= 4 )
-        return( 0 );
-
-    while( ( x & 0x8000 ) == 0 )
-        x <<= 1;
-    x &= 0x7fff;
-    return( x >= 0x1000 );
-}
-#endif
-
 /* END_HEADER */
 
 /* BEGIN_DEPENDENCIES
@@ -1323,113 +1287,6 @@
 }
 /* END_CASE */
 
-/* BEGIN_CASE depends_on:MBEDTLS_TEST_HOOKS:MBEDTLS_ECP_SHORT_WEIERSTRASS_ENABLED */
-void mpi_random_many( int min, data_t *bound_bytes, int iterations )
-{
-    /* Generate numbers in the range 1..bound-1. Do it iterations times.
-     * This function assumes that the value of bound is at least 2 and
-     * that iterations is large enough that a one-in-2^iterations chance
-     * effectively never occurs.
-     */
-
-    mbedtls_mpi upper_bound;
-    size_t n_bits;
-    mbedtls_mpi result;
-    size_t b;
-    /* If upper_bound is small, stats[b] is the number of times the value b
-     * has been generated. Otherwise stats[b] is the number of times a
-     * value with bit b set has been generated. */
-    size_t *stats = NULL;
-    size_t stats_len;
-    int full_stats;
-    size_t i;
-
-    mbedtls_mpi_init( &upper_bound );
-    mbedtls_mpi_init( &result );
-
-    TEST_EQUAL( 0, mbedtls_mpi_read_binary( &upper_bound,
-                                            bound_bytes->x, bound_bytes->len ) );
-    n_bits = mbedtls_mpi_bitlen( &upper_bound );
-    /* Consider a bound "small" if it's less than 2^5. This value is chosen
-     * to be small enough that the probability of missing one value is
-     * negligible given the number of iterations. It must be less than
-     * 256 because some of the code below assumes that "small" values
-     * fit in a byte. */
-    if( n_bits <= 5 )
-    {
-        full_stats = 1;
-        stats_len = bound_bytes->x[bound_bytes->len - 1];
-    }
-    else
-    {
-        full_stats = 0;
-        stats_len = n_bits;
-    }
-    ASSERT_ALLOC( stats, stats_len );
-
-    for( i = 0; i < (size_t) iterations; i++ )
-    {
-        mbedtls_test_set_step( i );
-        TEST_EQUAL( 0, mbedtls_mpi_random( &result, min, &upper_bound,
-                                           mbedtls_test_rnd_std_rand, NULL ) );
-
-        TEST_ASSERT( mbedtls_mpi_cmp_mpi( &result, &upper_bound ) < 0 );
-        TEST_ASSERT( mbedtls_mpi_cmp_int( &result, min ) >= 0 );
-        if( full_stats )
-        {
-            uint8_t value;
-            TEST_EQUAL( 0, mbedtls_mpi_write_binary( &result, &value, 1 ) );
-            TEST_ASSERT( value < stats_len );
-            ++stats[value];
-        }
-        else
-        {
-            for( b = 0; b < n_bits; b++ )
-                stats[b] += mbedtls_mpi_get_bit( &result, b );
-        }
-    }
-
-    if( full_stats )
-    {
-        for( b = 1; b < stats_len; b++ )
-        {
-            mbedtls_test_set_step( 1000000 + b );
-            /* Assert that each value has been reached at least once.
-             * This is almost guaranteed if the iteration count is large
-             * enough. This is a very crude way of checking the distribution.
-             */
-            TEST_ASSERT( stats[b] > 0 );
-        }
-    }
-    else
-    {
-        for( b = 0; b < n_bits; b++ )
-        {
-            mbedtls_test_set_step( 1000000 + b );
-            /* Assert that each bit has been set in at least one result and
-             * clear in at least one result. Provided that iterations is not
-             * too small, it would be extremely unlikely for this not to be
-             * the case if the results are uniformly distributed.
-             *
-             * As an exception, the top bit may legitimately never be set
-             * if bound is a power of 2 or only slightly above.
-             */
-            if( b != n_bits - 1 ||
-                is_significantly_above_a_power_of_2( bound_bytes ) )
-            {
-                TEST_ASSERT( stats[b] > 0 );
-            }
-            TEST_ASSERT( stats[b] < (size_t) iterations );
-        }
-    }
-
-exit:
-    mbedtls_mpi_free( &upper_bound );
-    mbedtls_mpi_free( &result );
-    mbedtls_free( stats );
-}
-/* END_CASE */
-
 /* BEGIN_CASE depends_on:MBEDTLS_SELF_TEST */
 void ecp_selftest(  )
 {
diff --git a/tests/suites/test_suite_mpi.data b/tests/suites/test_suite_mpi.data
index 59fd782..a057dcd 100644
--- a/tests/suites/test_suite_mpi.data
+++ b/tests/suites/test_suite_mpi.data
@@ -1033,6 +1033,105 @@
 Fill random: MAX_SIZE bytes, RNG failure after MAX_SIZE-1 bytes
 mpi_fill_random:MBEDTLS_MPI_MAX_SIZE:MBEDTLS_MPI_MAX_SIZE-1:MBEDTLS_ERR_ENTROPY_SOURCE_FAILED
 
+MPI random in range: 1..4
+mpi_random_many:1:"04":1000
+
+MPI random in range: 1..5
+mpi_random_many:1:"05":1000
+
+MPI random in range: 1..6
+mpi_random_many:1:"06":1000
+
+MPI random in range: 1..7
+mpi_random_many:1:"07":1000
+
+MPI random in range: 1..8
+mpi_random_many:1:"08":1000
+
+MPI random in range: 1..9
+mpi_random_many:1:"09":1000
+
+MPI random in range: 1..10
+mpi_random_many:1:"0a":1000
+
+MPI random in range: 1..11
+mpi_random_many:1:"0b":1000
+
+MPI random in range: 1..12
+mpi_random_many:1:"0c":1000
+
+MPI random in range: 1..255
+mpi_random_many:1:"ff":100
+
+MPI random in range: 1..256
+mpi_random_many:1:"0100":100
+
+MPI random in range: 1..257
+mpi_random_many:1:"0101":100
+
+MPI random in range: 1..272
+mpi_random_many:1:"0110":100
+
+MPI random in range: 1..2^64-1
+mpi_random_many:1:"ffffffffffffffff":100
+
+MPI random in range: 1..2^64
+mpi_random_many:1:"010000000000000000":100
+
+MPI random in range: 1..2^64+1
+mpi_random_many:1:"010000000000000001":100
+
+MPI random in range: 1..2^64+2^63
+mpi_random_many:1:"018000000000000000":100
+
+MPI random in range: 1..2^65-1
+mpi_random_many:1:"01ffffffffffffffff":100
+
+MPI random in range: 1..2^65
+mpi_random_many:1:"020000000000000000":100
+
+MPI random in range: 1..2^65+1
+mpi_random_many:1:"020000000000000001":100
+
+MPI random in range: 1..2^65+2^64
+mpi_random_many:1:"030000000000000000":100
+
+MPI random in range: 1..2^66+2^65
+mpi_random_many:1:"060000000000000000":100
+
+MPI random in range: 1..2^71-1
+mpi_random_many:1:"7fffffffffffffffff":100
+
+MPI random in range: 1..2^71
+mpi_random_many:1:"800000000000000000":100
+
+MPI random in range: 1..2^71+1
+mpi_random_many:1:"800000000000000001":100
+
+MPI random in range: 1..2^71+2^63
+mpi_random_many:1:"c00000000000000000":100
+
+MPI random in range: 1..2^72-1
+mpi_random_many:1:"ffffffffffffffffff":100
+
+MPI random in range: 1..2^72
+mpi_random_many:1:"01000000000000000000":100
+
+MPI random in range: 1..2^72+1
+mpi_random_many:1:"01000000000000000001":100
+
+MPI random in range: 1..2^72+2^63
+mpi_random_many:1:"01800000000000000000":100
+
+MPI random in range: 0..4
+mpi_random_many:0:"04":1000
+
+MPI random in range: 2..4
+mpi_random_many:1:"04":1000
+
+MPI random in range: 3..4
+mpi_random_many:1:"04":1000
+
 MPI Selftest
 depends_on:MBEDTLS_SELF_TEST
 mpi_selftest:
diff --git a/tests/suites/test_suite_mpi.function b/tests/suites/test_suite_mpi.function
index c5bb5a6..8cf53b3 100644
--- a/tests/suites/test_suite_mpi.function
+++ b/tests/suites/test_suite_mpi.function
@@ -64,6 +64,40 @@
     return( 0 );
 }
 
+/* Test whether bytes represents (in big-endian base 256) a number B that
+ * is "significantly" above a power of 2, which is defined as follows.
+ * Let n be the integer such that 2^n <= B < 2^{n+1}. B is significantly
+ * above a power of 2 if (B - 2^n) / 2^n is not negligible. "Negligible"
+ * is defined as having a negligible chance that if you draw an integer
+ * in the range [1, B-1] K times, the number will always be less than 2^n,
+ * where K is the iteration count passed to genkey_sw_many.
+ */
+static int is_significantly_above_a_power_of_2( data_t *bytes )
+{
+    const uint8_t *p = bytes->x;
+    size_t len = bytes->len;
+    unsigned x;
+    while( len > 0 && p[0] == 0 )
+    {
+        ++p;
+        --len;
+    }
+    if( len == 0 )
+        return( 0 );
+    else if( len == 1 )
+        x = p[0];
+    else
+        x = ( p[0] << 8 ) | p[1];
+
+    if( x <= 4 )
+        return( 0 );
+
+    while( ( x & 0x8000 ) == 0 )
+        x <<= 1;
+    x &= 0x7fff;
+    return( x >= 0x1000 );
+}
+
 /* END_HEADER */
 
 /* BEGIN_DEPENDENCIES
@@ -1396,6 +1430,113 @@
 }
 /* END_CASE */
 
+/* BEGIN_CASE */
+void mpi_random_many( int min, data_t *bound_bytes, int iterations )
+{
+    /* Generate numbers in the range 1..bound-1. Do it iterations times.
+     * This function assumes that the value of bound is at least 2 and
+     * that iterations is large enough that a one-in-2^iterations chance
+     * effectively never occurs.
+     */
+
+    mbedtls_mpi upper_bound;
+    size_t n_bits;
+    mbedtls_mpi result;
+    size_t b;
+    /* If upper_bound is small, stats[b] is the number of times the value b
+     * has been generated. Otherwise stats[b] is the number of times a
+     * value with bit b set has been generated. */
+    size_t *stats = NULL;
+    size_t stats_len;
+    int full_stats;
+    size_t i;
+
+    mbedtls_mpi_init( &upper_bound );
+    mbedtls_mpi_init( &result );
+
+    TEST_EQUAL( 0, mbedtls_mpi_read_binary( &upper_bound,
+                                            bound_bytes->x, bound_bytes->len ) );
+    n_bits = mbedtls_mpi_bitlen( &upper_bound );
+    /* Consider a bound "small" if it's less than 2^5. This value is chosen
+     * to be small enough that the probability of missing one value is
+     * negligible given the number of iterations. It must be less than
+     * 256 because some of the code below assumes that "small" values
+     * fit in a byte. */
+    if( n_bits <= 5 )
+    {
+        full_stats = 1;
+        stats_len = bound_bytes->x[bound_bytes->len - 1];
+    }
+    else
+    {
+        full_stats = 0;
+        stats_len = n_bits;
+    }
+    ASSERT_ALLOC( stats, stats_len );
+
+    for( i = 0; i < (size_t) iterations; i++ )
+    {
+        mbedtls_test_set_step( i );
+        TEST_EQUAL( 0, mbedtls_mpi_random( &result, min, &upper_bound,
+                                           mbedtls_test_rnd_std_rand, NULL ) );
+
+        TEST_ASSERT( mbedtls_mpi_cmp_mpi( &result, &upper_bound ) < 0 );
+        TEST_ASSERT( mbedtls_mpi_cmp_int( &result, min ) >= 0 );
+        if( full_stats )
+        {
+            uint8_t value;
+            TEST_EQUAL( 0, mbedtls_mpi_write_binary( &result, &value, 1 ) );
+            TEST_ASSERT( value < stats_len );
+            ++stats[value];
+        }
+        else
+        {
+            for( b = 0; b < n_bits; b++ )
+                stats[b] += mbedtls_mpi_get_bit( &result, b );
+        }
+    }
+
+    if( full_stats )
+    {
+        for( b = 1; b < stats_len; b++ )
+        {
+            mbedtls_test_set_step( 1000000 + b );
+            /* Assert that each value has been reached at least once.
+             * This is almost guaranteed if the iteration count is large
+             * enough. This is a very crude way of checking the distribution.
+             */
+            TEST_ASSERT( stats[b] > 0 );
+        }
+    }
+    else
+    {
+        for( b = 0; b < n_bits; b++ )
+        {
+            mbedtls_test_set_step( 1000000 + b );
+            /* Assert that each bit has been set in at least one result and
+             * clear in at least one result. Provided that iterations is not
+             * too small, it would be extremely unlikely for this not to be
+             * the case if the results are uniformly distributed.
+             *
+             * As an exception, the top bit may legitimately never be set
+             * if bound is a power of 2 or only slightly above.
+             */
+            if( b != n_bits - 1 ||
+                is_significantly_above_a_power_of_2( bound_bytes ) )
+            {
+                TEST_ASSERT( stats[b] > 0 );
+            }
+            TEST_ASSERT( stats[b] < (size_t) iterations );
+        }
+    }
+
+exit:
+    mbedtls_mpi_free( &upper_bound );
+    mbedtls_mpi_free( &result );
+    mbedtls_free( stats );
+}
+/* END_CASE */
+
 /* BEGIN_CASE depends_on:MBEDTLS_SELF_TEST */
 void mpi_selftest(  )
 {