1/* mpfr-mini-gmp.c -- Interface functions for mini-gmp. 2 3Copyright 2014-2023 Free Software Foundation, Inc. 4Contributed by the AriC and Caramba projects, INRIA. 5 6This file is part of the GNU MPFR Library. 7 8The GNU MPFR Library is free software; you can redistribute it and/or modify 9it under the terms of the GNU Lesser General Public License as published by 10the Free Software Foundation; either version 3 of the License, or (at your 11option) any later version. 12 13The GNU MPFR Library is distributed in the hope that it will be useful, but 14WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY 15or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public 16License for more details. 17 18You should have received a copy of the GNU Lesser General Public License 19along with the GNU MPFR Library; see the file COPYING.LESSER. If not, see 20https://www.gnu.org/licenses/ or write to the Free Software Foundation, Inc., 2151 Franklin St, Fifth Floor, Boston, MA 02110-1301, USA. */ 22 23/* The following include will do 2 things: include the config.h 24 if there is one (as it may define MPFR_USE_MINI_GMP), and avoid 25 an empty translation unit (see ISO C99, 6.9). */ 26#include "mpfr-impl.h" 27 28#ifdef MPFR_USE_MINI_GMP 29 30/************************ random generation functions ************************/ 31 32#ifdef WANT_gmp_randinit_default 33void 34gmp_randinit_default (gmp_randstate_t state) 35{ 36} 37#endif 38 39#ifdef WANT_gmp_randseed_ui 40void 41gmp_randseed_ui (gmp_randstate_t state, unsigned long int seed) 42{ 43 /* Note: We possibly ignore the high-order bits of seed. One should take 44 that into account when setting GMP_CHECK_RANDOMIZE for the tests. */ 45 srand ((unsigned int) seed); 46} 47#endif 48 49#ifdef WANT_gmp_randclear 50void 51gmp_randclear (gmp_randstate_t state) 52{ 53} 54#endif 55 56#ifdef WANT_gmp_randinit_set 57void 58gmp_randinit_set (gmp_randstate_t s1, gmp_randstate_t s2) 59{ 60} 61#endif 62 63static unsigned int 64rand15 (void) 65{ 66 /* With a good PRNG, we could use "rand () % 32768", but let's choose 67 the following from <https://c-faq.com/lib/randrange.html>. Note that 68 on most platforms, the compiler should generate a shift. */ 69 return rand () / (RAND_MAX / 32768 + 1); 70} 71 72static mp_limb_t 73random_limb (void) 74{ 75 mp_limb_t r = 0; 76 int i = GMP_NUMB_BITS; 77 78 while (i > 0) 79 { 80 r = (r << 15) | rand15 (); 81 i -= 15; 82 } 83 84 return r; 85} 86 87#ifdef WANT_mpz_urandomb 88void 89mpz_urandomb (mpz_t rop, gmp_randstate_t state, mp_bitcnt_t nbits) 90{ 91 unsigned long n, i; 92 93 mpz_realloc2 (rop, nbits); 94 n = (nbits - 1) / GMP_NUMB_BITS + 1; /* number of limbs */ 95 for (i = n; i-- > 0;) 96 rop->_mp_d[i] = random_limb (); 97 i = n * GMP_NUMB_BITS - nbits; 98 /* mask the upper i bits */ 99 if (i) 100 rop->_mp_d[n-1] = MPFR_LIMB_LSHIFT(rop->_mp_d[n-1], i) >> i; 101 while (n > 0 && (rop->_mp_d[n-1] == 0)) 102 n--; 103 rop->_mp_size = n; 104} 105#endif 106 107#ifdef WANT_gmp_urandomm_ui 108/* generates a random unsigned long */ 109static unsigned long 110random_ulong (void) 111{ 112#ifdef MPFR_LONG_WITHIN_LIMB 113 /* we assume a limb and an unsigned long have both a number of different 114 values that is a power of two, thus when we cast a random limb into 115 an unsigned long, we still get an uniform distribution */ 116 return random_limb (); 117#else 118 /* with the same assumption as above, we need to generate as many random 119 limbs needed to "fill" an unsigned long */ 120 unsigned long u, v; 121 122 v = MPFR_LIMB_MAX; 123 u = random_limb (); 124 while (v < ULONG_MAX) 125 { 126 v = (v << GMP_NUMB_BITS) + MPFR_LIMB_MAX; 127 u = (u << GMP_NUMB_BITS) + random_limb (); 128 } 129 return u; 130#endif 131} 132 133unsigned long 134gmp_urandomm_ui (gmp_randstate_t state, unsigned long n) 135{ 136 unsigned long p, q, r; 137 138 MPFR_ASSERTN (n > 0); 139 p = random_ulong (); /* p is in [0, ULONG_MAX], thus p is uniform among 140 ULONG_MAX+1 values */ 141 q = n * (ULONG_MAX / n); 142 r = ULONG_MAX % n; 143 if (r != n - 1) /* ULONG_MAX+1 is not multiple of n, will happen whenever 144 n is not a power of two */ 145 while (p >= q) 146 p = random_ulong (); 147 return p % n; 148} 149#endif 150 151#ifdef WANT_gmp_urandomb_ui 152unsigned long 153gmp_urandomb_ui (gmp_randstate_t state, unsigned long n) 154{ 155#ifdef MPFR_LONG_WITHIN_LIMB 156 /* Since n may be equal to the width of unsigned long, 157 we must not shift 1UL by n as this may be UB. */ 158 return n == 0 ? 0 : random_limb () & (((1UL << (n - 1)) << 1) - 1); 159#else 160 unsigned long res = 0; 161 int m = n; /* remaining bits to generate */ 162 while (m >= GMP_NUMB_BITS) 163 { 164 /* we can generate a full limb */ 165 res = (res << GMP_NUMB_BITS) | (unsigned long) random_limb (); 166 m -= GMP_NUMB_BITS; 167 } 168 MPFR_ASSERTD (m < GMP_NUMB_BITS); /* thus m < width(unsigned long) */ 169 if (m != 0) /* generate m extra bits */ 170 res = (res << m) | (unsigned long) (random_limb () % (1UL << m)); 171 return res; 172#endif 173} 174#endif 175 176/************************* division functions ********************************/ 177 178#ifdef WANT_mpn_divrem_1 179mp_limb_t 180mpn_divrem_1 (mp_limb_t *qp, mp_size_t qxn, mp_limb_t *np, mp_size_t nn, 181 mp_limb_t d0) 182{ 183 mpz_t q, r, n, d; 184 mp_limb_t ret, dd[1]; 185 186 d->_mp_d = dd; 187 d->_mp_d[0] = d0; 188 d->_mp_size = 1; 189 mpz_init (q); 190 mpz_init (r); 191 if (qxn == 0) 192 { 193 n->_mp_d = np; 194 n->_mp_size = nn; 195 } 196 else 197 { 198 mpz_init2 (n, (nn + qxn) * GMP_NUMB_BITS); 199 mpn_copyi (n->_mp_d + qxn, np, nn); 200 mpn_zero (n->_mp_d, qxn); 201 n->_mp_size = nn + qxn; 202 } 203 mpz_tdiv_qr (q, r, n, d); 204 if (q->_mp_size > 0) 205 mpn_copyi (qp, q->_mp_d, q->_mp_size); 206 if (q->_mp_size < nn + qxn) 207 mpn_zero (qp + q->_mp_size, nn + qxn - q->_mp_size); 208 ret = (r->_mp_size == 1) ? r->_mp_d[0] : 0; 209 mpz_clear (q); 210 mpz_clear (r); 211 if (qxn != 0) 212 mpz_clear (n); 213 return ret; 214} 215#endif 216 217#ifdef WANT_mpn_divrem 218mp_limb_t 219mpn_divrem (mp_limb_t *qp, mp_size_t qn, mp_limb_t *np, 220 mp_size_t nn, const mp_limb_t *dp, mp_size_t dn) 221{ 222 mpz_t q, r, n, d; 223 mp_limb_t ret; 224 225 MPFR_ASSERTN(qn == 0); 226 qn = nn - dn; 227 n->_mp_d = np; 228 n->_mp_size = nn; 229 d->_mp_d = (mp_limb_t*) dp; 230 d->_mp_size = dn; 231 mpz_init (q); 232 mpz_init (r); 233 mpz_tdiv_qr (q, r, n, d); 234 MPFR_ASSERTN(q->_mp_size == qn || q->_mp_size == qn + 1); 235 mpn_copyi (qp, q->_mp_d, qn); 236 ret = (q->_mp_size == qn) ? 0 : q->_mp_d[qn]; 237 if (r->_mp_size > 0) 238 mpn_copyi (np, r->_mp_d, r->_mp_size); 239 if (r->_mp_size < dn) 240 mpn_zero (np + r->_mp_size, dn - r->_mp_size); 241 mpz_clear (q); 242 mpz_clear (r); 243 return ret; 244} 245#endif 246 247#ifdef WANT_mpn_tdiv_qr 248void 249mpn_tdiv_qr (mp_limb_t *qp, mp_limb_t *rp, mp_size_t qxn, 250 const mp_limb_t *np, mp_size_t nn, 251 const mp_limb_t *dp, mp_size_t dn) 252{ 253 mpz_t q, r, n, d; 254 255 MPFR_ASSERTN(qxn == 0); 256 n->_mp_d = (mp_limb_t*) np; 257 n->_mp_size = nn; 258 d->_mp_d = (mp_limb_t*) dp; 259 d->_mp_size = dn; 260 mpz_init (q); 261 mpz_init (r); 262 mpz_tdiv_qr (q, r, n, d); 263 MPFR_ASSERTN(q->_mp_size > 0); 264 mpn_copyi (qp, q->_mp_d, q->_mp_size); 265 if (q->_mp_size < nn - dn + 1) 266 mpn_zero (qp + q->_mp_size, nn - dn + 1 - q->_mp_size); 267 if (r->_mp_size > 0) 268 mpn_copyi (rp, r->_mp_d, r->_mp_size); 269 if (r->_mp_size < dn) 270 mpn_zero (rp + r->_mp_size, dn - r->_mp_size); 271 mpz_clear (q); 272 mpz_clear (r); 273} 274#endif 275 276#if 0 /* this function is useful for debugging, thus please keep it here */ 277void 278mpz_dump (mpz_t z) 279{ 280 mp_size_t n = z->_mp_size; 281 282 MPFR_STAT_STATIC_ASSERT ((GMP_NUMB_BITS % 4) == 0); 283 284 if (n == 0) 285 printf ("0"); 286 else 287 { 288 int first = 1; 289 if (n < 0) 290 { 291 printf ("-"); 292 n = -n; 293 } 294 while (n > 0) 295 { 296 if (first) 297 { 298 printf ("%lx", (unsigned long) z->_mp_d[n-1]); 299 first = 0; 300 } 301 else 302 { 303 char s[17]; 304 int len; 305 sprintf (s, "%lx", (unsigned long) z->_mp_d[n-1]); 306 len = strlen (s); 307 /* one character should be printed for 4 bits */ 308 while (len++ < GMP_NUMB_BITS / 4) 309 printf ("0"); 310 printf ("%lx", (unsigned long) z->_mp_d[n-1]); 311 } 312 n--; 313 } 314 } 315 printf ("\n"); 316} 317#endif 318 319#endif /* MPFR_USE_MINI_GMP */ 320