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.1 - A Program for Computing Primitive Polynomials.
 23 |     Copyright (C) 1999-2021 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 
 62 using 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.h"         // Global functions.
 72 #include "ppArith.h"          // Basic arithmetic functions.
 73 #include "ppBigInt.h"         // Arbitrary precision integer arithmetic.
 74 #include "ppOperationCount.h" // OperationCount collection for factoring and poly finding.
 75 #include "ppFactor.h"         // Prime factorization and Euler Phi.
 76 #include "ppPolynomial.h"     // Polynomial operations and mod polynomial operations.
 77 #include "ppParser.h"         // Parsing of polynomials and I/O services.
 78 #include "ppUnitTest.h"       // 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 
 93 static const string legalNotice
 94 {
 95     "\n"
 96     "Primpoly Version 16.1 - A Program for Computing Primitive Polynomials.\n"
 97     "Copyright (C) 1999-2021 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 
105 static 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 
218 int 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             PolyOrder order( f ) ;
258             cout << f << " is " << (order.isPrimitive() ? "" : "NOT") << " primitive!" << endl ;
259 
260             if (parser.print_operation_count_)
261                 cout << order.statistics_ << endl ;
262 
263             // Also do a very slow maximal order test for primitivity, if asked to do so.
264             if (parser.slow_confirm_)
265             {
266                 cout << confirmWarning ;
267                 cout << " confirmed " << (order.maximal_order() ? "" : "NOT") << " primitive!" << endl ;
268             }
269         }
270         //  Find one primitive polynomial at random.  Optionally, find all primitive polynomials.
271         else
272         {
273             //   Generate and test all possible n th degree, monic, modulo p polynomials
274             //   f(x).  A polynomial is primitive if passes all the tests successfully.
275             //                       n
276             //   Initialize f(x) to x  + (-1).  Then, when f(x) passes through function
277             //                                                                        n
278             //   nextTrialPoly for the first time, it will have the correct value, x
279             Polynomial f ;
280             f.initialTrialPoly( parser.n, parser.p ) ;
281             PolyOrder order( f ) ;
282 
283             bool is_primitive_poly{ false } ;
284             bool tried_all_poly{ false } ;
285             bool stopTesting{ false } ;
286 
287             BigInt num_poly( 0u ) ;
288             BigInt num_primitive_poly( 0u ) ;
289 
290             if (parser.list_all_primitive_polynomials_)
291                 cout << "\n\nThere are " << order.getNumPrimPoly() << " primitive polynomials modulo " << f.modulus() << " of degree " << f.deg() << "\n\n" ;
292 
293             do {
294                 ++num_poly ;
295 
296                 #ifdef DEBUG_PP_POLYNOMIAL
297                 cout << "Testing polynomial # " << num_poly << ") p(x) = " << f << " for primitivity" << endl ;
298                 #endif // DEBUG_PP_POLYNOMIAL
299 
300                 order.resetPolynomial( f ) ;
301                 is_primitive_poly = order.isPrimitive() ;
302 
303                 if (is_primitive_poly)
304                 {
305                     ++num_primitive_poly ;
306                     cout << "\n\nPrimitive polynomial modulo " << f.modulus() << " of degree " << f.deg() << "\n\n" ;
307                     cout << f ;
308                     cout << endl << endl ;
309 
310                     // Do a very slow maximal order test for primitivity.
311                     if (parser.slow_confirm_)
312                     {
313                         cout << confirmWarning ;
314                         if (order.maximal_order())
315                             cout << f << " confirmed primitive!" << endl ;
316                         else
317                         {
318                             ostringstream os ;
319                             os << "Fast test says " << f << " is a primitive polynomial but slow test disagrees.\n" ;
320                             throw PolynomialError( os.str() ) ;
321                         }
322                     }
323 
324                     // Early out if we've found all the primitive polynomials.
325                     if (num_primitive_poly >= order.getNumPrimPoly())
326                         break ;
327                 }
328 
329                 tried_all_poly = (num_poly >= order.getMaxNumPoly()) ;
330                 stopTesting = tried_all_poly || (!parser.list_all_primitive_polynomials_ && is_primitive_poly) ;
331 
332                 f.nextTrialPoly() ;      // Try next polynomal in sequence.
333             } while( !stopTesting ) ;
334 
335             if (parser.print_operation_count_)
336                 cout << order.statistics_ << endl ;
337 
338             // Didn't find a primitive polynomial in the find-only-one-primitive-polynomial case, which is an error.
339             if (!parser.list_all_primitive_polynomials_ && !is_primitive_poly)
340             {
341                 ostringstream os ;
342                 os << "Tested all " << order.getMaxNumPoly() << " possible polynomials, but failed to find a primitive polynomial" ;
343                 throw PolynomialError( os.str() ) ;
344             }
345         }
346         return static_cast<int>( ReturnStatus::Success ) ;
347     }
348     // Catch all exceptions and report what happened to the user.
349     // First do the user-defined exceptions.
350     catch( PrimpolyError & e )
351     {
352         cerr << "\nTop Level Error: " << " [ " << e.what() << " ] " << endl << writeToAuthorMessage ;
353         return static_cast<int>( ReturnStatus::InternalError ) ;
354     }
355     catch( ParserError & e )
356     {
357         cerr << "Inputs are incorrect or out of range: " << " [ " << e.what() << " ] " << endl << helpText ;
358         return static_cast<int>( ReturnStatus::RangeError ) ;
359     }
360     catch( FactorError & e )
361     {
362         cerr << "Error in prime factorization:  " << " [ " << e.what() << " ] " << endl << writeToAuthorMessage ;
363         return static_cast<int>( ReturnStatus::InternalError ) ;
364     }
365     catch( BigIntRangeError & e )
366     {
367         cerr << "Internal range error in multiple precision arithmetic:  " << " [ " << e.what() << " ] " << endl << writeToAuthorMessage ;
368         return static_cast<int>( ReturnStatus::InternalError ) ;
369     }
370     catch( BigIntDomainError & e )
371     {
372         cerr << "Internal domain error in multiple precision arithmetic:  " << " [ " << e.what() << " ] " << endl << writeToAuthorMessage ;
373         return static_cast<int>( ReturnStatus::InternalError ) ;
374     }
375     catch( BigIntUnderflow & e )
376     {
377         cerr << "Internal underflow error in multiple precision arithmetic:  " << " [ " << e.what() << " ] " << endl << writeToAuthorMessage ;
378         return static_cast<int>( ReturnStatus::InternalError ) ;
379     }
380     catch( BigIntOverflow & e )
381     {
382         cerr << "Internal overflow error in multiple precision arithmetic:  " << " [ " << e.what() << " ] " << endl << writeToAuthorMessage ;
383         return static_cast<int>( ReturnStatus::InternalError ) ;
384     }
385     catch( BigIntZeroDivide & e )
386     {
387         cerr << "Internal zero divide error in multiple precision arithmetic:  " << " [ " << e.what() << " ] " << endl << writeToAuthorMessage ;
388         return static_cast<int>( ReturnStatus::InternalError ) ;
389     }
390     catch( BigIntMathError & e )
391     {
392         cerr << "Internal math error in multiple precision arithmetic:  " << " [ " << e.what() << " ] " << endl << writeToAuthorMessage ;
393         return static_cast<int>( ReturnStatus::InternalError ) ;
394     }
395     catch( ArithModPError & e )
396     {
397         cerr << "Internal modulo p arithmetic error:  " << " [ " << e.what() << " ] " << endl << writeToAuthorMessage ;
398         return static_cast<int>( ReturnStatus::InternalError ) ;
399     }
400     catch( PolynomialRangeError & e )
401     {
402         cout << "Error.  Polynomial has bad syntax or coefficients are out of range. " << " [ " << e.what() << " ] " << endl << helpText ;
403         return static_cast<int>( ReturnStatus::RangeError ) ;
404     }
405     catch( PolynomialError & e )
406     {
407         cerr << "Error in polynomial arithmetic:  " << " [ " << e.what() << " ] " << endl << writeToAuthorMessage ;
408         return static_cast<int>( ReturnStatus::InternalError ) ;
409     }
410     //  Standard library exceptions.
411     catch( bad_alloc & e )
412     {
413         cerr << "Error allocating memory:  " << " [ " << e.what() << " ] "<< endl ;
414         cerr << "Run on a different computer with more RAM or virtual memory." << writeToAuthorMessage ;
415         return static_cast<int>( ReturnStatus::InternalError ) ;
416     }
417     catch( exception & e )
418     {
419         cerr << "System error: " << e.what() << endl << writeToAuthorMessage ;
420         return static_cast<int>( ReturnStatus::InternalError ) ;
421     }
422     // Catch all other uncaught exceptions, which would otherwise call terminate()
423     // which in turn calls abort() and which would halt this program.
424     //
425     // Limitations:
426     //     We can't handle the case where terminate() gets called because the
427     //     exception mechanism finds a corrupt stack or catches a destructor
428     //     throwing an exception.
429     // 
430     //     Also we can't catch exceptions which are thrown by initializing or
431     //     destructing global variables.
432     catch( ... )
433     {
434         cerr << "Unexpected exception: " << endl << writeToAuthorMessage ;
435         return static_cast<int>( ReturnStatus::InternalError ) ;
436     }
437 
438     return static_cast<int>( ReturnStatus::Success ) ;
439 
440 } //========================== end of function main ========================