1 /* mpfr-mini-gmp.c -- Interface functions for mini-gmp.
2
3 Copyright 2014-2023 Free Software Foundation, Inc.
4 Contributed by the AriC and Caramba projects, INRIA.
5
6 This file is part of the GNU MPFR Library.
7
8 The GNU MPFR Library is free software; you can redistribute it and/or modify
9 it under the terms of the GNU Lesser General Public License as published by
10 the Free Software Foundation; either version 3 of the License, or (at your
11 option) any later version.
12
13 The GNU MPFR Library is distributed in the hope that it will be useful, but
14 WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
15 or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
16 License for more details.
17
18 You should have received a copy of the GNU Lesser General Public License
19 along with the GNU MPFR Library; see the file COPYING.LESSER. If not, see
20 https://www.gnu.org/licenses/ or write to the Free Software Foundation, Inc.,
21 51 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
33 void
34 gmp_randinit_default (gmp_randstate_t state)
35 {
36 }
37 #endif
38
39 #ifdef WANT_gmp_randseed_ui
40 void
41 gmp_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
50 void
51 gmp_randclear (gmp_randstate_t state)
52 {
53 }
54 #endif
55
56 #ifdef WANT_gmp_randinit_set
57 void
58 gmp_randinit_set (gmp_randstate_t s1, gmp_randstate_t s2)
59 {
60 }
61 #endif
62
63 static unsigned int
64 rand15 (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
72 static mp_limb_t
73 random_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
88 void
89 mpz_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 */
109 static unsigned long
110 random_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
133 unsigned long
134 gmp_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
152 unsigned long
153 gmp_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
179 mp_limb_t
180 mpn_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
218 mp_limb_t
219 mpn_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
248 void
249 mpn_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 */
277 void
278 mpz_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 */