IDA_new/gf-complete/tools/gf_poly.c

   1 /*
   2  * GF-Complete: A Comprehensive Open Source Library for Galois Field Arithmetic
   3  * James S. Plank, Ethan L. Miller, Kevin M. Greenan,
   4  * Benjamin A. Arnold, John A. Burnum, Adam W. Disney, Allen C. McBride.
   5  *
   6  * gf_poly.c - program to help find irreducible polynomials in composite fields,
   7  * using the Ben-Or algorithm.
   8  *
   9  * (This one was written by Jim)
  10  *
  11  * Please see the following paper for a description of the Ben-Or algorithm:
  12  *
  13  * author    S. Gao and D. Panario
  14  * title     Tests and Constructions of Irreducible Polynomials over Finite Fields
  15  * booktitle Foundations of Computational Mathematics
  16  * year      1997
  17  * publisher Springer Verlag
  18  * pages     346-361
  19  *
  20  * The basic technique is this.  You have a polynomial f(x) whose coefficients are
  21  * in a base field GF(2^w).  The polynomial is of degree n.  You need to do the
  22  * following for all i from 1 to n/2:
  23  *
  24  * Construct x^(2^w)^i modulo f.  That will be a polynomial of maximum degree n-1
  25  * with coefficients in GF(2^w).  You construct that polynomial by starting with x
  26  * and doubling it w times, each time taking the result modulo f.  Then you
  27  * multiply that by itself i times, again each time taking the result modulo f.
  28  *
  29  * When you're done, you need to "subtract" x -- since addition = subtraction =
  30  * XOR, that means XOR x.
  31  *
  32  * Now, find the GCD of that last polynomial and f, using Euclid's algorithm.  If
  33  * the GCD is not one, then f is reducible.  If it is not reducible for each of
  34  * those i, then it is irreducible.
  35  *
  36  * In this code, I am using a gf_general_t to represent elements of GF(2^w).  This
  37  * is so that I can use base fields that are GF(2^64) or GF(2^128).
  38  *
  39  * I have two main procedures.  The first is x_to_q_to_i_minus_x, which calculates
  40  * x^(2^w)^i - x, putting the result into a gf_general_t * called retval.
  41  *
  42  * The second is gcd_one, which takes a polynomial of degree n and a second one
  43  * of degree n-1, and uses Euclid's algorithm to decide if their GCD == 1.
  44  *
  45  * These can be made faster (e.g. calculate x^(2^w) once and store it).
  46  */
  47
  48 #include "gf_complete.h"
  49 #include "gf_method.h"
  50 #include "gf_general.h"
  51 #include "gf_int.h"
  52 #include <stdio.h>
  53 #include <stdlib.h>
  54 #include <string.h>
  55 #include <assert.h>
  56
  57 char *BM = "Bad Method: ";
  58
  59 void usage(char *s)
  60 {
  61   fprintf(stderr, "usage: gf_poly w(base-field) method power:coef [ power:coef .. ]\n");
  62   fprintf(stderr, "\n");
  63   fprintf(stderr, "       use - for the default method.\n");
  64   fprintf(stderr, "       use 0x in front of the coefficient if it's in hex\n");
  65   fprintf(stderr, "       \n");
  66   fprintf(stderr, "       For example, to test whether x^2 + 2x + 1 is irreducible\n");
  67   fprintf(stderr, "       in GF(2^16), the call is:\n");
  68   fprintf(stderr, "       \n");
  69   fprintf(stderr, "       gf_poly 16 - 2:1 1:2 0:1\n");
  70   fprintf(stderr, "       \n");
  71   fprintf(stderr, "       See the user's manual for more information.\n");
  72   if (s != NULL) {
  73     fprintf(stderr, "\n");
  74     if (s == BM) {
  75       fprintf(stderr, "%s", s);
  76       gf_error();
  77     } else {
  78       fprintf(stderr, "%s\n", s);
  79     }
  80   }
  81   exit(1);
  82 }
  83
  84 int gcd_one(gf_t *gf, int w, int n, gf_general_t *poly, gf_general_t *prod)
  85 {
  86   gf_general_t *a, *b, zero, factor, p;
  87   int i, j, da, db;
  88
  89   gf_general_set_zero(&zero, w);
  90
  91   a = (gf_general_t *) malloc(sizeof(gf_general_t) * n+1);
  92   b = (gf_general_t *) malloc(sizeof(gf_general_t) * n);
  93   for (i = 0; i <= n; i++) gf_general_add(gf, &zero, poly+i, a+i);
  94   for (i = 0; i < n; i++) gf_general_add(gf, &zero, prod+i, b+i);
  95
  96   da = n;
  97   while (1) {
  98     for (db = n-1; db >= 0 && gf_general_is_zero(b+db, w); db--) ;
  99     if (db < 0) return 0;
 100     if (db == 0) return 1;
 101     for (j = da; j >= db; j--) {
 102       if (!gf_general_is_zero(a+j, w)) {
 103         gf_general_divide(gf, a+j, b+db, &factor);
 104         for (i = 0; i <= db; i++) {
 105           gf_general_multiply(gf, b+i, &factor, &p);
 106           gf_general_add(gf, &p, a+(i+j-db), a+(i+j-db));
 107         }
 108       }
 109     }
 110     for (i = 0; i < n; i++) {
 111       gf_general_add(gf, a+i, &zero, &p);
 112       gf_general_add(gf, b+i, &zero, a+i);
 113       gf_general_add(gf, &p, &zero, b+i);
 114     }
 115   }
 116
 117 }
 118
 119 void x_to_q_to_i_minus_x(gf_t *gf, int w, int n, gf_general_t *poly, int logq, int i, gf_general_t *retval)
 120 {
 121   gf_general_t x;
 122   gf_general_t *x_to_q;
 123   gf_general_t *product;
 124   gf_general_t p, zero, factor;
 125   int j, k, lq;
 126
 127   gf_general_set_zero(&zero, w);
 128   product = (gf_general_t *) malloc(sizeof(gf_general_t) * n*2);
 129   x_to_q = (gf_general_t *) malloc(sizeof(gf_general_t) * n);
 130   for (j = 0; j < n; j++) gf_general_set_zero(x_to_q+j, w);
 131   gf_general_set_one(x_to_q+1, w);
 132
 133   for (lq = 0; lq < logq; lq++) {
 134     for (j = 0; j < n*2; j++) gf_general_set_zero(product+j, w);
 135     for (j = 0; j < n; j++) {
 136       for (k = 0; k < n; k++) {
 137         gf_general_multiply(gf, x_to_q+j, x_to_q+k, &p);
 138         gf_general_add(gf, product+(j+k), &p, product+(j+k));
 139       }
 140     }
 141     for (j = n*2-1; j >= n; j--) {
 142       if (!gf_general_is_zero(product+j, w)) {
 143         gf_general_add(gf, product+j, &zero, &factor);
 144         for (k = 0; k <= n; k++) {
 145           gf_general_multiply(gf, poly+k, &factor, &p);
 146           gf_general_add(gf, product+(j-n+k), &p, product+(j-n+k));
 147         }
 148       }
 149     }
 150     for (j = 0; j < n; j++) gf_general_add(gf, product+j, &zero, x_to_q+j);
 151   }
 152   for (j = 0; j < n; j++) gf_general_set_zero(retval+j, w);
 153   gf_general_set_one(retval, w);
 154
 155   while (i > 0) {
 156     for (j = 0; j < n*2; j++) gf_general_set_zero(product+j, w);
 157     for (j = 0; j < n; j++) {
 158       for (k = 0; k < n; k++) {
 159         gf_general_multiply(gf, x_to_q+j, retval+k, &p);
 160         gf_general_add(gf, product+(j+k), &p, product+(j+k));
 161       }
 162     }
 163     for (j = n*2-1; j >= n; j--) {
 164       if (!gf_general_is_zero(product+j, w)) {
 165         gf_general_add(gf, product+j, &zero, &factor);
 166         for (k = 0; k <= n; k++) {
 167           gf_general_multiply(gf, poly+k, &factor, &p);
 168           gf_general_add(gf, product+(j-n+k), &p, product+(j-n+k));
 169         }
 170       }
 171     }
 172     for (j = 0; j < n; j++) gf_general_add(gf, product+j, &zero, retval+j);
 173     i--;
 174   }
 175
 176   gf_general_set_one(&x, w);
 177   gf_general_add(gf, &x, retval+1, retval+1);
 178
 179   free(product);
 180   free(x_to_q);
 181 }
 182
 183 int main(int argc, char **argv)
 184 {
 185   int w, i, power, n, ap, success;
 186   gf_t gf;
 187   gf_general_t *poly, *prod;
 188   char *string, *ptr;
 189   char buf[100];
 190
 191   if (argc < 4) usage(NULL);
 192
 193   if (sscanf(argv[1], "%d", &w) != 1 || w <= 0) usage("Bad w.");
 194   ap = create_gf_from_argv(&gf, w, argc, argv, 2);
 195
 196   if (ap == 0) usage(BM);
 197
 198   if (ap == argc) usage("No powers/coefficients given.");
 199
 200   n = -1;
 201   for (i = ap; i < argc; i++) {
 202     if (strchr(argv[i], ':') == NULL || sscanf(argv[i], "%d:", &power) != 1) {
 203       string = (char *) malloc(sizeof(char)*(strlen(argv[i]+100)));
 204       sprintf(string, "Argument '%s' not in proper format of power:coefficient\n", argv[i]);
 205       usage(string);
 206     }
 207     if (power < 0) {
 208       usage("Can't have negative powers\n");
 209     } else {
 210       n = power;
 211     }
 212   }
 213   // in case the for-loop header fails
 214   assert (n >= 0);
 215
 216   poly = (gf_general_t *) malloc(sizeof(gf_general_t)*(n+1));
 217   for (i = 0; i <= n; i++) gf_general_set_zero(poly+i, w);
 218   prod = (gf_general_t *) malloc(sizeof(gf_general_t)*n);
 219
 220   for (i = ap; i < argc; i++) {
 221     sscanf(argv[i], "%d:", &power);
 222     ptr = strchr(argv[i], ':');
 223     ptr++;
 224     if (strncmp(ptr, "0x", 2) == 0) {
 225       success = gf_general_s_to_val(poly+power, w, ptr+2, 1);
 226     } else {
 227       success = gf_general_s_to_val(poly+power, w, ptr, 0);
 228     }
 229     if (success == 0) {
 230       string = (char *) malloc(sizeof(char)*(strlen(argv[i]+100)));
 231       sprintf(string, "Argument '%s' not in proper format of power:coefficient\n", argv[i]);
 232       usage(string);
 233     }
 234   }
 235
 236   printf("Poly:");
 237   for (power = n; power >= 0; power--) {
 238     if (!gf_general_is_zero(poly+power, w)) {
 239       printf("%s", (power == n) ? " " : " + ");
 240       if (!gf_general_is_one(poly+power, w)) {
 241         gf_general_val_to_s(poly+power, w, buf, 1);
 242         if (n > 0) {
 243           printf("(0x%s)", buf);
 244         } else {
 245           printf("0x%s", buf);
 246         }
 247       }
 248       if (power == 0) {
 249         if (gf_general_is_one(poly+power, w)) printf("1");
 250       } else if (power == 1) {
 251         printf("x");
 252       } else {
 253         printf("x^%d", power);
 254       }
 255     }
 256   }
 257   printf("\n");
 258
 259   if (!gf_general_is_one(poly+n, w)) {
 260     printf("\n");
 261     printf("Can't do Ben-Or, because the polynomial is not monic.\n");
 262     exit(0);
 263   }
 264
 265   for (i = 1; i <= n/2; i++) {
 266     x_to_q_to_i_minus_x(&gf, w, n, poly, w, i, prod);
 267     if (!gcd_one(&gf, w, n, poly, prod)) {
 268       printf("Reducible.\n");
 269       exit(0);
 270     }
 271   }
 272
 273   printf("Irreducible.\n");
 274   exit(0);
 275 }