1/* Test mpn_fib2m. 2 3Copyright 2018 Free Software Foundation, Inc. 4 5This file is part of the GNU MP Library test suite. 6 7The GNU MP Library test suite is free software; you can redistribute it 8and/or modify it under the terms of the GNU General Public License as 9published by the Free Software Foundation; either version 3 of the License, 10or (at your option) any later version. 11 12The GNU MP Library test suite is distributed in the hope that it will be 13useful, but WITHOUT ANY WARRANTY; without even the implied warranty of 14MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General 15Public License for more details. 16 17You should have received a copy of the GNU General Public License along with 18the 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 35static int 36test_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 170static int 171test_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 329int 330main (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} 345