```  1 /*==============================================================================
2 |
3 |  File Name:
4 |
5 |     ppOrder.c
6 |
7 |  Description:
8 |
9 |     Routines for testing the order of polynomials.
10 |
11 |  Functions:
12 |
13 |     order_m
14 |     order_r
15 |     maximal_order
16 |
17 |  LEGAL
18 |
19 |     Primpoly Version 16.1 - A Program for Computing Primitive Polynomials.
21 |
22 |     This program is free software: you can redistribute it and/or modify
24 |     the Free Software Foundation, either version 3 of the License, or
25 |     (at your option) any later version.
26 |
27 |     This program is distributed in the hope that it will be useful,
28 |     but WITHOUT ANY WARRANTY; without even the implied warranty of
29 |     MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
30 |     GNU General Public License for more details.
31 |
32 |     You should have received a copy of the GNU General Public License
33 |     along with this program.  If not, see <http://www.gnu.org/licenses/>.
34 |
35 |     The author's address is artificer!AT!seanerikoconnor!DOT!freeservers!DOT!com
36 |     with !DOT! replaced by . and the !AT! replaced by @
37 |
38 ==============================================================================*/
39
40 /*------------------------------------------------------------------------------
41 |                                Include Files                                 |
42 ------------------------------------------------------------------------------*/
43
44 #include <stdio.h>
45 #include <stdlib.h> /* for rand() function */
46
47 #include "Primpoly.h"
48
49
50 /*==============================================================================
51 |                                order_m                                       |
52 ================================================================================
53
54 DESCRIPTION
55                 m
56     Check that x  (mod f(x), p) is not an integer for m = r / p  but skip
57                                            n                   i
58                                           p  - 1
59     this test if p  | (p-1).  Recall r = -------.
60                   i                       p - 1
61
62 INPUT
63
64      power_table (int **)     x ^ k (mod f(x), p) for n <= k <= 2n-2, f monic.
65      n      (int, n >= 1)     Degree of f(x).
66      p      (int)             Modulo p coefficient arithmetic.
67      r (int)                  See above.
68      primes (bigint *)        Distinct prime factors of r.
69      prime_count              Number of primes.
70
71 RETURNS
72
73      YES    if all tests are passed.
74      NO     if any one test fails.
75
76 EXAMPLE
77                                            2
78      Let n = 4 and p = 5.  Then r = 156 = 2 * 3 * 13, and p = 2, p = 3,
79                                                             1      2
80      and p = 13.  m = { 156/2, 156/3, 156/13 } = { 78, 52, 12 }.  We can
81           3
82      skip the test for m = 78 because p = 2 divides p-1 = 4.  Exponentiation
83                                        1
84             52       3   2                                      12
85      gives x    = 2 x + x + 4 x, which is not an integer and x =
86         3       2
87      4 x  +  2 x  +  4 x  + 3 which is not an integer either, so we return
88
89      YES.
90
91 METHOD
92
93     Exponentiate x with x_to_power and test the result with is_integer.
94     Return right away if the result is not an integer.
95
96 BUGS
97
98     None.
99
100 --------------------------------------------------------------------------------
101 |                                Function Call                                 |
102 ------------------------------------------------------------------------------*/
103
104 int
105     order_m( int power_table[][ MAXDEGPOLY ], int n, int p, bigint r,
106              bigint * primes, int prime_count )
107 {
108
109 /*------------------------------------------------------------------------------
110 |                               Local Variables                                |
111 ------------------------------------------------------------------------------*/
112
113 int
114     i,                  /*  Loop counter.  */
115     g[ MAXDEGPOLY ] ;   /* g(x) = x ^ m (mod f(x), p) */
116
117 bigint
118     m ;                 /*  Exponent of m. */
119
120 /*------------------------------------------------------------------------------
121 |                                Function Body                                 |
122 ------------------------------------------------------------------------------*/
123
124 for (i = 0 ;  i <= prime_count ;  ++i)
125
126     if (!skip_test( i, primes, p ))
127     {
128         m = r / primes[ i ] ;
129
130         x_to_power( m, g, power_table, n, p ) ;
131
132         #ifdef DEBUG_PP_PRIMPOLY
133         printf( "    order m test for prime = %lld, x^ m = x ^ %lld = ", primes[i], m ) ;
134         write_poly( g, n-1 ) ;
135         printf( "\n\n" );
136         #endif
137
138         if (is_integer( g, n-1 ))
139
140             return( NO ) ;
141     }
142
143 return( YES ) ;
144
145 } /* ====================== end of function order_m ========================= */
146
147
148
149 /*==============================================================================
150 |                                order_r                                       |
151 ================================================================================
152
153 DESCRIPTION
154               r
155     Check if x  (mod f(x), p) = a, where a is an integer.
156
157 INPUT
158
159      power_table (int **)     x ^ k (mod f(x), p) for n <= k <= 2n-2, f monic.
160      n      (int, n >= 1)     Degree of f(x).
161      p      (int)             Modulo p coefficient arithmetic.
162      r (int)                  See above.
163      a (int *)                Pointer to value of a.
164
165 RETURNS
166
167      YES    if the formula above is true.
168      NO     if it isn't.
169
170 EXAMPLE
171                           4    2
172      Let r = 156, f(x) = x  + x  + 2 x + 3, n = 4 and p = 5.  It turns
173                156                           4
174      out that x    = 3 (mod f(x), 5) = (-1)  * 3, so we return YES.
175
176 METHOD
177                           r
178     First compute g(x) = x (mod f(x), p).
179     Then test if g(x) is a constant polynomial,
180
181 BUGS
182
183     None.
184
185 --------------------------------------------------------------------------------
186 |                                Function Call                                 |
187 ------------------------------------------------------------------------------*/
188
189 int
190     order_r( int power_table[][ MAXDEGPOLY ], int n, int p, bigint r, int * a )
191 {
192
193 /*------------------------------------------------------------------------------
194 |                               Local Variables                                |
195 ------------------------------------------------------------------------------*/
196
197 int
198     g[ MAXDEGPOLY ] ;   /* g(x) = x ^ m (mod f(x), p) */
199
200 /*------------------------------------------------------------------------------
201 |                                Function Body                                 |
202 ------------------------------------------------------------------------------*/
203
204 x_to_power( r, g, power_table, n, p ) ;
205
206 #ifdef DEBUG_PP_PRIMPOLY
207 printf( "    order r test for x^r = x ^ %lld = ", r ) ;
208 write_poly( g, n-1 ) ;
209 printf( "\n\n" );
210 #endif
211
212 /*  Return the value a = constant term of g(x) */
213 *a = g[ 0 ] ;
214
215 return( is_integer( g, n-1 ) ? YES : NO  ) ;
216
217 } /* ====================== end of function order_r ========================= */
218
219
220
221 /*==============================================================================
222 |                                  maximal_order                               |
223 ================================================================================
224
225 DESCRIPTION
226               k                                  n
227     Check if x  = 1 (mod f(x), p) only when k = p  - 1 and not for any smaller
228     power of k, i.e. that f(x) is a primitive polynomial.
229
230 INPUT
231
232      f (int *)                Monic polynomial f(x).
233      n      (int, n >= 1)     Degree of f(x).
234      p      (int)             Modulo p coefficient arithmetic.
235
236 RETURNS
237
238      YES    if f( x ) is primitive.
239      NO     if it isn't.
240
241 EXAMPLE
242
243                4
244      f( x ) = x  + x  +  1 is a primitive polynomial modulo p = 2,
245                                          4
246      because it generates the group GF( 2  ) with the exception of
247                             2         15
248      zero.  The powers {x, x , ... , x  } modulo f(x), mod 2 are
249                                         16       4
250      distinct and not equal to 1 until x   (mod x  + x + 1, 2) = 1.
251
252 METHOD
253
254     Confirm f(x) is primitive using the definition of primitive
255     polynomial as a generator of the Galois group
256          n                   n
257     GF( p ) by testing that p - 1 is the smallest power for which
258      k
259     x = 1 (mod f(x), p).
260
261 BUGS
262
263     None.
264
265 --------------------------------------------------------------------------------
266 |                                Function Call                                 |
267 ------------------------------------------------------------------------------*/
268
269 int maximal_order( int * f, int n, int p )
270 {
271     int g[ MAXDEGPOLY ] ;   /* g(x) = x ^ m (mod f(x), p) */
272     bigint maxOrder ;
273     bigint k ;
274     int power_table[ MAXDEGPOLY - 1 ] [ MAXDEGPOLY ] ;    /*  x ^ n , ... , x ^ 2n-2 (mod f(x), p) */
275
276     /*                         n         2n-2
277         Precompute the powers x ,  ..., x     (mod f(x), p)
278         for use in all later computations.
279     */
280     construct_power_table( power_table, f, n, p ) ;
281
282     /*  Highest possible order for x. */
283     maxOrder = power( p, n ) - 1 ;
284
285     for (k = 1 ;  k <= maxOrder ;  ++k)
286     {
287         x_to_power( k, g, power_table, n, p ) ;
288
289         if (is_integer( g, n-1 ) &&
290             g[0] == 1 &&
291             k < maxOrder)
292         {
293             return 0 ;
294         }
295
296     } /* end for k */
297
298     return 1 ;
299
300 } /* ================= end of function maximal_order ======================== */
```