1 /* mpfr_urandomb (rop, state, nbits) -- Generate a uniform pseudorandom
2 real number between 0 (inclusive) and 1 (exclusive) of size NBITS,
3 using STATE as the random state previously initialized by a call to
4 gmp_randinit_lc_2exp_size().
6 Copyright 2000, 2001, 2002, 2003, 2004, 2006, 2007, 2008, 2009, 2010, 2011 Free Software Foundation, Inc.
7 Contributed by the Arenaire and Caramel projects, INRIA.
9 This file is part of the GNU MPFR Library.
11 The GNU MPFR Library is free software; you can redistribute it and/or modify
12 it under the terms of the GNU Lesser General Public License as published by
13 the Free Software Foundation; either version 3 of the License, or (at your
14 option) any later version.
16 The GNU MPFR Library is distributed in the hope that it will be useful, but
17 WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
18 or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
19 License for more details.
21 You should have received a copy of the GNU Lesser General Public License
22 along with the GNU MPFR Library; see the file COPYING.LESSER. If not, see
23 http://www.gnu.org/licenses/ or write to the Free Software Foundation, Inc.,
24 51 Franklin St, Fifth Floor, Boston, MA 02110-1301, USA. */
27 #define MPFR_NEED_LONGLONG_H
28 #include "mpfr-impl.h"
30 /* generate nbits random bits into mp[], assuming mp was allocated to contain
31 a sufficient number of limbs */
33 mpfr_rand_raw (mpfr_limb_ptr mp, gmp_randstate_t rstate,
34 unsigned long int nbits)
38 /* To be sure to avoid the potential allocation of mpz_urandomb */
39 ALLOC(z) = SIZ(z) = ((nbits - 1) / GMP_NUMB_BITS) + 1;
41 mpz_urandomb (z, rstate, nbits);
45 mpfr_urandomb (mpfr_ptr rop, gmp_randstate_t rstate)
50 mp_size_t k; /* number of high zero limbs */
55 nbits = MPFR_PREC (rop);
56 nlimbs = MPFR_LIMB_SIZE (rop);
58 cnt = nlimbs * GMP_NUMB_BITS - nbits;
60 /* Uniform non-normalized significand */
61 /* generate exactly nbits so that the random generator stays in the same
62 state, independent of the machine word size GMP_NUMB_BITS */
63 mpfr_rand_raw (rp, rstate, nbits);
64 if (MPFR_LIKELY (cnt != 0)) /* this will put the low bits to zero */
65 mpn_lshift (rp, rp, nlimbs, cnt);
67 /* Count the null significant limbs and remaining limbs */
70 while (nlimbs != 0 && rp[nlimbs - 1] == 0)
77 if (MPFR_LIKELY (nlimbs != 0)) /* otherwise value is zero */
79 count_leading_zeros (cnt, rp[nlimbs - 1]);
81 if (mpfr_set_exp (rop, exp - cnt))
83 /* If the exponent is not in the current exponent range, we
84 choose to return a NaN as this is probably a user error.
85 Indeed this can happen only if the exponent range has been
86 reduced to a very small interval and/or the precision is
87 huge (very unlikely). */
89 __gmpfr_flags |= MPFR_FLAGS_NAN; /* Can't use MPFR_RET_NAN */
93 mpn_lshift (rp + k, rp, nlimbs, cnt);