CoCalc -- finfield.c

Real-time collaboration for Jupyter Notebooks, Linux Terminals, LaTeX, VS Code, R IDE, and more,
all in one place.

GAP 4.8.9 installation with standard packages -- copy to your CoCalc project to get it
Project: cocalc-sagemath-dev-slelievre
Views: ⁴¹⁸³⁴⁶
1
/****************************************************************************
2
**
3
*W  finfield.c                  GAP source                      Werner Nickel
4
*W                                                         & Martin Schönert
5
**
6
**
7
*Y  Copyright (C)  1996,  Lehrstuhl D für Mathematik,  RWTH Aachen,  Germany
8
*Y  (C) 1998 School Math and Comp. Sci., University of St Andrews, Scotland
9
*Y  Copyright (C) 2002 The GAP Group
10
**
11
**  This file contains the  functions  to compute  with elements  from  small
12
**  finite fields.
13
**
14
**  Finite fields  are an  important   domain in computational   group theory
15
**  because the classical matrix groups are defined over those finite fields.
16
**  In GAP we  support  small  finite  fields  with  up  to  65536  elements,
17
**  larger fields can be realized  as polynomial domains over smaller fields.
18
**
19
**  Elements in small finite fields are  represented  as  immediate  objects.
20
**
21
**      +----------------+-------------+---+
22
**      |     <value>    |   <field>   |010|
23
**      +----------------+-------------+---+
24
**
25
**  The least significant 3 bits of such an immediate object are always  010,
26
**  flagging the object as an object of a small finite field.
27
**
28
**  The next 13 bits represent the small finite field where the element lies.
29
**  They are simply an index into a global table of small finite fields.
30
**
31
**  The most significant 16 bits represent the value of the element.
32
**
33
**  If the  value is 0,  then the element is the  zero from the finite field.
34
**  Otherwise the integer is the logarithm of this  element with respect to a
35
**  fixed generator of the multiplicative group of the finite field plus one.
36
**  In the following desriptions we denote this generator always with $z$, it
37
**  is an element of order $o-1$, where $o$ is the order of the finite field.
38
**  Thus 1 corresponds to $z^{1-1} = z^0 = 1$, i.e., the  one from the field.
39
**  Likewise 2 corresponds to $z^{2-1} = z^1 = z$, i.e., the root itself.
40
**
41
**  This representation  makes multiplication very easy,  we only have to add
42
**  the values and subtract 1 , because  $z^{a-1} * z^{b-1} = z^{(a+b-1)-1}$.
43
**  Addition is reduced to * by the formula $z^a +  z^b = z^b * (z^{a-b}+1)$.
44
**  This makes it neccessary to know the successor $z^a + 1$ of every value.
45
**
46
**  The  finite field  bag contains  the  successor for  every nonzero value,
47
**  i.e., 'SUCC_FF(<ff>)[<a>]' is  the successor of  the element <a>, i.e, it
48
**  is the  logarithm  of $z^{a-1} +   1$.  This list  is  usually called the
49
**  Zech-Logarithm  table.  The zeroth  entry in the  finite field bag is the
50
**  order of the finite field minus one.
51
*/
52
#include        "system.h"              /* Ints, UInts                     */
53

54

55
#include        "gasman.h"              /* garbage collector               */
56
#include        "objects.h"             /* objects                         */
57
#include        "scanner.h"             /* scanner                         */
58

59
#include        "gap.h"                 /* error handling, initialisation  */
60

61
#include        "gvars.h"               /* global variables                */
62

63
#include        "calls.h"               /* generic call mechanism          */
64
#include        "opers.h"               /* generic operations              */
65

66
#include        "ariths.h"              /* basic arithmetic                */
67

68
#include        "bool.h"                /* booleans                        */
69

70
#include        "integer.h"             /* integers                        */
71

72
#include        "finfield.h"            /* finite fields and ff elements   */
73

74
#include        "records.h"             /* generic records                 */
75
#include        "precord.h"             /* plain records                   */
76

77
#include        "lists.h"               /* generic lists                   */
78
#include        "plist.h"               /* plain lists                     */
79
#include        "string.h"              /* strings                         */
80

81
#include	"code.h"		/* coder                           */
82
#include	"aobjects.h"		/* atomic access to plists	   */
83
#include	"thread.h"		/* threads			   */
84
#include	"tls.h"			/* thread-local storage		   */
85

86
#include	"ffdata.h"		/* pre-computed finite field data  */
87

88

89
/****************************************************************************
90
**
91

92
*T  FF  . . . . . . . . . . . . . . . . . . . . . type of small finite fields
93
**
94
**  'FF' is the type used to represent small finite fields.
95
**
96
**  Small finite fields are represented by an  index  into  a  global  table.
97
**
98
**  Since there are  only  6542 (prime) + 93 (nonprime)  small finite fields,
99
**  the index fits into a 'UInt2' (actually into 13 bits).
100
**
101
**  'FF' is defined in the declaration part of this package as follows
102
**
103
typedef UInt2       FF;
104
*/
105

106

107
/****************************************************************************
108
**
109
*F  CHAR_FF(<ff>) . . . . . . . . . . .  characteristic of small finite field
110
**
111
**  'CHAR_FF' returns the characteristic of the small finite field <ff>.
112
**
113
**  Note that  'CHAR_FF' is a macro,  so do not call  it  with arguments that
114
**  have side effects.
115
**
116
**  'CHAR_FF' is defined in the declaration part of this package as follows
117
**
118
#define CHAR_FF(ff)             INT_INTOBJ( ELM_PLIST( CharFF, ff ) )
119
*/
120
Obj             CharFF;
121

122

123
/****************************************************************************
124
**
125
*F  DEGR_FF(<ff>) . . . . . . . . . . . . . . .  degree of small finite field
126
**
127
**  'DEGR_FF' returns the degree of the small finite field <ff>.
128
**
129
**  Note that 'DEGR_FF' is  a macro, so do   not call it with  arguments that
130
**  have side effects.
131
**
132
**  'DEGR_FF' is defined in the declaration part of this package as follows
133
**
134
#define DEGR_FF(ff)             INT_INTOBJ( ELM_PLIST( DegrFF, ff ) )
135
*/
136
Obj             DegrFF;
137

138

139
/****************************************************************************
140
**
141
*F  SIZE_FF(<ff>) . . . . . . . . . . . . . . . .  size of small finite field
142
**
143
**  'SIZE_FF' returns the size of the small finite field <ff>.
144
**
145
**  Note that 'SIZE_FF' is a macro, so do not call  it  with  arguments  that
146
**  have side effects.
147
**
148
**  'SIZE_FF' is defined in the declaration part of this package as follows
149
**
150
#define SIZE_FF(ff)             (*SUCC_FF(ff)+1)
151
*/
152

153

154
/****************************************************************************
155
**
156
*F  SUCC_FF(<ff>) . . . . . . . . . . . successor table of small finite field
157
**
158
**  'SUCC_FF' returns a pointer to the successor table of  the  small  finite
159
**  field <ff>.
160
**
161
**  Note that 'SUCC_FF' is a macro, so do not call  it  with  arguments  that
162
**  side effects.
163
**
164
**  'SUCC_FF' is defined in the declaration part of this package as follows
165
**
166
#define SUCC_FF(ff)             ((FFV*)ADDR_OBJ( ELM_PLIST( SuccFF, ff ) ))
167
*/
168
Obj             SuccFF;
169

170

171
/****************************************************************************
172
**
173
*F  TYPE_FF(<ff>) . . . . . . . . . . . . . . .  type of a small finite field
174
**
175
**  'TYPE_FF' returns the type of elements of the small finite field <ff>.
176
**
177
**  Note that  'TYPE_FF' is a macro, so  do not call  it  with arguments that
178
**  have side effects.
179
**
180
**  'TYPE_FF' is defined in the declaration part of this package as follows
181
**
182
#define TYPE_FF(ff)             (ELM_PLIST( TypeFF, ff ))
183
*/
184
Obj             TypeFF;
185
Obj             TypeFF0;
186

187
Obj             TYPE_FFE;
188
Obj             TYPE_FFE0;
189

190

191
/****************************************************************************
192
**
193
*T  FFV . . . . . . . . type of the values of elements of small finite fields
194
**
195
**  'FFV'  is the  type used  to represent  the values  of elements  of small
196
**  finite fields.
197
**
198
**  Values of  elements  of  small  finite  fields  are  represented  by  the
199
**  logarithm of the element with respect to the root plus one.
200
**
201
**  Since small finite fields contain at most 65536 elements,  the value fits
202
**  into a 'UInt2'.
203
**
204
**  It may be possible to change this to 'UInt4' to allow small finite fields
205
**  with more than than  65536 elements.  The macros  and have been coded  in
206
**  such a  way that they work  without problems.  The exception is 'POW_FFV'
207
**  which  will only work if  the product of integers  of type 'FFV' does not
208
**  cause an overflow.  And of course the successor table stored for a finite
209
**  field will become quite large for fields with more than 65536 elements.
210
**
211
**  'FFV' is defined in the declaration part of this package as follows
212
**
213
typedef UInt2           FFV;
214
*/
215

216

217
/****************************************************************************
218
**
219
*F  SUM_FFV(<a>,<b>,<f>)  . . . . . . . . . . . .  sum of finite field values
220
**
221
**  'SUM_FFV' returns the sum of the two finite field values <a> and <b> from
222
**  the finite field pointed to by the pointer <f>.
223
**
224
**  Note that 'SUM_FFV' may only  be used if  the operands are represented in
225
**  the same finite field.  If you want to add two elements where one lies in
226
**  a subfield of the other use 'SumFFEFFE'.
227
**
228
**  Use  'SUM_FFV' only with arguments  that are variables or array elements,
229
**  because it is a macro and arguments with side effects will behave strange,
230
**  and because it  is a complex macro  so most C  compilers will be upset by
231
**  complex arguments.  Especially do not use 'SUM_FFV(a,NEG_FFV(b,f),f)'.
232
**
233
**  If either operand is 0, the sum is just the other operand.
234
**  If $a <= b$ we have
235
**  $a + b ~ z^{a-1}+z^{b-1} = z^{a-1} * (z^{(b-1)-(a-1)}+1) ~ a * f[b-a+1]$,
236
**  otherwise we have
237
**  $a + b ~ z^{b-1}+z^{a-1} = z^{b-1} * (z^{(a-1)-(b-1)}+1) ~ b * f[a-b+1]$.
238
**
239
**  'SUM_FFV' is defined in the declaration part of this package as follows
240
**
241
#define SUM2_FFV(a,b,f) PROD_FFV( a, (f)[(b)-(a)+1], f )
242
#define SUM1_FFV(a,b,f) ( (a)<=(b) ? SUM2_FFV(a,b,f) : SUM2_FFV(b,a,f) )
243
#define SUM_FFV(a,b,f)  ( (a)==0 || (b)==0 ? (a)+(b) : SUM1_FFV(a,b,f) )
244
*/
245

246

247
/****************************************************************************
248
**
249
*F  NEG_FFV(<a>,<f>)  . . . . . . . . . . . .  negative of finite field value
250
**
251
**  'NEG_FFV' returns  the negative of the   finite field value  <a> from the
252
**  finite field pointed to by the pointer <f>.
253
**
254
**  Use  'NEG_FFV' only with arguments  that are variables or array elements,
255
**  because it is a macro and arguments with side effects will behave strange,
256
**  and because it is  a complex macro so most  C compilers will be upset  by
257
**  complex arguments.  Especially do not use 'NEG_FFV(PROD_FFV(a,b,f),f)'.
258
**
259
**  If the characteristic is 2, every element is its  own  additive  inverse.
260
**  Otherwise note that $z^{o-1} = 1 = -1^2$ so $z^{(o-1)/2} = 1^{1/2} = -1$.
261
**  If $a <= (o-1)/2$ we have
262
**  $-a ~ -1 * z^{a-1} = z^{(o-1)/2} * z^{a-1} = z^{a+(o-1)/2-1} ~ a+(o-1)/2$
263
**  otherwise we have
264
**  $-a ~ -1 * z^{a-1} = z^{a+(o-1)/2-1} = z^{a+(o-1)/2-1-(o-1)} ~ a-(o-1)/2$
265
**
266
**  'NEG_FFV' is defined in the declaration part of this package as follows
267
**
268
#define NEG2_FFV(a,f)   ( (a)<=*(f)/2 ? (a)+*(f)/2 : (a)-*(f)/2 )
269
#define NEG1_FFV(a,f)   ( *(f)%2==1 ? (a) : NEG2_FFV(a,f) )
270
#define NEG_FFV(a,f)    ( (a)==0 ? 0 : NEG1_FFV(a,f) )
271
*/
272

273

274
/****************************************************************************
275
**
276
*F  PROD_FFV(<a>,<b>,<f>) . . . . . . . . . . . product of finite field value
277
**
278
**  'PROD_FFV' returns the product of the two finite field values <a> and <b>
279
**  from the finite field pointed to by the pointer <f>.
280
**
281
**  Note that 'PROD_FFV' may only be used if the  operands are represented in
282
**  the  same finite field.  If you  want to multiply  two elements where one
283
**  lies in a subfield of the other use 'ProdFFEFFE'.
284
**
285
**  Use 'PROD_FFV' only with arguments that are  variables or array elements,
286
**  because it is a macro and arguments with side effects will behave strange,
287
**  and  because it is  a complex macro so most  C compilers will be upset by
288
**  complex arguments.  Especially do not use 'NEG_FFV(PROD_FFV(a,b,f),f)'.
289
**
290
**  If one of the values is 0 the product is 0.
291
**  If $a+b <= o$ we have $a * b ~ z^{a-1} * z^{b-1} = z^{(a+b-1)-1} ~ a+b-1$
292
**  otherwise   we   have $a * b ~ z^{(a+b-2)-(o-1)} = z^{(a+b-o)-1} ~ a+b-o$
293
**
294
**  'PROD_FF' is defined in the declaration part of this package as follows
295
**
296
#define PROD1_FFV(a,b,f) ( (a)-1<=*(f)-(b) ? (a)-1+(b) : (a)-1-(*(f)-(b)) )
297
#define PROD_FFV(a,b,f) ( (a)==0 || (b)==0 ? 0 : PROD1_FFV(a,b,f) )
298
*/
299

300

301
/****************************************************************************
302
**
303
*F  QUO_FFV(<a>,<b>,<f>)  . . . . . . . . . . quotient of finite field values
304
**
305
**  'QUO_FFV' returns the quotient of the two finite field values <a> and <b>
306
**  from the finite field pointed to by the pointer <f>.
307
**
308
**  Note that 'QUO_FFV' may  only be used  if the operands are represented in
309
**  the same finite field.  If you want to divide two elements where one lies
310
**  in a subfield of the other use 'QuoFFEFFE'.
311
**
312
**  Use 'QUO_FFV' only with arguments  that are variables or array  elements,
313
**  because it is a macro and arguments with side effects will behave strange,
314
**  and  because it is  a complex macro so most  C compilers will be upset by
315
**  complex arguments.  Especially do not use 'NEG_FFV(PROD_FFV(a,b,f),f)'.
316
**
317
**  A division by 0 is an error,  and dividing 0 by a nonzero value gives  0.
318
**  If $0 <= a-b$ we have  $a / b ~ z^{a-1} / z^{b-1} = z^{a-b+1-1} ~ a-b+1$,
319
**  otherwise   we   have  $a / b ~ z^{a-b+1-1}  =  z^{a-b+(o-1)}   ~ a-b+o$.
320
**
321
**  'QUO_FFV' is defined in the declaration part of this package as follows
322
**
323
#define QUO1_FFV(a,b,f) ( (b)<=(a) ? (a)-(b)+1 : *(f)-(b)+1+(a) )
324
#define QUO_FFV(a,b,f)  ( (a)==0 ? 0 : QUO1_FFV(a,b,f) )
325
*/
326

327

328
/****************************************************************************
329
**
330
*F  POW_FFV(<a>,<n>,<f>)  . . . . . . . . . . . power of a finite field value
331
**
332
**  'POW_FFV' returns the <n>th power of the finite  field value <a> from the
333
**  the finite field pointed to by the pointer <f>.
334
**
335
**  Note that 'POW_FFV' may only be used  if the right  operand is an integer
336
**  in the range $0..order(f)-1$.
337
**
338
**  Finally 'POW_FFV' may only be used if the  product of two integers of the
339
**  size of 'FFV'   does  not cause an  overflow,   i.e.  only if  'FFV'   is
340
**  'unsigned short'.
341
**
342
**  Note  that 'POW_FFV' is a macro,  so do not call  it  with arguments that
343
**  have side effects.  For optimal performance  put the operands in registers
344
**  before calling 'POW_FFV'.
345
**
346
**  If the finite field element is 0 the power is also 0, otherwise  we  have
347
**  $a^n ~ (z^{a-1})^n = z^{(a-1)*n} = z^{(a-1)*n % (o-1)} ~ (a-1)*n % (o-1)$
348
**
349
**  'POW_FFV' is defined in the declaration part of this package as follows
350
**
351
#define POW1_FFV(a,n,f) ( (((a)-1) * (n)) % *(f) + 1 )
352
#define POW_FFV(a,n,f)  ( (n)==0 ? 1 : ( (a)==0 ? 0 : POW1_FFV(a,n,f) ) )
353
*/
354

355

356
/****************************************************************************
357
**
358
*F  FLD_FFE(<ffe>)  . . . . . . . field of an element of a small finite field
359
**
360
**  'FLD_FFE' returns the small finite field over which the element  <ffe> is
361
**  represented.
362
**
363
**  Note that 'FLD_FFE' is a macro, so do not call  it  with  arguments  that
364
**  have side effects.
365
**
366
**  'FLD_FFE' is defined in the declaration part of this package as follows
367
**
368
#define FLD_FFE(ffe)            ((((UInt)(ffe)) & 0xFFFF) >> 3)
369
*/
370

371

372
/****************************************************************************
373
**
374
*F  VAL_FFE(<ffe>)  . . . . . . . value of an element of a small finite field
375
**
376
**  'VAL_FFE' returns the value of the element <ffe> of a small finite field.
377
**  Thus,  if <ffe> is $0_F$, it returns 0;  if <ffe> is $1_F$, it returns 1;
378
**  and otherwise if <ffe> is $z^i$, it returns $i+1$.
379
**
380
**  Note that 'VAL_FFE' is a macro, so do not call  it  with  arguments  that
381
**  have side effects.
382
**
383
**  'VAL_FFE' is defined in the declaration part of this package as follows
384
**
385
#define VAL_FFE(ffe)            (((UInt)(ffe)) >> 16)
386
*/
387

388

389
/****************************************************************************
390
**
391
*F  NEW_FFE(<fld>,<val>)  . . . .  make a new element of a small finite field
392
**
393
**  'NEW_FFE' returns a new element  <ffe>  of the  small finite  field <fld>
394
**  with the value <val>.
395
**
396
**  Note that 'NEW_FFE' is a macro, so do not  call  it  with  arguments that
397
**  have side effects.
398
**
399
**  'NEW_FFE' is defined in the declaration part of this package as follows
400
**
401
#define NEW_FFE(fld,val)        ((Obj)(((val) << 16) + ((fld) << 3) + 0x02))
402
*/
403

404

405
/****************************************************************************
406
**
407
*V  PolsFF  . . . . . . . . . .  list of Conway polynomials for finite fields
408
**
409
**  'PolsFF' is a  list of  Conway  polynomials for finite fields.   The even
410
**  entries are the  proper prime powers,  odd entries are the  corresponding
411
**  conway polynomials.
412
*/
413
unsigned long   PolsFF [] = {
414
       4, 1+2,
415
       8, 1+2,
416
      16, 1+2,
417
      32, 1  +4,
418
      64, 1+2  +8+16,
419
     128, 1+2,
420
     256, 1  +4+8+16,
421
     512, 1      +16,
422
    1024, 1+2+4+8   +32+64,
423
    2048, 1  +4,
424
    4096, 1+2  +8   +32+64+128,
425
    8192, 1+2  +8+16,
426
   16384, 1    +8   +32   +128,
427
   32768, 1  +4  +16+32,
428
   65536, 1  +4+8   +32,
429
       9,  2 +2*3,
430
      27,  1 +2*3,
431
      81,  2           +2*27,
432
     243,  1 +2*3,
433
     729,  2 +2*3 +1*9       +2*81,
434
    2187,  1      +2*9,
435
    6561,  2 +2*3 +2*9       +1*81 +2*243,
436
   19683,  1 +1*3 +2*9 +2*27,
437
   59049,  2 +1*3            +2*81 +2*243 +2*729,
438
      25,  2 +4*5,
439
     125,  3 +3*5,
440
     625,  2 +4*5 +4*25,
441
    3125,  3 +4*5,
442
   15625,  2      +1*25 +4*125 +1*625,
443
      49,  3 +6*7,
444
     343,  4      +6*49,
445
    2401,  3 +4*7 +5*49,
446
   16807,  4 +1*7,
447
     121,  2 + 7*11,
448
    1331,  9 + 2*11,
449
   14641,  2 +10*11 +8*121,
450
     169,  2 +12*13,
451
    2197, 11 + 2*13,
452
   28561,  2 +12*13 +3*169,
453
     289,  3 +16*17,
454
    4913, 14 + 1*17,
455
     361,  2 +18*19,
456
    6859, 17 + 4*19,
457
     529,  5 +21*23,
458
   12167, 18 + 2*23,
459
     841,  2 +24*29,
460
   24389, 27 + 2*29,
461
     961,  3 +29*31,
462
   29791, 28 + 1*31,
463
    1369,  2 +33*37,
464
   50653, 35 + 6*37,
465
    1681,  6 + 38* 41,
466
    1849,  3 + 42* 43,
467
    2209,  5 + 45* 47,
468
    2809,  2 + 49* 53,
469
    3481,  2 + 58* 59,
470
    3721,  2 + 60* 61,
471
    4489,  2 + 63* 67,
472
    5041,  7 + 69* 71,
473
    5329,  5 + 70* 73,
474
    6241,  3 + 78* 79,
475
    6889,  2 + 82* 83,
476
    7921,  3 + 82* 89,
477
    9409,  5 + 96* 97,
478
   10201,  2 + 97*101,
479
   10609,  5 +102*103,
480
   11449,  2 +103*107,
481
   11881,  6 +108*109,
482
   12769,  3 +101*113,
483
   16129,  3 +126*127,
484
   17161,  2 +127*131,
485
   18769,  3 +131*137,
486
   19321,  2 +138*139,
487
   22201,  2 +145*149,
488
   22801,  6 +149*151,
489
   24649,  5 +152*157,
490
   26569,  2 +159*163,
491
   27889,  5 +166*167,
492
   29929,  2 +169*173,
493
   32041,  2 +172*179,
494
   32761,  2 +177*181,
495
   36481, 19 +190*191,
496
   37249,  5 +192*193,
497
   38809,  2 +192*197,
498
   39601,  3 +193*199,
499
   44521,  2 +207*211,
500
   49729,  3 +221*223,
501
   51529,  2 +220*227,
502
   52441,  6 +228*229,
503
   54289,  3 +232*233,
504
   57121,  7 +237*239,
505
   58081,  7 +238*241,
506
   63001,  6 +242*251,
507
};
508

509

510
/****************************************************************************
511
**
512
*F  FiniteField(<p>,<d>)  . . . make the small finite field with <q> elements
513
**
514
**  'FiniteField' returns the small finite field with <p>^<d> elements.
515
*/
516
FF              FiniteField (
517
    UInt                p,
518
    UInt                d )
519
{
520
    FF                  ff;             /* finite field, result            */
521
    Obj                 bag1, bag2;     /* successor table bags            */
522
    FFV *               succ;           /* successor table                 */
523
    FFV *               indx;           /* index table                     */
524
    UInt                q;              /* size of finite field            */
525
    UInt                poly;           /* Conway polynomial of extension  */
526
    UInt                i, l, f, n, e;  /* loop variables                  */
527

528
    /* search through the finite field table                               */
529
    for ( ff = 1; ff <= LEN_PLIST(SuccFF); ff++ ) {
530
        if ( CHAR_FF(ff) == p && DEGR_FF(ff) == d ) {
531
            return ff;
532
        }
533
    }
534

535
    /* check whether we can build such a finite field                      */
536
    if ( (  2 <= p && 17 <= d) || (  3 <= p && 11 <= d)
537
      || (  5 <= p &&  7 <= d) || (  7 <= p &&  6 <= d)
538
      || ( 11 <= p &&  5 <= d) || ( 17 <= p &&  4 <= d)
539
      || ( 41 <= p &&  3 <= d) || (257 <= p &&  2 <= d) ) {
540
        return 0;
541
    }
542

543
    /* allocate a bag for the finite field and one for a temporary         */
544
    q = 1;
545
    for ( i = 1; i <= d; i++ ) q *= p;
546
    bag1  = NewBag( T_PERM2, q * sizeof(FFV) );
547
    bag2  = NewBag( T_PERM2, q * sizeof(FFV) );
548
    indx = (FFV*)ADDR_OBJ( bag1 );
549
    succ = (FFV*)ADDR_OBJ( bag2 );
550

551
    /* if q is a prime find the smallest primitive root $e$, use $x - e$   */
552
    /*N 1990/02/04 mschoene this is likely to explode if 'FFV' is 'UInt4'  */
553
    /*N 1990/02/04 mschoene there are few dumber ways to find prim. roots  */
554
    if ( d == 1 ) {
555
        for ( e = 1, i = 1; i != p-1; ++e ) {
556
            for ( f = e, i = 1; f != 1; ++i )
557
                f = (f * e) % p;
558
        }
559
        poly = p-(e-1);
560
    }
561

562
    /* otherwise look up the polynomial used to construct this field       */
563
    else {
564
        for ( i = 0; PolsFF[i] != q; i += 2 ) ;
565
        poly = PolsFF[i+1];
566
    }
567

568
    /* construct 'indx' such that 'e = x^(indx[e]-1) % poly' for every e   */
569
    /*N 1990/02/04 mschoene this is likely to explode if 'FFV' is 'UInt4'  */
570
    indx[ 0 ] = 0;
571
    for ( e = 1, n = 0; n < q-1; ++n ) {
572
        indx[ e ] = n + 1;
573
        /* e =p*e mod poly =x*e mod poly =x*x^n mod poly =x^{n+1} mod poly */
574
        if ( p != 2 ) {
575
            f = p * (e % (q/p));  l = ((p-1) * (e / (q/p))) % p;  e = 0;
576
            for ( i = 1; i < q; i *= p )
577
                e = e + i * ((f/i + l * (poly/i)) % p);
578
        }
579
        else {
580
            if ( 2*e & q )  e = 2*e ^ poly ^ q;
581
            else            e = 2*e;
582
        }
583
    }
584

585
    /* construct 'succ' such that 'x^(n-1)+1 = x^(succ[n]-1)' for every n  */
586
    succ[ 0 ] = q-1;
587
    for ( e = 1, f = p-1; e < q; e++ ) {
588
        if ( e < f ) {
589
            succ[ indx[e] ] = indx[ e+1 ];
590
        }
591
        else {
592
            succ[ indx[e] ] = indx[ e+1-p ];
593
            f += p;
594
        }
595
    }
596

597
    /* enter the finite field in the tables                                */
598
    ASS_LIST( CharFF, ff, INTOBJ_INT(p) );
599
    ASS_LIST( DegrFF, ff, INTOBJ_INT(d) );
600
    ASS_LIST( SuccFF, ff, bag2 );
601
    CHANGED_BAG(SuccFF);
602
    bag1 = CALL_1ARGS( TYPE_FFE, INTOBJ_INT(p) );
603
    ASS_LIST( TypeFF, ff, bag1 );
604
    CHANGED_BAG(TypeFF);
605
    bag1 = CALL_1ARGS( TYPE_FFE0, INTOBJ_INT(p) );
606
    ASS_LIST( TypeFF0, ff, bag1 );
607
    CHANGED_BAG(TypeFF0);
608

609
    /* return the finite field                                             */
610
    return ff;
611
}
612

613

614
/****************************************************************************
615
**
616
*F  CommonFF(<f1>,<d1>,<f2>,<d2>) . . . . . . . . . . . . find a common field
617
**
618
**  'CommonFF' returns  a small finite field  that can represent  elements of
619
**  degree <d1> from the small finite field <f1> and  elements of degree <d2>
620
**  from the small finite field <f2>.  Note that this is not guaranteed to be
621
**  the smallest such field.  If  <f1> and <f2> have different characteristic
622
**  or the smallest common field, is too large, 'CommonFF' returns 0.
623
*/
624
FF              CommonFF (
625
    FF                  f1,
626
    UInt                d1,
627
    FF                  f2,
628
    UInt                d2 )
629
{
630
    UInt                p;              /* characteristic                  */
631
    UInt                d;              /* degree                          */
632

633
    /* trivial case first                                                  */
634
    if ( f1 == f2 ) {
635
        return f1;
636
    }
637

638
    /* get and check the characteristics                                   */
639
    p = CHAR_FF( f1 );
640
    if ( p != CHAR_FF( f2 ) ) {
641
        return 0;
642
    }
643

644
    /* check whether one of the fields will do                             */
645
    if ( DEGR_FF(f1) % d2 == 0 ) {
646
        return f1;
647
    }
648
    if ( DEGR_FF(f2) % d1 == 0 ) {
649
        return f2;
650
    }
651

652
    /* compute the neccessary degree                                       */
653
    d = d1;
654
    while ( d % d2 != 0 ) {
655
        d += d1;
656
    }
657

658
    /* try to build the field                                              */
659
    return FiniteField( p, d );
660
}
661

662

663
/****************************************************************************
664
**
665
*F  CharFFE(<ffe>)  . . . . . . . . .  characteristic of a small finite field
666
**
667
**  'CharFFE' returns the characteristic of the small finite field  in  which
668
**  the element <ffe> lies.
669
*/
670
UInt CharFFE (
671
    Obj                 ffe )
672
{
673
    return CHAR_FF( FLD_FFE(ffe) );
674
}
675

676
Obj FuncCHAR_FFE_DEFAULT (
677
    Obj                 self,
678
    Obj                 ffe )
679
{
680
    return INTOBJ_INT( CHAR_FF( FLD_FFE(ffe) ) );
681
}
682

683

684
/****************************************************************************
685
**
686
*F  DegreeFFE(<ffe>)  . . . . . . . . . . . .  degree of a small finite field
687
**
688
**  'DegreeFFE' returns the degree of the smallest finite field in which  the
689
**  element <ffe> lies.
690
*/
691
UInt DegreeFFE (
692
    Obj                 ffe )
693
{
694
    UInt                d;              /* degree, result                  */
695
    FFV                 val;            /* value of element                */
696
    FF                  fld;            /* field of element                */
697
    UInt                q;              /* size  of field                  */
698
    UInt                p;              /* char. of field                  */
699
    UInt                m;              /* size  of minimal field          */
700

701
    /* get the value, the field, the size, and the characteristic          */
702
    val = VAL_FFE( ffe );
703
    fld = FLD_FFE( ffe );
704
    q = SIZE_FF( fld );
705
    p = CHAR_FF( fld );
706

707
    /* the zero element has a degree of one                                */
708
    if ( val == 0 ) {
709
        return 1L;
710
    }
711

712
    /* compute the degree                                                  */
713
    m = p;
714
    d = 1;
715
    while ( (q-1) % (m-1) != 0 || (val-1) % ((q-1)/(m-1)) != 0 ) {
716
        m *= p;
717
        d += 1;
718
    }
719

720
    /* return the degree                                                   */
721
    return d;
722
}
723

724
Obj FunDEGREE_FFE_DEFAULT (
725
    Obj                 self,
726
    Obj                 ffe )
727
{
728
    return INTOBJ_INT( DegreeFFE( ffe ) );
729
}
730

731

732
/****************************************************************************
733
**
734
*F  TypeFFE(<ffe>)  . . . . . . . . . . type of element of small finite field
735
**
736
**  'TypeFFE' returns the type of the element <ffe> of a small finite field.
737
**
738
**  'TypeFFE' is the function in 'TypeObjFuncs' for  elements in small finite
739
**  fields.
740
*/
741
Obj             TypeFFE (
742
    Obj                 ffe )
743
{
744
  if (VAL_FFE(ffe) == 0)
745
    return TYPE_FF0( FLD_FFE( ffe ) );
746
  else
747
    return TYPE_FF( FLD_FFE( ffe ) );
748
}
749

750

751
/****************************************************************************
752
**
753
*F  EqFFE(<opL>,<opR>)  . . . . . . . test if finite field elements are equal
754
**
755
**  'EqFFE' returns  'True' if the two  finite field elements <opL> and <opR>
756
**  are equal and 'False' othwise.
757
**
758
**  This is complicated because it must account  for the following situation.
759
**  Suppose 'a' is 'Z(3)', 'b' is 'Z(3^2)^4' and  finally 'c' is 'Z(3^3)^13'.
760
**  Mathematically 'a' is equal to 'b', so we  want 'a =  b' to be 'true' and
761
**  since 'a' is represented over a  subfield of 'b'  this is no big problem.
762
**  Again  'a' is equal to  'c', and again we want  'a = c'  to be 'true' and
763
**  again this is no problem since 'a' is represented over a subfield of 'c'.
764
**  Since '=' ought  to be transitive we also  want 'b = c'  to be 'true' and
765
**  this is a problem, because they are represented over incompatible fields.
766
*/
767
Int             EqFFE (
768
    Obj                 opL,
769
    Obj                 opR )
770
{
771
    FFV                 vL, vR;         /* value of left and right         */
772
    FF                  fL, fR;         /* field of left and right         */
773
    UInt                pL, pR;         /* char. of left and right         */
774
    UInt                qL, qR;         /* size  of left and right         */
775
    UInt                mL, mR;         /* size  of minimal field          */
776

777
    /* get the values and the fields over which they are represented       */
778
    vL = VAL_FFE( opL );
779
    vR = VAL_FFE( opR );
780
    fL = FLD_FFE( opL );
781
    fR = FLD_FFE( opR );
782

783
    /* if the elements are represented over the same field, it is easy     */
784
    if ( fL == fR ) {
785
        return (vL == vR);
786
    }
787

788
    /* elements in fields of different characteristic are different too    */
789
    pL = CHAR_FF( fL );
790
    pR = CHAR_FF( fR );
791
    if ( pL != pR ) {
792
        return 0L;
793
    }
794

795
    /* the zero element is not equal to any other element                  */
796
    if ( vL == 0 || vR == 0 ) {
797
        return (vL == 0 && vR == 0);
798
    }
799

800
    /* compute the sizes of the minimal fields in which the elements lie   */
801
    qL = SIZE_FF( fL );
802
    mL = pL;
803
    while ( (qL-1) % (mL-1) != 0 || (vL-1) % ((qL-1)/(mL-1)) != 0 ) mL *= pL;
804
    qR = SIZE_FF( fR );
805
    mR = pR;
806
    while ( (qR-1) % (mR-1) != 0 || (vR-1) % ((qR-1)/(mR-1)) != 0 ) mR *= pR;
807

808
    /* elements in different fields are different too                      */
809
    if ( mL != mR ) {
810
        return 0L;
811
    }
812

813
    /* otherwise compare the elements in the common minimal field          */
814
    return ((vL-1)/((qL-1)/(mL-1)) == (vR-1)/((qR-1)/(mR-1)));
815
}
816

817

818
/****************************************************************************
819
**
820
*F  LtFFE(<opL>,<opR>)  . . . . . .  test if finite field elements is smaller
821
**
822
**  'LtFFEFFE' returns 'True' if the  finite field element <opL> is  strictly
823
**  less than the finite field element <opR> and 'False' otherwise.
824
*/
825
Int             LtFFE (
826
    Obj                 opL,
827
    Obj                 opR )
828
{
829
    FFV                 vL, vR;         /* value of left and right         */
830
    FF                  fL, fR;         /* field of left and right         */
831
    UInt                pL, pR;         /* char. of left and right         */
832
    UInt                qL, qR;         /* size  of left and right         */
833
    UInt                mL, mR;         /* size  of minimal field          */
834

835
    /* get the values and the fields over which they are represented       */
836
    vL = VAL_FFE( opL );
837
    vR = VAL_FFE( opR );
838
    fL = FLD_FFE( opL );
839
    fR = FLD_FFE( opR );
840

841
    /* elements in fields of different characteristic are not comparable   */
842
    pL = CHAR_FF( fL );
843
    pR = CHAR_FF( fR );
844
    if ( pL != pR ) {
845
        return (DoOperation2Args( LtOper, opL, opR ) == True);
846
    }
847

848
    /* the zero element is smaller than any other element                  */
849
    if ( vL == 0 || vR == 0 ) {
850
        return (vL == 0 && vR != 0);
851
    }
852

853
    /* compute the sizes of the minimal fields in which the elements lie   */
854
    qL = SIZE_FF( fL );
855
    mL = pL;
856
    while ( (qL-1) % (mL-1) != 0 || (vL-1) % ((qL-1)/(mL-1)) != 0 ) mL *= pL;
857
    qR = SIZE_FF( fR );
858
    mR = pR;
859
    while ( (qR-1) % (mR-1) != 0 || (vR-1) % ((qR-1)/(mR-1)) != 0 ) mR *= pR;
860

861
    /* elements in smaller fields are smaller too                          */
862
    if ( mL != mR ) {
863
        return (mL < mR);
864
    }
865

866
    /* otherwise compare the elements in the common minimal field          */
867
    return ((vL-1)/((qL-1)/(mL-1)) < (vR-1)/((qR-1)/(mR-1)));
868
}
869

870

871
/****************************************************************************
872
**
873
*F  PrFFV(<fld>,<val>)  . . . . . . . . . . . . .  print a finite field value
874
**
875
**  'PrFFV' prints the value <val> from the finite field <fld>.
876
**
877
*/
878
void            PrFFV (
879
    FF                  fld,
880
    FFV                 val )
881
{
882
    UInt                q;              /* size   of finite field          */
883
    UInt                p;              /* char.  of finite field          */
884
    UInt                m;              /* size   of minimal field         */
885
    UInt                d;              /* degree of minimal field         */
886

887
    /* get the characteristic, order of the minimal field and the degree   */
888
    q = SIZE_FF( fld );
889
    p = CHAR_FF( fld );
890

891
    /* print the zero                                                      */
892
    if ( val == 0 ) {
893
        Pr( "%>0*Z(%>%d%2<)", (Int)p, 0L );
894
    }
895

896
    /* print a nonzero element as power of the primitive root              */
897
    else {
898

899
        /* find the degree of the minimal field in that the element lies   */
900
        d = 1;  m = p;
901
        while ( (q-1) % (m-1) != 0 || (val-1) % ((q-1)/(m-1)) != 0 ) {
902
            d++;  m *= p;
903
        }
904
        val = (val-1) / ((q-1)/(m-1)) + 1;
905

906
        /* print the element                                               */
907
        Pr( "%>Z(%>%d%<", (Int)p, 0L );
908
        if ( d == 1 ) {
909
            Pr( "%<)", 0L, 0L );
910
        }
911
        else {
912
            Pr( "^%>%d%2<)", (Int)d, 0L );
913
        }
914
        if ( val != 2 ) {
915
            Pr( "^%>%d%<", (Int)(val-1), 0L );
916
        }
917
    }
918

919
}
920

921

922
/****************************************************************************
923
**
924
*F  PrFFE(<ffe>)  . . . . . . . . . . . . . . .  print a finite field element
925
**
926
**  'PrFFE' prints the finite field element <ffe>.
927
*/
928
void            PrFFE (
929
    Obj                 ffe )
930
{
931
    PrFFV( FLD_FFE(ffe), VAL_FFE(ffe) );
932
}
933

934

935
/****************************************************************************
936
**
937
*F  SumFFEFFE(<opL>,<opR>)  . . . . . . . . . .  sum of finite field elements
938
**
939
**  'SumFFEFFE' returns  the sum of the  two finite  field elements <opL> and
940
**  <opR>.  The sum is represented over the field over which the operands are
941
**  represented, even if it lies in a much smaller field.
942
**
943
**  If one of the  elements is represented over  a subfield of the field over
944
**  which the  other element  is  represented, it is  lifted  into the larger
945
**  field before the addition.
946
**
947
**  'SumFFEFFE' just does the conversions mentioned  above and then calls the
948
**  macro 'SUM_FFV' to do the actual addition.
949
*/
950
Obj             SUM_FFE_LARGE;
951

952
Obj             SumFFEFFE (
953
    Obj                 opL,
954
    Obj                 opR )
955
{
956
    FFV                 vL, vR, vX;     /* value of left, right, result    */
957
    FF                  fL, fR, fX;     /* field of left, right, result    */
958
    UInt                qL, qR, qX;     /* size  of left, right, result    */
959

960
    /* get the values, handle trivial cases                                */
961
    vL = VAL_FFE( opL );
962
    vR = VAL_FFE( opR );
963

964
    /* bring the two operands into a common field <fX>                     */
965
    fL = FLD_FFE( opL );
966
    qL = SIZE_FF( fL );
967
    fR = FLD_FFE( opR );
968
    qR = SIZE_FF( fR  );
969

970
    /*N 1997/01/04 mschoene this is likely to explode if 'FFV' is 'UInt4'  */
971
    if ( qL == qR ) {
972
        fX = fL;
973
    }
974
    else if ( qL % qR == 0 && (qL-1) % (qR-1) == 0 ) {
975
        fX = fL;
976
        if ( vR != 0 )  vR = (qL-1) / (qR-1) * (vR-1) + 1;
977
    }
978
    else if ( qR % qL == 0 && (qR-1) % (qL-1) == 0 ) {
979
        fX = fR;
980
        if ( vL != 0 )  vL = (qR-1) / (qL-1) * (vL-1) + 1;
981
    }
982
    else {
983
        fX = CommonFF( fL, DegreeFFE(opL), fR, DegreeFFE(opR) );
984
        if ( fX == 0 )  return CALL_2ARGS( SUM_FFE_LARGE, opL, opR );
985
        qX = SIZE_FF( fX );
986
        /* if ( vL != 0 )  vL = (qX-1) / (qL-1) * (vL-1) + 1; */
987
        if ( vL != 0 )  vL = ((qX-1) * (vL-1)) / (qL-1) + 1;
988
        /* if ( vR != 0 )  vR = (qX-1) / (qR-1) * (vR-1) + 1; */
989
        if ( vR != 0 )  vR = ((qX-1) * (vR-1)) / (qR-1) + 1;
990
    }
991

992
    /* compute and return the result                                       */
993
    vX = SUM_FFV( vL, vR, SUCC_FF(fX) );
994
    return NEW_FFE( fX, vX );
995
}
996

997
Obj             SumFFEInt (
998
    Obj                 opL,
999
    Obj                 opR )
1000
{
1001
    FFV                 vL, vR, vX;     /* value of left, right, result    */
1002
    FF                  fX;             /* field of result                 */
1003
    Int                 pX;             /* char. of result                 */
1004
    FF*                 sX;             /* successor table of result field */
1005

1006
    /* get the field for the result                                        */
1007
    fX = FLD_FFE( opL );
1008
    pX = CHAR_FF( fX );
1009
    sX = SUCC_FF( fX );
1010

1011
    /* get the right operand                                               */
1012
    vX = ((INT_INTOBJ( opR ) % pX) + pX) % pX;
1013
    if ( vX == 0 ) {
1014
        vR = 0;
1015
    }
1016
    else {
1017
        vR = 1;
1018
        for ( ; 1 < vX; vX-- )  vR = sX[vR];
1019
    }
1020

1021
    /* get the left operand                                                */
1022
    vL = VAL_FFE( opL );
1023

1024
    /* compute and return the result                                       */
1025
    vX = SUM_FFV( vL, vR, sX );
1026
    return NEW_FFE( fX, vX );
1027
}
1028

1029
Obj             SumIntFFE (
1030
    Obj                 opL,
1031
    Obj                 opR )
1032
{
1033
    FFV                 vL, vR, vX;     /* value of left, right, result    */
1034
    FF                  fX;             /* field of result                 */
1035
    Int                 pX;             /* char. of result                 */
1036
    FF*                 sX;             /* successor table of result field */
1037

1038
    /* get the field for the result                                        */
1039
    fX = FLD_FFE( opR );
1040
    pX = CHAR_FF( fX );
1041
    sX = SUCC_FF( fX );
1042

1043
    /* get the left operand                                                */
1044
    vX = ((INT_INTOBJ( opL ) % pX) + pX) % pX;
1045
    if ( vX == 0 ) {
1046
        vL = 0;
1047
    }
1048
    else {
1049
        vL = 1;
1050
        for ( ; 1 < vX; vX-- )  vL = sX[vL];
1051
    }
1052

1053
    /* get the right operand                                               */
1054
    vR = VAL_FFE( opR );
1055

1056
    /* compute and return the result                                       */
1057
    vX = SUM_FFV( vL, vR, sX );
1058
    return NEW_FFE( fX, vX );
1059
}
1060

1061

1062
/****************************************************************************
1063
**
1064
*F  ZeroFFE(<op>) . . . . . . . . . . . . . .  zero of a finite field element
1065
*/
1066
Obj             ZeroFFE (
1067
    Obj                 op )
1068
{
1069
    FF                  fX;             /* field of result                 */
1070

1071
    /* get the field for the result                                        */
1072
    fX = FLD_FFE( op );
1073

1074
    /* return the result                                                   */
1075
    return NEW_FFE( fX, 0 );
1076
}
1077

1078

1079
/****************************************************************************
1080
**
1081
*F  AInvFFE(<op>) . . . . . . . . . . additive inverse of finite field element
1082
*/
1083
Obj             AInvFFE (
1084
    Obj                 op )
1085
{
1086
    FFV                 v, vX;          /* value of operand, result        */
1087
    FF                  fX;             /* field of result                 */
1088
    FF*                 sX;             /* successor table of result field */
1089

1090
    /* get the field for the result                                        */
1091
    fX = FLD_FFE( op );
1092
    sX = SUCC_FF( fX );
1093

1094
    /* get the operand                                                     */
1095
    v = VAL_FFE( op );
1096

1097
    /* compute and return the result                                       */
1098
    vX = NEG_FFV( v, sX ); 
1099
    return NEW_FFE( fX, vX );
1100
}
1101

1102

1103
/****************************************************************************
1104
**
1105
*F  DiffFFEFFE(<opL>,<opR>) . . . . . . . difference of finite field elements
1106
**
1107
**  'DiffFFEFFE' returns  the difference  of  the two  finite  field elements
1108
**  <opL> and <opR>.  The difference is represented over the field over which
1109
**  the operands are represented, even if it lies in a much smaller field.
1110
**
1111
**  If one of the elements is  represented over a subfield  of the field over
1112
**  which the  other element is  represented,  it is  lifted into  the larger
1113
**  field before the subtraction.
1114
**
1115
**  'DiffFFEFFE' just does the conversions mentioned above and then calls the
1116
**  macros 'NEG_FFV' and 'SUM_FFV' to do the actual subtraction.
1117
*/
1118
Obj             DIFF_FFE_LARGE;
1119

1120
Obj             DiffFFEFFE (
1121
    Obj                 opL,
1122
    Obj                 opR )
1123
{
1124
    FFV                 vL, vR, vX;     /* value of left, right, result    */
1125
    FF                  fL, fR, fX;     /* field of left, right, result    */
1126
    UInt                qL, qR, qX;     /* size  of left, right, result    */
1127

1128
    /* get the values, handle trivial cases                                */
1129
    vL = VAL_FFE( opL );
1130
    vR = VAL_FFE( opR );
1131

1132
    /* bring the two operands into a common field <fX>                     */
1133
    fL = FLD_FFE( opL );
1134
    qL = SIZE_FF( fL );
1135
    fR = FLD_FFE( opR );
1136
    qR = SIZE_FF( fR  );
1137

1138
    /*N 1997/01/04 mschoene this is likely to explode if 'FFV' is 'UInt4'  */
1139
    if ( qL == qR ) {
1140
        fX = fL;
1141
    }
1142
    else if ( qL % qR == 0 && (qL-1) % (qR-1) == 0 ) {
1143
        fX = fL;
1144
        if ( vR != 0 )  vR = (qL-1) / (qR-1) * (vR-1) + 1;
1145
    }
1146
    else if ( qR % qL == 0 && (qR-1) % (qL-1) == 0 ) {
1147
        fX = fR;
1148
        if ( vL != 0 )  vL = (qR-1) / (qL-1) * (vL-1) + 1;
1149
    }
1150
    else {
1151
        fX = CommonFF( fL, DegreeFFE(opL), fR, DegreeFFE(opR) );
1152
        if ( fX == 0 )  return CALL_2ARGS( DIFF_FFE_LARGE, opL, opR );
1153
        qX = SIZE_FF( fX );
1154
        /* if ( vL != 0 )  vL = (qX-1) / (qL-1) * (vL-1) + 1; */
1155
        if ( vL != 0 )  vL = ((qX-1) * (vL-1)) / (qL-1) + 1;
1156
        /* if ( vR != 0 )  vR = (qX-1) / (qR-1) * (vR-1) + 1; */
1157
        if ( vR != 0 )  vR = ((qX-1) * (vR-1)) / (qR-1) + 1;
1158
    }
1159

1160
    /* compute and return the result                                       */
1161
    vR = NEG_FFV( vR, SUCC_FF(fX) );
1162
    vX = SUM_FFV( vL, vR, SUCC_FF(fX) );
1163
    return NEW_FFE( fX, vX );
1164
}
1165

1166
Obj             DiffFFEInt (
1167
    Obj                 opL,
1168
    Obj                 opR )
1169
{
1170
    FFV                 vL, vR, vX;     /* value of left, right, result    */
1171
    FF                  fX;             /* field of result                 */
1172
    Int                 pX;             /* char. of result                 */
1173
    FF*                 sX;             /* successor table of result field */
1174

1175
    /* get the field for the result                                        */
1176
    fX = FLD_FFE( opL );
1177
    pX = CHAR_FF( fX );
1178
    sX = SUCC_FF( fX );
1179

1180
    /* get the right operand                                               */
1181
    vX = ((INT_INTOBJ( opR ) % pX) + pX) % pX;
1182
    if ( vX == 0 ) {
1183
        vR = 0;
1184
    }
1185
    else {
1186
        vR = 1;
1187
        for ( ; 1 < vX; vX-- )  vR = sX[vR];
1188
    }
1189

1190
    /* get the left operand                                                */
1191
    vL = VAL_FFE( opL );
1192

1193
    /* compute and return the result                                       */
1194
    vR = NEG_FFV( vR, sX );
1195
    vX = SUM_FFV( vL, vR, sX );
1196
    return NEW_FFE( fX, vX );
1197
}
1198

1199
Obj             DiffIntFFE (
1200
    Obj                 opL,
1201
    Obj                 opR )
1202
{
1203
    FFV                 vL, vR, vX;     /* value of left, right, result    */
1204
    FF                  fX;             /* field of result                 */
1205
    Int                 pX;             /* char. of result                 */
1206
    FF*                 sX;             /* successor table of result field */
1207

1208
    /* get the field for the result                                        */
1209
    fX = FLD_FFE( opR );
1210
    pX = CHAR_FF( fX );
1211
    sX = SUCC_FF( fX );
1212

1213
    /* get the left operand                                                */
1214
    vX = ((INT_INTOBJ( opL ) % pX) + pX) % pX;
1215
    if ( vX == 0 ) {
1216
        vL = 0;
1217
    }
1218
    else {
1219
        vL = 1;
1220
        for ( ; 1 < vX; vX-- )  vL = sX[vL];
1221
    }
1222

1223
    /* get the right operand                                               */
1224
    vR = VAL_FFE( opR );
1225

1226
    /* compute and return the result                                       */
1227
    vR = NEG_FFV( vR, sX );
1228
    vX = SUM_FFV( vL, vR, sX );
1229
    return NEW_FFE( fX, vX );
1230
}
1231

1232

1233
/****************************************************************************
1234
**
1235
*F  ProdFFEFFE(<opL>,<opR>) . . . . . . . .  product of finite field elements
1236
**
1237
**  'ProdFFEFFE'  returns the product of  the two finite field elements <opL>
1238
**  and <opR>.   The  product is  represented  over the field  over which the
1239
**  operands are represented, even if it lies in a much smaller field.
1240
**
1241
**  If one of the elements  is represented over a  subfield of the field over
1242
**  which the other element  is  represented, it  is  lifted into the  larger
1243
**  field before the multiplication.
1244
**
1245
**  'ProdFFEFFE' just does the conversions mentioned above and then calls the
1246
**  macro 'PROD_FFV' to do the actual multiplication.
1247
*/
1248
Obj             PROD_FFE_LARGE;
1249

1250
Obj             ProdFFEFFE (
1251
    Obj                 opL,
1252
    Obj                 opR )
1253
{
1254
    FFV                 vL, vR, vX;     /* value of left, right, result    */
1255
    FF                  fL, fR, fX;     /* field of left, right, result    */
1256
    UInt                qL, qR, qX;     /* size  of left, right, result    */
1257

1258
    /* get the values, handle trivial cases                                */
1259
    vL = VAL_FFE( opL );
1260
    vR = VAL_FFE( opR );
1261

1262
    /* bring the two operands into a common field <fX>                     */
1263
    fL = FLD_FFE( opL );
1264
    qL = SIZE_FF( fL );
1265
    fR = FLD_FFE( opR );
1266
    qR = SIZE_FF( fR  );
1267

1268
    /*N 1997/01/04 mschoene this is likely to explode if 'FFV' is 'UInt4'  */
1269
    if ( qL == qR ) {
1270
        fX = fL;
1271
    }
1272
    else if ( qL % qR == 0 && (qL-1) % (qR-1) == 0 ) {
1273
        fX = fL;
1274
        if ( vR != 0 )  vR = (qL-1) / (qR-1) * (vR-1) + 1;
1275
    }
1276
    else if ( qR % qL == 0 && (qR-1) % (qL-1) == 0 ) {
1277
        fX = fR;
1278
        if ( vL != 0 )  vL = (qR-1) / (qL-1) * (vL-1) + 1;
1279
    }
1280
    else {
1281
        fX = CommonFF( fL, DegreeFFE(opL), fR, DegreeFFE(opR) );
1282
        if ( fX == 0 )  return CALL_2ARGS( PROD_FFE_LARGE, opL, opR );
1283
        qX = SIZE_FF( fX );
1284
        /* if ( vL != 0 )  vL = (qX-1) / (qL-1) * (vL-1) + 1; */
1285
        if ( vL != 0 )  vL = ((qX-1) * (vL-1)) / (qL-1) + 1;
1286
        /* if ( vR != 0 )  vR = (qX-1) / (qR-1) * (vR-1) + 1; */
1287
        if ( vR != 0 )  vR = ((qX-1) * (vR-1)) / (qR-1) + 1;
1288
    }
1289

1290
    /* compute and return the result                                       */
1291
    vX = PROD_FFV( vL, vR, SUCC_FF(fX) );
1292
    return NEW_FFE( fX, vX );
1293
}
1294

1295
Obj             ProdFFEInt (
1296
    Obj                 opL,
1297
    Obj                 opR )
1298
{
1299
    FFV                 vL, vR, vX;     /* value of left, right, result    */
1300
    FF                  fX;             /* field of result                 */
1301
    Int                 pX;             /* char. of result                 */
1302
    FF*                 sX;             /* successor table of result field */
1303

1304
    /* get the field for the result                                        */
1305
    fX = FLD_FFE( opL );
1306
    pX = CHAR_FF( fX );
1307
    sX = SUCC_FF( fX );
1308

1309
    /* get the right operand                                               */
1310
    vX = ((INT_INTOBJ( opR ) % pX) + pX) % pX;
1311
    if ( vX == 0 ) {
1312
        vR = 0;
1313
    }
1314
    else {
1315
        vR = 1;
1316
        for ( ; 1 < vX; vX-- )  vR = sX[vR];
1317
    }
1318

1319
    /* get the left operand                                                */
1320
    vL = VAL_FFE( opL );
1321

1322
    /* compute and return the result                                       */
1323
    vX = PROD_FFV( vL, vR, sX );
1324
    return NEW_FFE( fX, vX );
1325
}
1326

1327
Obj             ProdIntFFE (
1328
    Obj                 opL,
1329
    Obj                 opR )
1330
{
1331
    FFV                 vL, vR, vX;     /* value of left, right, result    */
1332
    FF                  fX;             /* field of result                 */
1333
    Int                 pX;             /* char. of result                 */
1334
    FF*                 sX;             /* successor table of result field */
1335

1336
    /* get the field for the result                                        */
1337
    fX = FLD_FFE( opR );
1338
    pX = CHAR_FF( fX );
1339
    sX = SUCC_FF( fX );
1340

1341
    /* get the left operand                                                */
1342
    vX = ((INT_INTOBJ( opL ) % pX) + pX) % pX;
1343
    if ( vX == 0 ) {
1344
        vL = 0;
1345
    }
1346
    else {
1347
        vL = 1;
1348
        for ( ; 1 < vX; vX-- )  vL = sX[vL];
1349
    }
1350

1351
    /* get the right operand                                               */
1352
    vR = VAL_FFE( opR );
1353

1354
    /* compute and return the result                                       */
1355
    vX = PROD_FFV( vL, vR, sX );
1356
    return NEW_FFE( fX, vX );
1357
}
1358

1359

1360
/****************************************************************************
1361
**
1362
*F  OneFFE(<op>)  . . . . . . . . . . . . . . . one of a finite field element
1363
*/
1364
Obj             OneFFE (
1365
    Obj                 op )
1366
{
1367
    FF                  fX;             /* field of result                 */
1368

1369
    /* get the field for the result                                        */
1370
    fX = FLD_FFE( op );
1371

1372
    /* return the result                                                   */
1373
    return NEW_FFE( fX, 1 );
1374
}
1375

1376

1377
/****************************************************************************
1378
**
1379
*F  InvFFE(<op>)  . . . . . . . . . . . . . . inverse of finite field element
1380
*/
1381
Obj             InvFFE (
1382
    Obj                 op )
1383
{
1384
    FFV                 v, vX;          /* value of operand, result        */
1385
    FF                  fX;             /* field of result                 */
1386
    FF*                 sX;             /* successor table of result field */
1387

1388
    /* get the field for the result                                        */
1389
    fX = FLD_FFE( op );
1390
    sX = SUCC_FF( fX );
1391

1392
    /* get the operand                                                     */
1393
    v = VAL_FFE( op );
1394
    if ( v == 0 ) return Fail;
1395

1396
    /* compute and return the result                                       */
1397
    vX = QUO_FFV( 1, v, sX ); 
1398
    return NEW_FFE( fX, vX );
1399
}
1400

1401

1402
/****************************************************************************
1403
**
1404
*F  QuoFFEFFE(<opL>,<opR>) . . . . . . . .  quotient of finite field elements
1405
**
1406
**  'QuoFFEFFE' returns  the quotient of  the two finite field elements <opL>
1407
**  and <opR>.  The  quotient is represented  over  the field  over which the
1408
**  operands are represented, even if it lies in a much smaller field.
1409
**
1410
**  If one of the  elements is represented over  a subfield of the field over
1411
**  which the  other  element is represented,  it is  lifted into the  larger
1412
**  field before the division.
1413
**
1414
**  'QuoFFEFFE' just does the conversions mentioned  above and then calls the
1415
**  macro 'QUO_FFV' to do the actual division.
1416
*/
1417
Obj             QUO_FFE_LARGE;
1418

1419
Obj             QuoFFEFFE (
1420
    Obj                 opL,
1421
    Obj                 opR )
1422
{
1423
    FFV                 vL, vR, vX;     /* value of left, right, result    */
1424
    FF                  fL, fR, fX;     /* field of left, right, result    */
1425
    UInt                qL, qR, qX;     /* size  of left, right, result    */
1426

1427
    /* get the values, handle trivial cases                                */
1428
    vL = VAL_FFE( opL );
1429
    vR = VAL_FFE( opR );
1430

1431
    /* bring the two operands into a common field <fX>                     */
1432
    fL = FLD_FFE( opL );
1433
    qL = SIZE_FF( fL );
1434
    fR = FLD_FFE( opR );
1435
    qR = SIZE_FF( fR  );
1436

1437
    /*N 1997/01/04 mschoene this is likely to explode if 'FFV' is 'UInt4'  */
1438
    if ( qL == qR ) {
1439
        fX = fL;
1440
    }
1441
    else if ( qL % qR == 0 && (qL-1) % (qR-1) == 0 ) {
1442
        fX = fL;
1443
        if ( vR != 0 )  vR = (qL-1) / (qR-1) * (vR-1) + 1;
1444
    }
1445
    else if ( qR % qL == 0 && (qR-1) % (qL-1) == 0 ) {
1446
        fX = fR;
1447
        if ( vL != 0 )  vL = (qR-1) / (qL-1) * (vL-1) + 1;
1448
    }
1449
    else {
1450
        fX = CommonFF( fL, DegreeFFE(opL), fR, DegreeFFE(opR) );
1451
        if ( fX == 0 )  return CALL_2ARGS( QUO_FFE_LARGE, opL, opR );
1452
        qX = SIZE_FF( fX );
1453
        /* if ( vL != 0 )  vL = (qX-1) / (qL-1) * (vL-1) + 1; */
1454
        if ( vL != 0 )  vL = ((qX-1) * (vL-1)) / (qL-1) + 1;
1455
        /* if ( vR != 0 )  vR = (qX-1) / (qR-1) * (vR-1) + 1; */
1456
        if ( vR != 0 )  vR = ((qX-1) * (vR-1)) / (qR-1) + 1;
1457
    }
1458

1459
    /* compute and return the result                                       */
1460
    if ( vR == 0 ) {
1461
        opR = ErrorReturnObj(
1462
            "FFE operations: <divisor> must not be zero",
1463
            0L, 0L,
1464
            "you can replace <divisor> via 'return <divisor>;'" );
1465
        return QUO( opL, opR );
1466
    }
1467
    vX = QUO_FFV( vL, vR, SUCC_FF(fX) );
1468
    return NEW_FFE( fX, vX );
1469
}
1470

1471
Obj             QuoFFEInt (
1472
    Obj                 opL,
1473
    Obj                 opR )
1474
{
1475
    FFV                 vL, vR, vX;     /* value of left, right, result    */
1476
    FF                  fX;             /* field of result                 */
1477
    Int                 pX;             /* char. of result                 */
1478
    FF*                 sX;             /* successor table of result field */
1479

1480
    /* get the field for the result                                        */
1481
    fX = FLD_FFE( opL );
1482
    pX = CHAR_FF( fX );
1483
    sX = SUCC_FF( fX );
1484

1485
    /* get the right operand                                               */
1486
    vX = ((INT_INTOBJ( opR ) % pX) + pX) % pX;
1487
    if ( vX == 0 ) {
1488
        vR = 0;
1489
    }
1490
    else {
1491
        vR = 1;
1492
        for ( ; 1 < vX; vX-- )  vR = sX[vR];
1493
    }
1494

1495
    /* get the left operand                                                */
1496
    vL = VAL_FFE( opL );
1497

1498
    /* compute and return the result                                       */
1499
    if ( vR == 0 ) {
1500
        opR = ErrorReturnObj(
1501
            "FFE operations: <divisor> must not be zero",
1502
            0L, 0L,
1503
            "you can replace <divisor> via 'return <divisor>;'" );
1504
        return QUO( opL, opR );
1505
    }
1506
    vX = QUO_FFV( vL, vR, sX );
1507
    return NEW_FFE( fX, vX );
1508
}
1509

1510
Obj             QuoIntFFE (
1511
    Obj                 opL,
1512
    Obj                 opR )
1513
{
1514
    FFV                 vL, vR, vX;     /* value of left, right, result    */
1515
    FF                  fX;             /* field of result                 */
1516
    Int                 pX;             /* char. of result                 */
1517
    FF*                 sX;             /* successor table of result field */
1518

1519
    /* get the field for the result                                        */
1520
    fX = FLD_FFE( opR );
1521
    pX = CHAR_FF( fX );
1522
    sX = SUCC_FF( fX );
1523

1524
    /* get the left operand                                                */
1525
    vX = ((INT_INTOBJ( opL ) % pX) + pX) % pX;
1526
    if ( vX == 0 ) {
1527
        vL = 0;
1528
    }
1529
    else {
1530
        vL = 1;
1531
        for ( ; 1 < vX; vX-- )  vL = sX[vL];
1532
    }
1533

1534
    /* get the right operand                                               */
1535
    vR = VAL_FFE( opR );
1536

1537
    /* compute and return the result                                       */
1538
    if ( vR == 0 ) {
1539
        opR = ErrorReturnObj(
1540
            "FFE operations: <divisor> must not be zero",
1541
            0L, 0L,
1542
            "you can replace <divisor> via 'return <divisor>;'" );
1543
        return QUO( opL, opR );
1544
    }
1545
    vX = QUO_FFV( vL, vR, sX );
1546
    return NEW_FFE( fX, vX );
1547
}
1548

1549

1550
/****************************************************************************
1551
**
1552
*F  PowFFEInt(<opL>,<opR>)  . . . . . . . . . power of a finite field element
1553
**
1554
**  'PowFFEInt' returns the  power of the finite  field element <opL> and the
1555
**  integer <opR>.  The power is  represented over the  field over which  the
1556
**  left operand is represented, even if it lies in a much smaller field.
1557
**
1558
**  'PowFFEInt' just does the conversions mentioned  above and then calls the
1559
**  macro 'POW_FFV' to do the actual exponentiation.
1560
*/
1561
Obj             PowFFEInt (
1562
    Obj                 opL,
1563
    Obj                 opR )
1564
{
1565
    FFV                 vL, vX;         /* value of left, result           */
1566
    Int                 vR;             /* value of right                  */
1567
    FF                  fX;             /* field of result                 */
1568
    FF*                 sX;             /* successor table of result field */
1569

1570
    /* get the field for the result                                        */
1571
    fX = FLD_FFE( opL );
1572
    sX = SUCC_FF( fX );
1573

1574
    /* get the right operand                                               */
1575
    vR = INT_INTOBJ( opR );
1576

1577
    /* get the left operand                                                */
1578
    vL = VAL_FFE( opL );
1579

1580
    /* if the exponent is negative, invert the left operand                */
1581
    if ( vR < 0 ) {
1582
        if ( vL == 0 ) {
1583
            opL = ErrorReturnObj(
1584
                "FFE operations: <divisor> must not be zero",
1585
                0L, 0L,
1586
                "you can replace <divisor> via 'return <divisor>;'" );
1587
            return POW( opL, opR );
1588
        }
1589
        vL = QUO_FFV( 1, vL, sX );
1590
        vR = -vR;
1591
    }
1592

1593
    /* catch the case when vL is zero.                                     */
1594
    if( vL == 0 ) return NEW_FFE( fX, (vR == 0 ? 1 : 0 ) );
1595

1596
    /* reduce vR modulo the order of the multiplicative group first.       */
1597
    vR %= *sX;
1598

1599
    /* compute and return the result                                       */
1600
    vX = POW_FFV( vL, vR, sX );
1601
    return NEW_FFE( fX, vX );
1602
}
1603

1604

1605
/****************************************************************************
1606
**
1607
*F  PowFFEFFE( <opL>, <opR> ) . . . . . . conjugate of a finite field element
1608
*/
1609
Obj PowFFEFFE (
1610
    Obj                 opL,
1611
    Obj                 opR )
1612
{
1613
    /* get the field for the result                                        */
1614
    if ( CHAR_FF( FLD_FFE(opL) ) != CHAR_FF( FLD_FFE(opR) ) ) {
1615
        opR = ErrorReturnObj(
1616
          "FFE operations: characteristic of conjugating element must be %d",
1617
          (Int)CHAR_FF(FLD_FFE(opL)), 0L,
1618
          "you can replace conjugating element <elt> via 'return <elt>;'" );
1619
        return POW( opL, opR );
1620
    }
1621

1622
    /* compute and return the result                                       */
1623
    return opL;
1624
}
1625

1626

1627
/****************************************************************************
1628
**
1629
*F  FuncIS_FFE( <self>, <obj> ) . . . . . . .  test for finite field elements
1630
**
1631
**  'FuncIsFFE' implements the internal function 'IsFFE( <obj> )'.
1632
**
1633
**  'IsFFE' returns  'true' if its argument  <obj> is a finite  field element
1634
**  and 'false' otherwise.   'IsFFE' will cause  an  error if  called with an
1635
**  unbound variable.
1636
*/
1637
Obj IsFFEFilt;
1638

1639
Obj FuncIS_FFE (
1640
    Obj                 self,
1641
    Obj                 obj )
1642
{
1643
    /* return 'true' if <obj> is a finite field element                    */
1644
    if ( TNUM_OBJ(obj) == T_FFE ) {
1645
        return True;
1646
    }
1647
    else if ( TNUM_OBJ(obj) < FIRST_EXTERNAL_TNUM ) {
1648
        return False;
1649
    }
1650
    else {
1651
        return DoFilter( self, obj );
1652
    }
1653
}
1654

1655

1656
/****************************************************************************
1657
**
1658
*F  FuncLOG_FFE_DEFAULT( <self>, <opZ>, <opR> ) .  logarithm of a ff constant
1659
**
1660
**  'FuncLOG_FFE_DEFAULT' implements the function 'LogFFE( <z>, <r> )'.
1661
**
1662
**  'LogFFE'  returns the logarithm of  the nonzero finite  field element <z>
1663
**  with respect to the root <r> which must lie in the same field like <z>.
1664
*/
1665
Obj LOG_FFE_LARGE;
1666
#include <stdio.h>
1667

1668
Obj FuncLOG_FFE_DEFAULT (
1669
    Obj                 self,
1670
    Obj                 opZ,
1671
    Obj                 opR )
1672
{
1673
    FFV                 vZ, vR;         /* value of left, right            */
1674
    FF                  fZ, fR, fX;     /* field of left, right, common    */
1675
    UInt                qZ, qR, qX;     /* size  of left, right, common    */
1676
    Int                 a, b, c, d, t;  /* temporaries                     */
1677

1678
    /* check the arguments                                                 */
1679
    if ( ! IS_FFE(opZ) || VAL_FFE(opZ) == 0 ) {
1680
        opZ = ErrorReturnObj(
1681
            "LogFFE: <z> must be a nonzero finite field element",
1682
             0L, 0L,
1683
             "you can replace <z> via 'return <z>;'" );
1684
        return FuncLOG_FFE_DEFAULT( self, opZ, opR );
1685
    }
1686
    if ( ! IS_FFE(opR) || VAL_FFE(opR) == 0 ) {
1687
        opR = ErrorReturnObj(
1688
            "LogFFE: <r> must be a nonzero finite field element",
1689
             0L, 0L,
1690
             "you can replace <r> via 'return <r>;'" );
1691
        return FuncLOG_FFE_DEFAULT( self, opZ, opR );
1692
    }
1693

1694
    /* get the values, handle trivial cases                                */
1695
    vZ = VAL_FFE( opZ );
1696
    vR = VAL_FFE( opR );
1697

1698
    /* bring the two operands into a common field <fX>                     */
1699
    fZ = FLD_FFE( opZ );
1700
    qZ = SIZE_FF( fZ );
1701
    fR = FLD_FFE( opR );
1702
    qR = SIZE_FF( fR  );
1703

1704
    /*N 1997/01/04 mschoene this is likely to explode if 'FFV' is 'UInt4'  */
1705
    if ( qZ == qR ) {
1706
        fX = fZ;
1707
        qX = qZ;
1708
    }
1709
    else if ( qZ % qR == 0 && (qZ-1) % (qR-1) == 0 ) {
1710
        fX = fZ;
1711
        qX = qZ;
1712
        if ( vR != 0 )  vR = (qZ-1) / (qR-1) * (vR-1) + 1;
1713
    }
1714
    else if ( qR % qZ == 0 && (qR-1) % (qZ-1) == 0 ) {
1715
        fX = fR;
1716
        qX = qR;
1717
        if ( vZ != 0 )  vZ = (qR-1) / (qZ-1) * (vZ-1) + 1;
1718
    }
1719
    else {
1720
        fX = CommonFF( fZ, DegreeFFE(opZ), fR, DegreeFFE(opR) );
1721
        if ( fX == 0 )  return CALL_2ARGS( LOG_FFE_LARGE, opZ, opR );
1722
        qX = SIZE_FF( fX );
1723
        /* if ( vZ != 0 )  vZ = (qX-1) / (qZ-1) * (vZ-1) + 1; */
1724
        if ( vZ != 0 )  vZ = ((qX-1) * (vZ-1)) / (qZ-1) + 1;
1725
        /* if ( vR != 0 )  vR = (qX-1) / (qR-1) * (vR-1) + 1; */
1726
        if ( vR != 0 )  vR = ((qX-1) * (vR-1)) / (qR-1) + 1;
1727
    }
1728

1729
    /* now solve <l> * (<vR>-1) = (<vZ>-1) % (<qX>-1)                      */
1730
    /*N 1990/02/04 mschoene this is likely to explode if 'FFV' is 'UInt4'  */
1731
    a = 1;             b = 0;
1732
    c = (Int) (vR-1);  d = (Int) (qX-1);
1733
    while ( d != 0 ) {
1734
        t = b;  b = a - (c/d) * b;  a = t;
1735
        t = d;  d = c - (c/d) * d;  c = t;
1736
    }
1737
    if ( ((Int) (vZ-1)) % c != 0 ) {
1738
      return Fail;
1739
    }
1740

1741

1742
    while (a < 0)
1743
      a+= (qX -1)/c;
1744
    /* return the logarithm                                                */
1745
    
1746

1747
    return INTOBJ_INT( (((Int) (vZ-1) / c) * a) % ((Int) (qX-1)) );
1748

1749
}
1750

1751

1752
/****************************************************************************
1753
**
1754
*F  FuncINT_FFE_DEFAULT( <self>, <z> )  . . . .   convert a ffe to an integer
1755
**
1756
**  'FuncINT_FFE_DEFAULT' implements the internal function 'IntFFE( <z> )'.
1757
**
1758
**  'IntFFE'  returns  the integer  that  corresponds  to  the  finite  field
1759
**  element <z>, which must of course be  an element  of a prime field, i.e.,
1760
**  the smallest integer <i> such that '<i> * <z>^0 = <z>'.
1761
*/
1762
Obj IntFF;
1763

1764
Obj INT_FF (
1765
    FF                  ff )
1766
{
1767
    Obj                 conv;           /* conversion table, result        */
1768
    Int                 q;              /* size of finite field            */
1769
    Int                 p;              /* char of finite field            */
1770
    FFV *               succ;           /* successor table of finite field */
1771
    FFV                 z;              /* one element of finite field     */
1772
    UInt                i;              /* loop variable                   */
1773

1774
    /* if the conversion table is not already known, construct it          */
1775
    if ( LEN_PLIST(IntFF) < ff || ELM_PLIST(IntFF,ff) == 0 ) {
1776
        q = SIZE_FF( ff );
1777
        p = CHAR_FF( ff );
1778
        conv = NEW_PLIST( T_PLIST, p-1 );
1779
        succ = SUCC_FF( ff );
1780
        SET_LEN_PLIST( conv, p-1 );
1781
        z = 1;
1782
        for ( i = 1; i < p; i++ ) {
1783
            SET_ELM_PLIST( conv, (z-1)/((q-1)/(p-1))+1, INTOBJ_INT(i) );
1784
            z = succ[ z ];
1785
        }
1786
        if ( LEN_PLIST(IntFF) < ff ) {
1787
            GROW_PLIST( IntFF, ff );
1788
            SET_LEN_PLIST( IntFF, (Int)ff );
1789
        }
1790
        SET_ELM_PLIST( IntFF, ff, conv );
1791
        CHANGED_BAG( IntFF );
1792
    }
1793

1794
    /* return the conversion table                                           */
1795
    return ELM_PLIST( IntFF, ff );
1796
}
1797

1798

1799

1800
Obj FuncINT_FFE_DEFAULT (
1801
    Obj                 self,
1802
    Obj                 z )
1803
{
1804
    FFV                 v;              /* value of finite field element   */
1805
    FF                  ff;             /* finite field                    */
1806
    Int                 q;              /* size of finite field            */
1807
    Int                 p;              /* char of finite field            */
1808
    Obj                 conv;           /* conversion table                */
1809

1810
    /* get the value                                                       */
1811
    v  = VAL_FFE( z );
1812

1813
    /* special case for 0                                                  */
1814
    if ( v == 0 ) {
1815
        return INTOBJ_INT( 0 );
1816
    }
1817

1818
    /* get the field, size, characteristic, and conversion table           */
1819
    ff   = FLD_FFE( z );
1820
    q    = SIZE_FF( ff );
1821
    p    = CHAR_FF( ff );
1822
    conv = INT_FF( ff );
1823

1824
    /* check the argument                                                  */
1825
    if ( (v-1) % ((q-1)/(p-1)) != 0 ) {
1826
        z = ErrorReturnObj(
1827
            "IntFFE: <z> must lie in prime field",
1828
            0L, 0L,
1829
            "you can replace <z> via 'return <z>;'" );
1830
        return FuncINT_FFE_DEFAULT( self, z );
1831
    }
1832

1833
    /* convert the value into the prime field                              */
1834
    v = (v-1) / ((q-1)/(p-1)) + 1;
1835

1836
    /* return the integer value                                            */
1837
    return ELM_PLIST( conv, v );
1838
}
1839

1840

1841
/****************************************************************************
1842
**
1843
*F  FuncZ( <self>, <q> )  . . .  return the generator of a small finite field
1844
**
1845
**  'FuncZ' implements the internal function 'Z( <q> )'.
1846
**
1847
**  'Z' returns the generators  of the small finite  field with <q> elements.
1848
**  <q> must be a positive prime power.
1849
*/
1850
static Obj ZOp;
1851

1852

1853

1854

1855
Obj FuncZ (
1856
    Obj                 self,
1857
    Obj                 q )
1858
{
1859
    FF                  ff;             /* the finite field                */
1860
    UInt                p;              /* characteristic                  */
1861
    UInt                d;              /* degree                          */
1862
    UInt                r;              /* temporary                       */
1863

1864
    /* check the argument                                                  */
1865
    if ( (TNUM_OBJ(q) == T_INT && (INT_INTOBJ(q) > 65536)) ||
1866
         (TNUM_OBJ(q) == T_INTPOS))
1867
      return CALL_1ARGS(ZOp, q);
1868
    
1869
    if ( TNUM_OBJ(q)!=T_INT || INT_INTOBJ(q)<=1 ) {
1870
        q = ErrorReturnObj(
1871
            "Z: <q> must be a positive prime power (not a %s)",
1872
            (Int)TNAM_OBJ(q), 0L,
1873
            "you can replace <q> via 'return <q>;'" );
1874
        return FuncZ( self, q );
1875
    }
1876

1877
    /* compute the prime and check that <q> is a prime power               */
1878
    if ( INT_INTOBJ(q) % 2 == 0 ) {
1879
        p = 2;
1880
    }
1881
    else {
1882
        p = 3;
1883
        while ( INT_INTOBJ(q) % p != 0 ) {
1884
            p += 2;
1885
        }
1886
    }
1887
    d = 1;
1888
    r = p;
1889
    while ( r < INT_INTOBJ(q) ) {
1890
        d = d + 1;
1891
        r = r * p;
1892
    }
1893
    if ( r != INT_INTOBJ(q) ) {
1894
        q = ErrorReturnObj(
1895
            "Z: <q> must be a positive prime power (not a %s)",
1896
            (Int)TNAM_OBJ(q), 0L,
1897
            "you can replace <q> via 'return <q>;'" );
1898
        return FuncZ( self, q );
1899
    }
1900

1901
    /* get the finite field                                                */
1902
    ff = FiniteField( p, d );
1903

1904
    /* make the root                                                       */
1905
    return NEW_FFE( ff, (p == 2 && d == 1 ? 1 : 2) );
1906
}
1907

1908
Obj FuncZ2 ( Obj self, Obj p, Obj d)
1909
{
1910
  FF ff;
1911
  Int ip,id,id1;
1912
  UInt q;
1913
  if (ARE_INTOBJS(p,d))
1914
    {
1915
      ip = INT_INTOBJ(p);
1916
      id = INT_INTOBJ(d);
1917
      if (ip > 1 && id > 0 && id <= 16 && ip <= 65536)
1918
        {
1919
          id1 = id;
1920
          q = ip;
1921
          while (--id1 > 0 && q <= 65536)
1922
            q *= ip;
1923
          if (q <= 65536)
1924
            {
1925
              /* get the finite field                                                */
1926
              ff = FiniteField( ip, id );
1927
              
1928
              /* make the root                                                       */
1929
              return NEW_FFE( ff, (ip == 2 && id == 1 ? 1 : 2) );
1930
            }
1931
        }
1932
    }
1933
  return CALL_2ARGS(ZOp, p, d);
1934
}
1935

1936

1937
/****************************************************************************
1938
**
1939

1940
*F * * * * * * * * * * * * * initialize package * * * * * * * * * * * * * * *
1941
*/
1942

1943

1944
/****************************************************************************
1945
**
1946

1947
*V  GVarFilts . . . . . . . . . . . . . . . . . . . list of filters to export
1948
*/
1949
static StructGVarFilt GVarFilts [] = {
1950

1951
    { "IS_FFE", "obj", &IsFFEFilt,
1952
      FuncIS_FFE, "src/finifield.c:IS_FFE" },
1953

1954
    { 0 }
1955

1956
};
1957

1958

1959
/****************************************************************************
1960
**
1961
*V  GVarFuncs . . . . . . . . . . . . . . . . . . list of functions to export
1962
*/
1963
static StructGVarFunc GVarFuncs [] = {
1964

1965
    { "CHAR_FFE_DEFAULT", 1, "z",
1966
      FuncCHAR_FFE_DEFAULT, "src/finifield.c:CHAR_FFE_DEFAULT" },
1967

1968
    { "DEGREE_FFE_DEFAULT", 1, "z",
1969
      FunDEGREE_FFE_DEFAULT, "src/finifield.c:DEGREE_FFE_DEFAULT" },
1970

1971
    { "LOG_FFE_DEFAULT", 2, "z, root",
1972
      FuncLOG_FFE_DEFAULT, "src/finifield.c:LOG_FFE_DEFAULT" },
1973

1974
    { "INT_FFE_DEFAULT", 1, "z",
1975
      FuncINT_FFE_DEFAULT, "src/finifield.c:INT_FFE_DEFAULT" },
1976

1977
    { "Z", 1, "q",
1978
      FuncZ, "src/finifield.c:Z" },
1979

1980
    { 0 }
1981

1982
};
1983

1984

1985
/****************************************************************************
1986
**
1987

1988
*F  InitKernel( <module> )  . . . . . . . . initialise kernel data structures
1989
*/
1990
static Int InitKernel (
1991
    StructInitInfo *    module )
1992
{
1993
    /* install the marking function                                        */
1994
    InfoBags[ T_FFE ].name = "ffe";
1995
    /* InitMarkFuncBags( T_FFE, MarkNoSubBags ); */
1996

1997
    /* install the type functions                                          */
1998
    ImportFuncFromLibrary( "TYPE_FFE", &TYPE_FFE );
1999
    ImportFuncFromLibrary( "TYPE_FFE0", &TYPE_FFE0 );
2000
    ImportFuncFromLibrary( "ZOp", &ZOp );
2001
    TypeObjFuncs[ T_FFE ] = TypeFFE;
2002

2003
    /* create the fields and integer conversion bags                       */
2004
    InitGlobalBag( &CharFF, "src/finfield.c:CharFF" );
2005
    InitGlobalBag( &DegrFF, "src/finfield.c:DegrFF" );
2006
    InitGlobalBag( &SuccFF, "src/finfield.c:SuccFF" );
2007
    InitGlobalBag( &TypeFF, "src/finfield.c:TypeFF" );
2008
    InitGlobalBag( &TypeFF0, "src/finfield.c:TypeFF0" );
2009
    InitGlobalBag( &IntFF, "src/finifield.c:IntFF" );
2010

2011
    /* install the functions that handle overflow                          */
2012
    ImportFuncFromLibrary( "SUM_FFE_LARGE",  &SUM_FFE_LARGE  );
2013
    ImportFuncFromLibrary( "DIFF_FFE_LARGE", &DIFF_FFE_LARGE );
2014
    ImportFuncFromLibrary( "PROD_FFE_LARGE", &PROD_FFE_LARGE );
2015
    ImportFuncFromLibrary( "QUO_FFE_LARGE",  &QUO_FFE_LARGE  );
2016
    ImportFuncFromLibrary( "LOG_FFE_LARGE",  &LOG_FFE_LARGE  );
2017

2018
    /* init filters and functions                                          */
2019
    InitHdlrFiltsFromTable( GVarFilts );
2020
    InitHdlrFuncsFromTable( GVarFuncs );
2021
    InitHandlerFunc( FuncZ2, "src/finfield.c: Z (2 args)");
2022

2023

2024
    /* install the printing method                                         */
2025
    PrintObjFuncs[ T_FFE ] = PrFFE;
2026

2027
    /* install the comparison methods                                      */
2028
    EqFuncs[   T_FFE ][ T_FFE ] = EqFFE;
2029
    LtFuncs[   T_FFE ][ T_FFE ] = LtFFE;
2030

2031
    /* install the arithmetic methods                                      */
2032
    ZeroFuncs[ T_FFE ] = ZeroFFE;
2033
    ZeroMutFuncs[ T_FFE ] = ZeroFFE;
2034
    AInvFuncs[ T_FFE ] = AInvFFE;
2035
    AInvMutFuncs[ T_FFE ] = AInvFFE;
2036
    OneFuncs [ T_FFE ] = OneFFE;
2037
    OneMutFuncs [ T_FFE ] = OneFFE;
2038
    InvFuncs [ T_FFE ] = InvFFE;
2039
    InvMutFuncs [ T_FFE ] = InvFFE;
2040
    SumFuncs[  T_FFE ][ T_FFE ] = SumFFEFFE;
2041
    SumFuncs[  T_FFE ][ T_INT ] = SumFFEInt;
2042
    SumFuncs[  T_INT ][ T_FFE ] = SumIntFFE;
2043
    DiffFuncs[ T_FFE ][ T_FFE ] = DiffFFEFFE;
2044
    DiffFuncs[ T_FFE ][ T_INT ] = DiffFFEInt;
2045
    DiffFuncs[ T_INT ][ T_FFE ] = DiffIntFFE;
2046
    ProdFuncs[ T_FFE ][ T_FFE ] = ProdFFEFFE;
2047
    ProdFuncs[ T_FFE ][ T_INT ] = ProdFFEInt;
2048
    ProdFuncs[ T_INT ][ T_FFE ] = ProdIntFFE;
2049
    QuoFuncs[  T_FFE ][ T_FFE ] = QuoFFEFFE;
2050
    QuoFuncs[  T_FFE ][ T_INT ] = QuoFFEInt;
2051
    QuoFuncs[  T_INT ][ T_FFE ] = QuoIntFFE;
2052
    PowFuncs[  T_FFE ][ T_INT ] = PowFFEInt;
2053
    PowFuncs[  T_FFE ][ T_FFE ] = PowFFEFFE;
2054

2055
    /* return success                                                      */
2056
    return 0;
2057
}
2058

2059

2060
/****************************************************************************
2061
**
2062
*F  InitLibrary( <module> ) . . . . . . .  initialise library data structures
2063
*/
2064
static Int InitLibrary (
2065
    StructInitInfo *    module )
2066
{
2067
    /* create the fields and integer conversion bags                       */
2068
    CharFF = NEW_PLIST( T_PLIST, 0 );
2069
    SET_LEN_PLIST( CharFF, 0 );
2070

2071
    DegrFF = NEW_PLIST( T_PLIST, 0 );
2072
    SET_LEN_PLIST( DegrFF, 0 );
2073

2074
    SuccFF = NEW_PLIST( T_PLIST, 0 );
2075
    SET_LEN_PLIST( SuccFF, 0 );
2076

2077
    TypeFF = NEW_PLIST( T_PLIST, 0 );
2078
    SET_LEN_PLIST( TypeFF, 0 );
2079

2080
    TypeFF0 = NEW_PLIST( T_PLIST, 0 );
2081
    SET_LEN_PLIST( TypeFF0, 0 );
2082

2083
    IntFF = NEW_PLIST( T_PLIST, 0 );
2084
    SET_LEN_PLIST( IntFF, 0 );
2085

2086
    /* init filters and functions                                          */
2087
    InitGVarFiltsFromTable( GVarFilts );
2088
    InitGVarFuncsFromTable( GVarFuncs );
2089
    HDLR_FUNC(VAL_GVAR(GVarName("Z")),2) = FuncZ2;
2090

2091
    /* return success                                                      */
2092
    return 0;
2093
}
2094

2095

2096
/****************************************************************************
2097
**
2098
*F  InitInfoFinfield()  . . . . . . . . . . . . . . . table of init functions
2099
*/
2100
static StructInitInfo module = {
2101
    MODULE_BUILTIN,                     /* type                           */
2102
    "finfield",                         /* name                           */
2103
    0,                                  /* revision entry of c file       */
2104
    0,                                  /* revision entry of h file       */
2105
    0,                                  /* version                        */
2106
    0,                                  /* crc                            */
2107
    InitKernel,                         /* initKernel                     */
2108
    InitLibrary,                        /* initLibrary                    */
2109
    0,                                  /* checkInit                      */
2110
    0,                                  /* preSave                        */
2111
    0,                                  /* postSave                       */
2112
    0                                   /* postRestore                    */
2113
};
2114

2115
StructInitInfo * InitInfoFinfield ( void )
2116
{
2117
    return &module;
2118
}
2119

2120

2121
/****************************************************************************
2122
**
2123

2124
*E  finfield.c  . . . . . . . . . . . . . . . . . . . . . . . . . . ends here
2125
*/
2126

2127
Real-time collaboration for Jupyter Notebooks, Linux Terminals, LaTeX, VS Code, R IDE, and more,
all in one place.

Product

Resources

Company

Real-time collaboration for Jupyter Notebooks, Linux Terminals, LaTeX, VS Code, R IDE, and more, all in one place.

Real-time collaboration for Jupyter Notebooks, Linux Terminals, LaTeX, VS Code, R IDE, and more,
all in one place.
