1 /* Test mpn_fib2m.
2
3 Copyright 2018 Free Software Foundation, Inc.
4
5 This file is part of the GNU MP Library test suite.
6
7 The GNU MP Library test suite is free software; you can redistribute it
8 and/or modify it under the terms of the GNU General Public License as
9 published by the Free Software Foundation; either version 3 of the License,
10 or (at your option) any later version.
11
12 The GNU MP Library test suite is distributed in the hope that it will be
13 useful, but WITHOUT ANY WARRANTY; without even the implied warranty of
14 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General
15 Public License for more details.
16
17 You should have received a copy of the GNU General Public License along with
18 the GNU MP Library test suite. If not, see https://www.gnu.org/licenses/. */
19
20 #include <stdio.h>
21 #include <stdlib.h>
22
23 #include "gmp-impl.h"
24 #include "tests.h"
25
26 #define MAX_K_BITS 16
27 #define MAX_K (1L << MAX_K_BITS)
28 #define MIN_K 1
29
30 #define MAX_MN 20
31 #define MAX_KN 30
32
33 #define COUNT 200
34
35 static int
36 test_fib2_fib2m (int count, gmp_randstate_ptr rands)
37 {
38 int test;
39 mp_ptr fk, fks1, fkm, fks1m, mp, qp;
40 mp_size_t mn, fn, size, max_mn;
41 TMP_DECL;
42
43 size = MPN_FIB2_SIZE (MAX_K);
44 max_mn = size / 4 + 10;
45 ASSERT (max_mn < size);
46
47 TMP_MARK;
48 fk = TMP_ALLOC_LIMBS (size);
49 fks1 = TMP_ALLOC_LIMBS (size);
50 qp = TMP_ALLOC_LIMBS (size);
51 mp = TMP_ALLOC_LIMBS (max_mn);
52 fkm = 1 + TMP_ALLOC_LIMBS (max_mn * 2 + 1 + 2);
53 fks1m = 1 + TMP_ALLOC_LIMBS (max_mn * 2 + 1 + 2);
54
55 for (test = 1; test <= count; ++test)
56 {
57 mp_limb_t fk_before, fk_after, fk1_before, fk1_after;
58 int signflip;
59 unsigned long k;
60
61 k = MIN_K +
62 gmp_urandomm_ui (rands, test < MAX_K_BITS ?
63 MAX_K >> test : (MAX_K - MIN_K));
64
65 fn = mpn_fib2_ui (fk, fks1, k);
66 do {
67 mn = gmp_urandomm_ui (rands, MAX_K) % (fn / 4 + 10);
68 } while (mn == 0);
69 ASSERT (mn <= max_mn);
70 mpn_random2 (mp, mn);
71 ASSERT (mp [mn - 1] != 0);
72
73 if (fn >= mn)
74 {
75 mpn_tdiv_qr (qp, fk, 0, fk, fn, mp, mn);
76 mpn_tdiv_qr (qp, fks1, 0, fks1, fn, mp, mn);
77 }
78 else
79 {
80 MPN_ZERO (fk + fn, mn - fn);
81 MPN_ZERO (fks1 + fn, mn - fn);
82 }
83
84 mpn_random2 (fkm - 1, 2*mn+1+2);
85 fk_before = fkm [-1];
86 fk_after = fkm [2 * mn + 1];
87
88 mpn_random2 (fks1m - 1, 2*mn+1+2);
89 fk1_before = fks1m [-1];
90 fk1_after = fks1m [2 * mn + 1];
91
92 qp [0] = k;
93 signflip = mpn_fib2m (fkm, fks1m, qp, 1, mp, mn);
94 if (fkm [-1] != fk_before || fkm [2 * mn + 1] != fk_after
95 || fks1m [-1] != fk1_before || fks1m [2 * mn + 1] != fk1_after)
96 {
97 printf ("REDZONE violation in test %d, k = %lu, mn = %u\n",
98 test, k, (unsigned) mn);
99 if (fkm[-1] != fk_before)
100 {
101 printf ("before fkm:"); mpn_dump (fkm - 1, 1);
102 printf ("keep: "); mpn_dump (&fk_before, 1);
103 }
104 if (fkm[2 * mn + 1] != fk_after)
105 {
106 printf ("after fkm:"); mpn_dump (fkm + 2 * mn + 1, 1);
107 printf ("keep: "); mpn_dump (&fk_after, 1);
108 }
109 if (fks1m[-1] != fk1_before)
110 {
111 printf ("before fks1m:"); mpn_dump (fks1m - 1, 1);
112 printf ("keep: "); mpn_dump (&fk1_before, 1);
113 }
114 if (fks1m[2 * mn + 1] != fk1_after)
115 {
116 printf ("after fks1m:"); mpn_dump (fks1m + 2 * mn + 1, 1);
117 printf ("keep: "); mpn_dump (&fk1_after, 1);
118 }
119 abort();
120 }
121
122 if (mpn_cmp (fkm, fk, mn) != 0)
123 {
124 if (mpn_sub_n (fk, mp, fk, mn) || mpn_cmp (fkm, fk, mn) != 0)
125 {
126 printf ("ERROR(k) in test %d, k = %lu, mn = %u\n",
127 test, k, (unsigned) mn);
128 mpn_dump (fk, mn);
129 mpn_dump (fkm, mn);
130 mpn_dump (mp, mn);
131 abort();
132 }
133 signflip ^= 1;
134 }
135
136 if (mpn_cmp (fks1m, fks1, mn) != 0)
137 {
138 if (mpn_sub_n (fks1, mp, fks1, mn) || mpn_cmp (fks1m, fks1, mn) != 0)
139 {
140 printf ("ERROR(k-1) in test %d, k = %lu, mn = %u\n",
141 test, k, (unsigned) mn);
142 mpn_dump (fks1, mn);
143 mpn_dump (fks1m, mn);
144 mpn_dump (mp, mn);
145 abort();
146 }
147 signflip ^= 1;
148 }
149
150 if (signflip != 0 && ! mpn_zero_p (fks1m, mn) && ! mpn_zero_p (fkm, mn))
151 {
152 if ((mp [0] & 1) == 0) /* Should we test only odd modulus-es? */
153 {
154 if (! mpn_lshift (fks1m, fks1m, mn, 1) &&
155 mpn_cmp (mp, fks1m, mn) == 0)
156 continue;
157 if (! mpn_lshift (fkm, fkm, mn, 1) &&
158 mpn_cmp (mp, fkm, mn) == 0)
159 continue;
160 }
161 printf ("ERROR(sign) in test %d, k = %lu, mn = %u\n",
162 test, k, (unsigned) mn);
163 abort();
164 }
165 }
166 TMP_FREE;
167 return 0;
168 }
169
170 static int
171 test_fib2m_2exp (int count, gmp_randstate_ptr rands)
172 {
173 int test;
174 mp_ptr fka, fks1a, fkb, fks1b, mp, kp;
175 TMP_DECL;
176
177 TMP_MARK;
178 kp = TMP_ALLOC_LIMBS (MAX_KN);
179 mp = TMP_ALLOC_LIMBS (MAX_MN);
180 fka = 1 + TMP_ALLOC_LIMBS (MAX_MN * 2 + 1 + 2);
181 fks1a = 1 + TMP_ALLOC_LIMBS (MAX_MN * 2 + 1 + 2);
182 fkb = 1 + TMP_ALLOC_LIMBS (MAX_MN * 2 + 1 + 2);
183 fks1b = 1 + TMP_ALLOC_LIMBS (MAX_MN * 2 + 1 + 2);
184
185 for (test = 1; test <= count; ++test)
186 {
187 mp_limb_t fka_before, fka_after, fk1a_before, fk1a_after;
188 mp_limb_t fkb_before, fkb_after, fk1b_before, fk1b_after;
189 mp_size_t mn, kn;
190 int signflip;
191 mp_bitcnt_t exp2;
192
193 mn = gmp_urandomm_ui (rands, MAX_MN - 1) + 1;
194 mpn_random2 (mp, mn);
195
196 exp2 = MIN_K + 1 + gmp_urandomm_ui (rands, MAX_KN * GMP_NUMB_BITS - MIN_K - 1);
197
198 kn = BITS_TO_LIMBS (exp2);
199 MPN_ZERO (kp, kn - 1);
200 kp [kn - 1] = CNST_LIMB (1) << ((exp2 - 1) % GMP_NUMB_BITS);
201
202 mpn_random2 (fka - 1, 2*mn+1+2);
203 fka_before = fka [-1];
204 fka_after = fka [2 * mn + 1];
205
206 mpn_random2 (fks1a - 1, 2*mn+1+2);
207 fk1a_before = fks1a [-1];
208 fk1a_after = fks1a [2 * mn + 1];
209
210 signflip = mpn_fib2m (fka, fks1a, kp, kn, mp, mn);
211 if (fka [-1] != fka_before || fka [2 * mn + 1] != fka_after
212 || fks1a [-1] != fk1a_before || fks1a [2 * mn + 1] != fk1a_after)
213 {
214 printf ("REDZONE(a) violation in test %d, exp2 = %lu\n", test, exp2);
215 if (fka[-1] != fka_before)
216 {
217 printf ("before fka:"); mpn_dump (fka - 1, 1);
218 printf ("keep: "); mpn_dump (&fka_before, 1);
219 }
220 if (fka[2 * mn + 1] != fka_after)
221 {
222 printf ("after fka:"); mpn_dump (fka + 2 * mn + 1, 1);
223 printf ("keep: "); mpn_dump (&fka_after, 1);
224 }
225 if (fks1a[-1] != fk1a_before)
226 {
227 printf ("before fks1a:"); mpn_dump (fks1a - 1, 1);
228 printf ("keep: "); mpn_dump (&fk1a_before, 1);
229 }
230 if (fks1a[2 * mn + 1] != fk1a_after)
231 {
232 printf ("after fks1a:"); mpn_dump (fks1a + 2 * mn + 1, 1);
233 printf ("keep: "); mpn_dump (&fk1a_after, 1);
234 }
235 abort();
236 }
237
238 if (signflip && ! mpn_zero_p (fks1a, mn))
239 mpn_sub_n (fks1a, mp, fks1a, mn);
240 if (mpn_sub_n (fka, fka, fks1a, mn))
241 ASSERT_CARRY (mpn_add_n (fka, fka, mp, mn));
242
243 mpn_sub_1 (kp, kp, kn, 1);
244 ASSERT (exp2 % GMP_NUMB_BITS == 1 || kp [kn - 1] != 0);
245 kn -= kp [kn - 1] == 0;
246
247 mpn_random2 (fkb - 1, 2*mn+1+2);
248 fkb_before = fkb [-1];
249 fkb_after = fkb [2 * mn + 1];
250
251 mpn_random2 (fks1b - 1, 2*mn+1+2);
252 fk1b_before = fks1b [-1];
253 fk1b_after = fks1b [2 * mn + 1];
254
255 signflip = mpn_fib2m (fkb, fks1b, kp, kn, mp, mn);
256 if (fkb [-1] != fkb_before || fkb [2 * mn + 1] != fkb_after
257 || fks1b [-1] != fk1b_before || fks1b [2 * mn + 1] != fk1b_after)
258 {
259 printf ("REDZONE(b) violation in test %d, exp2 = %lu\n", test, exp2);
260 if (fkb[-1] != fkb_before)
261 {
262 printf ("before fkb:"); mpn_dump (fkb - 1, 1);
263 printf ("keep: "); mpn_dump (&fkb_before, 1);
264 }
265 if (fkb[2 * mn + 1] != fkb_after)
266 {
267 printf ("after fkb:"); mpn_dump (fkb + 2 * mn + 1, 1);
268 printf ("keep: "); mpn_dump (&fkb_after, 1);
269 }
270 if (fks1b[-1] != fk1b_before)
271 {
272 printf ("before fks1b:"); mpn_dump (fks1b - 1, 1);
273 printf ("keep: "); mpn_dump (&fk1b_before, 1);
274 }
275 if (fks1b[2 * mn + 1] != fk1b_after)
276 {
277 printf ("after fks1b:"); mpn_dump (fks1b + 2 * mn + 1, 1);
278 printf ("keep: "); mpn_dump (&fk1b_after, 1);
279 }
280 abort();
281 }
282
283 if (mpn_cmp (fks1a, fkb, mn) != 0)
284 {
285 if (mpn_sub_n (fkb, mp, fkb, mn) || mpn_cmp (fks1a, fkb, mn) != 0)
286 {
287 printf ("ERROR(k) in test %d, exp2 = %lu\n", test, exp2);
288 mpn_dump (fks1a, mn);
289 mpn_dump (fkb, mn);
290 mpn_dump (mp, mn);
291 abort();
292 }
293 signflip ^= 1;
294 }
295
296 if (mpn_cmp (fka, fks1b, mn) != 0)
297 {
298 if (mpn_sub_n (fks1b, mp, fks1b, mn) || mpn_cmp (fka, fks1b, mn) != 0)
299 {
300 printf ("ERROR(k-1) in test %d, exp2 = %lu\n", test, exp2);
301 mpn_dump (fka, mn);
302 mpn_dump (fks1b, mn);
303 mpn_dump (mp, mn);
304 abort();
305 }
306 signflip ^= 1;
307 }
308
309 if (signflip != 0 && ! mpn_zero_p (fks1b, mn) && ! mpn_zero_p (fkb, mn))
310 {
311 if ((mp [0] & 1) == 0) /* Should we test only odd modulus-es? */
312 {
313 if (! mpn_lshift (fks1b, fks1b, mn, 1) &&
314 mpn_cmp (mp, fks1b, mn) == 0)
315 continue;
316 if (! mpn_lshift (fkb, fkb, mn, 1) &&
317 mpn_cmp (mp, fkb, mn) == 0)
318 continue;
319 }
320 printf ("ERROR(sign) in test %d, exp2 = %lu\n",
321 test, exp2);
322 abort();
323 }
324 }
325 TMP_FREE;
326 return 0;
327 }
328
329 int
330 main (int argc, char **argv)
331 {
332 int count = COUNT;
333 gmp_randstate_ptr rands;
334
335 tests_start ();
336 TESTS_REPS (count, argv, argc);
337 rands = RANDS;
338
339 test_fib2_fib2m (count / 2, rands);
340 test_fib2m_2exp (count / 2, rands);
341
342 tests_end ();
343 exit (0);
344 }