0
# Number Theory
1
2
gmpy2 provides extensive number theory functionality including GCD/LCM operations, modular arithmetic, primality testing, integer sequences, and advanced primality tests suitable for cryptographic applications.
3
4
## Capabilities
5
6
### GCD and LCM Operations
7
8
Greatest Common Divisor and Least Common Multiple operations with extended algorithms.
9
10
```python { .api }
11
def gcd(*args):
12
"""
13
Greatest Common Divisor of one or more integers.
14
15
Args:
16
*args: Integer values
17
18
Returns:
19
mpz: GCD of all arguments
20
21
Note:
22
For multiple arguments, computes GCD iteratively
23
"""
24
25
def gcdext(a, b):
26
"""
27
Extended Euclidean Algorithm.
28
29
Args:
30
a, b: Integer values
31
32
Returns:
33
tuple: (gcd, s, t) where gcd = s*a + t*b
34
35
Note:
36
Useful for computing modular inverses
37
"""
38
39
def lcm(*args):
40
"""
41
Least Common Multiple of one or more integers.
42
43
Args:
44
*args: Integer values
45
46
Returns:
47
mpz: LCM of all arguments
48
"""
49
```
50
51
### Modular Arithmetic
52
53
Modular inverse and division operations.
54
55
```python { .api }
56
def invert(x, y):
57
"""
58
Modular inverse of x modulo y.
59
60
Args:
61
x: Value to invert
62
y: Modulus
63
64
Returns:
65
mpz: z such that (x * z) ≡ 1 (mod y)
66
67
Raises:
68
ValueError: If gcd(x, y) != 1 (no inverse exists)
69
"""
70
71
def divm(a, b, m):
72
"""
73
Modular division: compute (a / b) mod m.
74
75
Args:
76
a: Dividend
77
b: Divisor
78
m: Modulus
79
80
Returns:
81
mpz: Result of (a * invert(b, m)) mod m
82
"""
83
```
84
85
### Primality Testing
86
87
Comprehensive primality testing with multiple algorithms.
88
89
```python { .api }
90
def is_prime(x, n=25):
91
"""
92
Miller-Rabin primality test.
93
94
Args:
95
x: Integer to test
96
n: Number of test rounds (default 25)
97
98
Returns:
99
bool: True if probably prime, False if composite
100
101
Note:
102
Probability of error is at most 4^(-n)
103
"""
104
105
def is_probab_prime(x, n=25):
106
"""
107
Probabilistic primality test with detailed result.
108
109
Args:
110
x: Integer to test
111
n: Number of test rounds
112
113
Returns:
114
int: 0 if composite, 1 if probably prime, 2 if definitely prime
115
"""
116
117
def next_prime(x):
118
"""
119
Find next prime number greater than x.
120
121
Args:
122
x: Starting value
123
124
Returns:
125
mpz: Next prime after x
126
"""
127
128
def prev_prime(x):
129
"""
130
Find previous prime number less than x (GMP 6.3.0+).
131
132
Args:
133
x: Starting value
134
135
Returns:
136
mpz: Previous prime before x
137
"""
138
```
139
140
### Advanced Primality Tests
141
142
Specialized primality tests for cryptographic and research applications.
143
144
```python { .api }
145
def is_bpsw_prp(x):
146
"""
147
Baillie-PSW probable prime test.
148
149
Args:
150
x: Integer to test
151
152
Returns:
153
bool: True if passes BPSW test
154
155
Note:
156
Very strong test with no known counterexamples < 2^64
157
"""
158
159
def is_strong_bpsw_prp(x):
160
"""
161
Strong Baillie-PSW probable prime test.
162
163
Args:
164
x: Integer to test
165
166
Returns:
167
bool: True if passes strong BPSW test
168
"""
169
170
def is_fermat_prp(x, a):
171
"""
172
Fermat probable prime test to base a.
173
174
Args:
175
x: Integer to test
176
a: Test base
177
178
Returns:
179
bool: True if x passes Fermat test to base a
180
"""
181
182
def is_euler_prp(x, a):
183
"""
184
Euler probable prime test to base a.
185
186
Args:
187
x: Integer to test
188
a: Test base
189
190
Returns:
191
bool: True if x passes Euler test to base a
192
"""
193
194
def is_strong_prp(x, a):
195
"""
196
Strong probable prime test (Miller-Rabin) to base a.
197
198
Args:
199
x: Integer to test
200
a: Test base
201
202
Returns:
203
bool: True if x passes strong test to base a
204
"""
205
206
def is_lucas_prp(x):
207
"""
208
Lucas probable prime test.
209
210
Args:
211
x: Integer to test
212
213
Returns:
214
bool: True if x passes Lucas test
215
"""
216
217
def is_strong_lucas_prp(x):
218
"""
219
Strong Lucas probable prime test.
220
221
Args:
222
x: Integer to test
223
224
Returns:
225
bool: True if x passes strong Lucas test
226
"""
227
228
def is_extra_strong_lucas_prp(x):
229
"""
230
Extra strong Lucas probable prime test.
231
232
Args:
233
x: Integer to test
234
235
Returns:
236
bool: True if x passes extra strong Lucas test
237
"""
238
239
def is_selfridge_prp(x):
240
"""
241
Selfridge probable prime test.
242
243
Args:
244
x: Integer to test
245
246
Returns:
247
bool: True if x passes Selfridge test
248
"""
249
250
def is_strong_selfridge_prp(x):
251
"""
252
Strong Selfridge probable prime test.
253
254
Args:
255
x: Integer to test
256
257
Returns:
258
bool: True if x passes strong Selfridge test
259
"""
260
261
def is_fibonacci_prp(x):
262
"""
263
Fibonacci probable prime test.
264
265
Args:
266
x: Integer to test
267
268
Returns:
269
bool: True if x passes Fibonacci test
270
"""
271
```
272
273
### Legendre, Jacobi, and Kronecker Symbols
274
275
Number theory symbols for quadratic residue testing.
276
277
```python { .api }
278
def legendre(a, p):
279
"""
280
Legendre symbol (a/p).
281
282
Args:
283
a: Integer
284
p: Odd prime
285
286
Returns:
287
int: -1, 0, or 1
288
289
Note:
290
Returns 1 if a is quadratic residue mod p, -1 if not, 0 if a ≡ 0 (mod p)
291
"""
292
293
def jacobi(a, b):
294
"""
295
Jacobi symbol (a/b).
296
297
Args:
298
a: Integer
299
b: Positive odd integer
300
301
Returns:
302
int: -1, 0, or 1
303
304
Note:
305
Generalization of Legendre symbol
306
"""
307
308
def kronecker(a, b):
309
"""
310
Kronecker symbol (a/b).
311
312
Args:
313
a: Integer
314
b: Integer
315
316
Returns:
317
int: -1, 0, or 1
318
319
Note:
320
Further generalization allowing even b
321
"""
322
```
323
324
### Combinatorial Functions
325
326
Factorial and binomial coefficient calculations.
327
328
```python { .api }
329
def fac(n):
330
"""
331
Factorial of n.
332
333
Args:
334
n: Non-negative integer
335
336
Returns:
337
mpz: n!
338
"""
339
340
def double_fac(n):
341
"""
342
Double factorial of n (n!!).
343
344
Args:
345
n: Non-negative integer
346
347
Returns:
348
mpz: n!! = n * (n-2) * (n-4) * ... * 1 or 2
349
"""
350
351
def multi_fac(n, m):
352
"""
353
Multi-factorial of n with step m.
354
355
Args:
356
n: Non-negative integer
357
m: Step size
358
359
Returns:
360
mpz: n * (n-m) * (n-2m) * ...
361
"""
362
363
def primorial(n):
364
"""
365
Primorial of n (product of primes <= n).
366
367
Args:
368
n: Non-negative integer
369
370
Returns:
371
mpz: Product of all primes <= n
372
"""
373
374
def bincoef(n, k):
375
"""
376
Binomial coefficient C(n, k).
377
378
Args:
379
n: Integer
380
k: Non-negative integer
381
382
Returns:
383
mpz: n! / (k! * (n-k)!)
384
385
Note:
386
Also available as comb(n, k)
387
"""
388
389
def comb(n, k):
390
"""
391
Binomial coefficient C(n, k) (alias for bincoef).
392
393
Args:
394
n: Integer
395
k: Non-negative integer
396
397
Returns:
398
mpz: n! / (k! * (n-k)!)
399
"""
400
```
401
402
### Integer Sequences
403
404
Fibonacci, Lucas, and related sequence generators.
405
406
```python { .api }
407
def fib(n):
408
"""
409
nth Fibonacci number.
410
411
Args:
412
n: Non-negative integer
413
414
Returns:
415
mpz: F(n) where F(0)=0, F(1)=1, F(n)=F(n-1)+F(n-2)
416
"""
417
418
def fib2(n):
419
"""
420
nth Fibonacci number and its predecessor.
421
422
Args:
423
n: Non-negative integer
424
425
Returns:
426
tuple: (F(n), F(n-1))
427
"""
428
429
def lucas(n):
430
"""
431
nth Lucas number.
432
433
Args:
434
n: Non-negative integer
435
436
Returns:
437
mpz: L(n) where L(0)=2, L(1)=1, L(n)=L(n-1)+L(n-2)
438
"""
439
440
def lucas2(n):
441
"""
442
nth Lucas number and its predecessor.
443
444
Args:
445
n: Non-negative integer
446
447
Returns:
448
tuple: (L(n), L(n-1))
449
"""
450
451
def lucasu(p, q, n):
452
"""
453
Lucas U sequence U_n(p, q).
454
455
Args:
456
p, q: Parameters
457
n: Index
458
459
Returns:
460
mpz: U_n(p, q)
461
"""
462
463
def lucasu_mod(p, q, n, m):
464
"""
465
Lucas U sequence U_n(p, q) mod m.
466
467
Args:
468
p, q: Parameters
469
n: Index
470
m: Modulus
471
472
Returns:
473
mpz: U_n(p, q) mod m
474
"""
475
476
def lucasv(p, q, n):
477
"""
478
Lucas V sequence V_n(p, q).
479
480
Args:
481
p, q: Parameters
482
n: Index
483
484
Returns:
485
mpz: V_n(p, q)
486
"""
487
488
def lucasv_mod(p, q, n, m):
489
"""
490
Lucas V sequence V_n(p, q) mod m.
491
492
Args:
493
p, q: Parameters
494
n: Index
495
m: Modulus
496
497
Returns:
498
mpz: V_n(p, q) mod m
499
"""
500
```
501
502
### Factor Manipulation
503
504
Functions for working with prime factors.
505
506
```python { .api }
507
def remove(x, f):
508
"""
509
Remove all occurrences of factor f from x.
510
511
Args:
512
x: Integer to factor
513
f: Factor to remove
514
515
Returns:
516
tuple: (result, count) where result = x / f^count
517
518
Note:
519
Efficiently removes highest power of f that divides x
520
"""
521
```
522
523
### Integer Property Tests
524
525
Tests for various integer properties.
526
527
```python { .api }
528
def is_even(x):
529
"""
530
Test if integer is even.
531
532
Args:
533
x: Integer value
534
535
Returns:
536
bool: True if x is even
537
"""
538
539
def is_odd(x):
540
"""
541
Test if integer is odd.
542
543
Args:
544
x: Integer value
545
546
Returns:
547
bool: True if x is odd
548
"""
549
550
def is_square(x):
551
"""
552
Test if integer is a perfect square.
553
554
Args:
555
x: Integer value
556
557
Returns:
558
bool: True if x = k^2 for some integer k
559
"""
560
561
def is_power(x):
562
"""
563
Test if integer is a perfect power.
564
565
Args:
566
x: Integer value
567
568
Returns:
569
bool: True if x = k^n for some integers k, n with n > 1
570
"""
571
572
def is_congruent(x, y, m):
573
"""
574
Test if x ≡ y (mod m).
575
576
Args:
577
x, y: Integer values
578
m: Modulus
579
580
Returns:
581
bool: True if x and y are congruent modulo m
582
"""
583
584
def is_divisible(x, y):
585
"""
586
Test if x is divisible by y.
587
588
Args:
589
x, y: Integer values
590
591
Returns:
592
bool: True if y divides x evenly
593
"""
594
```
595
596
## Usage Examples
597
598
### Basic Number Theory
599
600
```python
601
import gmpy2
602
603
# GCD and LCM
604
a = gmpy2.mpz(48)
605
b = gmpy2.mpz(18)
606
print(f"GCD: {gmpy2.gcd(a, b)}") # 6
607
print(f"LCM: {gmpy2.lcm(a, b)}") # 144
608
609
# Extended GCD
610
gcd, s, t = gmpy2.gcdext(a, b)
611
print(f"GCD: {gcd}, s: {s}, t: {t}")
612
print(f"Verification: {s * a + t * b}") # Should equal gcd
613
```
614
615
### Primality Testing
616
617
```python
618
import gmpy2
619
620
# Test various numbers for primality
621
candidates = [97, 98, 99, 101, 1009, 1013]
622
623
for n in candidates:
624
is_prime = gmpy2.is_prime(n)
625
prob_prime = gmpy2.is_probab_prime(n)
626
print(f"{n}: is_prime={is_prime}, probab_prime={prob_prime}")
627
628
# Find next prime
629
start = gmpy2.mpz(1000)
630
next_p = gmpy2.next_prime(start)
631
print(f"Next prime after {start}: {next_p}")
632
```
633
634
### Advanced Primality Testing
635
636
```python
637
import gmpy2
638
639
# Large number for testing
640
large_num = gmpy2.mpz("1000000000000000000000000000057")
641
642
# Various primality tests
643
tests = [
644
("BPSW", gmpy2.is_bpsw_prp),
645
("Strong BPSW", gmpy2.is_strong_bpsw_prp),
646
]
647
648
for name, test_func in tests:
649
result = test_func(large_num)
650
print(f"{name} test: {result}")
651
652
# Base-specific tests
653
base = 2
654
fermat_result = gmpy2.is_fermat_prp(large_num, base)
655
strong_result = gmpy2.is_strong_prp(large_num, base)
656
print(f"Fermat base {base}: {fermat_result}")
657
print(f"Strong base {base}: {strong_result}")
658
```
659
660
### Modular Arithmetic
661
662
```python
663
import gmpy2
664
665
# Modular inverse example
666
a = gmpy2.mpz(3)
667
m = gmpy2.mpz(11)
668
inv = gmpy2.invert(a, m)
669
print(f"Inverse of {a} modulo {m}: {inv}")
670
print(f"Verification: ({a} * {inv}) mod {m} = {(a * inv) % m}")
671
672
# Modular division
673
dividend = gmpy2.mpz(15)
674
divisor = gmpy2.mpz(7)
675
modulus = gmpy2.mpz(11)
676
result = gmpy2.divm(dividend, divisor, modulus)
677
print(f"({dividend} / {divisor}) mod {modulus} = {result}")
678
```
679
680
### Combinatorial Functions
681
682
```python
683
import gmpy2
684
685
# Factorials
686
n = 20
687
print(f"{n}! = {gmpy2.fac(n)}")
688
print(f"{n}!! = {gmpy2.double_fac(n)}")
689
690
# Binomial coefficients
691
n, k = 10, 3
692
coeff = gmpy2.bincoef(n, k)
693
print(f"C({n}, {k}) = {coeff}")
694
695
# Fibonacci numbers
696
fib_n = gmpy2.fib(50)
697
fib_n, fib_n_minus_1 = gmpy2.fib2(50)
698
print(f"F(50) = {fib_n}")
699
print(f"F(49) = {fib_n_minus_1}")
700
```
701
702
### Integer Sequences
703
704
```python
705
import gmpy2
706
707
# Generate Fibonacci sequence
708
print("First 10 Fibonacci numbers:")
709
for i in range(10):
710
print(f"F({i}) = {gmpy2.fib(i)}")
711
712
# Lucas numbers
713
print("\nFirst 10 Lucas numbers:")
714
for i in range(10):
715
print(f"L({i}) = {gmpy2.lucas(i)}")
716
717
# Custom Lucas sequences
718
p, q = 1, -1 # Parameters for Fibonacci-like sequence
719
for i in range(5):
720
u_val = gmpy2.lucasu(p, q, i)
721
v_val = gmpy2.lucasv(p, q, i)
722
print(f"U_{i}({p},{q}) = {u_val}, V_{i}({p},{q}) = {v_val}")
723
```