l4/pkg/ocaml/contrib/otherlibs/num/test/test_bng.c

   1 /***********************************************************************/
   2 /*                                                                     */
   3 /*                           Objective Caml                            */
   4 /*                                                                     */
   5 /*            Xavier Leroy, projet Cristal, INRIA Rocquencourt         */
   6 /*                                                                     */
   7 /*  Copyright 2003 Institut National de Recherche en Informatique et   */
   8 /*  en Automatique.  All rights reserved.  This file is distributed    */
   9 /*  under the terms of the GNU Library General Public License, with    */
  10 /*  the special exception on linking described in file ../../LICENSE.  */
  11 /*                                                                     */
  12 /***********************************************************************/
  13
  14 /* $Id: test_bng.c 5900 2003-11-07 07:59:10Z xleroy $ */
  15
  16 /* Test harness for the BNG primitives.  Use BigNum as a reference. */
  17
  18 #include <stdio.h>
  19 #include <stdlib.h>
  20 #include <string.h>
  21
  22 #include <BigNum.h>
  23
  24 #include "../../../config/m.h"
  25 #include "bng.h"
  26
  27 #if defined(__GNUC__) && BNG_ASM_LEVEL > 0
  28 #if defined(BNG_ARCH_ia32)
  29 #include "bng_ia32.c"
  30 #elif defined(BNG_ARCH_amd64)
  31 #include "bng_amd64.c"
  32 #elif defined(BNG_ARCH_ppc)
  33 #include "bng_ppc.c"
  34 #elif defined (BNG_ARCH_alpha)
  35 #include "bng_alpha.c"
  36 #elif defined (BNG_ARCH_sparc)
  37 #include "bng_sparc.c"
  38 #elif defined (BNG_ARCH_mips)
  39 #include "bng_mips.c"
  40 #endif
  41 #endif
  42
  43 #include "bng_digit.c"
  44
  45 /* Random generator for digits.  Can either generate "true" PRN numbers
  46    or numbers consisting of long sequences of 0 and 1 bits. */
  47
  48 static int rand_skewed = 0;
  49 static int rand_runlength = 0;
  50 static int rand_bit = 0;
  51 static bngdigit rand_seed = 0;
  52
  53 static bngdigit randdigit(void)
  54 {
  55   bngdigit res;
  56   int i;
  57
  58   if (rand_skewed) {
  59     for (i = 0, res = 0; i < BNG_BITS_PER_DIGIT; i++) {
  60       if (rand_runlength == 0) {
  61         rand_runlength = 1 + (rand() % (2 * BNG_BITS_PER_DIGIT));
  62         rand_bit ^= 1;
  63       }
  64       res = (res << 1) | rand_bit;
  65       rand_runlength--;
  66     }
  67     return res;
  68   } else {
  69     rand_seed = rand_seed * 69069 + 25173;
  70     return rand_seed;
  71   }
  72 }
  73
  74 /* Test the operations on digits.
  75    This uses double-width integer arithmetic as reference.
  76    This is only available on 32-bit platforms that support a 64-bit int type.
  77 */
  78
  79 #if defined(ARCH_UINT64_TYPE) && !defined(ARCH_SIXTYFOUR)
  80
  81 typedef ARCH_UINT64_TYPE dbldigit;
  82
  83 static int test_digit_ops(int i)
  84 {
  85   bngdigit a1, a2, a3, r1, r2;
  86   int ci, co, n;
  87
  88   a1 = randdigit();
  89   a2 = randdigit();
  90   a3 = randdigit();
  91   ci = randdigit() & 1;
  92
  93   BngAdd2(r1,co,a1,a2);
  94   if ((dbldigit) r1 + ((dbldigit) co << BNG_BITS_PER_DIGIT)
  95       != (dbldigit) a1 + (dbldigit) a2) {
  96     printf("Round %d, BngAdd2(%lx,%x,%lx, %lx)\n", i, r1, co, a1, a2);
  97     return 1;
  98   }
  99
 100   BngAdd2Carry(r1,co,a1,a2,ci);
 101   if ((dbldigit) r1 + ((dbldigit) co << BNG_BITS_PER_DIGIT)
 102       != (dbldigit) a1 + (dbldigit) a2 + (dbldigit) ci) {
 103     printf("Round %d, BngAdd2Carry(%lx,%x,%lx, %lx, %x)\n", i, r1, co, a1, a2, ci);
 104     return 1;
 105   }
 106
 107   r2 = 0;
 108   BngAdd3(r1,r2,a1,a2,a3);
 109   if ((dbldigit) r1 + ((dbldigit) r2 << BNG_BITS_PER_DIGIT)
 110       != (dbldigit) a1 + (dbldigit) a2 + (dbldigit) a3) {
 111     printf("Round %d, BngAdd3(%lx,%x,%lx, %lx, %lx)\n", i, r1, co, a1, a2, a3);
 112     return 1;
 113   }
 114
 115   BngSub2(r1,co,a1,a2);
 116   if ((dbldigit) r1 - ((dbldigit) co << BNG_BITS_PER_DIGIT)
 117       != (dbldigit) a1 - (dbldigit) a2) {
 118     printf("Round %d, BngSub2(%lx,%x,%lx, %lx)\n", i, r1, co, a1, a2);
 119     return 1;
 120   }
 121
 122   BngSub2Carry(r1,co,a1,a2,ci);
 123   if ((dbldigit) r1 - ((dbldigit) co << BNG_BITS_PER_DIGIT)
 124       != (dbldigit) a1 - (dbldigit) a2 - (dbldigit) ci) {
 125     printf("Round %d, BngSub2Carry(%lx,%x,%lx, %lx, %x)\n", i, r1, co, a1, a2, ci);
 126     return 1;
 127   }
 128
 129   r2 = 0;
 130   BngSub3(r1,r2,a1,a2,a3);
 131   if ((dbldigit) r1 - ((dbldigit) r2 << BNG_BITS_PER_DIGIT)
 132       != (dbldigit) a1 - (dbldigit) a2 - (dbldigit) a3) {
 133     printf("Round %d, BngSub3(%lx,%x,%lx, %lx, %lx)\n", i, r1, co, a1, a2, a3);
 134     return 1;
 135   }
 136
 137   BngMult(r1,r2,a1,a2);
 138   if ((((dbldigit) r1 << BNG_BITS_PER_DIGIT) | (dbldigit) r2)
 139       != (dbldigit) a1 * (dbldigit) a2) {
 140     printf("Round %d, BngMult(%lx,%lx,%lx, %lx)\n", i, r1, r2, a1, a2);
 141     return 1;
 142   }
 143
 144   /* Make sure a3 is normalized */
 145   a3 |= 1L << (BNG_BITS_PER_DIGIT - 1);
 146   if (a1 < a3) {
 147     BngDiv(r1,r2,a1,a2,a3);
 148     if (r1 != (((dbldigit) a1 << BNG_BITS_PER_DIGIT) | (dbldigit) a2) / a3
 149         ||
 150         r2 != (((dbldigit) a1 << BNG_BITS_PER_DIGIT) | (dbldigit) a2) % a3)
 151       {
 152         printf("Round %d, BngDiv(%lx,%lx,%lx, %lx, %lx)\n", i, r1, r2, a1, a2, a3);
 153         return 1;
 154       }
 155   }
 156
 157   n = bng_leading_zero_bits(a1);
 158   if (a1 == 0) {
 159     if (n != BNG_BITS_PER_DIGIT) {
 160       printf("Round %d, bng_leading_zero(bits(%lx) = %d", i, a1, n);
 161       return 1;
 162     }
 163   } else {
 164     if ((a1 << n) >> n != a1 ||
 165         ((a1 << n) & (1L << (BNG_BITS_PER_DIGIT - 1))) == 0) {
 166       printf("Round %d, bng_leading_zero(bits(%lx) = %d", i, a1, n);
 167       return 1;
 168     }
 169   }
 170   return 0;
 171 }
 172
 173 #endif
 174
 175 /* Test the bng operations.  Use BigNum as a reference. */
 176
 177 #define MAX_DIGITS 32
 178
 179 void randbng(bng a, bngsize n)
 180 {
 181   int i;
 182   for (i = 0; i < n; i++) a[i] = randdigit();
 183 }
 184
 185 char * bng2string(bng a, bngsize n)
 186 {
 187   char * buffer = malloc((BNG_BITS_PER_DIGIT / 4 + 1) * MAX_DIGITS);
 188   char temp[BNG_BITS_PER_DIGIT / 4 + 1];
 189   int i;
 190
 191   buffer[0] = 0;
 192   for (i = n - 1; i >= 0; i--) {
 193     sprintf(temp, "%lx", a[i]);
 194     strcat(buffer, temp);
 195     if (i > 0) strcat(buffer, "_");
 196   }
 197   return buffer;
 198 }
 199
 200 int bngsame(bng a, bng b, bngsize n)
 201 {
 202   int i;
 203   for (i = 0; i < n; i++)
 204     if (a[i] != b[i]) return 0;
 205   return 1;
 206 }
 207
 208 int test_bng_ops(int i)
 209 {
 210   bngsize p, q;
 211   bngdigit a[MAX_DIGITS], b[MAX_DIGITS], c[MAX_DIGITS], d[MAX_DIGITS];
 212   bngdigit f[2 * MAX_DIGITS], g[2 * MAX_DIGITS], h[2 * MAX_DIGITS];
 213   bngcarry ci, co, cp;
 214   bngdigit dg, do_, dp;
 215   int amount;
 216
 217   /* Determine random lengths p and q between 1 and MAX_DIGITS.
 218      Ensure p >= q. */
 219   p = 1 + (rand() % MAX_DIGITS);
 220   q = 1 + (rand() % MAX_DIGITS);
 221   if (q > p) { bngsize t = p; p = q; q = t; }
 222
 223   /* Randomly generate bignums a of size p, b of size q */
 224   randbng(a, p);
 225   randbng(b, q);
 226   ci = rand() & 1;
 227
 228   /* comparison */
 229   co = bng_compare(a, p, b, q);
 230   cp = BnnCompare(a, p, b, q);
 231   if (co != cp) {
 232     printf("Round %d, bng_compare(%s, %ld, %s, %ld) = %d\n",
 233            i, bng2string(a, p), p, bng2string(b, q), q, co);
 234     return 1;
 235   }
 236   co = bng_compare(b, q, a, p);
 237   cp = BnnCompare(b, q, a, p);
 238   if (co != cp) {
 239     printf("Round %d, bng_compare(%s, %ld, %s, %ld) = %d\n",
 240            i, bng2string(b, q), q, bng2string(a, p), p, co);
 241     return 1;
 242   }
 243   /* add carry */
 244   bng_assign(c, a, p);
 245   co = bng_add_carry(c, p, ci);
 246   BnnAssign(d, a, p);
 247   cp = BnnAddCarry(d, p, ci);
 248   if (co != cp || !bngsame(c, d, p)) {
 249     printf("Round %d, bng_add_carry(%s, %ld, %d) -> %s, %d\n",
 250            i, bng2string(a, p), p, ci, bng2string(c, p), co);
 251     return 1;
 252   }
 253   /* add */
 254   bng_assign(c, a, p);
 255   co = bng_add(c, p, b, q, ci);
 256   BnnAssign(d, a, p);
 257   cp = BnnAdd(d, p, b, q, ci);
 258   if (co != cp || !bngsame(c, d, p)) {
 259     printf("Round %d, bng_add(%s, %ld, %s, %ld, %d) -> %s, %d\n",
 260            i, bng2string(a, p), p, bng2string(b, q), q, ci,
 261            bng2string(c, p), co);
 262     return 1;
 263   }
 264   /* sub carry */
 265   bng_assign(c, a, p);
 266   co = bng_sub_carry(c, p, ci);
 267   BnnAssign(d, a, p);
 268   cp = BnnSubtractBorrow(d, p, ci ^ 1) ^ 1;
 269   if (co != cp || !bngsame(c, d, p)) {
 270     printf("Round %d, bng_sub_carry(%s, %ld, %d) -> %s, %d\n",
 271            i, bng2string(a, p), p, ci, bng2string(c, p), co);
 272     return 1;
 273   }
 274   /* sub */
 275   bng_assign(c, a, p);
 276   co = bng_sub(c, p, b, q, ci);
 277   BnnAssign(d, a, p);
 278   cp = BnnSubtract(d, p, b, q, ci ^ 1) ^ 1;
 279   if (co != cp || !bngsame(c, d, p)) {
 280     printf("Round %d, bng_sub(%s, %ld, %s, %ld, %d) -> %s, %d\n",
 281            i, bng2string(a, p), p, bng2string(b, q), q, ci,
 282            bng2string(c, p), co);
 283     return 1;
 284   }
 285   /* shift left */
 286   amount = rand() % BNG_BITS_PER_DIGIT;
 287   bng_assign(c, a, p);
 288   do_ = bng_shift_left(c, p, amount);
 289   BnnAssign(d, a, p);
 290   dp = BnnShiftLeft(d, p, amount);
 291   if (do_ != dp || !bngsame(c, d, p)) {
 292     printf("Round %d, bng_shift_left(%s, %ld, %d) -> %s, %ld\n",
 293            i, bng2string(a, p), p, amount, bng2string(c, p), do_);
 294     return 1;
 295   }
 296   /* shift right */
 297   amount = rand() % BNG_BITS_PER_DIGIT;
 298   bng_assign(c, a, p);
 299   do_ = bng_shift_right(c, p, amount);
 300   BnnAssign(d, a, p);
 301   dp = BnnShiftRight(d, p, amount);
 302   if (do_ != dp || !bngsame(c, d, p)) {
 303     printf("Round %d, bng_shift_right(%s, %ld, %d) -> %s, %ld\n",
 304            i, bng2string(a, p), p, amount, bng2string(c, p), do_);
 305     return 1;
 306   }
 307   /* mult_add_digit */
 308   dg = randdigit();
 309   if (p >= q + 1) {
 310     bng_assign(c, a, p);
 311     co = bng_mult_add_digit(c, p, b, q, dg);
 312     BnnAssign(d, a, p);
 313     cp = BnnMultiplyDigit(d, p, b, q, dg);
 314     if (co != cp || !bngsame(c, d, p)) {
 315       printf("Round %d, bng_mult_add_digit(%s, %ld, %s, %ld, %ld) -> %s, %d\n",
 316              i, bng2string(a, p), p, bng2string(b, q), q, dg,
 317              bng2string(c, p), co);
 318       return 1;
 319     }
 320   }
 321   /* mult_sub_digit */
 322   dg = randdigit();
 323   bng_assign(c, a, p);
 324   do_ = bng_mult_add_digit(c, p, b, q, dg);
 325   bng_assign(d, c, p);
 326   dp = bng_mult_sub_digit(d, p, b, q, dg);
 327   if (do_ != dp || !bngsame(a, d, p)) {
 328     printf("Round %d, bng_mult_sub_digit(%s, %ld, %s, %ld, %ld) -> %s, %ld\n",
 329            i, bng2string(c, p), p, bng2string(b, q), q, dg,
 330            bng2string(d, p), dp);
 331     return 1;
 332   }
 333   /* mult_add */
 334   randbng(f, 2*p);
 335   bng_assign(g, f, 2*p);
 336   co = bng_mult_add(g, 2*p, a, p, b, q);
 337   BnnAssign(h, f, 2*p);
 338   cp = BnnMultiply(h, 2*p, a, p, b, q);
 339   if (co != cp || !bngsame(g, h, 2*p)) {
 340     printf("Round %d, bng_mult_add(%s, %ld, %s, %ld, %s, %ld) -> %s, %d\n",
 341            i, bng2string(f, 2*p), 2*p,
 342            bng2string(a, p), p,
 343            bng2string(b, q), q,
 344            bng2string(g, 2*p), co);
 345     return 1;
 346   }
 347   /* square_add */
 348   randbng(f, 2*p);
 349   bng_assign(g, f, 2*p);
 350   co = bng_square_add(g, 2*p, b, q);
 351   BnnAssign(h, f, 2*p);
 352   cp = BnnAdd(h, 2*p, h, 2*p);
 353   cp += BnnMultiply(h, 2*p, b, q, b, q);
 354   if (co != cp || !bngsame(g, h, 2*p)) {
 355     printf("Round %d, bng_square_add(%s, %ld, %s, %ld) -> %s, %d\n",
 356            i, bng2string(f, 2*p), 2*p,
 357            bng2string(b, q), q,
 358            bng2string(g, 2*p), co);
 359     return 1;
 360   }
 361   /* div_rem_digit */
 362   if (a[p - 1] < dg) {
 363     do_ = bng_div_rem_digit(c, a, p, dg);
 364     dp = BnnDivideDigit(d, a, p, dg);
 365     if (do_ != dp || !bngsame(c, d, p-1)) {
 366       printf("Round %d, bng_div_rem_digit(%s, %s, %ld, %lx) -> %lx\n",
 367              i, bng2string(d, p-1), bng2string(a, p), p, dg, do_);
 368       return 1;
 369     }
 370   }
 371   /* div_rem */
 372   if (p > q && a[p - 1] < b[q - 1]) {
 373     bng_assign(c, a, p);
 374     bng_div_rem(c, p, b, q);
 375     BnnAssign(d, a, p);
 376     BnnDivide(d, p, b, q);
 377     if (!bngsame(c, d, p)) {
 378       printf("Round %d, bng_div_rem(%s, %ld, %s, %ld) -> %s, %s\n",
 379              i, bng2string(a, p), p, bng2string(b, q), q,
 380              bng2string(c + q, p - q),
 381              bng2string(c, q));
 382       return 1;
 383     }
 384   }
 385   return 0;
 386 }
 387
 388 int main(int argc, char ** argv)
 389 {
 390   int niter = 100000;
 391   int i, err;
 392
 393   bng_init();
 394   if (argc >= 2) niter = atoi(argv[1]);
 395 #if defined(ARCH_UINT64_TYPE) && !defined(ARCH_SIXTYFOUR)
 396   printf("Testing single-digit operations\n");
 397   for (err = 0, i = 1; i < niter; i++) err += test_digit_ops(i);
 398   printf("%d rounds performed, %d errors found\n", niter, err);
 399 #endif
 400   printf("Testing bignum operations\n");
 401   for (err = 0, i = 1; i < niter; i++) err += test_bng_ops(i);
 402   printf("%d rounds performed, %d errors found\n", niter, err);
 403   printf("Testing bignum operations with skewed PRNG\n");
 404   rand_skewed = 1;
 405   for (err = 0, i = 1; i < niter; i++) err += test_bng_ops(i);
 406   printf("%d rounds performed, %d errors found\n", niter, err);
 407   return 0;
 408 }