Annotation of OpenXM_contrib/gmp/mpz/rrandomb.c, Revision 1.1.1.2
1.1 maekawa 1: /* mpz_rrandomb -- Generate a positive random mpz_t of specified bit size, with
2: long runs of consecutive ones and zeros in the binary representation.
3: Meant for testing of other MP routines.
4:
1.1.1.2 ! ohara 5: Copyright 2000, 2001, 2002 Free Software Foundation, Inc.
1.1 maekawa 6:
7: This file is part of the GNU MP Library.
8:
9: The GNU MP Library is free software; you can redistribute it and/or modify
10: it under the terms of the GNU Lesser General Public License as published by
11: the Free Software Foundation; either version 2.1 of the License, or (at your
12: option) any later version.
13:
14: The GNU MP Library is distributed in the hope that it will be useful, but
15: WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
16: or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
17: License for more details.
18:
19: You should have received a copy of the GNU Lesser General Public License
20: along with the GNU MP Library; see the file COPYING.LIB. If not, write to
21: the Free Software Foundation, Inc., 59 Temple Place - Suite 330, Boston,
22: MA 02111-1307, USA. */
23:
24: #include "gmp.h"
25: #include "gmp-impl.h"
26:
27: static void gmp_rrandomb _PROTO ((mp_ptr rp, gmp_randstate_t rstate, unsigned long int nbits));
28:
29: void
30: mpz_rrandomb (mpz_ptr x, gmp_randstate_t rstate, unsigned long int nbits)
31: {
32: mp_size_t nl = 0;
33:
34: if (nbits != 0)
35: {
36: mp_ptr xp;
1.1.1.2 ! ohara 37: nl = (nbits + GMP_NUMB_BITS - 1) / GMP_NUMB_BITS;
! 38: MPZ_REALLOC (x, nl);
1.1 maekawa 39:
40: xp = PTR(x);
41: gmp_rrandomb (xp, rstate, nbits);
42: MPN_NORMALIZE (xp, nl);
43: }
44:
45: SIZ(x) = nl;
46: }
47:
1.1.1.2 ! ohara 48: /* It's a bit tricky to get this right, so please test the code well if you
! 49: hack with it. Some early versions of the function produced random numbers
! 50: with the leading limb == 0, and some versions never made the most
! 51: significant bit set.
! 52:
! 53: This code and mpn_random2 are almost identical, though the latter makes bit
! 54: runs of 1 to 32, and forces the first block to contain 1-bits.
! 55:
! 56: The random state produces some number of random bits per underlying lc
! 57: invocation (BITS_PER_RANDCALL). We should perhaps ask for that, instead of
! 58: asking for 32, presuming that limbs are at least 32 bits. FIXME: Handle
! 59: smaller limbs, such as 4-bit limbs useful for testing purposes, or limbs
! 60: truncated by nailing.
! 61:
! 62: For efficiency, we make sure to use most bits returned from _gmp_rand, since
! 63: the underlying random number generator is slow. We keep returned bits in
! 64: ranm/ran, and a count of how many bits remaining in ran_nbits. */
! 65:
! 66: #define LOGBITS_PER_BLOCK 4
! 67:
! 68: /* Ask _gmp_rand for 32 bits per call unless that's more than a limb can hold.
! 69: Thus, we get the same random number sequence in the common cases.
! 70: FIXME: We should always generate the same random number sequence! */
! 71: #if GMP_NUMB_BITS < 32
! 72: #define BITS_PER_RANDCALL GMP_NUMB_BITS
! 73: #else
! 74: #define BITS_PER_RANDCALL 32
! 75: #endif
1.1 maekawa 76:
77: static void
78: gmp_rrandomb (mp_ptr rp, gmp_randstate_t rstate, unsigned long int nbits)
79: {
80: int nb;
1.1.1.2 ! ohara 81: int bit_pos; /* bit number of least significant bit where
! 82: next bit field to be inserted */
! 83: mp_size_t ri; /* index in rp */
! 84: mp_limb_t ran, ranm; /* buffer for random bits */
! 85: mp_limb_t acc; /* accumulate output random data here */
! 86: int ran_nbits; /* number of valid bits in ran */
! 87:
! 88: ran_nbits = 0;
! 89: bit_pos = (nbits - 1) % GMP_NUMB_BITS;
! 90: ri = (nbits - 1) / GMP_NUMB_BITS;
1.1 maekawa 91:
92: acc = 0;
1.1.1.2 ! ohara 93: while (ri >= 0)
1.1 maekawa 94: {
1.1.1.2 ! ohara 95: if (ran_nbits < LOGBITS_PER_BLOCK + 1)
! 96: {
! 97: _gmp_rand (&ranm, rstate, BITS_PER_RANDCALL);
! 98: ran = ranm;
! 99: ran_nbits = BITS_PER_RANDCALL;
! 100: }
! 101:
! 102: nb = (ran >> 1) % (1 << LOGBITS_PER_BLOCK) + 1;
1.1 maekawa 103: if ((ran & 1) != 0)
104: {
1.1.1.2 ! ohara 105: /* Generate a string of nb ones. */
1.1 maekawa 106: if (nb > bit_pos)
107: {
1.1.1.2 ! ohara 108: rp[ri--] = acc | (((mp_limb_t) 2 << bit_pos) - 1);
! 109: bit_pos += GMP_NUMB_BITS;
1.1 maekawa 110: bit_pos -= nb;
1.1.1.2 ! ohara 111: acc = ((~(mp_limb_t) 1) << bit_pos) & GMP_NUMB_MASK;
1.1 maekawa 112: }
113: else
114: {
115: bit_pos -= nb;
1.1.1.2 ! ohara 116: acc |= (((mp_limb_t) 2 << nb) - 2) << bit_pos;
1.1 maekawa 117: }
118: }
119: else
120: {
1.1.1.2 ! ohara 121: /* Generate a string of nb zeroes. */
1.1 maekawa 122: if (nb > bit_pos)
123: {
1.1.1.2 ! ohara 124: rp[ri--] = acc;
1.1 maekawa 125: acc = 0;
1.1.1.2 ! ohara 126: bit_pos += GMP_NUMB_BITS;
1.1 maekawa 127: }
128: bit_pos -= nb;
129: }
1.1.1.2 ! ohara 130: ran_nbits -= LOGBITS_PER_BLOCK + 1;
! 131: ran >>= LOGBITS_PER_BLOCK + 1;
1.1 maekawa 132: }
133: }
FreeBSD-CVSweb <freebsd-cvsweb@FreeBSD.org>