class sympy.ntheory.generate.Sieve

In [ ]:

```
>>> from sympy import sieve
>>> from array import array # this line and next for doctest only
>>> sieve._list = array('l', [2, 3, 5, 7, 11, 13])
```

In [ ]:

```
>>> 25 in sieve
False
>>> sieve._list
array('l', [2, 3, 5, 7, 11, 13, 17, 19, 23])
```

extend(n)

Grow the sieve to cover all primes <= n (a real number).

In [ ]:

```
>>> from sympy import sieve
>>> from array import array # this line and next for doctest only
>>> sieve._list = array('l', [2, 3, 5, 7, 11, 13])
```

In [ ]:

```
>>> sieve.extend(30)
>>> sieve[10] == 29
True
```

extend_to_no(i)

Extend to include the ith prime number.

i must be an integer.

In [ ]:

```
>>> from sympy import sieve
>>> from array import array # this line and next for doctest only
>>> sieve._list = array('l', [2, 3, 5, 7, 11, 13])
```

In [ ]:

```
>>> sieve.extend_to_no(9)
>>> sieve._list
array('l', [2, 3, 5, 7, 11, 13, 17, 19, 23])
```

primerange(a, b)

Generate all prime numbers in the range [a, b).

In [ ]:

```
>>> from sympy import sieve
>>> print([i for i in sieve.primerange(7, 18)])
[7, 11, 13, 17]
```

search(n)

Return the indices i, j of the primes that bound n.

If n is prime then i == j.

In [ ]:

```
>>> from sympy import sieve
>>> sieve.search(25)
(9, 10)
>>> sieve.search(23)
(9, 9)
```

sympy.ntheory.generate.prime(nth)

In [ ]:

```
>>> from sympy import prime
>>> prime(10)
29
>>> prime(1)
2
```

sympy.ntheory.generate.primepi(n)

In [ ]:

```
>>> from sympy import primepi
>>> primepi(25)
9
```

sympy.ntheory.generate.nextprime(n, ith=1)

Return the ith prime greater than n.

i must be an integer.

Potential primes are located at 6*j +/- 1. This property is used during searching.

In [ ]:

```
>>> from sympy import nextprime
>>> [(i, nextprime(i)) for i in range(10, 15)]
[(10, 11), (11, 13), (12, 13), (13, 17), (14, 17)]
>>> nextprime(2, ith=2) # the 2nd prime after 2
5
```

sympy.ntheory.generate.prevprime(n)

Return the largest prime smaller than n.

Potential primes are located at 6*j +/- 1. This property is used during searching.

In [ ]:

```
>>> from sympy import prevprime
>>> [(i, prevprime(i)) for i in range(10, 15)]
[(10, 7), (11, 7), (12, 11), (13, 11), (14, 13)]
```

nextprime : Return the ith prime greater than n primerange : Generates all primes in a given range

sympy.ntheory.generate.primerange(a, b)

Generate a list of all prime numbers in the range [a, b).

Some famous conjectures about the occurence of primes in a given range are [1]:

Twin primes: though often not, the following will give 2 primes

an infinite number of times:

primerange(6*n - 1, 6*n + 2)

Legendre's: the following always yields at least one prime

primerange(n**2, (n+1)**2+1)

Bertrand's (proven): there is always a prime in the range

primerange(n, 2*n)

Brocard's: there are at least four primes in the range

primerange(prime(n)**2, prime(n+1)**2)

In [ ]:

```
>>> from sympy import primerange, sieve
>>> print([i for i in primerange(1, 30)])
[2, 3, 5, 7, 11, 13, 17, 19, 23, 29]
```

In [ ]:

```
>>> list(sieve.primerange(1, 30))
[2, 3, 5, 7, 11, 13, 17, 19, 23, 29]
```

or extend the sieve to contain the requested range.

sympy.ntheory.generate.randprime(a, b)

Return a random prime number in the range [a, b).

Bertrand's postulate assures that randprime(a, 2*a) will always succeed for a > 1.

In [ ]:

```
>>> from sympy import randprime, isprime
>>> randprime(1, 30) #doctest: +SKIP
13
>>> isprime(randprime(1, 30))
True
```

primerange : Generate all primes in a given range

sympy.ntheory.generate.primorial(n, nth=True)

In [ ]:

```
>>> from sympy.ntheory.generate import primorial, randprime, primerange
>>> from sympy import factorint, Mul, primefactors, sqrt
>>> primorial(4) # the first 4 primes are 2, 3, 5, 7
210
>>> primorial(4, nth=False) # primes <= 4 are 2 and 3
6
>>> primorial(1)
2
>>> primorial(1, nth=False)
1
>>> primorial(sqrt(101), nth=False)
210
```

In this case, the number itself is a new prime:

In [ ]:

```
>>> factorint(primorial(4) + 1)
{211: 1}
```

In this case two new primes are the factors:

In [ ]:

```
>>> factorint(primorial(4) - 1)
{11: 1, 19: 1}
```

Here, some primes smaller and larger than the primes multiplied together are obtained:

In [ ]:

```
>>> p = list(primerange(10, 20))
>>> sorted(set(primefactors(Mul(*p) + 1)).difference(set(p)))
[2, 5, 31, 149]
```

primerange : Generate all primes in a given range

sympy.ntheory.generate.cycle_length(f, x0, nmax=None, values=False)

In [ ]:

```
>>> from sympy.ntheory.generate import cycle_length
```

This will yield successive values of i <-- func(i):

In [ ]:

```
>>> def iter(func, i):
... while 1:
... ii = func(i)
... yield ii
... i = ii
...
```

A function is defined:

In [ ]:

```
>>> func = lambda i: (i**2 + 1) % 51
```

and given a seed of 4 and the mu and lambda terms calculated:

In [ ]:

```
>>> next(cycle_length(func, 4))
(6, 2)
```

We can see what is meant by looking at the output:

In [ ]:

```
>>> n = cycle_length(func, 4, values=True)
>>> list(ni for ni in n)
[17, 35, 2, 5, 26, 14, 44, 50, 2, 5, 26, 14]
```

There are 6 repeating values after the first 2.

In [ ]:

```
>>> next(cycle_length(func, 4, nmax = 4))
(4, None)
>>> [ni for ni in cycle_length(func, 4, nmax = 4, values=True)]
[17, 35, 2, 5]
```

Code modified from:

sympy.ntheory.factor_.smoothness(n)

Return the B-smooth and B-power smooth values of n.

In [ ]:

```
>>> from sympy.ntheory.factor_ import smoothness
>>> smoothness(2**7*3**2)
(3, 128)
>>> smoothness(2**4*13)
(13, 16)
>>> smoothness(2)
(2, 2)
```

factorint, smoothness_p

sympy.ntheory.factor_.smoothness_p(n, m=-1, power=0, visual=None)

Return a list of [m, (p, (M, sm(p + m), psm(p + m)))...] where:

p**M is the base-p divisor of n

sm(p + m) is the smoothness of p + m (m = -1 by default)

psm(p + m) is the power smoothness of p + m

The list is sorted according to smoothness (default) or by power smoothness if power=1.

In [ ]:

```
>>> from sympy.ntheory.factor_ import smoothness_p, factorint
>>> smoothness_p(10431, m=1)
(1, [(3, (2, 2, 4)), (19, (1, 5, 5)), (61, (1, 31, 31))])
>>> smoothness_p(10431)
(-1, [(3, (2, 2, 2)), (19, (1, 3, 9)), (61, (1, 5, 5))])
>>> smoothness_p(10431, power=1)
(-1, [(3, (2, 2, 2)), (61, (1, 5, 5)), (19, (1, 3, 9))])
```

If visual=True then an annotated string will be returned:

In [ ]:

```
>>> print(smoothness_p(21477639576571, visual=1))
p**i=4410317**1 has p-1 B=1787, B-pow=1787
p**i=4869863**1 has p-1 B=2434931, B-pow=2434931
```

This string can also be generated directly from a factorization dictionary and vice versa:

In [ ]:

```
>>> factorint(17*9)
{3: 2, 17: 1}
>>> smoothness_p(_)
'p**i=3**2 has p-1 B=2, B-pow=2\np**i=17**1 has p-1 B=2, B-pow=16'
>>> smoothness_p(_)
{3: 2, 17: 1}
```

The table of the output logic is:

In [ ]:

Visual

Input

True

False

other

dict

str

tuple

str

str

str

tuple

dict

tuple

str

tuple

str

n

str

tuple

tuple

mul

str

tuple

tuple

factorint, smoothness

sympy.ntheory.factor_.trailing(n)

In [ ]:

```
>>> from sympy import trailing
>>> trailing(128)
7
>>> trailing(63)
0
```

sympy.ntheory.factor_.multiplicity(p, n)

Find the greatest integer m such that p**m divides n.

In [ ]:

```
>>> from sympy.ntheory import multiplicity
>>> from sympy.core.numbers import Rational as R
>>> [multiplicity(5, n) for n in [8, 5, 25, 125, 250]]
[0, 1, 2, 3, 3]
>>> multiplicity(3, R(1, 9))
-2
```

sympy.ntheory.factor_.perfect_power(n, candidates=None, big=True, factor=True)

Return "(b, e)" such that "n" == "b**e" if "n" is a perfect power; otherwise return "False".

In [ ]:

```
>>> from sympy import perfect_power
>>> perfect_power(16)
(2, 4)
>>> perfect_power(16, big = False)
(4, 2)
```

sympy.ntheory.factor_.pollard_rho(n, s=2, a=1, retries=5, seed=1234, max_steps=None, F=None)

In [ ]:

```
>>> from sympy.ntheory.generate import cycle_length
>>> n = 16843009
>>> F = lambda x:(2048*pow(x, 2, n) + 32767) % n
>>> for s in range(5):
... print('loop length = %4i; leader length = %3i' % next(cycle_length(F, s)))
...
loop length = 2489; leader length = 42
loop length = 78; leader length = 120
loop length = 1482; leader length = 99
loop length = 1482; leader length = 285
loop length = 1482; leader length = 100
```

In [ ]:

```
>>> x=2
>>> for i in range(9):
... x=(x**2+12)%17
... print(x)
...
16
13
11
14
4
11
14
4
11
>>> next(cycle_length(lambda x: (x**2+12)%17, 2))
(3, 2)
>>> list(cycle_length(lambda x: (x**2+12)%17, 2, values=True))
[16, 13, 11, 14, 4]
```

In [ ]:

```
>>> from sympy import pollard_rho
>>> n=16843009
>>> F=lambda x:(2048*pow(x,2,n) + 32767) % n
>>> pollard_rho(n, F=F)
257
```

Use the default setting with a bad value of "a" and no retries:

In [ ]:

```
>>> pollard_rho(n, a=n-2, retries=0)
```

If retries is > 0 then perhaps the problem will correct itself when new values are generated for a:

In [ ]:

```
>>> pollard_rho(n, a=n-2, retries=1)
257
```

sympy.ntheory.factor_.pollard_pm1(n, B=10, a=2, retries=0, seed=1234)

Note: the value of M is lcm(1..B) = reduce(ilcm, range(2, B + 1)).

**M - 1, n) = N then a**M % q**r is 1
for every prime power divisor of N. But consider the following:

In [ ]:

```
>>> from sympy.ntheory.factor_ import smoothness_p, pollard_pm1
>>> n=257*1009
>>> smoothness_p(n)
(-1, [(257, (1, 2, 256)), (1009, (1, 7, 16))])
```

So we should (and can) find a root with B=16:

In [ ]:

```
>>> pollard_pm1(n, B=16, a=3)
1009
```

If we attempt to increase B to 256 we find that it doesn't work:

In [ ]:

```
>>> pollard_pm1(n, B=256)
>>>
```

But if the value of "a" is changed we find that only multiples of 257 work, e.g.:

In [ ]:

```
>>> pollard_pm1(n, B=256, a=257)
1009
```

In [ ]:

```
>>> from sympy.core.numbers import ilcm, igcd
>>> from sympy import factorint, Pow
>>> M = 1
>>> for i in range(2, 256):
... M = ilcm(M, i)
...
>>> set([igcd(pow(a, M, n) - 1, n) for a in range(2, 256) if
... igcd(pow(a, M, n) - 1, n) != n])
set([1009])
```

But does aM % d for every divisor of n give 1?

In [ ]:

```
>>> aM = pow(255, M, n)
>>> [(d, aM%Pow(*d.args)) for d in factorint(n, visual=True).args]
[(257**1, 1), (1009**1, 1)]
```

the power smoothness of the p - 1 value next to the root does not exceed B

a**M % p != 1 for any of the divisors of n.

By trying more than one "a" it is possible that one of them will yield a factor.

With the default smoothness bound, this number can't be cracked:

In [ ]:

```
>>> from sympy.ntheory import pollard_pm1, primefactors
>>> pollard_pm1(21477639576571)
```

Increasing the smoothness bound helps:

In [ ]:

```
>>> pollard_pm1(21477639576571, B=2000)
4410317
```

Looking at the smoothness of the factors of this number we find:

In [ ]:

```
>>> from sympy.utilities import flatten
>>> from sympy.ntheory.factor_ import smoothness_p, factorint
>>> print(smoothness_p(21477639576571, visual=1))
p**i=4410317**1 has p-1 B=1787, B-pow=1787
p**i=4869863**1 has p-1 B=2434931, B-pow=2434931
```

In [ ]:

```
>>> factorint(4410317 - 1)
{2: 2, 617: 1, 1787: 1}
>>> factorint(4869863-1)
{2: 1, 2434931: 1}
```

Note that until B reaches the B-pow value of 1787, the number is not cracked;

In [ ]:

```
>>> pollard_pm1(21477639576571, B=1786)
>>> pollard_pm1(21477639576571, B=1787)
4410317
```

**15 to 15**15 + 10**4 are 10**6 power smooth so a B of 10**6
will fail 85% of the time in that range. From 10**8 to 10**8 + 10**3 the
percentages are nearly reversed...but in that range the simple trial
division is quite fast.

In [ ]:

```
>>> from sympy.ntheory import factorint
>>> factorint(2000) # 2000 = (2**4) * (5**3)
{2: 4, 5: 3}
>>> factorint(65537) # This number is prime
{65537: 1}
```

For input less than 2, factorint behaves as follows:

"factorint(1)" returns the empty factorization, "{}"

"factorint(0)" returns "{0:1}"

"factorint(-n)" adds "-1:1" to the factors and then factors "n"

Partial Factorization:

In [ ]:

```
>>> from sympy.ntheory import isprime
>>> from sympy.core.compatibility import long
>>> a = 1407633717262338957430697921446883
>>> f = factorint(a, limit=10000)
>>> f == {991: 1, long(202916782076162456022877024859): 1, 7: 1}
True
>>> isprime(max(f))
False
```

This number has a small factor and a residual perfect power whose base is greater than the limit:

In [ ]:

```
>>> factorint(3*101**7, limit=5)
{3: 1, 101: 7}
```

Visual Factorization:

In [ ]:

```
>>> from sympy import pprint
>>> pprint(factorint(4200, visual=True))
3 1 2 1
2 *3 *5 *7
```

You can easily switch between the two forms by sending them back to factorint:

In [ ]:

```
>>> from sympy import Mul, Pow
>>> regular = factorint(1764); regular
{2: 2, 3: 2, 7: 2}
>>> pprint(factorint(regular))
2 2 2
2 *3 *7
```

In [ ]:

```
>>> visual = factorint(1764, visual=True); pprint(visual)
2 2 2
2 *3 *7
>>> print(factorint(visual))
{2: 2, 3: 2, 7: 2}
```

In [ ]:

```
>>> factorint(factorint({4: 2, 12: 3})) # twice to toggle to dict form
{2: 10, 3: 3}
>>> factorint(Mul(4, 12, evaluate=False))
{2: 4, 3: 1}
```

The table of the output logic is:

Input

True

False

other

dict

mul

dict

mul

n

mul

dict

dict

mul

mul

dict

dict

Algorithm:

In [ ]:

```
>>> factors = factorint(12345678910111213141516)
>>> for base, exp in sorted(factors.items()):
... print('%s %s' % (base, exp))
...
2 2
2507191691 1
1231026625769 1
```

Any of these methods can optionally be disabled with the following boolean parameters:

"use_trial": Toggle use of trial division

"use_rho": Toggle use of Pollard's rho method

"use_pm1": Toggle use of Pollard's p-1 method

If "verbose" is set to "True", detailed progress is printed.

smoothness, smoothness_p, divisors

sympy.ntheory.factor_.primefactors(n, limit=None, verbose=False)

In [ ]:

```
>>> from sympy.ntheory import primefactors, factorint, isprime
>>> primefactors(6)
[2, 3]
>>> primefactors(-5)
[5]
```

In [ ]:

```
>>> sorted(factorint(123456).items())
[(2, 6), (3, 1), (643, 1)]
>>> primefactors(123456)
[2, 3, 643]
```

In [ ]:

```
>>> sorted(factorint(10000000001, limit=200).items())
[(101, 1), (99009901, 1)]
>>> isprime(99009901)
False
>>> primefactors(10000000001, limit=300)
[101]
```

divisors

sympy.ntheory.factor_.divisors(n, generator=False)

In [ ]:

```
>>> from sympy import divisors, divisor_count
>>> divisors(24)
[1, 2, 3, 4, 6, 8, 12, 24]
>>> divisor_count(24)
8
```

In [ ]:

```
>>> list(divisors(120, generator=True))
[1, 2, 4, 8, 3, 6, 12, 24, 5, 10, 20, 40, 15, 30, 60, 120]
```

primefactors, factorint, divisor_count

sympy.ntheory.factor_.divisor_count(n, modulus=1)

In [ ]:

```
>>> from sympy import divisor_count
>>> divisor_count(6)
4
```

factorint, divisors, totient

sympy.ntheory.factor_.totient(*args, **kw_args)

Calculate the Euler totient function phi(n)

In [ ]:

```
>>> from sympy.ntheory import totient
>>> totient(1)
1
>>> totient(25)
20
```

divisor_count

sympy.ntheory.modular.symmetric_residue(a, m)

Return the residual mod m such that it is within half of the modulus.

In [ ]:

```
>>> from sympy.ntheory.modular import symmetric_residue
>>> symmetric_residue(1, 6)
1
>>> symmetric_residue(4, 6)
-2
```

sympy.ntheory.modular.crt(m, v, symmetric=False, check=True)

Chinese Remainder Theorem.

The keyword "check" can be set to False if it is known that the moduli are coprime.

In [ ]:

```
>>> from sympy.ntheory.modular import crt, solve_congruence
>>> crt([99, 97, 95], [49, 76, 65])
(639985, 912285)
```

This is the correct result because:

In [ ]:

```
>>> [639985 % m for m in [99, 97, 95]]
[49, 76, 65]
```

If the moduli are not co-prime, you may receive an incorrect result if you use "check=False":

In [ ]:

```
>>> crt([12, 6, 17], [3, 4, 2], check=False)
(954, 1224)
>>> [954 % m for m in [12, 6, 17]]
[6, 0, 2]
>>> crt([12, 6, 17], [3, 4, 2]) is None
True
>>> crt([3, 6], [2, 5])
(5, 6)
```

solve_congruence sympy.polys.galoistools.gf_crt : low level crt routine used by this routine

sympy.ntheory.modular.crt1(m)

First part of Chinese Remainder Theorem, for multiple application.

In [ ]:

```
>>> from sympy.ntheory.modular import crt1
>>> crt1([18, 42, 6])
(4536, [252, 108, 756], [0, 2, 0])
```

sympy.ntheory.modular.crt2(m, v, mm, e, s, symmetric=False)

Second part of Chinese Remainder Theorem, for multiple application.

In [ ]:

```
>>> from sympy.ntheory.modular import crt1, crt2
>>> mm, e, s = crt1([18, 42, 6])
>>> crt2([18, 42, 6], [0, 0, 0], mm, e, s)
(0, 4536)
```

sympy.ntheory.modular.solve_congruence(*remainder_modulus_pairs, **hint)

In [ ]:

```
>>> from sympy.ntheory.modular import solve_congruence
```

What number is 2 mod 3, 3 mod 5 and 2 mod 7?

In [ ]:

```
>>> solve_congruence((2, 3), (3, 5), (2, 7))
(23, 105)
>>> [23 % m for m in [3, 5, 7]]
[2, 3, 2]
```

In [ ]:

```
>>> solve_congruence(*zip((2, 3, 2), (3, 5, 7)))
(23, 105)
```

The moduli need not be co-prime; in this case there may or may not be a solution:

In [ ]:

```
>>> solve_congruence((2, 3), (4, 6)) is None
True
```

In [ ]:

```
>>> solve_congruence((2, 3), (5, 6))
(5, 6)
```

The symmetric flag will make the result be within 1/2 of the modulus:

In [ ]:

```
>>> solve_congruence((2, 3), (5, 6), symmetric=True)
(-1, 6)
```

crt : high level routine implementing the Chinese Remainder Theorem

sympy.ntheory.multinomial.binomial_coefficients(n)

Return a dictionary containing pairs where are binomial coefficients and .

In [ ]:

```
>>> from sympy.ntheory import binomial_coefficients
>>> binomial_coefficients(9)
{(0, 9): 1, (1, 8): 9, (2, 7): 36, (3, 6): 84,
(4, 5): 126, (5, 4): 126, (6, 3): 84, (7, 2): 36, (8, 1): 9, (9, 0): 1}
```

binomial_coefficients_list, multinomial_coefficients

sympy.ntheory.multinomial.binomial_coefficients_list(n)

Return a list of binomial coefficients as rows of the Pascal's triangle.

In [ ]:

```
>>> from sympy.ntheory import binomial_coefficients_list
>>> binomial_coefficients_list(9)
[1, 9, 36, 84, 126, 126, 84, 36, 9, 1]
```

binomial_coefficients, multinomial_coefficients

sympy.ntheory.multinomial.multinomial_coefficients(m, n)

For example:

In [ ]:

```
>>> from sympy.ntheory import multinomial_coefficients
>>> multinomial_coefficients(2, 5) # indirect doctest
{(0, 5): 1, (1, 4): 5, (2, 3): 10, (3, 2): 10, (4, 1): 5, (5, 0): 1}
```

The algorithm is based on the following result:

Code contributed to Sage by Yann Laigle-Chapuy, copied with permission of the author.

binomial_coefficients_list, binomial_coefficients

sympy.ntheory.multinomial.multinomial_coefficients_iterator(m, n, _tuple=

multinomial coefficient iterator

*m* large with respect to *n* by taking
advantage of the fact that when the monomial tuples *t* are stripped of
zeros, their coefficient is the same as that of the monomial tuples from
"multinomial_coefficients(n, n)". Therefore, the latter coefficients are
precomputed to save memory and time.

In [ ]:

```
>>> from sympy.ntheory.multinomial import multinomial_coefficients
>>> m53, m33 = multinomial_coefficients(5,3), multinomial_coefficients(3,3)
>>> m53[(0,0,0,1,2)] == m53[(0,0,1,0,2)] == m53[(1,0,2,0,0)] == m33[(0,1,2)]
True
```

In [ ]:

```
>>> from sympy.ntheory.multinomial import multinomial_coefficients_iterator
>>> it = multinomial_coefficients_iterator(20,3)
>>> next(it)
((3, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0), 1)
```

sympy.ntheory.partitions_.npartitions(n, verbose=False)

The correctness of this implementation has been tested for 10**n up to n = 8.

In [ ]:

```
>>> from sympy.ntheory import npartitions
>>> npartitions(25)
1958
```

sympy.ntheory.primetest.mr(n, bases)

Perform a Miller-Rabin strong pseudoprime test on n using a given list of bases/witnesses.

In [ ]:

```
>>> from sympy.ntheory.primetest import mr
>>> mr(1373651, [2, 3])
False
>>> mr(479001599, [31, 73])
True
```

sympy.ntheory.primetest.isprime(n)

Negative primes (e.g. -2) are not considered prime.

In [ ]:

```
>>> from sympy.ntheory import isprime
>>> isprime(13)
True
>>> isprime(15)
False
```

sympy.ntheory.residue_ntheory.n_order(a, n)

Returns the order of "a" modulo "n".

In [ ]:

```
>>> from sympy.ntheory import n_order
>>> n_order(3, 7)
6
>>> n_order(4, 7)
3
```

sympy.ntheory.residue_ntheory.is_primitive_root(a, p)

Returns True if "a" is a primitive root of "p"

a**totient(p) cong 1 mod(p)

In [ ]:

```
>>> from sympy.ntheory import is_primitive_root, n_order, totient
>>> is_primitive_root(3, 10)
True
>>> is_primitive_root(9, 10)
False
>>> n_order(3, 10) == totient(10)
True
>>> n_order(9, 10) == totient(10)
False
```

sympy.ntheory.residue_ntheory.primitive_root(p)

Returns the smallest primitive root or None

p : positive integer

In [ ]:

```
>>> from sympy.ntheory.residue_ntheory import primitive_root
>>> primitive_root(19)
2
```

sympy.ntheory.residue_ntheory.sqrt_mod(a, p, all_roots=False)

find a root of "x**2 = a mod p"

a : integer p : positive integer all_roots : if True the list of roots is returned or None

In [ ]:

```
>>> from sympy.ntheory import sqrt_mod
>>> sqrt_mod(11, 43)
21
>>> sqrt_mod(17, 32, True)
[7, 9, 23, 25]
```

sympy.ntheory.residue_ntheory.quadratic_residues(p)

Returns the list of quadratic residues.

In [ ]:

```
>>> from sympy.ntheory.residue_ntheory import quadratic_residues
>>> quadratic_residues(7)
[0, 1, 2, 4]
```

sympy.ntheory.residue_ntheory.nthroot_mod(a, n, p, all_roots=False)

find the solutions to "x**n = a mod p"

In [ ]:

```
>>> from sympy.ntheory.residue_ntheory import nthroot_mod
>>> nthroot_mod(11, 4, 19)
8
>>> nthroot_mod(11, 4, 19, True)
[8, 11]
>>> nthroot_mod(68, 3, 109)
23
```

sympy.ntheory.residue_ntheory.is_nthpow_residue(a, n, m)

Returns True if "x**n == a (mod m)" has solutions.

Hackman "Elementary Number Theory" (2009), page 76

sympy.ntheory.residue_ntheory.is_quad_residue(a, p)

In [ ]:

```
>>> from sympy.ntheory import is_quad_residue
>>> list(set([i**2 % 7 for i in range(7)]))
[0, 1, 2, 4]
>>> [j for j in range(7) if is_quad_residue(j, 7)]
[0, 1, 2, 4]
```

legendre_symbol, jacobi_symbol

sympy.ntheory.residue_ntheory.legendre_symbol(a, p)

0 if a is multiple of p

1 if a is a quadratic residue of p

-1 otherwise

p should be an odd prime by definition

In [ ]:

```
>>> from sympy.ntheory import legendre_symbol
>>> [legendre_symbol(i, 7) for i in range(7)]
[0, 1, 1, -1, 1, -1, -1]
>>> list(set([i**2 % 7 for i in range(7)]))
[0, 1, 2, 4]
```

is_quad_residue, jacobi_symbol

sympy.ntheory.residue_ntheory.jacobi_symbol(m, n)

Returns the product of the legendre_symbol(m, p) for all the prime factors, p, of n.

0 if m cong 0 mod(n)

1 if x**2 cong m mod(n) has a solution

-1 otherwise

In [ ]:

```
>>> from sympy.ntheory import jacobi_symbol, legendre_symbol
>>> from sympy import Mul, S
>>> jacobi_symbol(45, 77)
-1
>>> jacobi_symbol(60, 121)
1
```

The relationship between the jacobi_symbol and legendre_symbol can be demonstrated as follows:

In [ ]:

```
>>> L = legendre_symbol
>>> S(45).factors()
{3: 2, 5: 1}
>>> jacobi_symbol(7, 45) == L(7, 3)**2 * L(7, 5)**1
True
```

is_quad_residue, legendre_symbol