Move mbedtls_mpi_random to the bignum module

Since mbedtls_mpi_random() is not specific to ECC code, move it from
the ECP module to the bignum module.

This increases the code size in builds without short Weierstrass
curves (including builds without ECC at all) that do not optimize out
unused functions.

Signed-off-by: Gilles Peskine <Gilles.Peskine@arm.com>
diff --git a/include/mbedtls/bignum.h b/include/mbedtls/bignum.h
index 637360e..7f718e6 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 e00204b..1cde6b0 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 1447a5a..d67b78b 100644
--- a/library/ecp.c
+++ b/library/ecp.c
@@ -3074,62 +3074,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(  )
 {