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-2004, 2006-2015 Free Software Foundation, Inc.
7 Contributed by the AriC 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
,
38 MPFR_ASSERTN (nbits
>= 1);
39 /* To be sure to avoid the potential allocation of mpz_urandomb */
40 ALLOC(z
) = SIZ(z
) = MPFR_PREC2LIMBS (nbits
);
43 /* Check for integer overflow (unless mp_bitcnt_t is signed,
44 but according to the GMP manual, this shouldn't happen).
45 Note: mp_bitcnt_t has been introduced in GMP 5.0.0. */
46 MPFR_ASSERTN ((mp_bitcnt_t
) -1 < 0 || nbits
<= (mp_bitcnt_t
) -1);
48 mpz_urandomb (z
, rstate
, nbits
);
52 mpfr_urandomb (mpfr_ptr rop
, gmp_randstate_t rstate
)
57 mp_size_t k
; /* number of high zero limbs */
62 nbits
= MPFR_PREC (rop
);
63 nlimbs
= MPFR_LIMB_SIZE (rop
);
65 cnt
= nlimbs
* GMP_NUMB_BITS
- nbits
;
67 /* Uniform non-normalized significand */
68 /* generate exactly nbits so that the random generator stays in the same
69 state, independent of the machine word size GMP_NUMB_BITS */
70 mpfr_rand_raw (rp
, rstate
, nbits
);
71 if (MPFR_LIKELY (cnt
!= 0)) /* this will put the low bits to zero */
72 mpn_lshift (rp
, rp
, nlimbs
, cnt
);
74 /* Count the null significant limbs and remaining limbs */
77 while (nlimbs
!= 0 && rp
[nlimbs
- 1] == 0)
84 if (MPFR_LIKELY (nlimbs
!= 0)) /* otherwise value is zero */
86 count_leading_zeros (cnt
, rp
[nlimbs
- 1]);
88 if (mpfr_set_exp (rop
, exp
- cnt
))
90 /* If the exponent is not in the current exponent range, we
91 choose to return a NaN as this is probably a user error.
92 Indeed this can happen only if the exponent range has been
93 reduced to a very small interval and/or the precision is
94 huge (very unlikely). */
96 __gmpfr_flags
|= MPFR_FLAGS_NAN
; /* Can't use MPFR_RET_NAN */
100 mpn_lshift (rp
+ k
, rp
, nlimbs
, cnt
);
102 MPN_COPY (rp
+ k
, rp
, nlimbs
);