| 1 | <html> |
|---|
| 2 | <head> |
|---|
| 3 | <title> |
|---|
| 4 | A Tour of NTL: Examples: Big Integers </title> |
|---|
| 5 | </head> |
|---|
| 6 | |
|---|
| 7 | <body bgcolor="#fff9e6"> |
|---|
| 8 | |
|---|
| 9 | <center> |
|---|
| 10 | <img src="arrow1.gif" alt="[Previous]" align=bottom> |
|---|
| 11 | <a href="tour-examples.html"><img src="arrow2.gif" alt="[Up]" align=bottom></a> |
|---|
| 12 | <a href="tour-ex2.html"> <img src="arrow3.gif" alt="[Next]" align=bottom></a> |
|---|
| 13 | </center> |
|---|
| 14 | |
|---|
| 15 | <h1> |
|---|
| 16 | <p align=center> |
|---|
| 17 | A Tour of NTL: Examples: Big Integers |
|---|
| 18 | </p> |
|---|
| 19 | </h1> |
|---|
| 20 | |
|---|
| 21 | <p> <hr> <p> |
|---|
| 22 | |
|---|
| 23 | The first example makes use of the class |
|---|
| 24 | <tt>ZZ</tt>, |
|---|
| 25 | which |
|---|
| 26 | represents "big integers": signed, arbitrary length integers. |
|---|
| 27 | This program reads two big integers <tt>a</tt> and <tt>b</tt>, |
|---|
| 28 | and prints <tt>(a+1)*(b+1)</tt>. |
|---|
| 29 | |
|---|
| 30 | <p> |
|---|
| 31 | <pre> |
|---|
| 32 | #include <NTL/ZZ.h> |
|---|
| 33 | |
|---|
| 34 | NTL_CLIENT |
|---|
| 35 | |
|---|
| 36 | int main() |
|---|
| 37 | { |
|---|
| 38 | ZZ a, b, c; |
|---|
| 39 | |
|---|
| 40 | cin >> a; |
|---|
| 41 | cin >> b; |
|---|
| 42 | c = (a+1)*(b+1); |
|---|
| 43 | cout << c << "\n"; |
|---|
| 44 | } |
|---|
| 45 | </pre> |
|---|
| 46 | |
|---|
| 47 | <p> |
|---|
| 48 | |
|---|
| 49 | This program declares three variables <tt>a</tt>, <tt>b</tt>, |
|---|
| 50 | and <tt>c</tt> of type <tt>ZZ</tt>. |
|---|
| 51 | The values <tt>a</tt> and <tt>b</tt> are read from standard input. |
|---|
| 52 | The value <tt>c</tt> is then computed as <tt>(a+1)*(b+1)</tt>. |
|---|
| 53 | Finally, the value of <tt>c</tt> is printed to the standard output. |
|---|
| 54 | |
|---|
| 55 | <p> |
|---|
| 56 | Note that one can compute with <tt>ZZ</tt>s much as with ordinary |
|---|
| 57 | <tt>int</tt>s, in that most of the standard arithmetic and |
|---|
| 58 | assignment operators can be used in a direct and natural way. |
|---|
| 59 | The <tt>C++</tt> compiler and the NTL library routines |
|---|
| 60 | automatically take care |
|---|
| 61 | of all the bookkeeping involved |
|---|
| 62 | with memory management and temporary objects. |
|---|
| 63 | |
|---|
| 64 | <p> |
|---|
| 65 | Note the macro <tt>NTL_CLIENT</tt>. |
|---|
| 66 | When compiling NTL is ISO mode (the default), this |
|---|
| 67 | expands to |
|---|
| 68 | <pre> |
|---|
| 69 | using namespace std; |
|---|
| 70 | using namespace NTL; |
|---|
| 71 | </pre> |
|---|
| 72 | When compiling NTL in traditional mode, this expands |
|---|
| 73 | to the empty string. |
|---|
| 74 | More details <a href="tour-stdcxx.html">here</a>. |
|---|
| 75 | |
|---|
| 76 | <p> <hr> <p> |
|---|
| 77 | |
|---|
| 78 | Here's a program that reads a list of integers from standard |
|---|
| 79 | input and prints the sum of their squares. |
|---|
| 80 | |
|---|
| 81 | <p> |
|---|
| 82 | <pre> |
|---|
| 83 | #include <NTL/ZZ.h> |
|---|
| 84 | |
|---|
| 85 | NTL_CLIENT |
|---|
| 86 | |
|---|
| 87 | int main() |
|---|
| 88 | { |
|---|
| 89 | ZZ acc, val; |
|---|
| 90 | |
|---|
| 91 | acc = 0; |
|---|
| 92 | while (SkipWhiteSpace(cin)) { |
|---|
| 93 | cin >> val; |
|---|
| 94 | acc += val*val; |
|---|
| 95 | } |
|---|
| 96 | |
|---|
| 97 | cout << acc << "\n"; |
|---|
| 98 | } |
|---|
| 99 | </pre> |
|---|
| 100 | |
|---|
| 101 | <p> |
|---|
| 102 | |
|---|
| 103 | The function <tt>SkipWhiteSpace</tt> is defined by NTL. |
|---|
| 104 | It skips over white space, and returns 1 if there is something |
|---|
| 105 | following it. |
|---|
| 106 | This function is useful, because |
|---|
| 107 | NTL's input operators raise an error if an input |
|---|
| 108 | is missing or ill-formed. |
|---|
| 109 | This is different from the standard I/O library, |
|---|
| 110 | which does not raise an error. |
|---|
| 111 | Personally, I find that not raising an error, or at least |
|---|
| 112 | an exception, is a bad idea, since the caller of the I/O |
|---|
| 113 | routine must constantly check the status of the input |
|---|
| 114 | stream. |
|---|
| 115 | |
|---|
| 116 | |
|---|
| 117 | |
|---|
| 118 | |
|---|
| 119 | <p> |
|---|
| 120 | <hr> |
|---|
| 121 | <p> |
|---|
| 122 | |
|---|
| 123 | Here's a simple modular exponentiation routine for computing |
|---|
| 124 | <tt>a^e mod n</tt>. |
|---|
| 125 | NTL already provides a more sophisticated one, though. |
|---|
| 126 | |
|---|
| 127 | <p> |
|---|
| 128 | <pre> |
|---|
| 129 | ZZ PowerMod(const ZZ& a, const ZZ& e, const ZZ& n) |
|---|
| 130 | { |
|---|
| 131 | if (e == 0) return to_ZZ(1); |
|---|
| 132 | |
|---|
| 133 | long k = NumBits(e); |
|---|
| 134 | |
|---|
| 135 | ZZ res; |
|---|
| 136 | res = 1; |
|---|
| 137 | |
|---|
| 138 | for (long i = k-1; i >= 0; i--) { |
|---|
| 139 | res = (res*res) % n; |
|---|
| 140 | if (bit(e, i) == 1) res = (res*a) % n; |
|---|
| 141 | } |
|---|
| 142 | |
|---|
| 143 | if (e < 0) |
|---|
| 144 | return InvMod(res, n); |
|---|
| 145 | else |
|---|
| 146 | return res; |
|---|
| 147 | } |
|---|
| 148 | </pre> |
|---|
| 149 | <p> |
|---|
| 150 | |
|---|
| 151 | Note that as an alternative, we could implement the inner loop |
|---|
| 152 | as follows: |
|---|
| 153 | |
|---|
| 154 | <pre> |
|---|
| 155 | res = SqrMod(res, n); |
|---|
| 156 | if (bit(e, i) == 1) res = MulMod(res, a, n); |
|---|
| 157 | </pre> |
|---|
| 158 | |
|---|
| 159 | We could also write this as: |
|---|
| 160 | |
|---|
| 161 | <pre> |
|---|
| 162 | SqrMod(res, res, n); |
|---|
| 163 | if (bit(e, i) == 1) MulMod(res, res, a, n); |
|---|
| 164 | </pre> |
|---|
| 165 | |
|---|
| 166 | This illustrates an important point about NTL's programming interface. |
|---|
| 167 | For every function in NTL, there is a procedural version that |
|---|
| 168 | stores its result in its first argument. |
|---|
| 169 | The reason for using the procedural variant is efficieny: |
|---|
| 170 | on every iteration through the above loop, the functional form |
|---|
| 171 | of <tt>SqrMod</tt> will cause a temporary <tt>ZZ</tt> object to |
|---|
| 172 | be created and destroyed, whereas the procedural version |
|---|
| 173 | will not create any temporaries. |
|---|
| 174 | Where performance is critical, the procedural version |
|---|
| 175 | is to be preferred. |
|---|
| 176 | Although it is usually silly to get too worked up about performance, |
|---|
| 177 | it may be reasonable to argue that modular exponentiation |
|---|
| 178 | is an important enough routine that it should be as fast as possible. |
|---|
| 179 | |
|---|
| 180 | <p> |
|---|
| 181 | |
|---|
| 182 | Note that when the functional version of a function |
|---|
| 183 | can be naturally named with an operator, this is done. |
|---|
| 184 | So for example, NTL provides a 3-argument <tt>mul</tt> routine |
|---|
| 185 | for <tt>ZZ</tt> multiplication, and a functional version |
|---|
| 186 | whose name is <tt>operator *</tt>, and not <tt>mul</tt>. |
|---|
| 187 | |
|---|
| 188 | <p> |
|---|
| 189 | |
|---|
| 190 | While we are taking about temporaries, consider the first version |
|---|
| 191 | of the inner loop. |
|---|
| 192 | Execution of the statement |
|---|
| 193 | <pre> |
|---|
| 194 | res = (res*res) % n; |
|---|
| 195 | </pre> |
|---|
| 196 | will result in the creation of two temporary objects, |
|---|
| 197 | one for the product, and one for the result of the mod operation, |
|---|
| 198 | whose value is copied into <tt>res</tt>. |
|---|
| 199 | Of course, the compiler automatically generates the code for |
|---|
| 200 | cleaning up temporaries and other local objects at the right time. |
|---|
| 201 | The programmer does not have to worry about this. |
|---|
| 202 | |
|---|
| 203 | |
|---|
| 204 | <p> <hr> <p> |
|---|
| 205 | |
|---|
| 206 | This example is a bit more interesting. |
|---|
| 207 | The following program prompts the user for an input, |
|---|
| 208 | and applies a simple probabilistic primality test. |
|---|
| 209 | Note that NTL already provides a slightly more sophisticated |
|---|
| 210 | prime test. |
|---|
| 211 | |
|---|
| 212 | <p> |
|---|
| 213 | <pre> |
|---|
| 214 | #include <NTL/ZZ.h> |
|---|
| 215 | |
|---|
| 216 | NTL_CLIENT |
|---|
| 217 | |
|---|
| 218 | long witness(const ZZ& n, const ZZ& x) |
|---|
| 219 | { |
|---|
| 220 | ZZ m, y, z; |
|---|
| 221 | long j, k; |
|---|
| 222 | |
|---|
| 223 | if (x == 0) return 0; |
|---|
| 224 | |
|---|
| 225 | // compute m, k such that n-1 = 2^k * m, m odd: |
|---|
| 226 | |
|---|
| 227 | k = 1; |
|---|
| 228 | m = n/2; |
|---|
| 229 | while (m % 2 == 0) { |
|---|
| 230 | k++; |
|---|
| 231 | m /= 2; |
|---|
| 232 | } |
|---|
| 233 | |
|---|
| 234 | z = PowerMod(x, m, n); // z = x^m % n |
|---|
| 235 | if (z == 1) return 0; |
|---|
| 236 | |
|---|
| 237 | j = 0; |
|---|
| 238 | do { |
|---|
| 239 | y = z; |
|---|
| 240 | z = (y*y) % n; |
|---|
| 241 | j++; |
|---|
| 242 | } while (j < k && z != 1); |
|---|
| 243 | |
|---|
| 244 | return z != 1 || y != n-1; |
|---|
| 245 | } |
|---|
| 246 | |
|---|
| 247 | |
|---|
| 248 | long PrimeTest(const ZZ& n, long t) |
|---|
| 249 | { |
|---|
| 250 | if (n <= 1) return 0; |
|---|
| 251 | |
|---|
| 252 | // first, perform trial division by primes up to 2000 |
|---|
| 253 | |
|---|
| 254 | PrimeSeq s; // a class for quickly generating primes in sequence |
|---|
| 255 | long p; |
|---|
| 256 | |
|---|
| 257 | p = s.next(); // first prime is always 2 |
|---|
| 258 | while (p && p < 2000) { |
|---|
| 259 | if ((n % p) == 0) return (n == p); |
|---|
| 260 | p = s.next(); |
|---|
| 261 | } |
|---|
| 262 | |
|---|
| 263 | // second, perform t Miller-Rabin tests |
|---|
| 264 | |
|---|
| 265 | ZZ x; |
|---|
| 266 | long i; |
|---|
| 267 | |
|---|
| 268 | for (i = 0; i < t; i++) { |
|---|
| 269 | x = RandomBnd(n); // random number between 0 and n-1 |
|---|
| 270 | |
|---|
| 271 | if (witness(n, x)) |
|---|
| 272 | return 0; |
|---|
| 273 | } |
|---|
| 274 | |
|---|
| 275 | return 1; |
|---|
| 276 | } |
|---|
| 277 | |
|---|
| 278 | int main() |
|---|
| 279 | { |
|---|
| 280 | ZZ n; |
|---|
| 281 | |
|---|
| 282 | cout << "n: "; |
|---|
| 283 | cin >> n; |
|---|
| 284 | |
|---|
| 285 | if (PrimeTest(n, 10)) |
|---|
| 286 | cout << n << " is probably prime\n"; |
|---|
| 287 | else |
|---|
| 288 | cout << n << " is composite\n"; |
|---|
| 289 | } |
|---|
| 290 | </pre> |
|---|
| 291 | <p> |
|---|
| 292 | |
|---|
| 293 | Note that in NTL, there are typically a number of ways to |
|---|
| 294 | compute the same thing. |
|---|
| 295 | For example, consider the computation of <tt>m</tt> and <tt>k</tt> |
|---|
| 296 | in function <tt>witness</tt>. |
|---|
| 297 | We could have written it thusly: |
|---|
| 298 | |
|---|
| 299 | <pre> |
|---|
| 300 | k = 1; |
|---|
| 301 | m = n >> 1; |
|---|
| 302 | while (!IsOdd(m)) { |
|---|
| 303 | k++; |
|---|
| 304 | m >>= 1; |
|---|
| 305 | } |
|---|
| 306 | </pre> |
|---|
| 307 | |
|---|
| 308 | It turns out that this is actually not significantly more |
|---|
| 309 | efficient than the original version, because the implementation |
|---|
| 310 | optimizes multiplication and division by 2. |
|---|
| 311 | |
|---|
| 312 | <p> |
|---|
| 313 | |
|---|
| 314 | The following is more efficient: |
|---|
| 315 | |
|---|
| 316 | <pre> |
|---|
| 317 | k = 1; |
|---|
| 318 | while (bit(n, k) == 0) k++; |
|---|
| 319 | m = n >> k; |
|---|
| 320 | </pre> |
|---|
| 321 | |
|---|
| 322 | As it happens, there is a built-in NTL routine that does just what we want: |
|---|
| 323 | |
|---|
| 324 | <pre> |
|---|
| 325 | m = n-1; |
|---|
| 326 | k = MakeOdd(m); |
|---|
| 327 | </pre> |
|---|
| 328 | |
|---|
| 329 | |
|---|
| 330 | |
|---|
| 331 | <p> <hr> <p> |
|---|
| 332 | |
|---|
| 333 | Having seen a number of examples involving <tt>ZZ</tt>s, |
|---|
| 334 | let's look at the <tt>ZZ</tt> interface in a bit more detail. |
|---|
| 335 | |
|---|
| 336 | <p> |
|---|
| 337 | |
|---|
| 338 | <b> |
|---|
| 339 | Constructors, assignment, and conversions |
|---|
| 340 | </b> |
|---|
| 341 | |
|---|
| 342 | <p> |
|---|
| 343 | |
|---|
| 344 | When you declare an object of type <tt>ZZ</tt>, |
|---|
| 345 | the default constructor initializes to the value <tt>0</tt>. |
|---|
| 346 | As we have already seen, there is an assignment operator that |
|---|
| 347 | allows one to copy the value of one <tt>ZZ</tt> to another. |
|---|
| 348 | Note that these copies (like almost all copies in NTL) are "deep", |
|---|
| 349 | i.e., the actual data is copied, and not just a pointer. |
|---|
| 350 | Of course, if the amount of space allocated by the destination |
|---|
| 351 | of the assignment is insufficient to hold the value of the source, |
|---|
| 352 | space is automatically re-allocated. |
|---|
| 353 | |
|---|
| 354 | <p> |
|---|
| 355 | One can also assign a value of type <tt>long</tt> to a <tt>ZZ</tt>: |
|---|
| 356 | <pre> |
|---|
| 357 | ZZ x; |
|---|
| 358 | x = 1; |
|---|
| 359 | </pre> |
|---|
| 360 | |
|---|
| 361 | <p> |
|---|
| 362 | Note that one cannot write |
|---|
| 363 | <pre> |
|---|
| 364 | ZZ x = 1; // error |
|---|
| 365 | </pre> |
|---|
| 366 | to initialize a <tt>ZZ</tt>. |
|---|
| 367 | As a design principle, NTL avoids implicit conversions, and unfortunately, |
|---|
| 368 | the only way to allow initializations like this in <tt>C++</tt> |
|---|
| 369 | is to define an implicit conversion operator. |
|---|
| 370 | Instead, one could write |
|---|
| 371 | <pre> |
|---|
| 372 | ZZ x = to_ZZ(1); |
|---|
| 373 | </pre> |
|---|
| 374 | This is an example of one of NTL's conversion routines. |
|---|
| 375 | For very large constants, one can write: |
|---|
| 376 | <pre> |
|---|
| 377 | ZZ x = to_ZZ("99999999999999999999"); |
|---|
| 378 | </pre> |
|---|
| 379 | These examples illustrate conversion rountines in their |
|---|
| 380 | functional forms. |
|---|
| 381 | The corresponding procedural forms are all called <tt>conv</tt>, e.g., |
|---|
| 382 | <pre> |
|---|
| 383 | ZZ x; |
|---|
| 384 | conv(x, 1); |
|---|
| 385 | conv(x, "99999999999999999999"); |
|---|
| 386 | </pre> |
|---|
| 387 | |
|---|
| 388 | <p> |
|---|
| 389 | <b> |
|---|
| 390 | Functionality |
|---|
| 391 | </b> |
|---|
| 392 | <p> |
|---|
| 393 | |
|---|
| 394 | All of the basic arithmetic operators are supported, |
|---|
| 395 | including comparison, arithmetic, shift, and bit-wise logical operations. |
|---|
| 396 | One can mix <tt>ZZ</tt>s and <tt>long</tt>s in any expresion in |
|---|
| 397 | a natural way. |
|---|
| 398 | As was already mentioned, NTL does not support implicit type conversion; |
|---|
| 399 | rather, for basic operations, it simply overloads the operators |
|---|
| 400 | or functions in a way to achieve a kind of "promotion logic": |
|---|
| 401 | if one input is a <tt>ZZ</tt> and the other is a <tt>long</tt> |
|---|
| 402 | (or something that implicitly converts to a <tt>long</tt>, like |
|---|
| 403 | an <tt>int</tt>), the <tt>long</tt> input is effectively converted |
|---|
| 404 | to a <tt>ZZ</tt>. |
|---|
| 405 | Moreover, wherever possible, the implementation does this |
|---|
| 406 | as efficiently as possible, and usually avoids the creation |
|---|
| 407 | of a temporary <tt>ZZ</tt>. |
|---|
| 408 | |
|---|
| 409 | <p> |
|---|
| 410 | There are also procedural versions for all the basic arithmetic |
|---|
| 411 | operations: |
|---|
| 412 | <pre> |
|---|
| 413 | add, sub, negate, mul, sqr, div, rem, DivRem, |
|---|
| 414 | LeftShift, RightShift, |
|---|
| 415 | bit_and, bit_or, bit_xor |
|---|
| 416 | </pre> |
|---|
| 417 | |
|---|
| 418 | <p> |
|---|
| 419 | There are many other routines. |
|---|
| 420 | Here is a brief summary: |
|---|
| 421 | <ul> |
|---|
| 422 | <li> |
|---|
| 423 | <tt>GCD</tt> -- computes greatest common divisor of two integers |
|---|
| 424 | <li> |
|---|
| 425 | <tt>XGCD</tt> -- extended Euclidean algorithm |
|---|
| 426 | <li> |
|---|
| 427 | <tt>AddMod</tt>, <tt>SubMod</tt>, <tt>NegateMod</tt>, |
|---|
| 428 | <tt>MulMod</tt>, <tt>SqrMod</tt>, <tt>InvMod</tt>, |
|---|
| 429 | <tt>PowerMod</tt> -- routines for modular arithmetic, |
|---|
| 430 | including inversion and exponentiation |
|---|
| 431 | <li> |
|---|
| 432 | <tt>NumBits</tt> -- length of binary representation |
|---|
| 433 | <li> |
|---|
| 434 | <tt>bit</tt> -- extract a bit |
|---|
| 435 | <li> |
|---|
| 436 | <tt>ZZFromBytes</tt>, <tt>BytesFromZZ</tt> -- |
|---|
| 437 | convert between octet strings and <tt>ZZ</tt>s |
|---|
| 438 | <li> |
|---|
| 439 | <tt>RandomBnd</tt>, <tt>RandomBits</tt>, <tt>RandomLen</tt> -- |
|---|
| 440 | routines for generating pseudo-random numbers |
|---|
| 441 | <li> |
|---|
| 442 | <tt>GenPrime</tt>, <tt>ProbPrime</tt> -- routines for generating primes |
|---|
| 443 | and testing primality |
|---|
| 444 | <li> |
|---|
| 445 | <tt>power</tt> -- (non-modular) exponentiation |
|---|
| 446 | <li> |
|---|
| 447 | <tt>SqrRoot</tt> -- integer part of square root |
|---|
| 448 | <li> |
|---|
| 449 | <tt>Jacobi</tt>, <tt>SqrRootMod</tt> -- Jacobi symbol and modular |
|---|
| 450 | square root |
|---|
| 451 | </ul> |
|---|
| 452 | |
|---|
| 453 | <p> |
|---|
| 454 | Most of these functions also have pure <tt>long</tt> versions as |
|---|
| 455 | well, and as usual, there are both functional and procedural |
|---|
| 456 | variants. |
|---|
| 457 | |
|---|
| 458 | <p> |
|---|
| 459 | There are other functions as well. |
|---|
| 460 | See <a href="ZZ.txt"><tt>ZZ.txt</tt></a> for complete details. |
|---|
| 461 | Also see <a href="tools.txt"><tt>tools.txt</tt></a> for some basic |
|---|
| 462 | services provided by NTL. |
|---|
| 463 | |
|---|
| 464 | |
|---|
| 465 | <p> |
|---|
| 466 | <center> |
|---|
| 467 | <img src="arrow1.gif" alt="[Previous]" align=bottom> |
|---|
| 468 | <a href="tour-examples.html"><img src="arrow2.gif" alt="[Up]" align=bottom></a> |
|---|
| 469 | <a href="tour-ex2.html"> <img src="arrow3.gif" alt="[Next]" align=bottom></a> |
|---|
| 470 | </center> |
|---|
| 471 | |
|---|
| 472 | </body> |
|---|
| 473 | </html> |
|---|
