1/*==============================================================================
2|
3| NAME
4|
5| Primpoly.cpp
6|
7| DESCRIPTION
8|
9| Program for finding primitive polynomials of degree n modulo p for
10| any prime p >= 2 and any n >= 2.
11|
12| Useful for generating PN sequences and finite fields for error control coding.
13|
14| Please see the user manual and complete technical documentation at
15|
16| http://seanerikoconnor.freeservers.com/Mathematics/AbstractAlgebra/PrimitivePolynomials/overview.html
17|
18| This is a console app to be run in a terminal window. The OS calls main(), which returns a status to the OS.
19|
20| LEGAL
21|
22| Primpoly Version 16.3 - A Program for Computing Primitive Polynomials.
23| Copyright (C) 1999-2024 by Sean Erik O'Connor. All Rights Reserved.
24|
25| This program is free software: you can redistribute it and/or modify
26| it under the terms of the GNU General Public License as published by
27| the Free Software Foundation, either version 3 of the License, or
28| (at your option) any later version.
29|
30| This program is distributed in the hope that it will be useful,
31| but WITHOUT ANY WARRANTY; without even the implied warranty of
32| MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
33| GNU General Public License for more details.
34|
35| You should have received a copy of the GNU General Public License
36| along with this program. If not, see http://www.gnu.org/licenses/.
37|
38| The author's address is seanerikoconnor!AT!gmail!DOT!com
39| with the !DOT! replaced by . and the !AT! replaced by @
40|
41==============================================================================*/
42
43
44/*------------------------------------------------------------------------------
45| Include Files |
46------------------------------------------------------------------------------*/
47
48#include <cstdlib> // abort()
49#include <iostream> // Basic stream I/O.
50#include <new> // set_new_handler()
51#include <cmath> // Basic math functions e.g. sqrt()
52#include <limits> // Numeric limits.
53#include <complex> // Complex data type and operations.
54#include <fstream> // File stream I/O.
55#include <sstream> // String stream I/O.
56#include <vector> // STL vector class.
57#include <string> // STL string class.
58#include <algorithm> // Iterators.
59#include <stdexcept> // Exceptions.
60#include <cassert> // assert()
61
62using namespace std ; // I don't want to use the std:: prefix everywhere.
63
64#include "ctype.h" // C string functions.
65
66
67/*------------------------------------------------------------------------------
68| PP Include Files |
69------------------------------------------------------------------------------*/
70
71#include "Primpoly.hpp" // Global functions.
72#include "ppArith.hpp" // Basic arithmetic functions.
73#include "ppBigInt.hpp" // Arbitrary precision integer arithmetic.
74#include "ppOperationCount.hpp" // OperationCount collection for factoring and poly finding.
75#include "ppFactor.hpp" // Prime factorization and Euler Phi.
76#include "ppPolynomial.hpp" // Polynomial operations and mod polynomial operations.
77#include "ppParser.hpp" // Parsing of polynomials and I/O services.
78#include "ppUnitTest.hpp" // Complete unit test.
79
80
81/*=============================================================================
82|
83| NAME
84|
85| Message strings.
86|
87| DESCRIPTION
88|
89| Helpful information returned by top level main() function.
90|
91+============================================================================*/
92
93static const string legalNotice
94{
95 "\n"
96 "Primpoly Version 16.3 - A Program for Computing Primitive Polynomials.\n"
97 "Copyright (C) 1999-2024 by Sean Erik O'Connor. All Rights Reserved.\n"
98 "\n"
99 "Primpoly comes with ABSOLUTELY NO WARRANTY; for details see the\n"
100 "GNU General Public License. This is free software, and you are welcome\n"
101 "to redistribute it under certain conditions; see the GNU General Public License\n"
102 "for details.\n\n"
103} ;
104
105static const string helpText
106{
107 "This program generates primitive polynomials of degree n modulo p.\n"
108 "\n"
109 "Usage: Generate a single random polynomial of degree n modulo p where p is a prime >= 2 and n is an integer >= 2\n"
110 " Primpoly p n\n"
111 "Example:\n"
112 " Primpoly 2 4\n"
113 " Self-check passes...\n"
114 " Primitive polynomial modulo 2 of degree 4\n"
115 " x ^ 4 + x + 1, 2\n"
116 "Usage: Test whether a polynomial is primitive modulo p.\n"
117 " Primpoly -t ""<Polynomial to test>, p""\n"
118 " If you leave off the modulus p we default to p = 2\n"
119 "Examples:\n"
120 " Primpoly -t ""x^4 + x + 1, 2""\n"
121 " Self-check passes...\n"
122 " x ^ 4 + x + 1, 2 is primitive!\n"
123 "\n"
124 " Primpoly -t ""x^4 + x + 1""\n"
125 " Self-check passes...\n"
126 " x ^ 4 + x + 1, 2 is primitive!\n"
127 "Usage: Generate all primitive polynomial of degree n modulo p.\n"
128 " Primpoly -a p n\n"
129 "Example:\n"
130 " Primpoly -a 2 4\n"
131 " Self-check passes...\n"
132 " Primitive polynomial modulo 2 of degree 4\n"
133 " x ^ 4 + x + 1, 2\n"
134 " Primitive polynomial modulo 2 of degree 4\n"
135 " x ^ 4 + x ^ 3 + 1, 2\n"
136 "Usage: Same but show computation statistics.\n"
137 " Primpoly -s p n\n"
138 "Examples: \n"
139 "\n"
140 " Primpoly.exe -s 13 19\n"
141 " Self-check passes...\n"
142 " Primitive polynomial modulo 13 of degree 19\n"
143 " x ^ 19 + 9 x + 2, 13\n"
144 "\n"
145 " +--------- OperationCount --------------------------------\n"
146 " |\n"
147 " | Integer factorization: Table lookup + Trial division + Pollard Rho\n"
148 " |\n"
149 " | Number of trial divisions : 0\n"
150 " | Number of gcd's computed : 9027\n"
151 " | Number of primality tests : 2\n"
152 " | Number of squarings: 9026\n"
153 " |\n"
154 " | Polynomial Testing\n"
155 " |\n"
156 " | Total num. degree 19 poly mod 13 : 1461920290375446110677\n"
157 " | Number of possible primitive poly: 6411930599771980992\n"
158 " | Polynomials tested : 120\n"
159 " | Const. coeff. was primitive root : 46\n"
160 " | Free of linear factors : 11\n"
161 " | Irreducible to power >=1 : 1\n"
162 " | Had order r (x^r = integer) : 1\n"
163 " | Passed const. coeff. test : 1\n"
164 " | Had order m (x^m != integer) : 1\n"
165 " |\n"
166 " +-----------------------------------------------------\n"
167 "Usage: Print help message.\n"
168 " Primpoly -h\n"
169 " <Prints this help message.>\n"
170 "\n\n"
171 "Primitive polynomials find many uses in mathematics and communications\n"
172 "engineering:\n"
173 " * Generation of pseudonoise (PN) sequences for spread spectrum\n"
174 " communications and chip fault testing.\n"
175 " * Generating Sobol sequences for high dimensional numerical integration.\n"
176 " * Generation of CRC and Hamming codes.\n"
177 " * Generation of Galois (finite) fields for use in decoding Reed-Solomon\n"
178 " and BCH error correcting codes.\n"
179 "\n"
180 "For detailed technical information, see my web page\n"
181 " http://seanerikoconnor.freeservers.com/Mathematics/AbstractAlgebra/PrimitivePolynomials/overview.html\n"
182 "\n"
183} ;
184
185
186/*=============================================================================
187|
188| NAME
189|
190| main
191|
192| DESCRIPTION
193|
194| Program for finding primitive polynomials of degree n modulo p for
195| any prime p >= 2 and any n >= 2.
196|
197| Useful for generating PN sequences and finite fields for error control coding.
198|
199| Run the program by typing
200|
201| ./Primpoly p n
202|
203| You will get an nth degree primitive polynomial modulo p.
204|
205| To see all the options type
206|
207| ./Primpoly -h
208|
209| Please see the user manual and complete technical documentation at
210|
211| http://seanerikoconnor.freeservers.com/Mathematics/AbstractAlgebra/PrimitivePolynomials/overview.html
212|
213| The author's address is seanerikoconnor!AT!gmail!DOT!com
214| with the !DOT! replaced by . and the !AT! replaced by @
215|
216+============================================================================*/
217
218int main( int argc, const char * argv[] )
219{
220 try
221 {
222 // Make my lawyers happy.
223 cout << legalNotice ;
224
225 // Set up the full parser for both command line parsing and polynomial parsing.
226 PolyParser<PolySymbol, PolyValue> parser ;
227 parser.parseCommandLine( argc, argv ) ;
228
229 // Print the help message only and exit.
230 if (parser.print_help_)
231 {
232 cout << helpText ;
233 return static_cast<int>( ReturnStatus::AskForHelp ) ;
234 }
235
236 #ifdef SELF_CHECK
237 // Do a self check always. We might fail to pass one or more unit tests, or the unit test itself might fail.
238 try
239 {
240 UnitTest unitTest ;
241 if (!unitTest.run())
242 throw PrimpolyError( "Self-check failed!" ) ;
243 else
244 cout << "Self-check passes..." << endl ;
245 }
246 catch (UnitTestError & e)
247 {
248 throw PrimpolyError( static_cast<string>( "Could not run the self-check!\n" ) + " [ " + e.what() + " ] " ) ;
249 }
250 #endif
251
252 // The user input a polynomial. Test it for primitivity.
253 if (parser.test_polynomial_for_primitivity_)
254 {
255 // Test for primitivity with the quick test.
256 Polynomial f( parser.test_polynomial_ ) ;
257 #ifdef DEBUG_PP_PRIMITIVITY
258 cout << "Factoring into primes r = (p^n-1)/(p-1) = " << " for n = " << parser.n << " p = " << parser.p << endl ;
259 #endif // DEBUG_PP_PRIMITIVITY
260 PolyOrder order( f ) ;
261 cout << f << " is " << (order.isPrimitive() ? "" : "NOT") << " primitive!" << endl ;
262
263 if (parser.print_operation_count_)
264 cout << order.statistics_ << endl ;
265
266 // Also do a very slow maximal order test for primitivity, if asked to do so.
267 if (parser.slow_confirm_)
268 {
269 cout << confirmWarning ;
270 cout << " confirmed " << (order.maximal_order() ? "" : "NOT") << " primitive!" << endl ;
271 }
272 }
273 // Find one primitive polynomial at random. Optionally, find all primitive polynomials.
274 else
275 {
276 // Generate and test all possible n th degree, monic, modulo p polynomials
277 // f(x). A polynomial is primitive if passes all the tests successfully.
278 // n
279 // Initialize f(x) to x + (-1). Then, when f(x) passes through function
280 // n
281 // nextTrialPoly for the first time, it will have the correct value, x
282 TrialPolynomial f ;
283 f.initialTrialPoly( parser.n, parser.p ) ;
284 #ifdef DEBUG_PP_PRIMITIVITY
285 cout << "Factoring into primes r = (p^n-1)/(p-1) = " << " for n = " << parser.n << " p = " << parser.p << endl ;
286 #endif // DEBUG_PP_PRIMITIVITY
287 PolyOrder order( f ) ;
288
289 bool is_primitive_poly{ false } ;
290 bool tried_all_poly{ false } ;
291 bool stopTesting{ false } ;
292
293 BigInt num_poly( 0u ) ;
294 BigInt num_primitive_poly( 0u ) ;
295
296 if (parser.list_all_primitive_polynomials_)
297 cout << "\n\nThere are " << order.getNumPrimPoly() << " primitive polynomials modulo " << f.modulus() << " of degree " << f.deg() << "\n\n" ;
298
299 do {
300 ++num_poly ;
301
302 #ifdef DEBUG_PP_PRIMITIVITY
303 cout << "Testing polynomial # " << num_poly << ") p(x) = " << f << " for primitivity" << endl ;
304 #endif // DEBUG_PP_PRIMITIVITY
305
306 order.resetPolynomial( f ) ;
307 is_primitive_poly = order.isPrimitive() ;
308
309 if (is_primitive_poly)
310 {
311 ++num_primitive_poly ;
312 cout << "\n\nPrimitive polynomial modulo " << f.modulus() << " of degree " << f.deg() << "\n\n" ;
313 cout << f ;
314 cout << endl << endl ;
315
316 // Do a very slow maximal order test for primitivity.
317 if (parser.slow_confirm_)
318 {
319 cout << confirmWarning ;
320 if (order.maximal_order())
321 cout << f << " confirmed primitive!" << endl ;
322 else
323 {
324 ostringstream os ;
325 os << "Fast test says " << f << " is a primitive polynomial but slow test disagrees.\n" ;
326 throw PolynomialError( os.str() ) ;
327 }
328 }
329
330 // Early out if we've found all the primitive polynomials.
331 if (num_primitive_poly >= order.getNumPrimPoly())
332 break ;
333 }
334
335 tried_all_poly = (num_poly >= order.getMaxNumPoly()) ;
336 stopTesting = tried_all_poly || (!parser.list_all_primitive_polynomials_ && is_primitive_poly) ;
337
338 f.nextTrialPoly() ; // Try next polynomal in sequence.
339 } while( !stopTesting ) ;
340
341 if (parser.print_operation_count_)
342 cout << order.statistics_ << endl ;
343
344 // Didn't find a primitive polynomial in the find-only-one-primitive-polynomial case, which is an error.
345 if (!parser.list_all_primitive_polynomials_ && !is_primitive_poly)
346 {
347 ostringstream os ;
348 os << "Tested all " << order.getMaxNumPoly() << " possible polynomials, but failed to find a primitive polynomial" ;
349 throw PolynomialError( os.str() ) ;
350 }
351 }
352 return static_cast<int>( ReturnStatus::Success ) ;
353 }
354 // Catch all exceptions and report what happened to the user.
355 // First do the user-defined exceptions.
356 catch( PrimpolyError & e )
357 {
358 cerr << "\nTop Level Error: " << " [ " << e.what() << " ] " << endl << writeToAuthorMessage ;
359 return static_cast<int>( ReturnStatus::InternalError ) ;
360 }
361 catch( ParserError & e )
362 {
363 cerr << "Inputs are incorrect or out of range: " << " [ " << e.what() << " ] " << endl << helpText ;
364 return static_cast<int>( ReturnStatus::RangeError ) ;
365 }
366 catch( FactorError & e )
367 {
368 cerr << "Error in prime factorization: " << " [ " << e.what() << " ] " << endl << writeToAuthorMessage ;
369 return static_cast<int>( ReturnStatus::InternalError ) ;
370 }
371 catch( BigIntRangeError & e )
372 {
373 cerr << "Internal range error in multiple precision arithmetic: " << " [ " << e.what() << " ] " << endl << writeToAuthorMessage ;
374 return static_cast<int>( ReturnStatus::InternalError ) ;
375 }
376 catch( BigIntDomainError & e )
377 {
378 cerr << "Internal domain error in multiple precision arithmetic: " << " [ " << e.what() << " ] " << endl << writeToAuthorMessage ;
379 return static_cast<int>( ReturnStatus::InternalError ) ;
380 }
381 catch( BigIntUnderflow & e )
382 {
383 cerr << "Internal underflow error in multiple precision arithmetic: " << " [ " << e.what() << " ] " << endl << writeToAuthorMessage ;
384 return static_cast<int>( ReturnStatus::InternalError ) ;
385 }
386 catch( BigIntOverflow & e )
387 {
388 cerr << "Internal overflow error in multiple precision arithmetic: " << " [ " << e.what() << " ] " << endl << writeToAuthorMessage ;
389 return static_cast<int>( ReturnStatus::InternalError ) ;
390 }
391 catch( BigIntZeroDivide & e )
392 {
393 cerr << "Internal zero divide error in multiple precision arithmetic: " << " [ " << e.what() << " ] " << endl << writeToAuthorMessage ;
394 return static_cast<int>( ReturnStatus::InternalError ) ;
395 }
396 catch( BigIntMathError & e )
397 {
398 cerr << "Internal math error in multiple precision arithmetic: " << " [ " << e.what() << " ] " << endl << writeToAuthorMessage ;
399 return static_cast<int>( ReturnStatus::InternalError ) ;
400 }
401 catch( ArithModPError & e )
402 {
403 cerr << "Internal modulo p arithmetic error: " << " [ " << e.what() << " ] " << endl << writeToAuthorMessage ;
404 return static_cast<int>( ReturnStatus::InternalError ) ;
405 }
406 catch( PolynomialRangeError & e )
407 {
408 cout << "Error. Polynomial has bad syntax or coefficients are out of range. " << " [ " << e.what() << " ] " << endl << helpText ;
409 return static_cast<int>( ReturnStatus::RangeError ) ;
410 }
411 catch( PolynomialError & e )
412 {
413 cerr << "Error in polynomial arithmetic: " << " [ " << e.what() << " ] " << endl << writeToAuthorMessage ;
414 return static_cast<int>( ReturnStatus::InternalError ) ;
415 }
416 // Standard library exceptions.
417 catch( bad_alloc & e )
418 {
419 cerr << "Error allocating memory: " << " [ " << e.what() << " ] "<< endl ;
420 cerr << "Run on a different computer with more RAM or virtual memory." << writeToAuthorMessage ;
421 return static_cast<int>( ReturnStatus::InternalError ) ;
422 }
423 catch( exception & e )
424 {
425 cerr << "System error: " << e.what() << endl << writeToAuthorMessage ;
426 return static_cast<int>( ReturnStatus::InternalError ) ;
427 }
428 // Catch all other uncaught exceptions, which would otherwise call terminate()
429 // which in turn calls abort() and which would halt this program.
430 //
431 // Limitations:
432 // We can't handle the case where terminate() gets called because the
433 // exception mechanism finds a corrupt stack or catches a destructor
434 // throwing an exception.
435 //
436 // Also we can't catch exceptions which are thrown by initializing or
437 // destructing global variables.
438 catch( ... )
439 {
440 cerr << "Unexpected exception: " << endl << writeToAuthorMessage ;
441 return static_cast<int>( ReturnStatus::InternalError ) ;
442 }
443
444 return static_cast<int>( ReturnStatus::Success ) ;
445
446} //========================== end of function main ========================
