source: git/ntl/src/GF2X.c @ e6d2f67

fieker-DuValspielwiese
Last change on this file since e6d2f67 was de6a29, checked in by Hans Schönemann <hannes@…>, 19 years ago
* hannes: NTL-5.4 git-svn-id: file:///usr/local/Singular/svn/trunk@8693 2c84dea3-7e68-4137-9b89-c4e89433aadc
  • Property mode set to 100644
File size: 33.9 KB
Line 
1
2#include <NTL/GF2X.h>
3#include <NTL/vec_long.h>
4
5#include <NTL/new.h>
6
7NTL_START_IMPL
8
9long GF2X::HexOutput = 0;
10
11
12void GF2X::SetMaxLength(long n)
13{
14   if (n < 0) Error("GF2X::SetMaxLength: negative length");
15   if (NTL_OVERFLOW(n, 1, 0))
16      Error("GF2X::SetMaxLength: excessive length");
17   long w = (n + NTL_BITS_PER_LONG - 1)/NTL_BITS_PER_LONG;
18   xrep.SetMaxLength(w);
19}
20
21GF2X::GF2X(INIT_SIZE_TYPE, long n)
22{
23   SetMaxLength(n);
24}
25
26
27
28const GF2X& GF2X::zero()
29{
30   static GF2X z;
31   return z;
32}
33
34void GF2X::normalize()
35{
36   long n;
37   const _ntl_ulong *p;
38
39   n = xrep.length();
40   if (n == 0) return;
41   p = xrep.elts() + n;
42   while (n > 0 && (*--p) == 0) {
43      n--;
44   }
45   xrep.QuickSetLength(n);
46}
47
48long IsZero(const GF2X& a) 
49   { return a.xrep.length() == 0; }
50
51long IsOne(const GF2X& a)
52   { return a.xrep.length() == 1 && a.xrep[0] == 1; }
53
54
55
56
57
58
59
60long IsX(const GF2X& a)
61{
62   return a.xrep.length() == 1 && a.xrep[0] == 2;
63}
64
65GF2 coeff(const GF2X& a, long i)
66{
67   if (i < 0) return to_GF2(0);
68   long wi = i/NTL_BITS_PER_LONG;
69   if (wi >= a.xrep.length()) return to_GF2(0);
70   long bi = i - wi*NTL_BITS_PER_LONG;
71
72   return to_GF2((a.xrep[wi] & (1UL << bi)) != 0);
73}
74
75GF2 LeadCoeff(const GF2X& a)
76{
77   if (IsZero(a))
78      return to_GF2(0);
79   else
80      return to_GF2(1);
81}
82
83GF2 ConstTerm(const GF2X& a)
84{
85   if (IsZero(a))
86      return to_GF2(0);
87   else
88      return to_GF2((a.xrep[0] & 1) != 0);
89}
90
91
92void set(GF2X& x)
93{
94   x.xrep.SetLength(1);
95   x.xrep[0] = 1;
96}
97
98void SetX(GF2X& x)
99{
100   x.xrep.SetLength(1);
101   x.xrep[0] = 2;
102}
103
104void SetCoeff(GF2X& x, long i)
105{
106   if (i < 0) {
107      Error("SetCoeff: negative index");
108      return;
109   }
110
111   long n, j;
112
113   n = x.xrep.length();
114   long wi = i/NTL_BITS_PER_LONG;
115
116   if (wi >= n) {
117      x.xrep.SetLength(wi+1);
118      for (j = n; j <= wi; j++)
119         x.xrep[j] = 0;
120   }
121
122   long bi = i - wi*NTL_BITS_PER_LONG;
123
124   x.xrep[wi] |= (1UL << bi);
125}
126   
127
128
129void SetCoeff(GF2X& x, long i, long val)
130{
131   if (i < 0) {
132      Error("SetCoeff: negative index");
133      return;
134   }
135
136   val = val & 1;
137
138   if (val) {
139      SetCoeff(x, i);
140      return;
141   }
142
143   // we want to clear position i
144
145   long n;
146
147   n = x.xrep.length();
148   long wi = i/NTL_BITS_PER_LONG;
149
150   if (wi >= n) 
151      return;
152
153   long bi = i - wi*NTL_BITS_PER_LONG;
154
155   x.xrep[wi] &= ~(1UL << bi);
156   if (wi == n-1) x.normalize();
157}
158
159void SetCoeff(GF2X& x, long i, GF2 a)
160{
161   SetCoeff(x, i, rep(a));
162}
163
164
165void swap(GF2X& a, GF2X& b)
166{
167   swap(a.xrep, b.xrep);
168}
169
170
171
172long deg(const GF2X& aa)
173{
174   long n = aa.xrep.length();
175
176   if (n == 0)
177      return -1;
178
179   _ntl_ulong a = aa.xrep[n-1];
180   long i = 0;
181
182   if (a == 0) Error("GF2X: unnormalized polynomial detected in deg");
183
184   while (a>=256)
185      i += 8, a >>= 8;
186   if (a >=16)
187      i += 4, a >>= 4;
188   if (a >= 4)
189      i += 2, a >>= 2;
190   if (a >= 2)
191      i += 2;
192   else if (a >= 1)
193      i++;
194
195   return NTL_BITS_PER_LONG*(n-1) + i - 1;
196}
197   
198
199long operator==(const GF2X& a, const GF2X& b)
200{
201   return a.xrep == b.xrep;
202}
203
204long operator==(const GF2X& a, long b)
205{
206   if (b & 1) 
207      return IsOne(a);
208   else
209      return IsZero(a);
210}
211
212long operator==(const GF2X& a, GF2 b)
213{
214   if (b == 1) 
215      return IsOne(a);
216   else
217      return IsZero(a);
218}
219
220void random(GF2X& x, long n)
221{
222   if (n < 0) Error("GF2X random: negative length");
223
224   if (NTL_OVERFLOW(n, 1, 0))
225      Error("GF2X random: excessive length");
226
227   long wl = (n+NTL_BITS_PER_LONG-1)/NTL_BITS_PER_LONG;
228
229   x.xrep.SetLength(wl);
230
231   long i;
232   for (i = 0; i < wl-1; i++) {
233      x.xrep[i] = RandomWord();
234   }
235
236   if (n > 0) {
237      long pos = n % NTL_BITS_PER_LONG;
238      if (pos == 0) pos = NTL_BITS_PER_LONG;
239      x.xrep[wl-1] = RandomBits_ulong(pos);
240   }
241
242   x.normalize();
243}
244
245void add(GF2X& x, const GF2X& a, const GF2X& b)
246{
247   long sa = a.xrep.length();
248   long sb = b.xrep.length();
249
250   long i;
251
252   if (sa == sb) {
253      x.xrep.SetLength(sa);
254      if (sa == 0) return;
255
256      _ntl_ulong *xp = x.xrep.elts();
257      const _ntl_ulong *ap = a.xrep.elts();
258      const _ntl_ulong *bp = b.xrep.elts();
259
260      for (i = 0; i < sa; i++)
261         xp[i] = ap[i] ^ bp[i];
262
263      i = sa-1;
264      while (i >= 0 && !xp[i]) i--;
265      x.xrep.QuickSetLength(i+1);
266   }
267   
268   else if (sa < sb) {
269      x.xrep.SetLength(sb);
270      _ntl_ulong *xp = x.xrep.elts();
271      const _ntl_ulong *ap = a.xrep.elts();
272      const _ntl_ulong *bp = b.xrep.elts();
273
274      for (i = 0; i < sa; i++)
275         xp[i] = ap[i] ^ bp[i];
276
277      for (; i < sb; i++)
278         xp[i] = bp[i];
279   }
280   else { // sa > sb
281      x.xrep.SetLength(sa);
282      _ntl_ulong *xp = x.xrep.elts();
283      const _ntl_ulong *ap = a.xrep.elts();
284      const _ntl_ulong *bp = b.xrep.elts();
285
286      for (i = 0; i < sb; i++)
287         xp[i] = ap[i] ^ bp[i];
288
289      for (; i < sa; i++)
290         xp[i] = ap[i];
291   }
292}
293
294
295
296
297
298/*
299 * The bodies of mul1, Mul1, AddMul1, and mul_half
300 * are generated by the MakeDesc program, and the
301 * macros NTL_BB_MUL_CODE... are defined in mach_desc.h.
302 * Thanks to Paul Zimmermann for providing improvements
303 * to this approach.
304 */
305
306
307#if (defined(NTL_GF2X_ALTCODE1))
308
309#define NTL_EFF_BB_MUL_CODE0 NTL_ALT1_BB_MUL_CODE0
310#define NTL_EFF_BB_MUL_CODE1 NTL_ALT1_BB_MUL_CODE1
311#define NTL_EFF_BB_MUL_CODE2 NTL_ALT1_BB_MUL_CODE2
312#define NTL_EFF_SHORT_BB_MUL_CODE1 NTL_ALT1_SHORT_BB_MUL_CODE1
313#define NTL_EFF_HALF_BB_MUL_CODE0 NTL_ALT1_HALF_BB_MUL_CODE0
314
315#elif (defined(NTL_GF2X_ALTCODE))
316
317#define NTL_EFF_BB_MUL_CODE0 NTL_ALT_BB_MUL_CODE0
318#define NTL_EFF_BB_MUL_CODE1 NTL_ALT_BB_MUL_CODE1
319#define NTL_EFF_BB_MUL_CODE2 NTL_ALT_BB_MUL_CODE2
320#define NTL_EFF_SHORT_BB_MUL_CODE1 NTL_ALT_SHORT_BB_MUL_CODE1
321#define NTL_EFF_HALF_BB_MUL_CODE0 NTL_ALT_HALF_BB_MUL_CODE0
322
323#else
324
325#define NTL_EFF_BB_MUL_CODE0 NTL_BB_MUL_CODE0
326#define NTL_EFF_BB_MUL_CODE1 NTL_BB_MUL_CODE1
327#define NTL_EFF_BB_MUL_CODE2 NTL_BB_MUL_CODE2
328#define NTL_EFF_SHORT_BB_MUL_CODE1 NTL_SHORT_BB_MUL_CODE1
329#define NTL_EFF_HALF_BB_MUL_CODE0 NTL_HALF_BB_MUL_CODE0
330
331#endif
332
333
334
335static 
336void mul1(_ntl_ulong *c, _ntl_ulong a, _ntl_ulong b)
337{
338
339NTL_EFF_BB_MUL_CODE0
340
341
342}
343
344
345#ifdef NTL_GF2X_NOINLINE
346
347#define mul1_IL mul1
348
349#else
350
351static inline
352void mul1_inline(_ntl_ulong *c, _ntl_ulong a, _ntl_ulong b)
353{
354
355
356NTL_EFF_BB_MUL_CODE0
357
358
359}
360
361#define mul1_IL mul1_inline
362#endif
363
364
365static 
366void Mul1(_ntl_ulong *cp, const _ntl_ulong *bp, long sb, _ntl_ulong a)
367{
368 
369
370NTL_EFF_BB_MUL_CODE1
371
372
373}
374
375static 
376void AddMul1(_ntl_ulong *cp, const _ntl_ulong* bp, long sb, _ntl_ulong a)
377{
378
379
380NTL_EFF_BB_MUL_CODE2
381
382
383}
384
385
386static 
387void Mul1_short(_ntl_ulong *cp, const _ntl_ulong *bp, long sb, _ntl_ulong a)
388{
389 
390
391NTL_EFF_SHORT_BB_MUL_CODE1
392
393
394}
395
396
397
398
399
400static 
401void mul_half(_ntl_ulong *c, _ntl_ulong a, _ntl_ulong b)
402{
403
404
405NTL_EFF_HALF_BB_MUL_CODE0
406
407
408}
409
410
411// mul2...mul8 hard-code 2x2...8x8 word multiplies.
412// I adapted these routines from LiDIA (except mul3, see below).
413// Inlining these seems to hurt, not help.
414
415static
416void mul2(_ntl_ulong *c, const _ntl_ulong *a, const _ntl_ulong *b)
417{
418   _ntl_ulong hs0, hs1;
419   _ntl_ulong hl2[2];
420
421   hs0 = a[0] ^ a[1];
422   hs1 = b[0] ^ b[1];
423
424   mul1_IL(c, a[0], b[0]);
425   mul1_IL(c+2, a[1], b[1]);
426   mul1_IL(hl2, hs0, hs1);
427
428
429   hl2[0] = hl2[0] ^ c[0] ^ c[2];
430   hl2[1] = hl2[1] ^ c[1] ^ c[3];
431
432   c[1] ^= hl2[0];
433   c[2] ^= hl2[1];
434}
435
436
437/*
438 * This version of mul3 I got from Weimerskirch, Stebila,
439 * and Shantz, "Generic GF(2^m) arithmetic in software
440 * an its application to ECC" (ACISP 2003).
441 */
442
443static
444void mul3 (_ntl_ulong *c, const _ntl_ulong *a, const _ntl_ulong *b)
445{
446   _ntl_ulong d0[2], d1[2], d2[2], d01[2], d02[2], d12[2];
447
448   mul1_IL(d0, a[0], b[0]);
449   mul1_IL(d1, a[1], b[1]);
450   mul1_IL(d2, a[2], b[2]);
451   mul1_IL(d01, a[0]^a[1], b[0]^b[1]);
452   mul1_IL(d02, a[0]^a[2], b[0]^b[2]);
453   mul1_IL(d12, a[1]^a[2], b[1]^b[2]);
454
455   
456   c[0] = d0[0];
457   c[1] = d0[1] ^ d01[0] ^ d1[0] ^ d0[0];
458   c[2] = d01[1] ^ d1[1] ^ d0[1] ^ d02[0] ^ d2[0] ^ d0[0] ^ d1[0];
459   c[3] = d02[1] ^ d2[1] ^ d0[1] ^ d1[1] ^ d12[0] ^ d1[0] ^ d2[0];
460   c[4] = d12[1] ^ d1[1] ^ d2[1] ^ d2[0];
461   c[5] = d2[1];
462
463}
464
465static
466void mul4(_ntl_ulong *c, const _ntl_ulong *a, const _ntl_ulong *b)
467{
468   _ntl_ulong hs0[2], hs1[2];
469   _ntl_ulong hl2[4];
470
471   hs0[0] = a[0] ^ a[2];
472   hs0[1] = a[1] ^ a[3];
473   hs1[0] = b[0] ^ b[2];
474   hs1[1] = b[1] ^ b[3];
475
476   mul2(c, a, b);
477   mul2(c+4, a+2, b+2); 
478   mul2(hl2, hs0, hs1);
479
480   hl2[0] = hl2[0] ^ c[0] ^ c[4];
481   hl2[1] = hl2[1] ^ c[1] ^ c[5];
482   hl2[2] = hl2[2] ^ c[2] ^ c[6];
483   hl2[3] = hl2[3] ^ c[3] ^ c[7];
484   
485   c[2] ^= hl2[0];
486   c[3] ^= hl2[1];
487   c[4] ^= hl2[2];
488   c[5] ^= hl2[3];
489}
490
491static
492void mul5 (_ntl_ulong *c, const _ntl_ulong *a, const _ntl_ulong *b)
493{
494   _ntl_ulong hs0[3], hs1[3];
495   _ntl_ulong hl2[6];
496
497   hs0[0] = a[0] ^ a[3]; 
498   hs0[1] = a[1] ^ a[4];
499   hs0[2] = a[2];
500   hs1[0] = b[0] ^ b[3]; 
501   hs1[1] = b[1] ^ b[4];
502   hs1[2] = b[2];
503
504   mul3(c, a, b); 
505   mul3(hl2, hs0, hs1);
506   mul2(c+6, a+3, b+3);
507
508   hl2[0] = hl2[0] ^ c[0] ^ c[6]; 
509   hl2[1] = hl2[1] ^ c[1] ^ c[7];
510   hl2[2] = hl2[2] ^ c[2] ^ c[8];
511   hl2[3] = hl2[3] ^ c[3] ^ c[9];
512   hl2[4] = hl2[4] ^ c[4];
513   hl2[5] = hl2[5] ^ c[5];
514
515 
516   c[3] ^= hl2[0];
517   c[4] ^= hl2[1];
518   c[5] ^= hl2[2];
519   c[6] ^= hl2[3];
520   c[7] ^= hl2[4];
521   c[8] ^= hl2[5];
522}
523
524static
525void mul6(_ntl_ulong *c, const _ntl_ulong *a, const _ntl_ulong *b)
526{
527   _ntl_ulong hs0[3], hs1[3];
528   _ntl_ulong hl2[6];
529
530   hs0[0] = a[0] ^ a[3];   
531   hs0[1] = a[1] ^ a[4];
532   hs0[2] = a[2] ^ a[5];
533   hs1[0] = b[0] ^ b[3];   
534   hs1[1] = b[1] ^ b[4];
535   hs1[2] = b[2] ^ b[5];
536
537   mul3(c, a, b);   
538   mul3(c+6, a+3, b+3); 
539   mul3(hl2, hs0, hs1); 
540
541   hl2[0] = hl2[0] ^ c[0] ^ c[6];
542   hl2[1] = hl2[1] ^ c[1] ^ c[7];
543   hl2[2] = hl2[2] ^ c[2] ^ c[8];
544   hl2[3] = hl2[3] ^ c[3] ^ c[9];
545   hl2[4] = hl2[4] ^ c[4] ^ c[10];
546   hl2[5] = hl2[5] ^ c[5] ^ c[11];
547   
548   c[3] ^= hl2[0];
549   c[4] ^= hl2[1];
550   c[5] ^= hl2[2];
551   c[6] ^= hl2[3];
552   c[7] ^= hl2[4];
553   c[8] ^= hl2[5];
554}
555
556static
557void mul7(_ntl_ulong *c, const _ntl_ulong *a, const _ntl_ulong *b)
558{
559   _ntl_ulong hs0[4], hs1[4];
560   _ntl_ulong hl2[8];
561
562   hs0[0] = a[0] ^ a[4]; 
563   hs0[1] = a[1] ^ a[5];
564   hs0[2] = a[2] ^ a[6];
565   hs0[3] = a[3];
566   hs1[0] = b[0] ^ b[4]; 
567   hs1[1] = b[1] ^ b[5];
568   hs1[2] = b[2] ^ b[6];
569   hs1[3] = b[3];
570
571   mul4(c, a, b); 
572   mul4(hl2, hs0, hs1);
573   mul3(c+8, a+4, b+4); 
574
575   hl2[0] = hl2[0] ^ c[0] ^ c[8];
576   hl2[1] = hl2[1] ^ c[1] ^ c[9];
577   hl2[2] = hl2[2] ^ c[2] ^ c[10];
578   hl2[3] = hl2[3] ^ c[3] ^ c[11];
579   hl2[4] = hl2[4] ^ c[4] ^ c[12];
580   hl2[5] = hl2[5] ^ c[5] ^ c[13];
581   hl2[6] = hl2[6] ^ c[6];
582   hl2[7] = hl2[7] ^ c[7];
583   
584   c[4]  ^= hl2[0];
585   c[5]  ^= hl2[1];
586   c[6]  ^= hl2[2];
587   c[7]  ^= hl2[3];
588   c[8]  ^= hl2[4];
589   c[9]  ^= hl2[5];
590   c[10] ^= hl2[6];
591   c[11] ^= hl2[7];
592}
593
594static
595void mul8(_ntl_ulong *c, const _ntl_ulong *a, const _ntl_ulong *b)
596{
597   _ntl_ulong hs0[4], hs1[4];
598   _ntl_ulong hl2[8];
599
600   hs0[0] = a[0] ^ a[4]; 
601   hs0[1] = a[1] ^ a[5];
602   hs0[2] = a[2] ^ a[6];
603   hs0[3] = a[3] ^ a[7];
604   hs1[0] = b[0] ^ b[4]; 
605   hs1[1] = b[1] ^ b[5];
606   hs1[2] = b[2] ^ b[6];
607   hs1[3] = b[3] ^ b[7];
608
609   mul4(c, a, b); 
610   mul4(c+8, a+4, b+4);
611   mul4(hl2, hs0, hs1); 
612
613   hl2[0] = hl2[0] ^ c[0] ^ c[8];
614   hl2[1] = hl2[1] ^ c[1] ^ c[9];
615   hl2[2] = hl2[2] ^ c[2] ^ c[10];
616   hl2[3] = hl2[3] ^ c[3] ^ c[11];
617   hl2[4] = hl2[4] ^ c[4] ^ c[12];
618   hl2[5] = hl2[5] ^ c[5] ^ c[13];
619   hl2[6] = hl2[6] ^ c[6] ^ c[14];
620   hl2[7] = hl2[7] ^ c[7] ^ c[15];
621   
622   c[4]  ^= hl2[0];
623   c[5]  ^= hl2[1];
624   c[6]  ^= hl2[2];
625   c[7]  ^= hl2[3];
626   c[8]  ^= hl2[4];
627   c[9]  ^= hl2[5];
628   c[10] ^= hl2[6];
629   c[11] ^= hl2[7];
630}
631
632static
633void KarMul(_ntl_ulong *c, const _ntl_ulong *a, const _ntl_ulong *b, 
634            long len, _ntl_ulong *stk)
635{
636   if (len <= 8) {
637      switch (len) {
638         case 1: mul1(c, a[0], b[0]); break;
639         case 2: mul2(c, a, b); break;
640         case 3: mul3(c, a, b); break;
641         case 4: mul4(c, a, b); break;
642         case 5: mul5(c, a, b); break;
643         case 6: mul6(c, a, b); break;
644         case 7: mul7(c, a, b); break;
645         case 8: mul8(c, a, b); break;
646      }
647
648      return;
649   }
650
651   long ll, lh, i, ll2, lh2;
652   const _ntl_ulong *a0, *a1, *b0, *b1;
653   _ntl_ulong *a01, *b01, *h;
654
655
656   lh = len >> 1;
657   ll = (len+1) >> 1;
658
659   ll2 = ll << 1;
660   lh2 = lh << 1;
661
662   a01 = stk;  stk += ll+1;
663   b01 = stk;  stk += ll+1;
664   h = stk; stk += ll2+1;
665
666   a0 = a;
667   a1 = a+ll;
668   b0 = b;
669   b1 = b+ll;
670
671   KarMul(c, a0, b0, ll, stk);
672   KarMul(c+ll2, a1, b1, lh, stk);
673
674   for (i = 0; i < lh; i++) {
675      a01[i] = a[i] ^ a[i+ll];
676      b01[i] = b[i] ^ b[i+ll];
677   }
678
679   if (lh < ll) {
680      a01[lh] = a[lh];
681      b01[lh] = b[lh];
682   }
683
684   KarMul(h, a01, b01, ll, stk);
685
686   for (i = 0; i < ll2; i++)
687      h[i] ^= c[i];
688
689   for (i = 0; i < lh2; i++)
690      h[i] ^= c[i+ll2];
691
692   for (i = 0; i < ll2; i++)
693      c[i+ll] ^= h[i];
694}
695
696
697
698void mul(GF2X& c, const GF2X& a, const GF2X& b)
699{
700   long sa = a.xrep.length();
701   long sb = b.xrep.length();
702
703   if (sa <= 0 || sb <= 0) {
704      clear(c);
705      return;
706   }
707 
708   _ntl_ulong a0 = a.xrep[0];
709   _ntl_ulong b0 = b.xrep[0];
710
711   if (sb == 1 && b0 == 1) {
712      c = a;
713      return;
714   }
715
716   if (sa == 1 && a0 == 1) {
717      c = b;
718      return;
719   }
720
721   if (&a == &b) {
722      sqr(c, a);
723      return;
724   }
725
726   if (sa == sb && sa <= 8) {
727      // we treat these cases specially for efficiency reasons
728
729      switch (sa) {
730         case 1: {
731            _ntl_ulong v[2];
732            if (!(a0 >> NTL_BITS_PER_LONG/2))
733               mul_half(v, b0, a0);
734            else if (!(b0 >> NTL_BITS_PER_LONG/2))
735               mul_half(v, a0, b0);
736            else
737               mul1(v, a0, b0);
738
739            if (v[1]) {
740               c.xrep.SetLength(2);
741               _ntl_ulong *cp = &c.xrep[0];
742               cp[0] = v[0];
743               cp[1] = v[1];
744            }
745            else {
746               c.xrep.SetLength(1);
747               _ntl_ulong *cp = &c.xrep[0];
748               cp[0] = v[0];
749            }
750         }
751         return;
752
753         case 2: {
754            _ntl_ulong v[4];
755            mul2(v, &a.xrep[0], &b.xrep[0]);
756            if (v[3]) {
757               c.xrep.SetLength(4);
758               _ntl_ulong *cp = &c.xrep[0];
759               cp[0] = v[0];
760               cp[1] = v[1];
761               cp[2] = v[2];
762               cp[3] = v[3];
763            }
764            else {
765               c.xrep.SetLength(3);
766               _ntl_ulong *cp = &c.xrep[0];
767               cp[0] = v[0];
768               cp[1] = v[1];
769               cp[2] = v[2];
770            }
771         }
772         return;
773
774         case 3: {
775            _ntl_ulong v[6];
776            mul3(v, &a.xrep[0], &b.xrep[0]);
777            if (v[5]) {
778               c.xrep.SetLength(6);
779               _ntl_ulong *cp = &c.xrep[0];
780               cp[0] = v[0];
781               cp[1] = v[1];
782               cp[2] = v[2];
783               cp[3] = v[3];
784               cp[4] = v[4];
785               cp[5] = v[5];
786            }
787            else {
788               c.xrep.SetLength(5);
789               _ntl_ulong *cp = &c.xrep[0];
790               cp[0] = v[0];
791               cp[1] = v[1];
792               cp[2] = v[2];
793               cp[3] = v[3];
794               cp[4] = v[4];
795            }
796         }
797         return;
798
799         case 4: {
800            _ntl_ulong v[8];
801            mul4(v, &a.xrep[0], &b.xrep[0]);
802            if (v[7]) {
803               c.xrep.SetLength(8);
804               _ntl_ulong *cp = &c.xrep[0];
805               cp[0] = v[0];
806               cp[1] = v[1];
807               cp[2] = v[2];
808               cp[3] = v[3];
809               cp[4] = v[4];
810               cp[5] = v[5];
811               cp[6] = v[6];
812               cp[7] = v[7];
813            }
814            else {
815               c.xrep.SetLength(7);
816               _ntl_ulong *cp = &c.xrep[0];
817               cp[0] = v[0];
818               cp[1] = v[1];
819               cp[2] = v[2];
820               cp[3] = v[3];
821               cp[4] = v[4];
822               cp[5] = v[5];
823               cp[6] = v[6];
824            }
825         }
826         return;
827
828         case 5: {
829            _ntl_ulong v[10];
830            mul5(v, &a.xrep[0], &b.xrep[0]);
831            if (v[9]) {
832               c.xrep.SetLength(10);
833               _ntl_ulong *cp = &c.xrep[0];
834               cp[0] = v[0];
835               cp[1] = v[1];
836               cp[2] = v[2];
837               cp[3] = v[3];
838               cp[4] = v[4];
839               cp[5] = v[5];
840               cp[6] = v[6];
841               cp[7] = v[7];
842               cp[8] = v[8];
843               cp[9] = v[9];
844            }
845            else {
846               c.xrep.SetLength(9);
847               _ntl_ulong *cp = &c.xrep[0];
848               cp[0] = v[0];
849               cp[1] = v[1];
850               cp[2] = v[2];
851               cp[3] = v[3];
852               cp[4] = v[4];
853               cp[5] = v[5];
854               cp[6] = v[6];
855               cp[7] = v[7];
856               cp[8] = v[8];
857            }
858         }
859         return;
860
861         case 6: {
862            _ntl_ulong v[12];
863            mul6(v, &a.xrep[0], &b.xrep[0]);
864            if (v[11]) {
865               c.xrep.SetLength(12);
866               _ntl_ulong *cp = &c.xrep[0];
867               cp[0] = v[0];
868               cp[1] = v[1];
869               cp[2] = v[2];
870               cp[3] = v[3];
871               cp[4] = v[4];
872               cp[5] = v[5];
873               cp[6] = v[6];
874               cp[7] = v[7];
875               cp[8] = v[8];
876               cp[9] = v[9];
877               cp[10] = v[10];
878               cp[11] = v[11];
879            }
880            else {
881               c.xrep.SetLength(11);
882               _ntl_ulong *cp = &c.xrep[0];
883               cp[0] = v[0];
884               cp[1] = v[1];
885               cp[2] = v[2];
886               cp[3] = v[3];
887               cp[4] = v[4];
888               cp[5] = v[5];
889               cp[6] = v[6];
890               cp[7] = v[7];
891               cp[8] = v[8];
892               cp[9] = v[9];
893               cp[10] = v[10];
894            }
895         }
896         return; 
897
898         case 7: {
899            _ntl_ulong v[14];
900            mul7(v, &a.xrep[0], &b.xrep[0]);
901            if (v[13]) {
902               c.xrep.SetLength(14);
903               _ntl_ulong *cp = &c.xrep[0];
904               cp[0] = v[0];
905               cp[1] = v[1];
906               cp[2] = v[2];
907               cp[3] = v[3];
908               cp[4] = v[4];
909               cp[5] = v[5];
910               cp[6] = v[6];
911               cp[7] = v[7];
912               cp[8] = v[8];
913               cp[9] = v[9];
914               cp[10] = v[10];
915               cp[11] = v[11];
916               cp[12] = v[12];
917               cp[13] = v[13];
918            }
919            else {
920               c.xrep.SetLength(13);
921               _ntl_ulong *cp = &c.xrep[0];
922               cp[0] = v[0];
923               cp[1] = v[1];
924               cp[2] = v[2];
925               cp[3] = v[3];
926               cp[4] = v[4];
927               cp[5] = v[5];
928               cp[6] = v[6];
929               cp[7] = v[7];
930               cp[8] = v[8];
931               cp[9] = v[9];
932               cp[10] = v[10];
933               cp[11] = v[11];
934               cp[12] = v[12];
935            }
936         }
937         return; 
938
939         case 8: {
940            _ntl_ulong v[16];
941            mul8(v, &a.xrep[0], &b.xrep[0]);
942            if (v[15]) {
943               c.xrep.SetLength(16);
944               _ntl_ulong *cp = &c.xrep[0];
945               cp[0] = v[0];
946               cp[1] = v[1];
947               cp[2] = v[2];
948               cp[3] = v[3];
949               cp[4] = v[4];
950               cp[5] = v[5];
951               cp[6] = v[6];
952               cp[7] = v[7];
953               cp[8] = v[8];
954               cp[9] = v[9];
955               cp[10] = v[10];
956               cp[11] = v[11];
957               cp[12] = v[12];
958               cp[13] = v[13];
959               cp[14] = v[14];
960               cp[15] = v[15];
961            }
962            else {
963               c.xrep.SetLength(15);
964               _ntl_ulong *cp = &c.xrep[0];
965               cp[0] = v[0];
966               cp[1] = v[1];
967               cp[2] = v[2];
968               cp[3] = v[3];
969               cp[4] = v[4];
970               cp[5] = v[5];
971               cp[6] = v[6];
972               cp[7] = v[7];
973               cp[8] = v[8];
974               cp[9] = v[9];
975               cp[10] = v[10];
976               cp[11] = v[11];
977               cp[12] = v[12];
978               cp[13] = v[13];
979               cp[14] = v[14];
980            }
981         }
982         return; 
983
984      }
985   }
986
987   // another special case:  one of the two inputs
988   // has length 1 (or less).
989
990   if (sa == 1) {
991      c.xrep.SetLength(sb + 1);
992      _ntl_ulong *cp = c.xrep.elts();
993      const _ntl_ulong *bp = b.xrep.elts();
994
995      if (a0 >> (NTL_BITS_PER_LONG-NTL_BB_MUL1_BITS+1))
996         Mul1(cp, bp, sb, a0);
997      else
998         Mul1_short(cp, bp, sb, a0);
999
1000
1001      c.normalize();
1002      return;
1003   }
1004
1005   if (sb == 1) {
1006      c.xrep.SetLength(sa + 1);
1007      _ntl_ulong *cp = c.xrep.elts();
1008      const _ntl_ulong *ap = a.xrep.elts();
1009
1010
1011      if (b0 >> (NTL_BITS_PER_LONG-NTL_BB_MUL1_BITS+1))
1012         Mul1(cp, ap, sa, b0);
1013      else
1014         Mul1_short(cp, ap, sa, b0);
1015
1016      c.normalize();
1017      return;
1018   }
1019
1020   // finally: the general case
1021
1022   
1023   static WordVector mem;
1024   static WordVector stk;
1025   static WordVector vec;
1026
1027   const _ntl_ulong *ap, *bp;
1028   _ntl_ulong *cp;
1029
1030
1031   long sc = sa + sb;
1032   long in_mem = 0;
1033
1034   if (&a == &c || &b == &c) {
1035      mem.SetLength(sc);
1036      cp = mem.elts();
1037      in_mem = 1;
1038   }
1039   else {
1040      c.xrep.SetLength(sc);
1041      cp = c.xrep.elts();
1042   }
1043
1044
1045   long n, hn, sp;
1046
1047   n = min(sa, sb);
1048   sp = 0;
1049   while (n > 8) {
1050      hn = (n+1) >> 1;
1051      sp += (hn << 2) + 3;
1052      n = hn;
1053   }
1054
1055   stk.SetLength(sp);
1056   _ntl_ulong *stk_p = stk.elts();
1057
1058   if (sa > sb) {
1059      { long t; t = sa; sa = sb; sb = t; }
1060      ap = b.xrep.elts();
1061      bp = a.xrep.elts();
1062   }
1063   else {
1064      ap = a.xrep.elts();
1065      bp = b.xrep.elts();
1066   }
1067
1068
1069   vec.SetLength(2*sa);
1070
1071   _ntl_ulong *v = vec.elts();
1072
1073   long i, j;
1074
1075   for (i = 0; i < sc; i++)
1076      cp[i] = 0;
1077
1078   do {
1079      if (sa == 0) break;
1080
1081      if (sa == 1) {
1082         AddMul1(cp, bp, sb, ap[0]);
1083         
1084         break;
1085      }
1086
1087      // general case
1088
1089      for (i = 0; i+sa <= sb; i += sa) {
1090         KarMul(v, ap, bp + i, sa, stk_p);
1091         for (j = 0; j < 2*sa; j++)
1092            cp[i+j] ^= v[j]; 
1093      }
1094
1095      { const _ntl_ulong *t; t = ap; ap = bp + i; bp = t; }
1096      { long t; t = sa; sa = sb - i; sb = t; }
1097      cp = cp + i;
1098   } while (1);
1099
1100   if (in_mem)
1101      c.xrep = mem;
1102
1103   c.normalize();
1104}
1105
1106
1107
1108
1109void mul(GF2X& c, const GF2X& a, long b)
1110{
1111   if (b & 1)
1112      c = a;
1113   else
1114      clear(c);
1115}
1116
1117void mul(GF2X& c, const GF2X& a, GF2 b)
1118{
1119   if (b == 1)
1120      c = a;
1121   else
1122      clear(c);
1123}
1124
1125
1126void trunc(GF2X& x, const GF2X& a, long m)
1127{
1128   if (m < 0) Error("trunc: bad args");
1129
1130   long n = a.xrep.length();
1131   if (n == 0 || m == 0) {
1132      clear(x);
1133      return;
1134   }
1135
1136   if (&x == &a) {
1137      if (n*NTL_BITS_PER_LONG > m) {
1138         long wm = (m-1)/NTL_BITS_PER_LONG;
1139         long bm = m - NTL_BITS_PER_LONG*wm;
1140         _ntl_ulong msk;
1141         if (bm == NTL_BITS_PER_LONG)
1142            msk = ~(0UL);
1143         else
1144            msk = ((1UL << bm) - 1UL);
1145         x.xrep[wm] &= msk;
1146         x.xrep.QuickSetLength(wm+1);
1147         x.normalize();
1148      }
1149   }
1150   else if (n*NTL_BITS_PER_LONG <= m) 
1151      x = a;
1152   else {
1153      long wm = (m-1)/NTL_BITS_PER_LONG;
1154      long bm = m - NTL_BITS_PER_LONG*wm;
1155      x.xrep.SetLength(wm+1);
1156      _ntl_ulong *xp = &x.xrep[0];
1157      const _ntl_ulong *ap = &a.xrep[0];
1158      long i;
1159      for (i = 0; i < wm; i++)
1160         xp[i] = ap[i];
1161      _ntl_ulong msk;
1162      if (bm == NTL_BITS_PER_LONG)
1163         msk = ~(0UL);
1164      else
1165         msk = ((1UL << bm) - 1UL);
1166      xp[wm] = ap[wm] & msk;
1167      x.normalize();
1168   }
1169}
1170     
1171
1172/****** implementation of vec_GF2X ******/
1173
1174
1175NTL_vector_impl(GF2X,vec_GF2X)
1176
1177NTL_eq_vector_impl(GF2X,vec_GF2X)
1178
1179
1180
1181void MulByX(GF2X& x, const GF2X& a)
1182{
1183   long n = a.xrep.length();
1184   if (n == 0) {
1185      clear(x);
1186      return;
1187   }
1188
1189   if (a.xrep[n-1] & (1UL << (NTL_BITS_PER_LONG-1))) {
1190      x.xrep.SetLength(n+1);
1191      x.xrep[n] = 1;
1192   }
1193   else if (&x != &a)
1194      x.xrep.SetLength(n);
1195
1196   _ntl_ulong *xp = &x.xrep[0];
1197   const _ntl_ulong *ap = &a.xrep[0];
1198
1199   long i;
1200   for (i = n-1; i > 0; i--)
1201      xp[i] = (ap[i] << 1) | (ap[i-1] >> (NTL_BITS_PER_LONG-1));
1202
1203   xp[0] = ap[0] << 1;
1204
1205   // no need to normalize
1206}
1207
1208
1209
1210static _ntl_ulong sqrtab[256] = {
1211
12120UL, 1UL, 4UL, 5UL, 16UL, 17UL, 20UL, 21UL, 64UL, 
121365UL, 68UL, 69UL, 80UL, 81UL, 84UL, 85UL, 256UL, 
1214257UL, 260UL, 261UL, 272UL, 273UL, 276UL, 277UL, 320UL, 
1215321UL, 324UL, 325UL, 336UL, 337UL, 340UL, 341UL, 1024UL, 
12161025UL, 1028UL, 1029UL, 1040UL, 1041UL, 1044UL, 1045UL, 1088UL, 
12171089UL, 1092UL, 1093UL, 1104UL, 1105UL, 1108UL, 1109UL, 1280UL, 
12181281UL, 1284UL, 1285UL, 1296UL, 1297UL, 1300UL, 1301UL, 1344UL, 
12191345UL, 1348UL, 1349UL, 1360UL, 1361UL, 1364UL, 1365UL, 4096UL, 
12204097UL, 4100UL, 4101UL, 4112UL, 4113UL, 4116UL, 4117UL, 4160UL, 
12214161UL, 4164UL, 4165UL, 4176UL, 4177UL, 4180UL, 4181UL, 4352UL, 
12224353UL, 4356UL, 4357UL, 4368UL, 4369UL, 4372UL, 4373UL, 4416UL, 
12234417UL, 4420UL, 4421UL, 4432UL, 4433UL, 4436UL, 4437UL, 5120UL, 
12245121UL, 5124UL, 5125UL, 5136UL, 5137UL, 5140UL, 5141UL, 5184UL, 
12255185UL, 5188UL, 5189UL, 5200UL, 5201UL, 5204UL, 5205UL, 5376UL, 
12265377UL, 5380UL, 5381UL, 5392UL, 5393UL, 5396UL, 5397UL, 5440UL, 
12275441UL, 5444UL, 5445UL, 5456UL, 5457UL, 5460UL, 5461UL, 16384UL, 
122816385UL, 16388UL, 16389UL, 16400UL, 16401UL, 16404UL, 16405UL, 16448UL, 
122916449UL, 16452UL, 16453UL, 16464UL, 16465UL, 16468UL, 16469UL, 16640UL, 
123016641UL, 16644UL, 16645UL, 16656UL, 16657UL, 16660UL, 16661UL, 16704UL, 
123116705UL, 16708UL, 16709UL, 16720UL, 16721UL, 16724UL, 16725UL, 17408UL, 
123217409UL, 17412UL, 17413UL, 17424UL, 17425UL, 17428UL, 17429UL, 17472UL, 
123317473UL, 17476UL, 17477UL, 17488UL, 17489UL, 17492UL, 17493UL, 17664UL, 
123417665UL, 17668UL, 17669UL, 17680UL, 17681UL, 17684UL, 17685UL, 17728UL, 
123517729UL, 17732UL, 17733UL, 17744UL, 17745UL, 17748UL, 17749UL, 20480UL, 
123620481UL, 20484UL, 20485UL, 20496UL, 20497UL, 20500UL, 20501UL, 20544UL, 
123720545UL, 20548UL, 20549UL, 20560UL, 20561UL, 20564UL, 20565UL, 20736UL, 
123820737UL, 20740UL, 20741UL, 20752UL, 20753UL, 20756UL, 20757UL, 20800UL, 
123920801UL, 20804UL, 20805UL, 20816UL, 20817UL, 20820UL, 20821UL, 21504UL, 
124021505UL, 21508UL, 21509UL, 21520UL, 21521UL, 21524UL, 21525UL, 21568UL, 
124121569UL, 21572UL, 21573UL, 21584UL, 21585UL, 21588UL, 21589UL, 21760UL, 
124221761UL, 21764UL, 21765UL, 21776UL, 21777UL, 21780UL, 21781UL, 21824UL, 
124321825UL, 21828UL, 21829UL, 21840UL, 21841UL, 21844UL, 21845UL };
1244
1245
1246
1247
1248static inline 
1249void sqr1(_ntl_ulong *c, _ntl_ulong a)
1250{
1251   _ntl_ulong hi, lo;
1252
1253   NTL_BB_SQR_CODE
1254
1255   c[0] = lo;
1256   c[1] = hi;
1257}
1258
1259
1260
1261
1262void sqr(GF2X& c, const GF2X& a)
1263{
1264   long sa = a.xrep.length();
1265   if (sa <= 0) {
1266      clear(c);
1267      return;
1268   }
1269
1270   c.xrep.SetLength(sa << 1);
1271   _ntl_ulong *cp = c.xrep.elts();
1272   const _ntl_ulong *ap = a.xrep.elts();
1273   long i;
1274
1275   for (i = sa-1; i >= 0; i--)
1276      sqr1(cp + (i << 1), ap[i]);
1277
1278   c.normalize();
1279   return;
1280}
1281
1282
1283
1284void LeftShift(GF2X& c, const GF2X& a, long n)
1285{
1286   if (IsZero(a)) {
1287      clear(c);
1288      return;
1289   }
1290
1291   if (n == 1) {
1292      MulByX(c, a);
1293      return;
1294   }
1295
1296   if (n < 0) {
1297      if (n < -NTL_MAX_LONG) 
1298         clear(c);
1299      else
1300         RightShift(c, a, -n);
1301      return;
1302   }
1303
1304   if (NTL_OVERFLOW(n, 1, 0))
1305      Error("overflow in LeftShift");
1306
1307   if (n == 0) {
1308      c = a;
1309      return;
1310   }
1311
1312   long sa = a.xrep.length();
1313
1314   long wn = n/NTL_BITS_PER_LONG;
1315   long bn = n - wn*NTL_BITS_PER_LONG;
1316
1317   long sc = sa + wn;
1318   if (bn) sc++;
1319
1320   c.xrep.SetLength(sc);
1321
1322   _ntl_ulong *cp = c.xrep.elts();
1323   const _ntl_ulong *ap = a.xrep.elts();
1324
1325   long i;
1326
1327   if (bn == 0) {
1328      for (i = sa+wn-1; i >= wn; i--)
1329         cp[i] = ap[i-wn];
1330      for (i = wn-1; i >= 0; i--)
1331         cp[i] = 0;
1332   }
1333   else {
1334      cp[sa+wn] = ap[sa-1] >> (NTL_BITS_PER_LONG-bn);
1335      for (i = sa+wn-1; i >= wn+1; i--)
1336         cp[i] = (ap[i-wn] << bn) | (ap[i-wn-1] >> (NTL_BITS_PER_LONG-bn));
1337      cp[wn] = ap[0] << bn;
1338      for (i = wn-1; i >= 0; i--)
1339         cp[i] = 0;
1340   }
1341
1342   c.normalize();
1343}
1344
1345void ShiftAdd(GF2X& c, const GF2X& a, long n)
1346// c = c + a*X^n
1347{
1348   if (n < 0) Error("ShiftAdd: negative argument");
1349
1350   if (n == 0) {
1351      add(c, c, a);
1352      return;
1353   }
1354
1355   if (NTL_OVERFLOW(n, 1, 0))
1356      Error("overflow in ShiftAdd");
1357
1358   long sa = a.xrep.length();
1359   if (sa <= 0) {
1360      return;
1361   }
1362
1363   long sc = c.xrep.length();
1364
1365   long wn = n/NTL_BITS_PER_LONG;
1366   long bn = n - wn*NTL_BITS_PER_LONG;
1367
1368   long ss = sa + wn;
1369   if (bn) ss++;
1370
1371   if (ss > sc) 
1372      c.xrep.SetLength(ss);
1373
1374   _ntl_ulong *cp = c.xrep.elts();
1375   const _ntl_ulong *ap = a.xrep.elts();
1376
1377   long i;
1378
1379   for (i = sc; i < ss; i++)
1380      cp[i] = 0;
1381
1382   if (bn == 0) {
1383      for (i = sa+wn-1; i >= wn; i--)
1384         cp[i] ^= ap[i-wn];
1385   }
1386   else {
1387      cp[sa+wn] ^= ap[sa-1] >> (NTL_BITS_PER_LONG-bn);
1388      for (i = sa+wn-1; i >= wn+1; i--)
1389         cp[i] ^= (ap[i-wn] << bn) | (ap[i-wn-1] >> (NTL_BITS_PER_LONG-bn));
1390      cp[wn] ^= ap[0] << bn;
1391   }
1392
1393   c.normalize();
1394}
1395
1396
1397
1398
1399void RightShift(GF2X& c, const GF2X& a, long n)
1400{
1401   if (IsZero(a)) {
1402      clear(c);
1403      return;
1404   }
1405
1406   if (n < 0) {
1407      if (n < -NTL_MAX_LONG) Error("overflow in RightShift");
1408      LeftShift(c, a, -n);
1409      return;
1410   }
1411
1412   if (n == 0) {
1413      c = a;
1414      return;
1415   }
1416
1417   long sa = a.xrep.length();
1418   long wn = n/NTL_BITS_PER_LONG;
1419   long bn = n - wn*NTL_BITS_PER_LONG;
1420
1421   if (wn >= sa) {
1422      clear(c);
1423      return;
1424   }
1425
1426   c.xrep.SetLength(sa-wn);
1427   _ntl_ulong *cp = c.xrep.elts();
1428   const _ntl_ulong *ap = a.xrep.elts();
1429
1430   long i;
1431
1432   if (bn == 0) {
1433      for (i = 0; i < sa-wn; i++)
1434         cp[i] = ap[i+wn];
1435   }
1436   else {
1437      for (i = 0; i < sa-wn-1; i++)
1438         cp[i] = (ap[i+wn] >> bn) | (ap[i+wn+1] << (NTL_BITS_PER_LONG - bn));
1439
1440      cp[sa-wn-1] = ap[sa-1] >> bn;
1441   }
1442
1443   c.normalize();
1444}
1445
1446
1447
1448static _ntl_ulong revtab[256] = {
1449
14500UL, 128UL, 64UL, 192UL, 32UL, 160UL, 96UL, 224UL, 16UL, 144UL, 
145180UL, 208UL, 48UL, 176UL, 112UL, 240UL, 8UL, 136UL, 72UL, 200UL, 
145240UL, 168UL, 104UL, 232UL, 24UL, 152UL, 88UL, 216UL, 56UL, 184UL, 
1453120UL, 248UL, 4UL, 132UL, 68UL, 196UL, 36UL, 164UL, 100UL, 228UL, 
145420UL, 148UL, 84UL, 212UL, 52UL, 180UL, 116UL, 244UL, 12UL, 140UL, 
145576UL, 204UL, 44UL, 172UL, 108UL, 236UL, 28UL, 156UL, 92UL, 220UL, 
145660UL, 188UL, 124UL, 252UL, 2UL, 130UL, 66UL, 194UL, 34UL, 162UL, 
145798UL, 226UL, 18UL, 146UL, 82UL, 210UL, 50UL, 178UL, 114UL, 242UL, 
145810UL, 138UL, 74UL, 202UL, 42UL, 170UL, 106UL, 234UL, 26UL, 154UL, 
145990UL, 218UL, 58UL, 186UL, 122UL, 250UL, 6UL, 134UL, 70UL, 198UL, 
146038UL, 166UL, 102UL, 230UL, 22UL, 150UL, 86UL, 214UL, 54UL, 182UL, 
1461118UL, 246UL, 14UL, 142UL, 78UL, 206UL, 46UL, 174UL, 110UL, 238UL, 
146230UL, 158UL, 94UL, 222UL, 62UL, 190UL, 126UL, 254UL, 1UL, 129UL, 
146365UL, 193UL, 33UL, 161UL, 97UL, 225UL, 17UL, 145UL, 81UL, 209UL, 
146449UL, 177UL, 113UL, 241UL, 9UL, 137UL, 73UL, 201UL, 41UL, 169UL, 
1465105UL, 233UL, 25UL, 153UL, 89UL, 217UL, 57UL, 185UL, 121UL, 249UL, 
14665UL, 133UL, 69UL, 197UL, 37UL, 165UL, 101UL, 229UL, 21UL, 149UL, 
146785UL, 213UL, 53UL, 181UL, 117UL, 245UL, 13UL, 141UL, 77UL, 205UL, 
146845UL, 173UL, 109UL, 237UL, 29UL, 157UL, 93UL, 221UL, 61UL, 189UL, 
1469125UL, 253UL, 3UL, 131UL, 67UL, 195UL, 35UL, 163UL, 99UL, 227UL, 
147019UL, 147UL, 83UL, 211UL, 51UL, 179UL, 115UL, 243UL, 11UL, 139UL, 
147175UL, 203UL, 43UL, 171UL, 107UL, 235UL, 27UL, 155UL, 91UL, 219UL, 
147259UL, 187UL, 123UL, 251UL, 7UL, 135UL, 71UL, 199UL, 39UL, 167UL, 
1473103UL, 231UL, 23UL, 151UL, 87UL, 215UL, 55UL, 183UL, 119UL, 247UL, 
147415UL, 143UL, 79UL, 207UL, 47UL, 175UL, 111UL, 239UL, 31UL, 159UL, 
147595UL, 223UL, 63UL, 191UL, 127UL, 255UL  }; 
1476
1477static inline 
1478_ntl_ulong rev1(_ntl_ulong a)
1479{
1480   return NTL_BB_REV_CODE;
1481}
1482
1483
1484
1485void CopyReverse(GF2X& c, const GF2X& a, long hi)
1486// c[0..hi] = reverse(a[0..hi]), with zero fill as necessary
1487// input may alias output
1488
1489{
1490   if (hi < 0) { clear(c); return; }
1491
1492   if (NTL_OVERFLOW(hi, 1, 0))
1493      Error("overflow in CopyReverse");
1494
1495   long n = hi+1;
1496   long sa = a.xrep.length();
1497   if (n <= 0 || sa <= 0) {
1498      clear(c);
1499      return;
1500   }
1501
1502   long wn = n/NTL_BITS_PER_LONG;
1503   long bn = n - wn*NTL_BITS_PER_LONG;
1504
1505   if (bn != 0) {
1506      wn++;
1507      bn = NTL_BITS_PER_LONG - bn;
1508   }
1509
1510   c.xrep.SetLength(wn);
1511 
1512   _ntl_ulong *cp = c.xrep.elts();
1513   const _ntl_ulong *ap = a.xrep.elts();
1514
1515   long mm = min(sa, wn);
1516   long i;
1517
1518   for (i = 0; i < mm; i++)
1519      cp[i] = ap[i];
1520
1521   for (i = mm; i < wn; i++)
1522      cp[i] = 0;
1523
1524   if (bn != 0) {
1525      for (i = wn-1; i >= 1; i--)
1526         cp[i] = (cp[i] << bn) | (cp[i-1] >> (NTL_BITS_PER_LONG-bn));
1527      cp[0] = cp[0] << bn;
1528   }
1529
1530   for (i = 0; i < wn/2; i++) {
1531      _ntl_ulong t; t = cp[i]; cp[i] = cp[wn-1-i]; cp[wn-1-i] = t;
1532   }
1533
1534   for (i = 0; i < wn; i++)
1535      cp[i] = rev1(cp[i]);
1536
1537   c.normalize();
1538}
1539
1540
1541void div(GF2X& q, const GF2X& a, GF2 b)
1542{
1543   if (b == 0)
1544      Error("div: division by zero");
1545
1546   q = a;
1547}
1548
1549void div(GF2X& q, const GF2X& a, long b)
1550{
1551   if ((b & 1) == 0)
1552      Error("div: division by zero");
1553
1554   q = a;
1555}
1556
1557
1558
1559void GF2XFromBytes(GF2X& x, const unsigned char *p, long n)
1560{
1561   if (n <= 0) {
1562      x = 0;
1563      return;
1564   }
1565
1566   const long BytesPerLong = NTL_BITS_PER_LONG/8;
1567
1568   long lw, r, i, j;
1569
1570   lw = n/BytesPerLong;
1571   r = n - lw*BytesPerLong;
1572
1573   if (r != 0) 
1574      lw++;
1575   else
1576      r = BytesPerLong;
1577
1578   x.xrep.SetLength(lw);
1579   unsigned long *xp = x.xrep.elts();
1580
1581   for (i = 0; i < lw-1; i++) {
1582      unsigned long t = 0;
1583      for (j = 0; j < BytesPerLong; j++) {
1584         t >>= 8;
1585         t += (((unsigned long)(*p)) & 255UL) << ((BytesPerLong-1)*8);
1586         p++;
1587      }
1588      xp[i] = t;
1589   }
1590
1591   unsigned long t = 0;
1592   for (j = 0; j < r; j++) {
1593      t >>= 8;
1594      t += (((unsigned long)(*p)) & 255UL) << ((BytesPerLong-1)*8);
1595      p++;
1596   }
1597
1598   t >>= (BytesPerLong-r)*8;
1599   xp[lw-1] = t;
1600
1601   x.normalize();
1602}
1603
1604void BytesFromGF2X(unsigned char *p, const GF2X& a, long n)
1605{
1606   if (n < 0) n = 0;
1607
1608   const long BytesPerLong = NTL_BITS_PER_LONG/8;
1609
1610   long lbits = deg(a) + 1;
1611   long lbytes = (lbits+7)/8;
1612
1613   long min_bytes = min(lbytes, n);
1614
1615   long min_words = min_bytes/BytesPerLong;
1616   long r = min_bytes - min_words*BytesPerLong;
1617   if (r != 0)
1618      min_words++;
1619   else
1620      r = BytesPerLong;
1621
1622   const unsigned long *ap = a.xrep.elts();
1623
1624   long i, j;
1625
1626   for (i = 0; i < min_words-1; i++) {
1627      unsigned long t = ap[i];
1628      for (j = 0; j < BytesPerLong; j++) {
1629         *p = t & 255UL;
1630         t >>= 8;
1631         p++;
1632      }
1633   }
1634
1635   if (min_words > 0) {
1636      unsigned long t = ap[min_words-1];
1637      for (j = 0; j < r; j++) {
1638         *p = t & 255UL;
1639         t >>= 8;
1640         p++;
1641      }
1642   }
1643
1644   for (j = min_bytes; j < n; j++) {
1645      *p = 0;
1646      p++;
1647   }
1648}
1649
1650NTL_END_IMPL
Note: See TracBrowser for help on using the repository browser.