Miscellaneous arithmetic functions¶
AUTHORS:
Kevin Stueve (2010-01-17): in
is_prime(n), delegated calculation ton.is_prime()
- sage.arith.misc.CRT(a, b, m=None, n=None)[source]¶
Return a solution to a Chinese Remainder Theorem problem.
INPUT:
a,b– two residues (elements of some ring for which extended gcd is available), or two lists, one of residues and one of modulim,n– (default:None) two moduli, orNone
OUTPUT:
If
m,nare notNone, returns a solution \(x\) to the simultaneous congruences \(x\equiv a \bmod m\) and \(x\equiv b \bmod n\), if one exists. By the Chinese Remainder Theorem, a solution to the simultaneous congruences exists if and only if \(a\equiv b\pmod{\gcd(m,n)}\). The solution \(x\) is only well-defined modulo \(\text{lcm}(m,n)\).If
aandbare lists, returns a simultaneous solution to the congruences \(x\equiv a_i\pmod{b_i}\), if one exists.See also
EXAMPLES:
Using
crtby giving it pairs of residues and moduli:sage: crt(2, 1, 3, 5) 11 sage: crt(13, 20, 100, 301) 28013 sage: crt([2, 1], [3, 5]) 11 sage: crt([13, 20], [100, 301]) 28013
>>> from sage.all import * >>> crt(Integer(2), Integer(1), Integer(3), Integer(5)) 11 >>> crt(Integer(13), Integer(20), Integer(100), Integer(301)) 28013 >>> crt([Integer(2), Integer(1)], [Integer(3), Integer(5)]) 11 >>> crt([Integer(13), Integer(20)], [Integer(100), Integer(301)]) 28013
crt(2, 1, 3, 5) crt(13, 20, 100, 301) crt([2, 1], [3, 5]) crt([13, 20], [100, 301])
You can also use upper case:
sage: c = CRT(2,3, 3, 5); c 8 sage: c % 3 == 2 True sage: c % 5 == 3 True
>>> from sage.all import * >>> c = CRT(Integer(2),Integer(3), Integer(3), Integer(5)); c 8 >>> c % Integer(3) == Integer(2) True >>> c % Integer(5) == Integer(3) True
c = CRT(2,3, 3, 5); c c % 3 == 2 c % 5 == 3
Note that this also works for polynomial rings:
sage: # needs sage.rings.number_field sage: x = polygen(ZZ, 'x') sage: K.<a> = NumberField(x^3 - 7) sage: R.<y> = K[] sage: f = y^2 + 3 sage: g = y^3 - 5 sage: CRT(1, 3, f, g) -3/26*y^4 + 5/26*y^3 + 15/26*y + 53/26 sage: CRT(1, a, f, g) (-3/52*a + 3/52)*y^4 + (5/52*a - 5/52)*y^3 + (15/52*a - 15/52)*y + 27/52*a + 25/52
>>> from sage.all import * >>> # needs sage.rings.number_field >>> x = polygen(ZZ, 'x') >>> K = NumberField(x**Integer(3) - Integer(7), names=('a',)); (a,) = K._first_ngens(1) >>> R = K['y']; (y,) = R._first_ngens(1) >>> f = y**Integer(2) + Integer(3) >>> g = y**Integer(3) - Integer(5) >>> CRT(Integer(1), Integer(3), f, g) -3/26*y^4 + 5/26*y^3 + 15/26*y + 53/26 >>> CRT(Integer(1), a, f, g) (-3/52*a + 3/52)*y^4 + (5/52*a - 5/52)*y^3 + (15/52*a - 15/52)*y + 27/52*a + 25/52
# needs sage.rings.number_field x = polygen(ZZ, 'x') K.<a> = NumberField(x^3 - 7) R.<y> = K[] f = y^2 + 3 g = y^3 - 5 CRT(1, 3, f, g) CRT(1, a, f, g)
You can also do this for any number of moduli:
sage: # needs sage.rings.number_field sage: K.<a> = NumberField(x^3 - 7) sage: R.<x> = K[] sage: CRT([], []) 0 sage: CRT([a], [x]) a sage: f = x^2 + 3 sage: g = x^3 - 5 sage: h = x^5 + x^2 - 9 sage: k = CRT([1, a, 3], [f, g, h]); k (127/26988*a - 5807/386828)*x^9 + (45/8996*a - 33677/1160484)*x^8 + (2/173*a - 6/173)*x^7 + (133/6747*a - 5373/96707)*x^6 + (-6/2249*a + 18584/290121)*x^5 + (-277/8996*a + 38847/386828)*x^4 + (-135/4498*a + 42673/193414)*x^3 + (-1005/8996*a + 470245/1160484)*x^2 + (-1215/8996*a + 141165/386828)*x + 621/8996*a + 836445/386828 sage: k.mod(f) 1 sage: k.mod(g) a sage: k.mod(h) 3
>>> from sage.all import * >>> # needs sage.rings.number_field >>> K = NumberField(x**Integer(3) - Integer(7), names=('a',)); (a,) = K._first_ngens(1) >>> R = K['x']; (x,) = R._first_ngens(1) >>> CRT([], []) 0 >>> CRT([a], [x]) a >>> f = x**Integer(2) + Integer(3) >>> g = x**Integer(3) - Integer(5) >>> h = x**Integer(5) + x**Integer(2) - Integer(9) >>> k = CRT([Integer(1), a, Integer(3)], [f, g, h]); k (127/26988*a - 5807/386828)*x^9 + (45/8996*a - 33677/1160484)*x^8 + (2/173*a - 6/173)*x^7 + (133/6747*a - 5373/96707)*x^6 + (-6/2249*a + 18584/290121)*x^5 + (-277/8996*a + 38847/386828)*x^4 + (-135/4498*a + 42673/193414)*x^3 + (-1005/8996*a + 470245/1160484)*x^2 + (-1215/8996*a + 141165/386828)*x + 621/8996*a + 836445/386828 >>> k.mod(f) 1 >>> k.mod(g) a >>> k.mod(h) 3
# needs sage.rings.number_field K.<a> = NumberField(x^3 - 7) R.<x> = K[] CRT([], []) CRT([a], [x]) f = x^2 + 3 g = x^3 - 5 h = x^5 + x^2 - 9 k = CRT([1, a, 3], [f, g, h]); k k.mod(f) k.mod(g) k.mod(h)
If the moduli are not coprime, a solution may not exist:
sage: crt(4, 8, 8, 12) 20 sage: crt(4, 6, 8, 12) Traceback (most recent call last): ... ValueError: no solution to crt problem since gcd(8,12) does not divide 4-6 sage: x = polygen(QQ) sage: crt(2, 3, x - 1, x + 1) -1/2*x + 5/2 sage: crt(2, x, x^2 - 1, x^2 + 1) -1/2*x^3 + x^2 + 1/2*x + 1 sage: crt(2, x, x^2 - 1, x^3 - 1) Traceback (most recent call last): ... ValueError: no solution to crt problem since gcd(x^2 - 1,x^3 - 1) does not divide 2-x sage: crt(int(2), int(3), int(7), int(11)) 58
>>> from sage.all import * >>> crt(Integer(4), Integer(8), Integer(8), Integer(12)) 20 >>> crt(Integer(4), Integer(6), Integer(8), Integer(12)) Traceback (most recent call last): ... ValueError: no solution to crt problem since gcd(8,12) does not divide 4-6 >>> x = polygen(QQ) >>> crt(Integer(2), Integer(3), x - Integer(1), x + Integer(1)) -1/2*x + 5/2 >>> crt(Integer(2), x, x**Integer(2) - Integer(1), x**Integer(2) + Integer(1)) -1/2*x^3 + x^2 + 1/2*x + 1 >>> crt(Integer(2), x, x**Integer(2) - Integer(1), x**Integer(3) - Integer(1)) Traceback (most recent call last): ... ValueError: no solution to crt problem since gcd(x^2 - 1,x^3 - 1) does not divide 2-x >>> crt(int(Integer(2)), int(Integer(3)), int(Integer(7)), int(Integer(11))) 58
crt(4, 8, 8, 12) crt(4, 6, 8, 12) x = polygen(QQ) crt(2, 3, x - 1, x + 1) crt(2, x, x^2 - 1, x^2 + 1) crt(2, x, x^2 - 1, x^3 - 1) crt(int(2), int(3), int(7), int(11))
crt also work with numpy and gmpy2 numbers:
sage: import numpy # needs numpy sage: crt(numpy.int8(2), numpy.int8(3), numpy.int8(7), numpy.int8(11)) # needs numpy 58 sage: from gmpy2 import mpz sage: crt(mpz(2), mpz(3), mpz(7), mpz(11)) 58 sage: crt(mpz(2), 3, mpz(7), numpy.int8(11)) # needs numpy 58
>>> from sage.all import * >>> import numpy # needs numpy >>> crt(numpy.int8(Integer(2)), numpy.int8(Integer(3)), numpy.int8(Integer(7)), numpy.int8(Integer(11))) # needs numpy 58 >>> from gmpy2 import mpz >>> crt(mpz(Integer(2)), mpz(Integer(3)), mpz(Integer(7)), mpz(Integer(11))) 58 >>> crt(mpz(Integer(2)), Integer(3), mpz(Integer(7)), numpy.int8(Integer(11))) # needs numpy 58
import numpy # needs numpy crt(numpy.int8(2), numpy.int8(3), numpy.int8(7), numpy.int8(11)) # needs numpy from gmpy2 import mpz crt(mpz(2), mpz(3), mpz(7), mpz(11)) crt(mpz(2), 3, mpz(7), numpy.int8(11)) # needs numpy
- sage.arith.misc.CRT_basis(moduli)[source]¶
Return a CRT basis for the given moduli.
INPUT:
moduli– list of pairwise coprime moduli \(m\) which admit an extended Euclidean algorithm
OUTPUT:
a list of elements \(a_i\) of the same length as \(m\) such that \(a_i\) is congruent to 1 modulo \(m_i\) and to 0 modulo \(m_j\) for \(j\not=i\).
EXAMPLES:
sage: a1 = ZZ(mod(42,5)) sage: a2 = ZZ(mod(42,13)) sage: c1,c2 = CRT_basis([5,13]) sage: mod(a1*c1+a2*c2,5*13) 42
>>> from sage.all import * >>> a1 = ZZ(mod(Integer(42),Integer(5))) >>> a2 = ZZ(mod(Integer(42),Integer(13))) >>> c1,c2 = CRT_basis([Integer(5),Integer(13)]) >>> mod(a1*c1+a2*c2,Integer(5)*Integer(13)) 42
a1 = ZZ(mod(42,5)) a2 = ZZ(mod(42,13)) c1,c2 = CRT_basis([5,13]) mod(a1*c1+a2*c2,5*13)
A polynomial example:
sage: x=polygen(QQ) sage: mods = [x,x^2+1,2*x-3] sage: b = CRT_basis(mods) sage: b [-2/3*x^3 + x^2 - 2/3*x + 1, 6/13*x^3 - x^2 + 6/13*x, 8/39*x^3 + 8/39*x] sage: [[bi % mj for mj in mods] for bi in b] [[1, 0, 0], [0, 1, 0], [0, 0, 1]]
>>> from sage.all import * >>> x=polygen(QQ) >>> mods = [x,x**Integer(2)+Integer(1),Integer(2)*x-Integer(3)] >>> b = CRT_basis(mods) >>> b [-2/3*x^3 + x^2 - 2/3*x + 1, 6/13*x^3 - x^2 + 6/13*x, 8/39*x^3 + 8/39*x] >>> [[bi % mj for mj in mods] for bi in b] [[1, 0, 0], [0, 1, 0], [0, 0, 1]]
x=polygen(QQ) mods = [x,x^2+1,2*x-3] b = CRT_basis(mods) b [[bi % mj for mj in mods] for bi in b]
- sage.arith.misc.CRT_list(values, moduli)[source]¶
Given a list
valuesof elements and a list of correspondingmoduli, find a single element that reduces to each element ofvaluesmodulo the corresponding moduli.See also
EXAMPLES:
sage: CRT_list([2,3,2], [3,5,7]) 23 sage: x = polygen(QQ) sage: c = CRT_list([3], [x]); c 3 sage: c.parent() Univariate Polynomial Ring in x over Rational Field
>>> from sage.all import * >>> CRT_list([Integer(2),Integer(3),Integer(2)], [Integer(3),Integer(5),Integer(7)]) 23 >>> x = polygen(QQ) >>> c = CRT_list([Integer(3)], [x]); c 3 >>> c.parent() Univariate Polynomial Ring in x over Rational Field
CRT_list([2,3,2], [3,5,7]) x = polygen(QQ) c = CRT_list([3], [x]); c c.parent()
It also works if the moduli are not coprime:
sage: CRT_list([32,2,2],[60,90,150]) 452
>>> from sage.all import * >>> CRT_list([Integer(32),Integer(2),Integer(2)],[Integer(60),Integer(90),Integer(150)]) 452
CRT_list([32,2,2],[60,90,150])
But with non coprime moduli there is not always a solution:
sage: CRT_list([32,2,1],[60,90,150]) Traceback (most recent call last): ... ValueError: no solution to crt problem since gcd(180,150) does not divide 92-1
>>> from sage.all import * >>> CRT_list([Integer(32),Integer(2),Integer(1)],[Integer(60),Integer(90),Integer(150)]) Traceback (most recent call last): ... ValueError: no solution to crt problem since gcd(180,150) does not divide 92-1
CRT_list([32,2,1],[60,90,150])
The arguments must be lists:
sage: CRT_list([1,2,3],"not a list") Traceback (most recent call last): ... ValueError: arguments to CRT_list should be lists sage: CRT_list("not a list",[2,3]) Traceback (most recent call last): ... ValueError: arguments to CRT_list should be lists
>>> from sage.all import * >>> CRT_list([Integer(1),Integer(2),Integer(3)],"not a list") Traceback (most recent call last): ... ValueError: arguments to CRT_list should be lists >>> CRT_list("not a list",[Integer(2),Integer(3)]) Traceback (most recent call last): ... ValueError: arguments to CRT_list should be lists
CRT_list([1,2,3],"not a list") CRT_list("not a list",[2,3])The list of moduli must have the same length as the list of elements:
sage: CRT_list([1,2,3],[2,3,5]) 23 sage: CRT_list([1,2,3],[2,3]) Traceback (most recent call last): ... ValueError: arguments to CRT_list should be lists of the same length sage: CRT_list([1,2,3],[2,3,5,7]) Traceback (most recent call last): ... ValueError: arguments to CRT_list should be lists of the same length
>>> from sage.all import * >>> CRT_list([Integer(1),Integer(2),Integer(3)],[Integer(2),Integer(3),Integer(5)]) 23 >>> CRT_list([Integer(1),Integer(2),Integer(3)],[Integer(2),Integer(3)]) Traceback (most recent call last): ... ValueError: arguments to CRT_list should be lists of the same length >>> CRT_list([Integer(1),Integer(2),Integer(3)],[Integer(2),Integer(3),Integer(5),Integer(7)]) Traceback (most recent call last): ... ValueError: arguments to CRT_list should be lists of the same length
CRT_list([1,2,3],[2,3,5]) CRT_list([1,2,3],[2,3]) CRT_list([1,2,3],[2,3,5,7])
- sage.arith.misc.CRT_vectors(X, moduli)[source]¶
Vector form of the Chinese Remainder Theorem: given a list of integer vectors \(v_i\) and a list of coprime moduli \(m_i\), find a vector \(w\) such that \(w = v_i \pmod m_i\) for all \(i\). This is more efficient than applying
CRT()to each entry.INPUT:
X– list or tuple, consisting of lists/tuples/vectors/etc of integers of the same lengthmoduli– list of len(X) moduli
OUTPUT: list; application of CRT componentwise
EXAMPLES:
sage: CRT_vectors([[3,5,7],[3,5,11]], [2,3]) [3, 5, 5] sage: CRT_vectors([vector(ZZ, [2,3,1]), Sequence([1,7,8], ZZ)], [8,9]) # needs sage.modules [10, 43, 17]
>>> from sage.all import * >>> CRT_vectors([[Integer(3),Integer(5),Integer(7)],[Integer(3),Integer(5),Integer(11)]], [Integer(2),Integer(3)]) [3, 5, 5] >>> CRT_vectors([vector(ZZ, [Integer(2),Integer(3),Integer(1)]), Sequence([Integer(1),Integer(7),Integer(8)], ZZ)], [Integer(8),Integer(9)]) # needs sage.modules [10, 43, 17]
CRT_vectors([[3,5,7],[3,5,11]], [2,3]) CRT_vectors([vector(ZZ, [2,3,1]), Sequence([1,7,8], ZZ)], [8,9]) # needs sage.modules
- class sage.arith.misc.Euler_Phi[source]¶
Bases:
objectReturn the value of the Euler phi function on the integer \(n\). We defined this to be the number of positive integers <= n that are relatively prime to \(n\). Thus if \(n \leq 0\) then
euler_phi(n)is defined and equals 0.INPUT:
n– integer
EXAMPLES:
sage: euler_phi(1) 1 sage: euler_phi(2) 1 sage: euler_phi(3) # needs sage.libs.pari 2 sage: euler_phi(12) # needs sage.libs.pari 4 sage: euler_phi(37) # needs sage.libs.pari 36
>>> from sage.all import * >>> euler_phi(Integer(1)) 1 >>> euler_phi(Integer(2)) 1 >>> euler_phi(Integer(3)) # needs sage.libs.pari 2 >>> euler_phi(Integer(12)) # needs sage.libs.pari 4 >>> euler_phi(Integer(37)) # needs sage.libs.pari 36
euler_phi(1) euler_phi(2) euler_phi(3) # needs sage.libs.pari euler_phi(12) # needs sage.libs.pari euler_phi(37) # needs sage.libs.pari
Notice that
euler_phiis defined to be 0 on negative numbers and 0:sage: euler_phi(-1) 0 sage: euler_phi(0) 0 sage: type(euler_phi(0)) <class 'sage.rings.integer.Integer'>
>>> from sage.all import * >>> euler_phi(-Integer(1)) 0 >>> euler_phi(Integer(0)) 0 >>> type(euler_phi(Integer(0))) <class 'sage.rings.integer.Integer'>
euler_phi(-1) euler_phi(0) type(euler_phi(0))
We verify directly that the phi function is correct for 21:
sage: euler_phi(21) # needs sage.libs.pari 12 sage: [i for i in range(21) if gcd(21,i) == 1] [1, 2, 4, 5, 8, 10, 11, 13, 16, 17, 19, 20]
>>> from sage.all import * >>> euler_phi(Integer(21)) # needs sage.libs.pari 12 >>> [i for i in range(Integer(21)) if gcd(Integer(21),i) == Integer(1)] [1, 2, 4, 5, 8, 10, 11, 13, 16, 17, 19, 20]
euler_phi(21) # needs sage.libs.pari [i for i in range(21) if gcd(21,i) == 1]
The length of the list of integers ‘i’ in
range(n)such that thegcd(i,n) == 1equalseuler_phi(n):sage: len([i for i in range(21) if gcd(21,i) == 1]) == euler_phi(21) # needs sage.libs.pari True
>>> from sage.all import * >>> len([i for i in range(Integer(21)) if gcd(Integer(21),i) == Integer(1)]) == euler_phi(Integer(21)) # needs sage.libs.pari True
len([i for i in range(21) if gcd(21,i) == 1]) == euler_phi(21) # needs sage.libs.pari
The phi function also has a special plotting method:
sage: P = plot(euler_phi, -3, 71) # needs sage.libs.pari sage.plot
>>> from sage.all import * >>> P = plot(euler_phi, -Integer(3), Integer(71)) # needs sage.libs.pari sage.plot
P = plot(euler_phi, -3, 71) # needs sage.libs.pari sage.plot
Tests with numpy and gmpy2 numbers:
sage: from numpy import int8 # needs numpy sage: euler_phi(int8(37)) # needs numpy sage.libs.pari 36 sage: from gmpy2 import mpz sage: euler_phi(mpz(37)) # needs sage.libs.pari 36
>>> from sage.all import * >>> from numpy import int8 # needs numpy >>> euler_phi(int8(Integer(37))) # needs numpy sage.libs.pari 36 >>> from gmpy2 import mpz >>> euler_phi(mpz(Integer(37))) # needs sage.libs.pari 36
from numpy import int8 # needs numpy euler_phi(int8(37)) # needs numpy sage.libs.pari from gmpy2 import mpz euler_phi(mpz(37)) # needs sage.libs.pari
AUTHORS:
William Stein
Alex Clemesha (2006-01-10): some examples
- plot(xmin=1, xmax=50, pointsize=30, rgbcolor=(0, 0, 1), join=True, **kwds)[source]¶
Plot the Euler phi function.
INPUT:
xmin– (default: 1)xmax– (default: 50)pointsize– (default: 30)rgbcolor– (default: (0,0,1))join– boolean (default:True); whether to join the points**kwds– passed on
EXAMPLES:
sage: from sage.arith.misc import Euler_Phi sage: p = Euler_Phi().plot() # needs sage.libs.pari sage.plot sage: p.ymax() # needs sage.libs.pari sage.plot 46.0
>>> from sage.all import * >>> from sage.arith.misc import Euler_Phi >>> p = Euler_Phi().plot() # needs sage.libs.pari sage.plot >>> p.ymax() # needs sage.libs.pari sage.plot 46.0
from sage.arith.misc import Euler_Phi p = Euler_Phi().plot() # needs sage.libs.pari sage.plot p.ymax() # needs sage.libs.pari sage.plot
- sage.arith.misc.GCD(a, b=None, **kwargs)[source]¶
Return the greatest common divisor of
aandb.If
ais a list andbis omitted, return instead the greatest common divisor of all elements ofa.INPUT:
a,b– two elements of a ring with gcd ora– list or tuple of elements of a ring with gcd
Additional keyword arguments are passed to the respectively called methods.
OUTPUT:
The given elements are first coerced into a common parent. Then, their greatest common divisor in that common parent is returned.
EXAMPLES:
sage: GCD(97,100) 1 sage: GCD(97*10^15, 19^20*97^2) 97 sage: GCD(2/3, 4/5) 2/15 sage: GCD([2,4,6,8]) 2 sage: GCD(srange(0,10000,10)) # fast !! 10
>>> from sage.all import * >>> GCD(Integer(97),Integer(100)) 1 >>> GCD(Integer(97)*Integer(10)**Integer(15), Integer(19)**Integer(20)*Integer(97)**Integer(2)) 97 >>> GCD(Integer(2)/Integer(3), Integer(4)/Integer(5)) 2/15 >>> GCD([Integer(2),Integer(4),Integer(6),Integer(8)]) 2 >>> GCD(srange(Integer(0),Integer(10000),Integer(10))) # fast !! 10
GCD(97,100) GCD(97*10^15, 19^20*97^2) GCD(2/3, 4/5) GCD([2,4,6,8]) GCD(srange(0,10000,10)) # fast !!
Note that to take the gcd of \(n\) elements for \(n \not= 2\) you must put the elements into a list by enclosing them in
[..]. Before Issue #4988 the following wrongly returned 3 since the third parameter was just ignored:sage: gcd(3, 6, 2) Traceback (most recent call last): ... TypeError: ...gcd() takes ... sage: gcd([3, 6, 2]) 1
>>> from sage.all import * >>> gcd(Integer(3), Integer(6), Integer(2)) Traceback (most recent call last): ... TypeError: ...gcd() takes ... >>> gcd([Integer(3), Integer(6), Integer(2)]) 1
gcd(3, 6, 2) gcd([3, 6, 2])
Similarly, giving just one element (which is not a list) gives an error:
sage: gcd(3) Traceback (most recent call last): ... TypeError: 'sage.rings.integer.Integer' object is not iterable
>>> from sage.all import * >>> gcd(Integer(3)) Traceback (most recent call last): ... TypeError: 'sage.rings.integer.Integer' object is not iterable
gcd(3)
By convention, the gcd of the empty list is (the integer) 0:
sage: gcd([]) 0 sage: type(gcd([])) <class 'sage.rings.integer.Integer'>
>>> from sage.all import * >>> gcd([]) 0 >>> type(gcd([])) <class 'sage.rings.integer.Integer'>
gcd([]) type(gcd([]))
- class sage.arith.misc.Moebius[source]¶
Bases:
objectReturn the value of the Möbius function of abs(n), where n is an integer.
DEFINITION: \(\mu(n)\) is 0 if \(n\) is not square free, and otherwise equals \((-1)^r\), where \(n\) has \(r\) distinct prime factors.
For simplicity, if \(n=0\) we define \(\mu(n) = 0\).
IMPLEMENTATION: Factors or - for integers - uses the PARI C library.
INPUT:
n– anything that can be factored
OUTPUT: 0, 1, or -1
EXAMPLES:
sage: # needs sage.libs.pari sage: moebius(-5) -1 sage: moebius(9) 0 sage: moebius(12) 0 sage: moebius(-35) 1 sage: moebius(-1) 1 sage: moebius(7) -1
>>> from sage.all import * >>> # needs sage.libs.pari >>> moebius(-Integer(5)) -1 >>> moebius(Integer(9)) 0 >>> moebius(Integer(12)) 0 >>> moebius(-Integer(35)) 1 >>> moebius(-Integer(1)) 1 >>> moebius(Integer(7)) -1
# needs sage.libs.pari moebius(-5) moebius(9) moebius(12) moebius(-35) moebius(-1) moebius(7)
sage: moebius(0) # potentially nonstandard! 0
>>> from sage.all import * >>> moebius(Integer(0)) # potentially nonstandard! 0
moebius(0) # potentially nonstandard!
The moebius function even makes sense for non-integer inputs.
sage: x = GF(7)['x'].0 sage: moebius(x + 2) # needs sage.libs.pari -1
>>> from sage.all import * >>> x = GF(Integer(7))['x'].gen(0) >>> moebius(x + Integer(2)) # needs sage.libs.pari -1
x = GF(7)['x'].0 moebius(x + 2) # needs sage.libs.pari
Tests with numpy and gmpy2 numbers:
sage: from numpy import int8 # needs numpy sage: moebius(int8(-5)) # needs numpy sage.libs.pari -1 sage: from gmpy2 import mpz sage: moebius(mpz(-5)) # needs sage.libs.pari -1
>>> from sage.all import * >>> from numpy import int8 # needs numpy >>> moebius(int8(-Integer(5))) # needs numpy sage.libs.pari -1 >>> from gmpy2 import mpz >>> moebius(mpz(-Integer(5))) # needs sage.libs.pari -1
from numpy import int8 # needs numpy moebius(int8(-5)) # needs numpy sage.libs.pari from gmpy2 import mpz moebius(mpz(-5)) # needs sage.libs.pari
- plot(xmin=0, xmax=50, pointsize=30, rgbcolor=(0, 0, 1), join=True, **kwds)[source]¶
Plot the Möbius function.
INPUT:
xmin– (default: 0)xmax– (default: 50)pointsize– (default: 30)rgbcolor– (default: (0,0,1))join– (default:True) whether to join the points(very helpful in seeing their order)
**kwds– passed on
EXAMPLES:
sage: from sage.arith.misc import Moebius sage: p = Moebius().plot() # needs sage.libs.pari sage.plot sage: p.ymax() # needs sage.libs.pari sage.plot 1.0
>>> from sage.all import * >>> from sage.arith.misc import Moebius >>> p = Moebius().plot() # needs sage.libs.pari sage.plot >>> p.ymax() # needs sage.libs.pari sage.plot 1.0
from sage.arith.misc import Moebius p = Moebius().plot() # needs sage.libs.pari sage.plot p.ymax() # needs sage.libs.pari sage.plot
- range(start, stop=None, step=None)[source]¶
Return the Möbius function evaluated at the given range of values, i.e., the image of the list range(start, stop, step) under the Möbius function.
This is much faster than directly computing all these values with a list comprehension.
EXAMPLES:
sage: # needs sage.libs.pari sage: v = moebius.range(-10, 10); v [1, 0, 0, -1, 1, -1, 0, -1, -1, 1, 0, 1, -1, -1, 0, -1, 1, -1, 0, 0] sage: v == [moebius(n) for n in range(-10, 10)] True sage: v = moebius.range(-1000, 2000, 4) sage: v == [moebius(n) for n in range(-1000, 2000, 4)] True
>>> from sage.all import * >>> # needs sage.libs.pari >>> v = moebius.range(-Integer(10), Integer(10)); v [1, 0, 0, -1, 1, -1, 0, -1, -1, 1, 0, 1, -1, -1, 0, -1, 1, -1, 0, 0] >>> v == [moebius(n) for n in range(-Integer(10), Integer(10))] True >>> v = moebius.range(-Integer(1000), Integer(2000), Integer(4)) >>> v == [moebius(n) for n in range(-Integer(1000), Integer(2000), Integer(4))] True
# needs sage.libs.pari v = moebius.range(-10, 10); v v == [moebius(n) for n in range(-10, 10)] v = moebius.range(-1000, 2000, 4) v == [moebius(n) for n in range(-1000, 2000, 4)]
- class sage.arith.misc.Sigma[source]¶
Bases:
objectReturn the sum of the \(k\)-th powers of the divisors of \(n\).
INPUT:
n– integerk– integer (default: 1)
OUTPUT: integer
EXAMPLES:
sage: sigma(5) 6 sage: sigma(5,2) 26
>>> from sage.all import * >>> sigma(Integer(5)) 6 >>> sigma(Integer(5),Integer(2)) 26
sigma(5) sigma(5,2)
The sigma function also has a special plotting method.
sage: P = plot(sigma, 1, 100) # needs sage.plot
>>> from sage.all import * >>> P = plot(sigma, Integer(1), Integer(100)) # needs sage.plot
P = plot(sigma, 1, 100) # needs sage.plot
This method also works with \(k\)-th powers.
sage: P = plot(sigma, 1, 100, k=2) # needs sage.plot
>>> from sage.all import * >>> P = plot(sigma, Integer(1), Integer(100), k=Integer(2)) # needs sage.plot
P = plot(sigma, 1, 100, k=2) # needs sage.plot
AUTHORS:
William Stein: original implementation
Craig Citro (2007-06-01): rewrote for huge speedup
- plot(xmin=1, xmax=50, k=1, pointsize=30, rgbcolor=(0, 0, 1), join=True, **kwds)[source]¶
Plot the sigma (sum of \(k\)-th powers of divisors) function.
INPUT:
xmin– (default: 1)xmax– (default: 50)k– (default: 1)pointsize– (default: 30)rgbcolor– (default: (0,0,1))join– (default:True) whether to join the points**kwds– passed on
EXAMPLES:
sage: from sage.arith.misc import Sigma sage: p = Sigma().plot() # needs sage.libs.pari sage.plot sage: p.ymax() # needs sage.libs.pari sage.plot 124.0
>>> from sage.all import * >>> from sage.arith.misc import Sigma >>> p = Sigma().plot() # needs sage.libs.pari sage.plot >>> p.ymax() # needs sage.libs.pari sage.plot 124.0
from sage.arith.misc import Sigma p = Sigma().plot() # needs sage.libs.pari sage.plot p.ymax() # needs sage.libs.pari sage.plot
- sage.arith.misc.XGCD(a, b=None)[source]¶
Return the greatest common divisor and the Bézout coefficients of the input arguments.
When both
aandbare given, then return a triple(g,s,t)such that \(g = s\cdot a+t\cdot b = \gcd(a,b)\). When onlyais given, then return a tuplerof lengthlen(a) + 1such that \(r_0 = \sum_{i = 0}^{len(a) - 1} r_{i + 1}a_i = gcd(a_0, \dots, a_{len(a) - 1})\)Note
One exception is if the elements are not in a principal ideal domain (see Wikipedia article Principal_ideal_domain), e.g., they are both polynomials over the integers. Then this function can’t in general return
(g,s,t)orras above, since they need not exist. Instead, over the integers, whenaandbare given, we first multiply \(g\) by a divisor of the resultant of \(a/g\) and \(b/g\), up to sign.INPUT:
One of the following:
a, b– integers or more generally, element of a ring for which the xgcd make sense (e.g. a field or univariate polynomials).a– a list or tuple of at least two integers or more generally, elements of a ring which the xgcd make sense.
OUTPUT:
One of the following:
g, s, t– when two inputsa, bare given. They satisfy \(g = s\cdot a + t\cdot b\).r– a tuple, when onlyais given (andb = None). Its first entryr[0]is the gcd of the inputs, and has length one longer than the length ofa. Its entries satisfy \(r_0 = \sum_{i = 0}^{len(a) - 1} r_{i + 1}a_i\).
Note
There is no guarantee that the returned cofactors (s and t) are minimal.
EXAMPLES:
sage: xgcd(56, 44) (4, 4, -5) sage: 4*56 + (-5)*44 4 sage: xgcd([56, 44]) (4, 4, -5) sage: r = xgcd([30, 105, 70, 42]); r (1, -255, 85, -17, -2) sage: (-255)*30 + 85*105 + (-17)*70 + (-2)*42 1 sage: xgcd([]) (0,) sage: xgcd([42]) (42, 1) sage: g, a, b = xgcd(5/1, 7/1); g, a, b (1, 3, -2) sage: a*(5/1) + b*(7/1) == g True sage: x = polygen(QQ) sage: xgcd(x^3 - 1, x^2 - 1) (x - 1, 1, -x) sage: g, a, b, c = xgcd([x^4 - x, x^6 - 1, x^4 - 1]); g, a, b, c (x - 1, x^3, -x, 1) sage: a*(x^4 - x) + b*(x^6 - 1) + c*(x^4 - 1) == g True sage: K.<g> = NumberField(x^2 - 3) # needs sage.rings.number_field sage: g.xgcd(g + 2) # needs sage.rings.number_field (1, 1/3*g, 0) sage: # needs sage.rings.number_field sage: R.<a,b> = K[] sage: S.<y> = R.fraction_field()[] sage: xgcd(y^2, a*y + b) (1, a^2/b^2, ((-a)/b^2)*y + 1/b) sage: xgcd((b+g)*y^2, (a+g)*y + b) (1, (a^2 + (2*g)*a + 3)/(b^3 + g*b^2), ((-a + (-g))/b^2)*y + 1/b)
>>> from sage.all import * >>> xgcd(Integer(56), Integer(44)) (4, 4, -5) >>> Integer(4)*Integer(56) + (-Integer(5))*Integer(44) 4 >>> xgcd([Integer(56), Integer(44)]) (4, 4, -5) >>> r = xgcd([Integer(30), Integer(105), Integer(70), Integer(42)]); r (1, -255, 85, -17, -2) >>> (-Integer(255))*Integer(30) + Integer(85)*Integer(105) + (-Integer(17))*Integer(70) + (-Integer(2))*Integer(42) 1 >>> xgcd([]) (0,) >>> xgcd([Integer(42)]) (42, 1) >>> g, a, b = xgcd(Integer(5)/Integer(1), Integer(7)/Integer(1)); g, a, b (1, 3, -2) >>> a*(Integer(5)/Integer(1)) + b*(Integer(7)/Integer(1)) == g True >>> x = polygen(QQ) >>> xgcd(x**Integer(3) - Integer(1), x**Integer(2) - Integer(1)) (x - 1, 1, -x) >>> g, a, b, c = xgcd([x**Integer(4) - x, x**Integer(6) - Integer(1), x**Integer(4) - Integer(1)]); g, a, b, c (x - 1, x^3, -x, 1) >>> a*(x**Integer(4) - x) + b*(x**Integer(6) - Integer(1)) + c*(x**Integer(4) - Integer(1)) == g True >>> K = NumberField(x**Integer(2) - Integer(3), names=('g',)); (g,) = K._first_ngens(1)# needs sage.rings.number_field >>> g.xgcd(g + Integer(2)) # needs sage.rings.number_field (1, 1/3*g, 0) >>> # needs sage.rings.number_field >>> R = K['a, b']; (a, b,) = R._first_ngens(2) >>> S = R.fraction_field()['y']; (y,) = S._first_ngens(1) >>> xgcd(y**Integer(2), a*y + b) (1, a^2/b^2, ((-a)/b^2)*y + 1/b) >>> xgcd((b+g)*y**Integer(2), (a+g)*y + b) (1, (a^2 + (2*g)*a + 3)/(b^3 + g*b^2), ((-a + (-g))/b^2)*y + 1/b)
xgcd(56, 44) 4*56 + (-5)*44 xgcd([56, 44]) r = xgcd([30, 105, 70, 42]); r (-255)*30 + 85*105 + (-17)*70 + (-2)*42 xgcd([]) xgcd([42]) g, a, b = xgcd(5/1, 7/1); g, a, b a*(5/1) + b*(7/1) == g x = polygen(QQ) xgcd(x^3 - 1, x^2 - 1) g, a, b, c = xgcd([x^4 - x, x^6 - 1, x^4 - 1]); g, a, b, c a*(x^4 - x) + b*(x^6 - 1) + c*(x^4 - 1) == g K.<g> = NumberField(x^2 - 3) # needs sage.rings.number_field g.xgcd(g + 2) # needs sage.rings.number_field # needs sage.rings.number_field R.<a,b> = K[] S.<y> = R.fraction_field()[] xgcd(y^2, a*y + b) xgcd((b+g)*y^2, (a+g)*y + b)
Here is an example of a xgcd for two polynomials over the integers, where the linear combination is not the gcd but the gcd multiplied by the resultant:
sage: R.<x> = ZZ[] sage: gcd(2*x*(x-1), x^2) x sage: xgcd(2*x*(x-1), x^2) (2*x, -1, 2) sage: (2*(x-1)).resultant(x) # needs sage.libs.pari 2
>>> from sage.all import * >>> R = ZZ['x']; (x,) = R._first_ngens(1) >>> gcd(Integer(2)*x*(x-Integer(1)), x**Integer(2)) x >>> xgcd(Integer(2)*x*(x-Integer(1)), x**Integer(2)) (2*x, -1, 2) >>> (Integer(2)*(x-Integer(1))).resultant(x) # needs sage.libs.pari 2
R.<x> = ZZ[] gcd(2*x*(x-1), x^2) xgcd(2*x*(x-1), x^2) (2*(x-1)).resultant(x) # needs sage.libs.pari
Tests with numpy and gmpy2 types:
sage: from numpy import int8 # needs numpy sage: xgcd(4, int8(8)) # needs numpy (4, 1, 0) sage: xgcd(int8(4), int8(8)) # needs numpy (4, 1, 0) sage: xgcd([int8(4), int8(8), int(10)]) # needs numpy (2, -2, 0, 1) sage: from gmpy2 import mpz sage: xgcd(mpz(4), mpz(8)) (4, 1, 0) sage: xgcd(4, mpz(8)) (4, 1, 0) sage: xgcd([4, mpz(8), mpz(10)]) (2, -2, 0, 1)
>>> from sage.all import * >>> from numpy import int8 # needs numpy >>> xgcd(Integer(4), int8(Integer(8))) # needs numpy (4, 1, 0) >>> xgcd(int8(Integer(4)), int8(Integer(8))) # needs numpy (4, 1, 0) >>> xgcd([int8(Integer(4)), int8(Integer(8)), int(Integer(10))]) # needs numpy (2, -2, 0, 1) >>> from gmpy2 import mpz >>> xgcd(mpz(Integer(4)), mpz(Integer(8))) (4, 1, 0) >>> xgcd(Integer(4), mpz(Integer(8))) (4, 1, 0) >>> xgcd([Integer(4), mpz(Integer(8)), mpz(Integer(10))]) (2, -2, 0, 1)
from numpy import int8 # needs numpy xgcd(4, int8(8)) # needs numpy xgcd(int8(4), int8(8)) # needs numpy xgcd([int8(4), int8(8), int(10)]) # needs numpy from gmpy2 import mpz xgcd(mpz(4), mpz(8)) xgcd(4, mpz(8)) xgcd([4, mpz(8), mpz(10)])
- sage.arith.misc.algdep(z, degree, known_bits=None, use_bits=None, known_digits=None, use_digits=None, height_bound=None, proof=False)[source]¶
Return an irreducible polynomial of degree at most \(degree\) which is approximately satisfied by the number \(z\).
You can specify the number of known bits or digits of \(z\) with
known_bits=korknown_digits=k. PARI is then told to compute the result using \(0.8k\) of these bits/digits. Or, you can specify the precision to use directly withuse_bits=koruse_digits=k. If none of these are specified, then the precision is taken from the input value.A height bound may be specified to indicate the maximum coefficient size of the returned polynomial; if a sufficiently small polynomial is not found, then
Nonewill be returned. Ifproof=Truethen the result is returned only if it can be proved correct (i.e. the only possible minimal polynomial satisfying the height bound, or no such polynomial exists). Otherwise aValueErroris raised indicating that higher precision is required.ALGORITHM: Uses LLL for real/complex inputs, PARI C-library
algdepcommand otherwise.Note that
algebraic_dependencyis a synonym foralgdep.INPUT:
z– real, complex, or \(p\)-adic numberdegree– integerheight_bound– integer (default:None); specifying the maximum coefficient size for the returned polynomialproof– boolean (default:False); requires height_bound to be set
EXAMPLES:
sage: algdep(1.888888888888888, 1) # needs sage.libs.pari 9*x - 17 sage: algdep(0.12121212121212, 1) # needs sage.libs.pari 33*x - 4 sage: algdep(sqrt(2), 2) # needs sage.libs.pari sage.symbolic x^2 - 2
>>> from sage.all import * >>> algdep(RealNumber('1.888888888888888'), Integer(1)) # needs sage.libs.pari 9*x - 17 >>> algdep(RealNumber('0.12121212121212'), Integer(1)) # needs sage.libs.pari 33*x - 4 >>> algdep(sqrt(Integer(2)), Integer(2)) # needs sage.libs.pari sage.symbolic x^2 - 2
algdep(1.888888888888888, 1) # needs sage.libs.pari algdep(0.12121212121212, 1) # needs sage.libs.pari algdep(sqrt(2), 2) # needs sage.libs.pari sage.symbolic
This example involves a complex number:
sage: z = (1/2) * (1 + RDF(sqrt(3)) * CC.0); z # needs sage.symbolic 0.500000000000000 + 0.866025403784439*I sage: algdep(z, 6) # needs sage.symbolic x^2 - x + 1
>>> from sage.all import * >>> z = (Integer(1)/Integer(2)) * (Integer(1) + RDF(sqrt(Integer(3))) * CC.gen(0)); z # needs sage.symbolic 0.500000000000000 + 0.866025403784439*I >>> algdep(z, Integer(6)) # needs sage.symbolic x^2 - x + 1
z = (1/2) * (1 + RDF(sqrt(3)) * CC.0); z # needs sage.symbolic algdep(z, 6) # needs sage.symbolic
This example involves a \(p\)-adic number:
sage: K = Qp(3, print_mode='series') # needs sage.rings.padics sage: a = K(7/19); a # needs sage.rings.padics 1 + 2*3 + 3^2 + 3^3 + 2*3^4 + 2*3^5 + 3^8 + 2*3^9 + 3^11 + 3^12 + 2*3^15 + 2*3^16 + 3^17 + 2*3^19 + O(3^20) sage: algdep(a, 1) # needs sage.rings.padics 19*x - 7
>>> from sage.all import * >>> K = Qp(Integer(3), print_mode='series') # needs sage.rings.padics >>> a = K(Integer(7)/Integer(19)); a # needs sage.rings.padics 1 + 2*3 + 3^2 + 3^3 + 2*3^4 + 2*3^5 + 3^8 + 2*3^9 + 3^11 + 3^12 + 2*3^15 + 2*3^16 + 3^17 + 2*3^19 + O(3^20) >>> algdep(a, Integer(1)) # needs sage.rings.padics 19*x - 7
K = Qp(3, print_mode='series') # needs sage.rings.padics a = K(7/19); a # needs sage.rings.padics algdep(a, 1) # needs sage.rings.padics
These examples show the importance of proper precision control. We compute a 200-bit approximation to \(sqrt(2)\) which is wrong in the 33’rd bit:
sage: # needs sage.libs.pari sage.rings.real_mpfr sage: z = sqrt(RealField(200)(2)) + (1/2)^33 sage: p = algdep(z, 4); p 227004321085*x^4 - 216947902586*x^3 - 99411220986*x^2 + 82234881648*x - 211871195088 sage: factor(p) 227004321085*x^4 - 216947902586*x^3 - 99411220986*x^2 + 82234881648*x - 211871195088 sage: algdep(z, 4, known_bits=32) x^2 - 2 sage: algdep(z, 4, known_digits=10) x^2 - 2 sage: algdep(z, 4, use_bits=25) x^2 - 2 sage: algdep(z, 4, use_digits=8) x^2 - 2
>>> from sage.all import * >>> # needs sage.libs.pari sage.rings.real_mpfr >>> z = sqrt(RealField(Integer(200))(Integer(2))) + (Integer(1)/Integer(2))**Integer(33) >>> p = algdep(z, Integer(4)); p 227004321085*x^4 - 216947902586*x^3 - 99411220986*x^2 + 82234881648*x - 211871195088 >>> factor(p) 227004321085*x^4 - 216947902586*x^3 - 99411220986*x^2 + 82234881648*x - 211871195088 >>> algdep(z, Integer(4), known_bits=Integer(32)) x^2 - 2 >>> algdep(z, Integer(4), known_digits=Integer(10)) x^2 - 2 >>> algdep(z, Integer(4), use_bits=Integer(25)) x^2 - 2 >>> algdep(z, Integer(4), use_digits=Integer(8)) x^2 - 2
# needs sage.libs.pari sage.rings.real_mpfr z = sqrt(RealField(200)(2)) + (1/2)^33 p = algdep(z, 4); p factor(p) algdep(z, 4, known_bits=32) algdep(z, 4, known_digits=10) algdep(z, 4, use_bits=25) algdep(z, 4, use_digits=8)
Using the
height_boundandproofparameters, we can see that \(pi\) is not the root of an integer polynomial of degree at most 5 and coefficients bounded above by 10:sage: algdep(pi.n(), 5, height_bound=10, proof=True) is None # needs sage.libs.pari sage.symbolic True
>>> from sage.all import * >>> algdep(pi.n(), Integer(5), height_bound=Integer(10), proof=True) is None # needs sage.libs.pari sage.symbolic True
algdep(pi.n(), 5, height_bound=10, proof=True) is None # needs sage.libs.pari sage.symbolic
For stronger results, we need more precision:
sage: # needs sage.libs.pari sage.symbolic sage: algdep(pi.n(), 5, height_bound=100, proof=True) is None Traceback (most recent call last): ... ValueError: insufficient precision for non-existence proof sage: algdep(pi.n(200), 5, height_bound=100, proof=True) is None True sage: algdep(pi.n(), 10, height_bound=10, proof=True) is None Traceback (most recent call last): ... ValueError: insufficient precision for non-existence proof sage: algdep(pi.n(200), 10, height_bound=10, proof=True) is None True
>>> from sage.all import * >>> # needs sage.libs.pari sage.symbolic >>> algdep(pi.n(), Integer(5), height_bound=Integer(100), proof=True) is None Traceback (most recent call last): ... ValueError: insufficient precision for non-existence proof >>> algdep(pi.n(Integer(200)), Integer(5), height_bound=Integer(100), proof=True) is None True >>> algdep(pi.n(), Integer(10), height_bound=Integer(10), proof=True) is None Traceback (most recent call last): ... ValueError: insufficient precision for non-existence proof >>> algdep(pi.n(Integer(200)), Integer(10), height_bound=Integer(10), proof=True) is None True
# needs sage.libs.pari sage.symbolic algdep(pi.n(), 5, height_bound=100, proof=True) is None algdep(pi.n(200), 5, height_bound=100, proof=True) is None algdep(pi.n(), 10, height_bound=10, proof=True) is None algdep(pi.n(200), 10, height_bound=10, proof=True) is None
We can also use
proof=Trueto get positive results:sage: # needs sage.libs.pari sage.symbolic sage: a = sqrt(2) + sqrt(3) + sqrt(5) sage: algdep(a.n(), 8, height_bound=1000, proof=True) Traceback (most recent call last): ... ValueError: insufficient precision for uniqueness proof sage: f = algdep(a.n(1000), 8, height_bound=1000, proof=True); f x^8 - 40*x^6 + 352*x^4 - 960*x^2 + 576 sage: f(a).expand() 0
>>> from sage.all import * >>> # needs sage.libs.pari sage.symbolic >>> a = sqrt(Integer(2)) + sqrt(Integer(3)) + sqrt(Integer(5)) >>> algdep(a.n(), Integer(8), height_bound=Integer(1000), proof=True) Traceback (most recent call last): ... ValueError: insufficient precision for uniqueness proof >>> f = algdep(a.n(Integer(1000)), Integer(8), height_bound=Integer(1000), proof=True); f x^8 - 40*x^6 + 352*x^4 - 960*x^2 + 576 >>> f(a).expand() 0
# needs sage.libs.pari sage.symbolic a = sqrt(2) + sqrt(3) + sqrt(5) algdep(a.n(), 8, height_bound=1000, proof=True) f = algdep(a.n(1000), 8, height_bound=1000, proof=True); f f(a).expand()
- sage.arith.misc.algebraic_dependency(z, degree, known_bits=None, use_bits=None, known_digits=None, use_digits=None, height_bound=None, proof=False)[source]¶
Return an irreducible polynomial of degree at most \(degree\) which is approximately satisfied by the number \(z\).
You can specify the number of known bits or digits of \(z\) with
known_bits=korknown_digits=k. PARI is then told to compute the result using \(0.8k\) of these bits/digits. Or, you can specify the precision to use directly withuse_bits=koruse_digits=k. If none of these are specified, then the precision is taken from the input value.A height bound may be specified to indicate the maximum coefficient size of the returned polynomial; if a sufficiently small polynomial is not found, then
Nonewill be returned. Ifproof=Truethen the result is returned only if it can be proved correct (i.e. the only possible minimal polynomial satisfying the height bound, or no such polynomial exists). Otherwise aValueErroris raised indicating that higher precision is required.ALGORITHM: Uses LLL for real/complex inputs, PARI C-library
algdepcommand otherwise.Note that
algebraic_dependencyis a synonym foralgdep.INPUT:
z– real, complex, or \(p\)-adic numberdegree– integerheight_bound– integer (default:None); specifying the maximum coefficient size for the returned polynomialproof– boolean (default:False); requires height_bound to be set
EXAMPLES:
sage: algdep(1.888888888888888, 1) # needs sage.libs.pari 9*x - 17 sage: algdep(0.12121212121212, 1) # needs sage.libs.pari 33*x - 4 sage: algdep(sqrt(2), 2) # needs sage.libs.pari sage.symbolic x^2 - 2
>>> from sage.all import * >>> algdep(RealNumber('1.888888888888888'), Integer(1)) # needs sage.libs.pari 9*x - 17 >>> algdep(RealNumber('0.12121212121212'), Integer(1)) # needs sage.libs.pari 33*x - 4 >>> algdep(sqrt(Integer(2)), Integer(2)) # needs sage.libs.pari sage.symbolic x^2 - 2
algdep(1.888888888888888, 1) # needs sage.libs.pari algdep(0.12121212121212, 1) # needs sage.libs.pari algdep(sqrt(2), 2) # needs sage.libs.pari sage.symbolic
This example involves a complex number:
sage: z = (1/2) * (1 + RDF(sqrt(3)) * CC.0); z # needs sage.symbolic 0.500000000000000 + 0.866025403784439*I sage: algdep(z, 6) # needs sage.symbolic x^2 - x + 1
>>> from sage.all import * >>> z = (Integer(1)/Integer(2)) * (Integer(1) + RDF(sqrt(Integer(3))) * CC.gen(0)); z # needs sage.symbolic 0.500000000000000 + 0.866025403784439*I >>> algdep(z, Integer(6)) # needs sage.symbolic x^2 - x + 1
z = (1/2) * (1 + RDF(sqrt(3)) * CC.0); z # needs sage.symbolic algdep(z, 6) # needs sage.symbolic
This example involves a \(p\)-adic number:
sage: K = Qp(3, print_mode='series') # needs sage.rings.padics sage: a = K(7/19); a # needs sage.rings.padics 1 + 2*3 + 3^2 + 3^3 + 2*3^4 + 2*3^5 + 3^8 + 2*3^9 + 3^11 + 3^12 + 2*3^15 + 2*3^16 + 3^17 + 2*3^19 + O(3^20) sage: algdep(a, 1) # needs sage.rings.padics 19*x - 7
>>> from sage.all import * >>> K = Qp(Integer(3), print_mode='series') # needs sage.rings.padics >>> a = K(Integer(7)/Integer(19)); a # needs sage.rings.padics 1 + 2*3 + 3^2 + 3^3 + 2*3^4 + 2*3^5 + 3^8 + 2*3^9 + 3^11 + 3^12 + 2*3^15 + 2*3^16 + 3^17 + 2*3^19 + O(3^20) >>> algdep(a, Integer(1)) # needs sage.rings.padics 19*x - 7
K = Qp(3, print_mode='series') # needs sage.rings.padics a = K(7/19); a # needs sage.rings.padics algdep(a, 1) # needs sage.rings.padics
These examples show the importance of proper precision control. We compute a 200-bit approximation to \(sqrt(2)\) which is wrong in the 33’rd bit:
sage: # needs sage.libs.pari sage.rings.real_mpfr sage: z = sqrt(RealField(200)(2)) + (1/2)^33 sage: p = algdep(z, 4); p 227004321085*x^4 - 216947902586*x^3 - 99411220986*x^2 + 82234881648*x - 211871195088 sage: factor(p) 227004321085*x^4 - 216947902586*x^3 - 99411220986*x^2 + 82234881648*x - 211871195088 sage: algdep(z, 4, known_bits=32) x^2 - 2 sage: algdep(z, 4, known_digits=10) x^2 - 2 sage: algdep(z, 4, use_bits=25) x^2 - 2 sage: algdep(z, 4, use_digits=8) x^2 - 2
>>> from sage.all import * >>> # needs sage.libs.pari sage.rings.real_mpfr >>> z = sqrt(RealField(Integer(200))(Integer(2))) + (Integer(1)/Integer(2))**Integer(33) >>> p = algdep(z, Integer(4)); p 227004321085*x^4 - 216947902586*x^3 - 99411220986*x^2 + 82234881648*x - 211871195088 >>> factor(p) 227004321085*x^4 - 216947902586*x^3 - 99411220986*x^2 + 82234881648*x - 211871195088 >>> algdep(z, Integer(4), known_bits=Integer(32)) x^2 - 2 >>> algdep(z, Integer(4), known_digits=Integer(10)) x^2 - 2 >>> algdep(z, Integer(4), use_bits=Integer(25)) x^2 - 2 >>> algdep(z, Integer(4), use_digits=Integer(8)) x^2 - 2
# needs sage.libs.pari sage.rings.real_mpfr z = sqrt(RealField(200)(2)) + (1/2)^33 p = algdep(z, 4); p factor(p) algdep(z, 4, known_bits=32) algdep(z, 4, known_digits=10) algdep(z, 4, use_bits=25) algdep(z, 4, use_digits=8)
Using the
height_boundandproofparameters, we can see that \(pi\) is not the root of an integer polynomial of degree at most 5 and coefficients bounded above by 10:sage: algdep(pi.n(), 5, height_bound=10, proof=True) is None # needs sage.libs.pari sage.symbolic True
>>> from sage.all import * >>> algdep(pi.n(), Integer(5), height_bound=Integer(10), proof=True) is None # needs sage.libs.pari sage.symbolic True
algdep(pi.n(), 5, height_bound=10, proof=True) is None # needs sage.libs.pari sage.symbolic
For stronger results, we need more precision:
sage: # needs sage.libs.pari sage.symbolic sage: algdep(pi.n(), 5, height_bound=100, proof=True) is None Traceback (most recent call last): ... ValueError: insufficient precision for non-existence proof sage: algdep(pi.n(200), 5, height_bound=100, proof=True) is None True sage: algdep(pi.n(), 10, height_bound=10, proof=True) is None Traceback (most recent call last): ... ValueError: insufficient precision for non-existence proof sage: algdep(pi.n(200), 10, height_bound=10, proof=True) is None True
>>> from sage.all import * >>> # needs sage.libs.pari sage.symbolic >>> algdep(pi.n(), Integer(5), height_bound=Integer(100), proof=True) is None Traceback (most recent call last): ... ValueError: insufficient precision for non-existence proof >>> algdep(pi.n(Integer(200)), Integer(5), height_bound=Integer(100), proof=True) is None True >>> algdep(pi.n(), Integer(10), height_bound=Integer(10), proof=True) is None Traceback (most recent call last): ... ValueError: insufficient precision for non-existence proof >>> algdep(pi.n(Integer(200)), Integer(10), height_bound=Integer(10), proof=True) is None True
# needs sage.libs.pari sage.symbolic algdep(pi.n(), 5, height_bound=100, proof=True) is None algdep(pi.n(200), 5, height_bound=100, proof=True) is None algdep(pi.n(), 10, height_bound=10, proof=True) is None algdep(pi.n(200), 10, height_bound=10, proof=True) is None
We can also use
proof=Trueto get positive results:sage: # needs sage.libs.pari sage.symbolic sage: a = sqrt(2) + sqrt(3) + sqrt(5) sage: algdep(a.n(), 8, height_bound=1000, proof=True) Traceback (most recent call last): ... ValueError: insufficient precision for uniqueness proof sage: f = algdep(a.n(1000), 8, height_bound=1000, proof=True); f x^8 - 40*x^6 + 352*x^4 - 960*x^2 + 576 sage: f(a).expand() 0
>>> from sage.all import * >>> # needs sage.libs.pari sage.symbolic >>> a = sqrt(Integer(2)) + sqrt(Integer(3)) + sqrt(Integer(5)) >>> algdep(a.n(), Integer(8), height_bound=Integer(1000), proof=True) Traceback (most recent call last): ... ValueError: insufficient precision for uniqueness proof >>> f = algdep(a.n(Integer(1000)), Integer(8), height_bound=Integer(1000), proof=True); f x^8 - 40*x^6 + 352*x^4 - 960*x^2 + 576 >>> f(a).expand() 0
# needs sage.libs.pari sage.symbolic a = sqrt(2) + sqrt(3) + sqrt(5) algdep(a.n(), 8, height_bound=1000, proof=True) f = algdep(a.n(1000), 8, height_bound=1000, proof=True); f f(a).expand()
- sage.arith.misc.bernoulli(n, algorithm='default', num_threads=1)[source]¶
Return the \(n\)-th Bernoulli number, as a rational number.
INPUT:
n– integeralgorithm:'default'– use ‘flint’ for n <= 20000, then ‘arb’ for n <= 300000 and ‘bernmm’ for larger values (this is just a heuristic, and not guaranteed to be optimal on all hardware)'arb'– use thebernoulli_fmpq_uifunction (formerly part of Arb) of the FLINT library'flint'– use thearith_bernoulli_numberfunction of the FLINT library'pari'– use the PARI C library'gap'– use GAP'gp'– use PARI/GP interpreter'magma'– use MAGMA (optional)'bernmm'– use bernmm package (a multimodular algorithm)
num_threads– positive integer, number of threads to use (only used for bernmm algorithm)
EXAMPLES:
sage: bernoulli(12) # needs sage.libs.flint -691/2730 sage: bernoulli(50) # needs sage.libs.flint 495057205241079648212477525/66
>>> from sage.all import * >>> bernoulli(Integer(12)) # needs sage.libs.flint -691/2730 >>> bernoulli(Integer(50)) # needs sage.libs.flint 495057205241079648212477525/66
bernoulli(12) # needs sage.libs.flint bernoulli(50) # needs sage.libs.flint
We demonstrate each of the alternative algorithms:
sage: bernoulli(12, algorithm='arb') # needs sage.libs.flint -691/2730 sage: bernoulli(12, algorithm='flint') # needs sage.libs.flint -691/2730 sage: bernoulli(12, algorithm='gap') # needs sage.libs.gap -691/2730 sage: bernoulli(12, algorithm='gp') # needs sage.libs.pari -691/2730 sage: bernoulli(12, algorithm='magma') # optional - magma -691/2730 sage: bernoulli(12, algorithm='pari') # needs sage.libs.pari -691/2730 sage: bernoulli(12, algorithm='bernmm') # needs sage.libs.ntl -691/2730 sage: bernoulli(12, algorithm='bernmm', num_threads=4) # needs sage.libs.ntl -691/2730
>>> from sage.all import * >>> bernoulli(Integer(12), algorithm='arb') # needs sage.libs.flint -691/2730 >>> bernoulli(Integer(12), algorithm='flint') # needs sage.libs.flint -691/2730 >>> bernoulli(Integer(12), algorithm='gap') # needs sage.libs.gap -691/2730 >>> bernoulli(Integer(12), algorithm='gp') # needs sage.libs.pari -691/2730 >>> bernoulli(Integer(12), algorithm='magma') # optional - magma -691/2730 >>> bernoulli(Integer(12), algorithm='pari') # needs sage.libs.pari -691/2730 >>> bernoulli(Integer(12), algorithm='bernmm') # needs sage.libs.ntl -691/2730 >>> bernoulli(Integer(12), algorithm='bernmm', num_threads=Integer(4)) # needs sage.libs.ntl -691/2730
bernoulli(12, algorithm='arb') # needs sage.libs.flint bernoulli(12, algorithm='flint') # needs sage.libs.flint bernoulli(12, algorithm='gap') # needs sage.libs.gap bernoulli(12, algorithm='gp') # needs sage.libs.pari bernoulli(12, algorithm='magma') # optional - magma bernoulli(12, algorithm='pari') # needs sage.libs.pari bernoulli(12, algorithm='bernmm') # needs sage.libs.ntl bernoulli(12, algorithm='bernmm', num_threads=4) # needs sage.libs.ntl
AUTHOR:
David Joyner and William Stein
- sage.arith.misc.binomial(x, m, **kwds)[source]¶
Return the binomial coefficient.
\[\binom{x}{m} = x (x-1) \cdots (x-m+1) / m!\]which is defined for \(m \in \ZZ\) and any \(x\). We extend this definition to include cases when \(x-m\) is an integer but \(m\) is not by
\[\binom{x}{m} = \binom{x}{x-m}\]If \(m < 0\), return \(0\).
INPUT:
x,m– numbers or symbolic expressions; eithermorx-mmust be an integer
OUTPUT: number or symbolic expression (if input is symbolic)
EXAMPLES:
sage: from sage.arith.misc import binomial sage: binomial(5, 2) 10 sage: binomial(2, 0) 1 sage: binomial(1/2, 0) # needs sage.libs.pari 1 sage: binomial(3, -1) 0 sage: binomial(20, 10) 184756 sage: binomial(-2, 5) -6 sage: binomial(-5, -2) 0 sage: binomial(RealField()('2.5'), 2) # needs sage.rings.real_mpfr 1.87500000000000 sage: binomial(Zp(5)(99),50) 3 + 4*5^3 + 2*5^4 + 4*5^5 + 4*5^6 + 4*5^7 + 4*5^8 + 5^9 + 2*5^10 + 3*5^11 + 4*5^12 + 4*5^13 + 2*5^14 + 3*5^15 + 3*5^16 + 4*5^17 + 4*5^18 + 2*5^19 + O(5^20) sage: binomial(Qp(3)(2/3),2) 2*3^-2 + 2*3^-1 + 2 + 2*3 + 2*3^2 + 2*3^3 + 2*3^4 + 2*3^5 + 2*3^6 + 2*3^7 + 2*3^8 + 2*3^9 + 2*3^10 + 2*3^11 + 2*3^12 + 2*3^13 + 2*3^14 + 2*3^15 + 2*3^16 + 2*3^17 + O(3^18) sage: n = var('n'); binomial(n, 2) # needs sage.symbolic 1/2*(n - 1)*n sage: n = var('n'); binomial(n, n) # needs sage.symbolic 1 sage: n = var('n'); binomial(n, n - 1) # needs sage.symbolic n sage: binomial(2^100, 2^100) 1 sage: x = polygen(ZZ) sage: binomial(x, 3) 1/6*x^3 - 1/2*x^2 + 1/3*x sage: binomial(x, x - 3) 1/6*x^3 - 1/2*x^2 + 1/3*x
>>> from sage.all import * >>> from sage.arith.misc import binomial >>> binomial(Integer(5), Integer(2)) 10 >>> binomial(Integer(2), Integer(0)) 1 >>> binomial(Integer(1)/Integer(2), Integer(0)) # needs sage.libs.pari 1 >>> binomial(Integer(3), -Integer(1)) 0 >>> binomial(Integer(20), Integer(10)) 184756 >>> binomial(-Integer(2), Integer(5)) -6 >>> binomial(-Integer(5), -Integer(2)) 0 >>> binomial(RealField()('2.5'), Integer(2)) # needs sage.rings.real_mpfr 1.87500000000000 >>> binomial(Zp(Integer(5))(Integer(99)),Integer(50)) 3 + 4*5^3 + 2*5^4 + 4*5^5 + 4*5^6 + 4*5^7 + 4*5^8 + 5^9 + 2*5^10 + 3*5^11 + 4*5^12 + 4*5^13 + 2*5^14 + 3*5^15 + 3*5^16 + 4*5^17 + 4*5^18 + 2*5^19 + O(5^20) >>> binomial(Qp(Integer(3))(Integer(2)/Integer(3)),Integer(2)) 2*3^-2 + 2*3^-1 + 2 + 2*3 + 2*3^2 + 2*3^3 + 2*3^4 + 2*3^5 + 2*3^6 + 2*3^7 + 2*3^8 + 2*3^9 + 2*3^10 + 2*3^11 + 2*3^12 + 2*3^13 + 2*3^14 + 2*3^15 + 2*3^16 + 2*3^17 + O(3^18) >>> n = var('n'); binomial(n, Integer(2)) # needs sage.symbolic 1/2*(n - 1)*n >>> n = var('n'); binomial(n, n) # needs sage.symbolic 1 >>> n = var('n'); binomial(n, n - Integer(1)) # needs sage.symbolic n >>> binomial(Integer(2)**Integer(100), Integer(2)**Integer(100)) 1 >>> x = polygen(ZZ) >>> binomial(x, Integer(3)) 1/6*x^3 - 1/2*x^2 + 1/3*x >>> binomial(x, x - Integer(3)) 1/6*x^3 - 1/2*x^2 + 1/3*x
from sage.arith.misc import binomial binomial(5, 2) binomial(2, 0) binomial(1/2, 0) # needs sage.libs.pari binomial(3, -1) binomial(20, 10) binomial(-2, 5) binomial(-5, -2) binomial(RealField()('2.5'), 2) # needs sage.rings.real_mpfr binomial(Zp(5)(99),50) binomial(Qp(3)(2/3),2) n = var('n'); binomial(n, 2) # needs sage.symbolic n = var('n'); binomial(n, n) # needs sage.symbolic n = var('n'); binomial(n, n - 1) # needs sage.symbolic binomial(2^100, 2^100) x = polygen(ZZ) binomial(x, 3) binomial(x, x - 3)If \(x \in \ZZ\), there is an optional ‘algorithm’ parameter, which can be ‘gmp’ (faster for small values; alias: ‘mpir’) or ‘pari’ (faster for large values):
sage: a = binomial(100, 45, algorithm='gmp') sage: b = binomial(100, 45, algorithm='pari') # needs sage.libs.pari sage: a == b # needs sage.libs.pari True
>>> from sage.all import * >>> a = binomial(Integer(100), Integer(45), algorithm='gmp') >>> b = binomial(Integer(100), Integer(45), algorithm='pari') # needs sage.libs.pari >>> a == b # needs sage.libs.pari True
a = binomial(100, 45, algorithm='gmp') b = binomial(100, 45, algorithm='pari') # needs sage.libs.pari a == b # needs sage.libs.pari
- sage.arith.misc.binomial_coefficients(n)[source]¶
Return a dictionary containing pairs \(\{(k_1,k_2) : C_{k,n}\}\) where \(C_{k_n}\) are binomial coefficients and \(n = k_1 + k_2\).
INPUT:
n– integer
OUTPUT: dictionary
EXAMPLES:
sage: sorted(binomial_coefficients(3).items()) [((0, 3), 1), ((1, 2), 3), ((2, 1), 3), ((3, 0), 1)]
>>> from sage.all import * >>> sorted(binomial_coefficients(Integer(3)).items()) [((0, 3), 1), ((1, 2), 3), ((2, 1), 3), ((3, 0), 1)]
sorted(binomial_coefficients(3).items())
Notice the coefficients above are the same as below:
sage: R.<x,y> = QQ[] sage: (x+y)^3 x^3 + 3*x^2*y + 3*x*y^2 + y^3
>>> from sage.all import * >>> R = QQ['x, y']; (x, y,) = R._first_ngens(2) >>> (x+y)**Integer(3) x^3 + 3*x^2*y + 3*x*y^2 + y^3
R.<x,y> = QQ[] (x+y)^3
Tests with numpy and gmpy2 numbers:
sage: from numpy import int8 # needs numpy sage: sorted(binomial_coefficients(int8(3)).items()) # needs numpy [((0, 3), 1), ((1, 2), 3), ((2, 1), 3), ((3, 0), 1)] sage: from gmpy2 import mpz sage: sorted(binomial_coefficients(mpz(3)).items()) [((0, 3), 1), ((1, 2), 3), ((2, 1), 3), ((3, 0), 1)]
>>> from sage.all import * >>> from numpy import int8 # needs numpy >>> sorted(binomial_coefficients(int8(Integer(3))).items()) # needs numpy [((0, 3), 1), ((1, 2), 3), ((2, 1), 3), ((3, 0), 1)] >>> from gmpy2 import mpz >>> sorted(binomial_coefficients(mpz(Integer(3))).items()) [((0, 3), 1), ((1, 2), 3), ((2, 1), 3), ((3, 0), 1)]
from numpy import int8 # needs numpy sorted(binomial_coefficients(int8(3)).items()) # needs numpy from gmpy2 import mpz sorted(binomial_coefficients(mpz(3)).items())
AUTHORS:
Fredrik Johansson
- sage.arith.misc.carmichael_lambda(n)[source]¶
Return the Carmichael function of a positive integer
n.The Carmichael function of \(n\), denoted \(\lambda(n)\), is the smallest positive integer \(k\) such that \(a^k \equiv 1 \pmod{n}\) for all \(a \in \ZZ/n\ZZ\) satisfying \(\gcd(a, n) = 1\). Thus, \(\lambda(n) = k\) is the exponent of the multiplicative group \((\ZZ/n\ZZ)^{\ast}\).
INPUT:
n– positive integer
OUTPUT: the Carmichael function of
nALGORITHM:
If \(n = 2, 4\) then \(\lambda(n) = \varphi(n)\). Let \(p \geq 3\) be an odd prime and let \(k\) be a positive integer. Then \(\lambda(p^k) = p^{k - 1}(p - 1) = \varphi(p^k)\). If \(k \geq 3\), then \(\lambda(2^k) = 2^{k - 2}\). Now consider the case where \(n > 3\) is composite and let \(n = p_1^{k_1} p_2^{k_2} \cdots p_t^{k_t}\) be the prime factorization of \(n\). Then
\[\lambda(n) = \lambda(p_1^{k_1} p_2^{k_2} \cdots p_t^{k_t}) = \text{lcm}(\lambda(p_1^{k_1}), \lambda(p_2^{k_2}), \dots, \lambda(p_t^{k_t}))\]EXAMPLES:
The Carmichael function of all positive integers up to and including 10:
sage: from sage.arith.misc import carmichael_lambda sage: list(map(carmichael_lambda, [1..10])) [1, 1, 2, 2, 4, 2, 6, 2, 6, 4]
>>> from sage.all import * >>> from sage.arith.misc import carmichael_lambda >>> list(map(carmichael_lambda, (ellipsis_range(Integer(1),Ellipsis,Integer(10))))) [1, 1, 2, 2, 4, 2, 6, 2, 6, 4]
from sage.arith.misc import carmichael_lambda list(map(carmichael_lambda, [1..10]))
The Carmichael function of the first ten primes:
sage: list(map(carmichael_lambda, primes_first_n(10))) # needs sage.libs.pari [1, 2, 4, 6, 10, 12, 16, 18, 22, 28]
>>> from sage.all import * >>> list(map(carmichael_lambda, primes_first_n(Integer(10)))) # needs sage.libs.pari [1, 2, 4, 6, 10, 12, 16, 18, 22, 28]
list(map(carmichael_lambda, primes_first_n(10))) # needs sage.libs.pari
Cases where the Carmichael function is equivalent to the Euler phi function:
sage: carmichael_lambda(2) == euler_phi(2) True sage: carmichael_lambda(4) == euler_phi(4) # needs sage.libs.pari True sage: p = random_prime(1000, lbound=3, proof=True) # needs sage.libs.pari sage: k = randint(1, 1000) sage: carmichael_lambda(p^k) == euler_phi(p^k) # needs sage.libs.pari True
>>> from sage.all import * >>> carmichael_lambda(Integer(2)) == euler_phi(Integer(2)) True >>> carmichael_lambda(Integer(4)) == euler_phi(Integer(4)) # needs sage.libs.pari True >>> p = random_prime(Integer(1000), lbound=Integer(3), proof=True) # needs sage.libs.pari >>> k = randint(Integer(1), Integer(1000)) >>> carmichael_lambda(p**k) == euler_phi(p**k) # needs sage.libs.pari True
carmichael_lambda(2) == euler_phi(2) carmichael_lambda(4) == euler_phi(4) # needs sage.libs.pari p = random_prime(1000, lbound=3, proof=True) # needs sage.libs.pari k = randint(1, 1000) carmichael_lambda(p^k) == euler_phi(p^k) # needs sage.libs.pari
A case where \(\lambda(n) \neq \varphi(n)\):
sage: k = randint(3, 1000) sage: carmichael_lambda(2^k) == 2^(k - 2) # needs sage.libs.pari True sage: carmichael_lambda(2^k) == 2^(k - 2) == euler_phi(2^k) # needs sage.libs.pari False
>>> from sage.all import * >>> k = randint(Integer(3), Integer(1000)) >>> carmichael_lambda(Integer(2)**k) == Integer(2)**(k - Integer(2)) # needs sage.libs.pari True >>> carmichael_lambda(Integer(2)**k) == Integer(2)**(k - Integer(2)) == euler_phi(Integer(2)**k) # needs sage.libs.pari False
k = randint(3, 1000) carmichael_lambda(2^k) == 2^(k - 2) # needs sage.libs.pari carmichael_lambda(2^k) == 2^(k - 2) == euler_phi(2^k) # needs sage.libs.pari
Verifying the current implementation of the Carmichael function using another implementation. The other implementation that we use for verification is an exhaustive search for the exponent of the multiplicative group \((\ZZ/n\ZZ)^{\ast}\).
sage: from sage.arith.misc import carmichael_lambda sage: n = randint(1, 500) sage: c = carmichael_lambda(n) sage: def coprime(n): ....: return [i for i in range(n) if gcd(i, n) == 1] sage: def znpower(n, k): ....: L = coprime(n) ....: return list(map(power_mod, L, [k]*len(L), [n]*len(L))) sage: def my_carmichael(n): ....: if n == 1: ....: return 1 ....: for k in range(1, n): ....: L = znpower(n, k) ....: ones = [1] * len(L) ....: T = [L[i] == ones[i] for i in range(len(L))] ....: if all(T): ....: return k sage: c == my_carmichael(n) True
>>> from sage.all import * >>> from sage.arith.misc import carmichael_lambda >>> n = randint(Integer(1), Integer(500)) >>> c = carmichael_lambda(n) >>> def coprime(n): ... return [i for i in range(n) if gcd(i, n) == Integer(1)] >>> def znpower(n, k): ... L = coprime(n) ... return list(map(power_mod, L, [k]*len(L), [n]*len(L))) >>> def my_carmichael(n): ... if n == Integer(1): ... return Integer(1) ... for k in range(Integer(1), n): ... L = znpower(n, k) ... ones = [Integer(1)] * len(L) ... T = [L[i] == ones[i] for i in range(len(L))] ... if all(T): ... return k >>> c == my_carmichael(n) True
from sage.arith.misc import carmichael_lambda n = randint(1, 500) c = carmichael_lambda(n) def coprime(n): return [i for i in range(n) if gcd(i, n) == 1] def znpower(n, k): L = coprime(n) return list(map(power_mod, L, [k]*len(L), [n]*len(L))) def my_carmichael(n): if n == 1: return 1 for k in range(1, n): L = znpower(n, k) ones = [1] * len(L) T = [L[i] == ones[i] for i in range(len(L))] if all(T): return k c == my_carmichael(n)Carmichael’s theorem states that \(a^{\lambda(n)} \equiv 1 \pmod{n}\) for all elements \(a\) of the multiplicative group \((\ZZ/n\ZZ)^{\ast}\). Here, we verify Carmichael’s theorem.
sage: from sage.arith.misc import carmichael_lambda sage: n = randint(2, 1000) sage: c = carmichael_lambda(n) sage: ZnZ = IntegerModRing(n) sage: M = ZnZ.list_of_elements_of_multiplicative_group() sage: ones = [1] * len(M) sage: P = [power_mod(a, c, n) for a in M] sage: P == ones True
>>> from sage.all import * >>> from sage.arith.misc import carmichael_lambda >>> n = randint(Integer(2), Integer(1000)) >>> c = carmichael_lambda(n) >>> ZnZ = IntegerModRing(n) >>> M = ZnZ.list_of_elements_of_multiplicative_group() >>> ones = [Integer(1)] * len(M) >>> P = [power_mod(a, c, n) for a in M] >>> P == ones True
from sage.arith.misc import carmichael_lambda n = randint(2, 1000) c = carmichael_lambda(n) ZnZ = IntegerModRing(n) M = ZnZ.list_of_elements_of_multiplicative_group() ones = [1] * len(M) P = [power_mod(a, c, n) for a in M] P == ones
REFERENCES:
- sage.arith.misc.continuant(v, n=None)[source]¶
Function returns the continuant of the sequence \(v\) (list or tuple).
Definition: see Graham, Knuth and Patashnik, Concrete Mathematics, section 6.7: Continuants. The continuant is defined by
\(K_0() = 1\)
\(K_1(x_1) = x_1\)
\(K_n(x_1, \cdots, x_n) = K_{n-1}(x_n, \cdots x_{n-1})x_n + K_{n-2}(x_1, \cdots, x_{n-2})\)
If
n = Noneorn > len(v)the defaultn = len(v)is used.INPUT:
v– list or tuple of elements of a ringn– (optional) integer
OUTPUT: element of ring (integer, polynomial, etcetera).
EXAMPLES:
sage: continuant([1,2,3]) 10 sage: p = continuant([2, 1, 2, 1, 1, 4, 1, 1, 6, 1, 1, 8, 1, 1, 10]) sage: q = continuant([1, 2, 1, 1, 4, 1, 1, 6, 1, 1, 8, 1, 1, 10]) sage: p/q 517656/190435 sage: F = continued_fraction([2, 1, 2, 1, 1, 4, 1, 1, 6, 1, 1, 8, 1, 1, 10]) sage: F.convergent(14) 517656/190435 sage: x = PolynomialRing(RationalField(), 'x', 5).gens() sage: continuant(x) x0*x1*x2*x3*x4 + x0*x1*x2 + x0*x1*x4 + x0*x3*x4 + x2*x3*x4 + x0 + x2 + x4 sage: continuant(x, 3) x0*x1*x2 + x0 + x2 sage: continuant(x, 2) x0*x1 + 1
>>> from sage.all import * >>> continuant([Integer(1),Integer(2),Integer(3)]) 10 >>> p = continuant([Integer(2), Integer(1), Integer(2), Integer(1), Integer(1), Integer(4), Integer(1), Integer(1), Integer(6), Integer(1), Integer(1), Integer(8), Integer(1), Integer(1), Integer(10)]) >>> q = continuant([Integer(1), Integer(2), Integer(1), Integer(1), Integer(4), Integer(1), Integer(1), Integer(6), Integer(1), Integer(1), Integer(8), Integer(1), Integer(1), Integer(10)]) >>> p/q 517656/190435 >>> F = continued_fraction([Integer(2), Integer(1), Integer(2), Integer(1), Integer(1), Integer(4), Integer(1), Integer(1), Integer(6), Integer(1), Integer(1), Integer(8), Integer(1), Integer(1), Integer(10)]) >>> F.convergent(Integer(14)) 517656/190435 >>> x = PolynomialRing(RationalField(), 'x', Integer(5)).gens() >>> continuant(x) x0*x1*x2*x3*x4 + x0*x1*x2 + x0*x1*x4 + x0*x3*x4 + x2*x3*x4 + x0 + x2 + x4 >>> continuant(x, Integer(3)) x0*x1*x2 + x0 + x2 >>> continuant(x, Integer(2)) x0*x1 + 1
continuant([1,2,3]) p = continuant([2, 1, 2, 1, 1, 4, 1, 1, 6, 1, 1, 8, 1, 1, 10]) q = continuant([1, 2, 1, 1, 4, 1, 1, 6, 1, 1, 8, 1, 1, 10]) p/q F = continued_fraction([2, 1, 2, 1, 1, 4, 1, 1, 6, 1, 1, 8, 1, 1, 10]) F.convergent(14) x = PolynomialRing(RationalField(), 'x', 5).gens() continuant(x) continuant(x, 3) continuant(x, 2)
We verify the identity
\[K_n(z,z,\cdots,z) = \sum_{k=0}^n \binom{n-k}{k} z^{n-2k}\]for \(n = 6\) using polynomial arithmetic:
sage: z = QQ['z'].0 sage: continuant((z,z,z,z,z,z,z,z,z,z,z,z,z,z,z), 6) z^6 + 5*z^4 + 6*z^2 + 1 sage: continuant(9) Traceback (most recent call last): ... TypeError: object of type 'sage.rings.integer.Integer' has no len()
>>> from sage.all import * >>> z = QQ['z'].gen(0) >>> continuant((z,z,z,z,z,z,z,z,z,z,z,z,z,z,z), Integer(6)) z^6 + 5*z^4 + 6*z^2 + 1 >>> continuant(Integer(9)) Traceback (most recent call last): ... TypeError: object of type 'sage.rings.integer.Integer' has no len()
z = QQ['z'].0 continuant((z,z,z,z,z,z,z,z,z,z,z,z,z,z,z), 6) continuant(9)
Tests with numpy and gmpy2 numbers:
sage: from numpy import int8 # needs numpy sage: continuant([int8(1), int8(2), int8(3)]) # needs numpy 10 sage: from gmpy2 import mpz sage: continuant([mpz(1), mpz(2), mpz(3)]) mpz(10)
>>> from sage.all import * >>> from numpy import int8 # needs numpy >>> continuant([int8(Integer(1)), int8(Integer(2)), int8(Integer(3))]) # needs numpy 10 >>> from gmpy2 import mpz >>> continuant([mpz(Integer(1)), mpz(Integer(2)), mpz(Integer(3))]) mpz(10)
from numpy import int8 # needs numpy continuant([int8(1), int8(2), int8(3)]) # needs numpy from gmpy2 import mpz continuant([mpz(1), mpz(2), mpz(3)])
AUTHORS:
Jaap Spies (2007-02-06)
- sage.arith.misc.coprime_part(x, base)[source]¶
Given an element
xof a Euclidean domain and a factor basebase, return the largest divisor ofxthat is not divisible by any element ofbase.ALGORITHM: Divide \(x\) by the
smooth_part().EXAMPLES:
sage: from sage.arith.misc import coprime_part, smooth_part sage: from sage.rings.generic import ProductTree sage: coprime_part(10^77+1, primes(10000)) 2159827213801295896328509719222460043196544298056155507343412527 sage: tree = ProductTree(primes(10000)) sage: coprime_part(10^55+1, tree) 6426667196963538873896485804232411 sage: coprime_part(10^55+1, tree).factor() 20163494891 * 318727841165674579776721 sage: prod(smooth_part(10^55+1, tree)) * coprime_part(10^55+1, tree) 10000000000000000000000000000000000000000000000000000001
>>> from sage.all import * >>> from sage.arith.misc import coprime_part, smooth_part >>> from sage.rings.generic import ProductTree >>> coprime_part(Integer(10)**Integer(77)+Integer(1), primes(Integer(10000))) 2159827213801295896328509719222460043196544298056155507343412527 >>> tree = ProductTree(primes(Integer(10000))) >>> coprime_part(Integer(10)**Integer(55)+Integer(1), tree) 6426667196963538873896485804232411 >>> coprime_part(Integer(10)**Integer(55)+Integer(1), tree).factor() 20163494891 * 318727841165674579776721 >>> prod(smooth_part(Integer(10)**Integer(55)+Integer(1), tree)) * coprime_part(Integer(10)**Integer(55)+Integer(1), tree) 10000000000000000000000000000000000000000000000000000001
from sage.arith.misc import coprime_part, smooth_part from sage.rings.generic import ProductTree coprime_part(10^77+1, primes(10000)) tree = ProductTree(primes(10000)) coprime_part(10^55+1, tree) coprime_part(10^55+1, tree).factor() prod(smooth_part(10^55+1, tree)) * coprime_part(10^55+1, tree)
- sage.arith.misc.crt(a, b, m=None, n=None)[source]¶
Return a solution to a Chinese Remainder Theorem problem.
INPUT:
a,b– two residues (elements of some ring for which extended gcd is available), or two lists, one of residues and one of modulim,n– (default:None) two moduli, orNone
OUTPUT:
If
m,nare notNone, returns a solution \(x\) to the simultaneous congruences \(x\equiv a \bmod m\) and \(x\equiv b \bmod n\), if one exists. By the Chinese Remainder Theorem, a solution to the simultaneous congruences exists if and only if \(a\equiv b\pmod{\gcd(m,n)}\). The solution \(x\) is only well-defined modulo \(\text{lcm}(m,n)\).If
aandbare lists, returns a simultaneous solution to the congruences \(x\equiv a_i\pmod{b_i}\), if one exists.See also
EXAMPLES:
Using
crtby giving it pairs of residues and moduli:sage: crt(2, 1, 3, 5) 11 sage: crt(13, 20, 100, 301) 28013 sage: crt([2, 1], [3, 5]) 11 sage: crt([13, 20], [100, 301]) 28013
>>> from sage.all import * >>> crt(Integer(2), Integer(1), Integer(3), Integer(5)) 11 >>> crt(Integer(13), Integer(20), Integer(100), Integer(301)) 28013 >>> crt([Integer(2), Integer(1)], [Integer(3), Integer(5)]) 11 >>> crt([Integer(13), Integer(20)], [Integer(100), Integer(301)]) 28013
crt(2, 1, 3, 5) crt(13, 20, 100, 301) crt([2, 1], [3, 5]) crt([13, 20], [100, 301])
You can also use upper case:
sage: c = CRT(2,3, 3, 5); c 8 sage: c % 3 == 2 True sage: c % 5 == 3 True
>>> from sage.all import * >>> c = CRT(Integer(2),Integer(3), Integer(3), Integer(5)); c 8 >>> c % Integer(3) == Integer(2) True >>> c % Integer(5) == Integer(3) True
c = CRT(2,3, 3, 5); c c % 3 == 2 c % 5 == 3
Note that this also works for polynomial rings:
sage: # needs sage.rings.number_field sage: x = polygen(ZZ, 'x') sage: K.<a> = NumberField(x^3 - 7) sage: R.<y> = K[] sage: f = y^2 + 3 sage: g = y^3 - 5 sage: CRT(1, 3, f, g) -3/26*y^4 + 5/26*y^3 + 15/26*y + 53/26 sage: CRT(1, a, f, g) (-3/52*a + 3/52)*y^4 + (5/52*a - 5/52)*y^3 + (15/52*a - 15/52)*y + 27/52*a + 25/52
>>> from sage.all import * >>> # needs sage.rings.number_field >>> x = polygen(ZZ, 'x') >>> K = NumberField(x**Integer(3) - Integer(7), names=('a',)); (a,) = K._first_ngens(1) >>> R = K['y']; (y,) = R._first_ngens(1) >>> f = y**Integer(2) + Integer(3) >>> g = y**Integer(3) - Integer(5) >>> CRT(Integer(1), Integer(3), f, g) -3/26*y^4 + 5/26*y^3 + 15/26*y + 53/26 >>> CRT(Integer(1), a, f, g) (-3/52*a + 3/52)*y^4 + (5/52*a - 5/52)*y^3 + (15/52*a - 15/52)*y + 27/52*a + 25/52
# needs sage.rings.number_field x = polygen(ZZ, 'x') K.<a> = NumberField(x^3 - 7) R.<y> = K[] f = y^2 + 3 g = y^3 - 5 CRT(1, 3, f, g) CRT(1, a, f, g)
You can also do this for any number of moduli:
sage: # needs sage.rings.number_field sage: K.<a> = NumberField(x^3 - 7) sage: R.<x> = K[] sage: CRT([], []) 0 sage: CRT([a], [x]) a sage: f = x^2 + 3 sage: g = x^3 - 5 sage: h = x^5 + x^2 - 9 sage: k = CRT([1, a, 3], [f, g, h]); k (127/26988*a - 5807/386828)*x^9 + (45/8996*a - 33677/1160484)*x^8 + (2/173*a - 6/173)*x^7 + (133/6747*a - 5373/96707)*x^6 + (-6/2249*a + 18584/290121)*x^5 + (-277/8996*a + 38847/386828)*x^4 + (-135/4498*a + 42673/193414)*x^3 + (-1005/8996*a + 470245/1160484)*x^2 + (-1215/8996*a + 141165/386828)*x + 621/8996*a + 836445/386828 sage: k.mod(f) 1 sage: k.mod(g) a sage: k.mod(h) 3
>>> from sage.all import * >>> # needs sage.rings.number_field >>> K = NumberField(x**Integer(3) - Integer(7), names=('a',)); (a,) = K._first_ngens(1) >>> R = K['x']; (x,) = R._first_ngens(1) >>> CRT([], []) 0 >>> CRT([a], [x]) a >>> f = x**Integer(2) + Integer(3) >>> g = x**Integer(3) - Integer(5) >>> h = x**Integer(5) + x**Integer(2) - Integer(9) >>> k = CRT([Integer(1), a, Integer(3)], [f, g, h]); k (127/26988*a - 5807/386828)*x^9 + (45/8996*a - 33677/1160484)*x^8 + (2/173*a - 6/173)*x^7 + (133/6747*a - 5373/96707)*x^6 + (-6/2249*a + 18584/290121)*x^5 + (-277/8996*a + 38847/386828)*x^4 + (-135/4498*a + 42673/193414)*x^3 + (-1005/8996*a + 470245/1160484)*x^2 + (-1215/8996*a + 141165/386828)*x + 621/8996*a + 836445/386828 >>> k.mod(f) 1 >>> k.mod(g) a >>> k.mod(h) 3
# needs sage.rings.number_field K.<a> = NumberField(x^3 - 7) R.<x> = K[] CRT([], []) CRT([a], [x]) f = x^2 + 3 g = x^3 - 5 h = x^5 + x^2 - 9 k = CRT([1, a, 3], [f, g, h]); k k.mod(f) k.mod(g) k.mod(h)
If the moduli are not coprime, a solution may not exist:
sage: crt(4, 8, 8, 12) 20 sage: crt(4, 6, 8, 12) Traceback (most recent call last): ... ValueError: no solution to crt problem since gcd(8,12) does not divide 4-6 sage: x = polygen(QQ) sage: crt(2, 3, x - 1, x + 1) -1/2*x + 5/2 sage: crt(2, x, x^2 - 1, x^2 + 1) -1/2*x^3 + x^2 + 1/2*x + 1 sage: crt(2, x, x^2 - 1, x^3 - 1) Traceback (most recent call last): ... ValueError: no solution to crt problem since gcd(x^2 - 1,x^3 - 1) does not divide 2-x sage: crt(int(2), int(3), int(7), int(11)) 58
>>> from sage.all import * >>> crt(Integer(4), Integer(8), Integer(8), Integer(12)) 20 >>> crt(Integer(4), Integer(6), Integer(8), Integer(12)) Traceback (most recent call last): ... ValueError: no solution to crt problem since gcd(8,12) does not divide 4-6 >>> x = polygen(QQ) >>> crt(Integer(2), Integer(3), x - Integer(1), x + Integer(1)) -1/2*x + 5/2 >>> crt(Integer(2), x, x**Integer(2) - Integer(1), x**Integer(2) + Integer(1)) -1/2*x^3 + x^2 + 1/2*x + 1 >>> crt(Integer(2), x, x**Integer(2) - Integer(1), x**Integer(3) - Integer(1)) Traceback (most recent call last): ... ValueError: no solution to crt problem since gcd(x^2 - 1,x^3 - 1) does not divide 2-x >>> crt(int(Integer(2)), int(Integer(3)), int(Integer(7)), int(Integer(11))) 58
crt(4, 8, 8, 12) crt(4, 6, 8, 12) x = polygen(QQ) crt(2, 3, x - 1, x + 1) crt(2, x, x^2 - 1, x^2 + 1) crt(2, x, x^2 - 1, x^3 - 1) crt(int(2), int(3), int(7), int(11))
crt also work with numpy and gmpy2 numbers:
sage: import numpy # needs numpy sage: crt(numpy.int8(2), numpy.int8(3), numpy.int8(7), numpy.int8(11)) # needs numpy 58 sage: from gmpy2 import mpz sage: crt(mpz(2), mpz(3), mpz(7), mpz(11)) 58 sage: crt(mpz(2), 3, mpz(7), numpy.int8(11)) # needs numpy 58
>>> from sage.all import * >>> import numpy # needs numpy >>> crt(numpy.int8(Integer(2)), numpy.int8(Integer(3)), numpy.int8(Integer(7)), numpy.int8(Integer(11))) # needs numpy 58 >>> from gmpy2 import mpz >>> crt(mpz(Integer(2)), mpz(Integer(3)), mpz(Integer(7)), mpz(Integer(11))) 58 >>> crt(mpz(Integer(2)), Integer(3), mpz(Integer(7)), numpy.int8(Integer(11))) # needs numpy 58
import numpy # needs numpy crt(numpy.int8(2), numpy.int8(3), numpy.int8(7), numpy.int8(11)) # needs numpy from gmpy2 import mpz crt(mpz(2), mpz(3), mpz(7), mpz(11)) crt(mpz(2), 3, mpz(7), numpy.int8(11)) # needs numpy
- sage.arith.misc.dedekind_psi(N)[source]¶
Return the value of the Dedekind psi function at
N.INPUT:
N– positive integer
OUTPUT: integer
The Dedekind psi function is the multiplicative function defined by
\[\psi(n) = n \prod_{p|n, p prime} (1 + 1/p).\]See Wikipedia article Dedekind_psi_function and OEIS sequence A001615.
EXAMPLES:
sage: from sage.arith.misc import dedekind_psi sage: [dedekind_psi(d) for d in range(1, 12)] [1, 3, 4, 6, 6, 12, 8, 12, 12, 18, 12]
>>> from sage.all import * >>> from sage.arith.misc import dedekind_psi >>> [dedekind_psi(d) for d in range(Integer(1), Integer(12))] [1, 3, 4, 6, 6, 12, 8, 12, 12, 18, 12]
from sage.arith.misc import dedekind_psi [dedekind_psi(d) for d in range(1, 12)]
- sage.arith.misc.dedekind_sum(p, q, algorithm='default')[source]¶
Return the Dedekind sum \(s(p,q)\) defined for integers \(p\), \(q\) as
\[s(p,q) = \sum_{i=0}^{q-1} \left(\!\left(\frac{i}{q}\right)\!\right) \left(\!\left(\frac{pi}{q}\right)\!\right)\]where
\[\begin{split}((x))=\begin{cases} x-\lfloor x \rfloor - \frac{1}{2} &\mbox{if } x \in \QQ \setminus \ZZ \\ 0 & \mbox{if } x \in \ZZ. \end{cases}\end{split}\]Warning
Caution is required as the Dedekind sum sometimes depends on the algorithm or is left undefined when \(p\) and \(q\) are not coprime.
INPUT:
p,q– integersalgorithm– must be one of the following'default'– (default) use FLINT'flint'– use FLINT'pari'– use PARI (gives different results if \(p\) and \(q\)are not coprime)
OUTPUT: a rational number
EXAMPLES:
Several small values:
sage: for q in range(10): print([dedekind_sum(p,q) for p in range(q+1)]) # needs sage.libs.flint [0] [0, 0] [0, 0, 0] [0, 1/18, -1/18, 0] [0, 1/8, 0, -1/8, 0] [0, 1/5, 0, 0, -1/5, 0] [0, 5/18, 1/18, 0, -1/18, -5/18, 0] [0, 5/14, 1/14, -1/14, 1/14, -1/14, -5/14, 0] [0, 7/16, 1/8, 1/16, 0, -1/16, -1/8, -7/16, 0] [0, 14/27, 4/27, 1/18, -4/27, 4/27, -1/18, -4/27, -14/27, 0]
>>> from sage.all import * >>> for q in range(Integer(10)): print([dedekind_sum(p,q) for p in range(q+Integer(1))]) # needs sage.libs.flint [0] [0, 0] [0, 0, 0] [0, 1/18, -1/18, 0] [0, 1/8, 0, -1/8, 0] [0, 1/5, 0, 0, -1/5, 0] [0, 5/18, 1/18, 0, -1/18, -5/18, 0] [0, 5/14, 1/14, -1/14, 1/14, -1/14, -5/14, 0] [0, 7/16, 1/8, 1/16, 0, -1/16, -1/8, -7/16, 0] [0, 14/27, 4/27, 1/18, -4/27, 4/27, -1/18, -4/27, -14/27, 0]
for q in range(10): print([dedekind_sum(p,q) for p in range(q+1)]) # needs sage.libs.flint
Check relations for restricted arguments:
sage: q = 23; dedekind_sum(1, q); (q-1)*(q-2)/(12*q) # needs sage.libs.flint 77/46 77/46 sage: p, q = 100, 723 # must be coprime sage: dedekind_sum(p, q) + dedekind_sum(q, p) # needs sage.libs.flint 31583/86760 sage: -1/4 + (p/q + q/p + 1/(p*q))/12 31583/86760
>>> from sage.all import * >>> q = Integer(23); dedekind_sum(Integer(1), q); (q-Integer(1))*(q-Integer(2))/(Integer(12)*q) # needs sage.libs.flint 77/46 77/46 >>> p, q = Integer(100), Integer(723) # must be coprime >>> dedekind_sum(p, q) + dedekind_sum(q, p) # needs sage.libs.flint 31583/86760 >>> -Integer(1)/Integer(4) + (p/q + q/p + Integer(1)/(p*q))/Integer(12) 31583/86760
q = 23; dedekind_sum(1, q); (q-1)*(q-2)/(12*q) # needs sage.libs.flint p, q = 100, 723 # must be coprime dedekind_sum(p, q) + dedekind_sum(q, p) # needs sage.libs.flint -1/4 + (p/q + q/p + 1/(p*q))/12
We check that evaluation works with large input:
sage: dedekind_sum(3^54 - 1, 2^93 + 1) # needs sage.libs.flint 459340694971839990630374299870/29710560942849126597578981379 sage: dedekind_sum(3^54 - 1, 2^93 + 1, algorithm='pari') # needs sage.libs.pari 459340694971839990630374299870/29710560942849126597578981379
>>> from sage.all import * >>> dedekind_sum(Integer(3)**Integer(54) - Integer(1), Integer(2)**Integer(93) + Integer(1)) # needs sage.libs.flint 459340694971839990630374299870/29710560942849126597578981379 >>> dedekind_sum(Integer(3)**Integer(54) - Integer(1), Integer(2)**Integer(93) + Integer(1), algorithm='pari') # needs sage.libs.pari 459340694971839990630374299870/29710560942849126597578981379
dedekind_sum(3^54 - 1, 2^93 + 1) # needs sage.libs.flint dedekind_sum(3^54 - 1, 2^93 + 1, algorithm='pari') # needs sage.libs.pari
We check consistency of the results:
sage: dedekind_sum(5, 7, algorithm='default') # needs sage.libs.flint -1/14 sage: dedekind_sum(5, 7, algorithm='flint') # needs sage.libs.flint -1/14 sage: dedekind_sum(5, 7, algorithm='pari') # needs sage.libs.pari -1/14 sage: dedekind_sum(6, 8, algorithm='default') # needs sage.libs.flint -1/8 sage: dedekind_sum(6, 8, algorithm='flint') # needs sage.libs.flint -1/8 sage: dedekind_sum(6, 8, algorithm='pari') # needs sage.libs.pari -1/8
>>> from sage.all import * >>> dedekind_sum(Integer(5), Integer(7), algorithm='default') # needs sage.libs.flint -1/14 >>> dedekind_sum(Integer(5), Integer(7), algorithm='flint') # needs sage.libs.flint -1/14 >>> dedekind_sum(Integer(5), Integer(7), algorithm='pari') # needs sage.libs.pari -1/14 >>> dedekind_sum(Integer(6), Integer(8), algorithm='default') # needs sage.libs.flint -1/8 >>> dedekind_sum(Integer(6), Integer(8), algorithm='flint') # needs sage.libs.flint -1/8 >>> dedekind_sum(Integer(6), Integer(8), algorithm='pari') # needs sage.libs.pari -1/8
dedekind_sum(5, 7, algorithm='default') # needs sage.libs.flint dedekind_sum(5, 7, algorithm='flint') # needs sage.libs.flint dedekind_sum(5, 7, algorithm='pari') # needs sage.libs.pari dedekind_sum(6, 8, algorithm='default') # needs sage.libs.flint dedekind_sum(6, 8, algorithm='flint') # needs sage.libs.flint dedekind_sum(6, 8, algorithm='pari') # needs sage.libs.pari
Tests with numpy and gmpy2 numbers:
sage: from numpy import int8 # needs numpy sage: dedekind_sum(int8(5), int8(7), algorithm='default') # needs numpy sage.libs.flint -1/14 sage: from gmpy2 import mpz sage: dedekind_sum(mpz(5), mpz(7), algorithm='default') # needs sage.libs.flint -1/14
>>> from sage.all import * >>> from numpy import int8 # needs numpy >>> dedekind_sum(int8(Integer(5)), int8(Integer(7)), algorithm='default') # needs numpy sage.libs.flint -1/14 >>> from gmpy2 import mpz >>> dedekind_sum(mpz(Integer(5)), mpz(Integer(7)), algorithm='default') # needs sage.libs.flint -1/14
from numpy import int8 # needs numpy dedekind_sum(int8(5), int8(7), algorithm='default') # needs numpy sage.libs.flint from gmpy2 import mpz dedekind_sum(mpz(5), mpz(7), algorithm='default') # needs sage.libs.flint
REFERENCES:
- sage.arith.misc.differences(lis, n=1)[source]¶
Return the \(n\) successive differences of the elements in
lis.EXAMPLES:
sage: differences(prime_range(50)) # needs sage.libs.pari [1, 2, 2, 4, 2, 4, 2, 4, 6, 2, 6, 4, 2, 4] sage: differences([i^2 for i in range(1,11)]) [3, 5, 7, 9, 11, 13, 15, 17, 19] sage: differences([i^3 + 3*i for i in range(1,21)]) [10, 22, 40, 64, 94, 130, 172, 220, 274, 334, 400, 472, 550, 634, 724, 820, 922, 1030, 1144] sage: differences([i^3 - i^2 for i in range(1,21)], 2) [10, 16, 22, 28, 34, 40, 46, 52, 58, 64, 70, 76, 82, 88, 94, 100, 106, 112] sage: differences([p - i^2 for i, p in enumerate(prime_range(50))], 3) # needs sage.libs.pari [-1, 2, -4, 4, -4, 4, 0, -6, 8, -6, 0, 4]
>>> from sage.all import * >>> differences(prime_range(Integer(50))) # needs sage.libs.pari [1, 2, 2, 4, 2, 4, 2, 4, 6, 2, 6, 4, 2, 4] >>> differences([i**Integer(2) for i in range(Integer(1),Integer(11))]) [3, 5, 7, 9, 11, 13, 15, 17, 19] >>> differences([i**Integer(3) + Integer(3)*i for i in range(Integer(1),Integer(21))]) [10, 22, 40, 64, 94, 130, 172, 220, 274, 334, 400, 472, 550, 634, 724, 820, 922, 1030, 1144] >>> differences([i**Integer(3) - i**Integer(2) for i in range(Integer(1),Integer(21))], Integer(2)) [10, 16, 22, 28, 34, 40, 46, 52, 58, 64, 70, 76, 82, 88, 94, 100, 106, 112] >>> differences([p - i**Integer(2) for i, p in enumerate(prime_range(Integer(50)))], Integer(3)) # needs sage.libs.pari [-1, 2, -4, 4, -4, 4, 0, -6, 8, -6, 0, 4]
differences(prime_range(50)) # needs sage.libs.pari differences([i^2 for i in range(1,11)]) differences([i^3 + 3*i for i in range(1,21)]) differences([i^3 - i^2 for i in range(1,21)], 2) differences([p - i^2 for i, p in enumerate(prime_range(50))], 3) # needs sage.libs.pari
Tests with numpy and gmpy2 numbers:
sage: from numpy import int8 # needs numpy sage: differences([int8(1), int8(4), int8(6), int8(19)]) # needs numpy [3, 2, 13] sage: from gmpy2 import mpz sage: differences([mpz(1), mpz(4), mpz(6), mpz(19)]) [mpz(3), mpz(2), mpz(13)]
>>> from sage.all import * >>> from numpy import int8 # needs numpy >>> differences([int8(Integer(1)), int8(Integer(4)), int8(Integer(6)), int8(Integer(19))]) # needs numpy [3, 2, 13] >>> from gmpy2 import mpz >>> differences([mpz(Integer(1)), mpz(Integer(4)), mpz(Integer(6)), mpz(Integer(19))]) [mpz(3), mpz(2), mpz(13)]
from numpy import int8 # needs numpy differences([int8(1), int8(4), int8(6), int8(19)]) # needs numpy from gmpy2 import mpz differences([mpz(1), mpz(4), mpz(6), mpz(19)])
AUTHORS:
Timothy Clemans (2008-03-09)
- sage.arith.misc.divisors(n)[source]¶
Return the list of all divisors (up to units) of this element of a unique factorization domain.
For an integer, the list of all positive integer divisors of this integer, sorted in increasing order, is returned.
INPUT:
n– the element
EXAMPLES:
Divisors of integers:
sage: divisors(-3) [1, 3] sage: divisors(6) [1, 2, 3, 6] sage: divisors(28) [1, 2, 4, 7, 14, 28] sage: divisors(2^5) [1, 2, 4, 8, 16, 32] sage: divisors(100) [1, 2, 4, 5, 10, 20, 25, 50, 100] sage: divisors(1) [1] sage: divisors(0) Traceback (most recent call last): ... ValueError: n must be nonzero sage: divisors(2^3 * 3^2 * 17) [1, 2, 3, 4, 6, 8, 9, 12, 17, 18, 24, 34, 36, 51, 68, 72, 102, 136, 153, 204, 306, 408, 612, 1224]
>>> from sage.all import * >>> divisors(-Integer(3)) [1, 3] >>> divisors(Integer(6)) [1, 2, 3, 6] >>> divisors(Integer(28)) [1, 2, 4, 7, 14, 28] >>> divisors(Integer(2)**Integer(5)) [1, 2, 4, 8, 16, 32] >>> divisors(Integer(100)) [1, 2, 4, 5, 10, 20, 25, 50, 100] >>> divisors(Integer(1)) [1] >>> divisors(Integer(0)) Traceback (most recent call last): ... ValueError: n must be nonzero >>> divisors(Integer(2)**Integer(3) * Integer(3)**Integer(2) * Integer(17)) [1, 2, 3, 4, 6, 8, 9, 12, 17, 18, 24, 34, 36, 51, 68, 72, 102, 136, 153, 204, 306, 408, 612, 1224]
divisors(-3) divisors(6) divisors(28) divisors(2^5) divisors(100) divisors(1) divisors(0) divisors(2^3 * 3^2 * 17)
This function works whenever one has unique factorization:
sage: # needs sage.rings.number_field sage: K.<a> = QuadraticField(7) sage: divisors(K.ideal(7)) [Fractional ideal (1), Fractional ideal (a), Fractional ideal (7)] sage: divisors(K.ideal(3)) [Fractional ideal (1), Fractional ideal (3), Fractional ideal (a - 2), Fractional ideal (a + 2)] sage: divisors(K.ideal(35)) [Fractional ideal (1), Fractional ideal (5), Fractional ideal (a), Fractional ideal (7), Fractional ideal (5*a), Fractional ideal (35)]
>>> from sage.all import * >>> # needs sage.rings.number_field >>> K = QuadraticField(Integer(7), names=('a',)); (a,) = K._first_ngens(1) >>> divisors(K.ideal(Integer(7))) [Fractional ideal (1), Fractional ideal (a), Fractional ideal (7)] >>> divisors(K.ideal(Integer(3))) [Fractional ideal (1), Fractional ideal (3), Fractional ideal (a - 2), Fractional ideal (a + 2)] >>> divisors(K.ideal(Integer(35))) [Fractional ideal (1), Fractional ideal (5), Fractional ideal (a), Fractional ideal (7), Fractional ideal (5*a), Fractional ideal (35)]
# needs sage.rings.number_field K.<a> = QuadraticField(7) divisors(K.ideal(7)) divisors(K.ideal(3)) divisors(K.ideal(35))
- sage.arith.misc.eratosthenes(n)[source]¶
Return a list of the primes \(\leq n\).
This is extremely slow and is for educational purposes only.
INPUT:
n– positive integer
OUTPUT: list of primes less than or equal to \(n\)
EXAMPLES:
sage: eratosthenes(3) [2, 3] sage: eratosthenes(50) [2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31, 37, 41, 43, 47] sage: len(eratosthenes(100)) 25 sage: eratosthenes(213) == prime_range(213) # needs sage.libs.pari True
>>> from sage.all import * >>> eratosthenes(Integer(3)) [2, 3] >>> eratosthenes(Integer(50)) [2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31, 37, 41, 43, 47] >>> len(eratosthenes(Integer(100))) 25 >>> eratosthenes(Integer(213)) == prime_range(Integer(213)) # needs sage.libs.pari True
eratosthenes(3) eratosthenes(50) len(eratosthenes(100)) eratosthenes(213) == prime_range(213) # needs sage.libs.pari
- sage.arith.misc.factor(n, proof=None, int_=False, algorithm='pari', verbose=0, **kwds)[source]¶
Return the factorization of
n. The result depends on the type ofn.If
nis an integer, returns the factorization as an object of typeFactorization.If
nis not an integer,n.factor(proof=proof, **kwds)gets called. Seen.factor??for more documentation in this case.Warning
This means that applying
factor()to an integer result of a symbolic computation will not factor the integer, because it is considered as an element of a larger symbolic ring.EXAMPLES:
sage: f(n) = n^2 # needs sage.symbolic sage: is_prime(f(3)) # needs sage.symbolic False sage: factor(f(3)) # needs sage.symbolic 9
>>> from sage.all import * >>> __tmp__=var("n"); f = symbolic_expression(n**Integer(2) ).function(n)# needs sage.symbolic >>> is_prime(f(Integer(3))) # needs sage.symbolic False >>> factor(f(Integer(3))) # needs sage.symbolic 9
f(n) = n^2 # needs sage.symbolic is_prime(f(3)) # needs sage.symbolic factor(f(3)) # needs sage.symbolic
INPUT:
n– nonzero integerproof– boolean orNone(default:None)int_– boolean (default:False); whether to return answers as Python integersalgorithm– string'pari'– (default) use the PARI C library'kash'– use KASH computer algebra system (requires that kash be installed)'magma'– use Magma (requires magma be installed)
verbose– integer (default: 0); PARI’s debug variable is set to this. E.g., set to 4 or 8 to see lots of output during factorization.
OUTPUT: factorization of \(n\)
The qsieve and ecm commands give access to highly optimized implementations of algorithms for doing certain integer factorization problems. These implementations are not used by the generic
factor()command, which currently just calls PARI (note that PARI also implements sieve and ecm algorithms, but they are not as optimized). Thus you might consider using them instead for certain numbers.The factorization returned is an element of the class
Factorization; useFactorization??to see more details, and examples below for usage. AFactorizationcontains both the unit factor (\(+1\) or \(-1\)) and a sorted list of(prime, exponent)pairs.The factorization displays in pretty-print format but it is easy to obtain access to the
(prime, exponent)pairs and the unit, to recover the number from its factorization, and even to multiply two factorizations. See examples below.EXAMPLES:
sage: factor(500) 2^2 * 5^3 sage: factor(-20) -1 * 2^2 * 5 sage: f=factor(-20) sage: list(f) [(2, 2), (5, 1)] sage: f.unit() -1 sage: f.value() -20 sage: factor(-next_prime(10^2) * next_prime(10^7)) # needs sage.libs.pari -1 * 101 * 10000019
>>> from sage.all import * >>> factor(Integer(500)) 2^2 * 5^3 >>> factor(-Integer(20)) -1 * 2^2 * 5 >>> f=factor(-Integer(20)) >>> list(f) [(2, 2), (5, 1)] >>> f.unit() -1 >>> f.value() -20 >>> factor(-next_prime(Integer(10)**Integer(2)) * next_prime(Integer(10)**Integer(7))) # needs sage.libs.pari -1 * 101 * 10000019
factor(500) factor(-20) f=factor(-20) list(f) f.unit() f.value() factor(-next_prime(10^2) * next_prime(10^7)) # needs sage.libs.pari
sage: factor(293292629867846432923017396246429, algorithm='flint') # needs sage.libs.flint 3 * 4852301647696687 * 20148007492971089
>>> from sage.all import * >>> factor(Integer(293292629867846432923017396246429), algorithm='flint') # needs sage.libs.flint 3 * 4852301647696687 * 20148007492971089
factor(293292629867846432923017396246429, algorithm='flint') # needs sage.libs.flint
sage: factor(-500, algorithm='kash') -1 * 2^2 * 5^3
>>> from sage.all import * >>> factor(-Integer(500), algorithm='kash') -1 * 2^2 * 5^3
factor(-500, algorithm='kash')
sage: factor(-500, algorithm='magma') # optional - magma -1 * 2^2 * 5^3
>>> from sage.all import * >>> factor(-Integer(500), algorithm='magma') # optional - magma -1 * 2^2 * 5^3
factor(-500, algorithm='magma') # optional - magma
sage: factor(0) Traceback (most recent call last): ... ArithmeticError: factorization of 0 is not defined sage: factor(1) 1 sage: factor(-1) -1 sage: factor(2^(2^7) + 1) # needs sage.libs.pari 59649589127497217 * 5704689200685129054721
>>> from sage.all import * >>> factor(Integer(0)) Traceback (most recent call last): ... ArithmeticError: factorization of 0 is not defined >>> factor(Integer(1)) 1 >>> factor(-Integer(1)) -1 >>> factor(Integer(2)**(Integer(2)**Integer(7)) + Integer(1)) # needs sage.libs.pari 59649589127497217 * 5704689200685129054721
factor(0) factor(1) factor(-1) factor(2^(2^7) + 1) # needs sage.libs.pari
Sage calls PARI’s pari:factor, which has
proof=Falseby default. Sage has a global proof flag, set toTrueby default (seesage.structure.proof.proof, or useproof.[tab]). To override the default, call this function withproof=False.sage: factor(3^89 - 1, proof=False) # needs sage.libs.pari 2 * 179 * 1611479891519807 * 5042939439565996049162197
>>> from sage.all import * >>> factor(Integer(3)**Integer(89) - Integer(1), proof=False) # needs sage.libs.pari 2 * 179 * 1611479891519807 * 5042939439565996049162197
factor(3^89 - 1, proof=False) # needs sage.libs.pari
sage: factor(2^197 + 1) # long time (2s) # needs sage.libs.pari 3 * 197002597249 * 1348959352853811313 * 251951573867253012259144010843
>>> from sage.all import * >>> factor(Integer(2)**Integer(197) + Integer(1)) # long time (2s) # needs sage.libs.pari 3 * 197002597249 * 1348959352853811313 * 251951573867253012259144010843
factor(2^197 + 1) # long time (2s) # needs sage.libs.pari
Any object which has a factor method can be factored like this:
sage: K.<i> = QuadraticField(-1) # needs sage.rings.number_field sage: factor(122 - 454*i) # needs sage.rings.number_field (-i) * (-i - 2)^3 * (i + 1)^3 * (-2*i + 3) * (i + 4)
>>> from sage.all import * >>> K = QuadraticField(-Integer(1), names=('i',)); (i,) = K._first_ngens(1)# needs sage.rings.number_field >>> factor(Integer(122) - Integer(454)*i) # needs sage.rings.number_field (-i) * (-i - 2)^3 * (i + 1)^3 * (-2*i + 3) * (i + 4)
K.<i> = QuadraticField(-1) # needs sage.rings.number_field factor(122 - 454*i) # needs sage.rings.number_field
To access the data in a factorization:
sage: # needs sage.libs.pari sage: f = factor(420); f 2^2 * 3 * 5 * 7 sage: [x for x in f] [(2, 2), (3, 1), (5, 1), (7, 1)] sage: [p for p,e in f] [2, 3, 5, 7] sage: [e for p,e in f] [2, 1, 1, 1] sage: [p^e for p,e in f] [4, 3, 5, 7]
>>> from sage.all import * >>> # needs sage.libs.pari >>> f = factor(Integer(420)); f 2^2 * 3 * 5 * 7 >>> [x for x in f] [(2, 2), (3, 1), (5, 1), (7, 1)] >>> [p for p,e in f] [2, 3, 5, 7] >>> [e for p,e in f] [2, 1, 1, 1] >>> [p**e for p,e in f] [4, 3, 5, 7]
# needs sage.libs.pari f = factor(420); f [x for x in f] [p for p,e in f] [e for p,e in f] [p^e for p,e in f]
We can factor Python, numpy and gmpy2 numbers:
sage: factor(math.pi) 3.141592653589793 sage: import numpy # needs numpy sage: factor(numpy.int8(30)) # needs numpy sage.libs.pari 2 * 3 * 5 sage: import gmpy2 sage: factor(gmpy2.mpz(30)) 2 * 3 * 5
>>> from sage.all import * >>> factor(math.pi) 3.141592653589793 >>> import numpy # needs numpy >>> factor(numpy.int8(Integer(30))) # needs numpy sage.libs.pari 2 * 3 * 5 >>> import gmpy2 >>> factor(gmpy2.mpz(Integer(30))) 2 * 3 * 5
factor(math.pi) import numpy # needs numpy factor(numpy.int8(30)) # needs numpy sage.libs.pari import gmpy2 factor(gmpy2.mpz(30))
- sage.arith.misc.factorial(n, algorithm='gmp')[source]¶
Compute the factorial of \(n\), which is the product \(1\cdot 2\cdot 3 \cdots (n-1)\cdot n\).
INPUT:
n– integeralgorithm– string (default:'gmp'):'gmp'– use the GMP C-library factorial function'pari'– use PARI’s factorial function
OUTPUT: integer
EXAMPLES:
sage: from sage.arith.misc import factorial sage: factorial(0) 1 sage: factorial(4) 24 sage: factorial(10) 3628800 sage: factorial(1) == factorial(0) True sage: factorial(6) == 6*5*4*3*2 True sage: factorial(1) == factorial(0) True sage: factorial(71) == 71* factorial(70) True sage: factorial(-32) Traceback (most recent call last): ... ValueError: factorial -- must be nonnegative
>>> from sage.all import * >>> from sage.arith.misc import factorial >>> factorial(Integer(0)) 1 >>> factorial(Integer(4)) 24 >>> factorial(Integer(10)) 3628800 >>> factorial(Integer(1)) == factorial(Integer(0)) True >>> factorial(Integer(6)) == Integer(6)*Integer(5)*Integer(4)*Integer(3)*Integer(2) True >>> factorial(Integer(1)) == factorial(Integer(0)) True >>> factorial(Integer(71)) == Integer(71)* factorial(Integer(70)) True >>> factorial(-Integer(32)) Traceback (most recent call last): ... ValueError: factorial -- must be nonnegative
from sage.arith.misc import factorial factorial(0) factorial(4) factorial(10) factorial(1) == factorial(0) factorial(6) == 6*5*4*3*2 factorial(1) == factorial(0) factorial(71) == 71* factorial(70) factorial(-32)
Tests with numpy and gmpy2 numbers:
sage: from numpy import int8 # needs numpy sage: factorial(int8(4)) # needs numpy 24 sage: from gmpy2 import mpz sage: factorial(mpz(4)) 24
>>> from sage.all import * >>> from numpy import int8 # needs numpy >>> factorial(int8(Integer(4))) # needs numpy 24 >>> from gmpy2 import mpz >>> factorial(mpz(Integer(4))) 24
from numpy import int8 # needs numpy factorial(int8(4)) # needs numpy from gmpy2 import mpz factorial(mpz(4))
PERFORMANCE: This discussion is valid as of April 2006. All timings below are on a Pentium Core Duo 2Ghz MacBook Pro running Linux with a 2.6.16.1 kernel.
It takes less than a minute to compute the factorial of \(10^7\) using the GMP algorithm, and the factorial of \(10^6\) takes less than 4 seconds.
The GMP algorithm is faster and more memory efficient than the PARI algorithm. E.g., PARI computes \(10^7\) factorial in 100 seconds on the core duo 2Ghz.
For comparison, computation in Magma \(\leq\) 2.12-10 of \(n!\) is best done using
*[1..n]. It takes 113 seconds to compute the factorial of \(10^7\) and 6 seconds to compute the factorial of \(10^6\). Mathematica V5.2 compute the factorial of \(10^7\) in 136 seconds and the factorial of \(10^6\) in 7 seconds. (Mathematica is notably very efficient at memory usage when doing factorial calculations.)
- sage.arith.misc.falling_factorial(x, a)[source]¶
Return the falling factorial \((x)_a\).
The notation in the literature is a mess: often \((x)_a\), but there are many other notations: GKP: Concrete Mathematics uses \(x^{\underline{a}}\).
Definition: for integer \(a \ge 0\) we have \(x(x-1) \cdots (x-a+1)\). In all other cases we use the GAMMA-function: \(\frac {\Gamma(x+1)} {\Gamma(x-a+1)}\).
INPUT:
x– element of a ringa– nonnegative integer orx,a– any numbers
OUTPUT: the falling factorial
See also
EXAMPLES:
sage: falling_factorial(10, 3) 720 sage: falling_factorial(10, 10) 3628800 sage: factorial(10) 3628800 sage: # needs sage.symbolic sage: falling_factorial(10, RR('3.0')) 720.000000000000 sage: falling_factorial(10, RR('3.3')) 1310.11633396601 sage: a = falling_factorial(1 + I, I); a gamma(I + 2) sage: CC(a) 0.652965496420167 + 0.343065839816545*I sage: falling_factorial(1 + I, 4) 4*I + 2 sage: falling_factorial(I, 4) -10 sage: M = MatrixSpace(ZZ, 4, 4) # needs sage.modules sage: A = M([1,0,1,0, 1,0,1,0, 1,0,10,10, 1,0,1,1]) # needs sage.modules sage: falling_factorial(A, 2) # A(A - I) # needs sage.modules [ 1 0 10 10] [ 1 0 10 10] [ 20 0 101 100] [ 2 0 11 10] sage: x = ZZ['x'].0 sage: falling_factorial(x, 4) x^4 - 6*x^3 + 11*x^2 - 6*x
>>> from sage.all import * >>> falling_factorial(Integer(10), Integer(3)) 720 >>> falling_factorial(Integer(10), Integer(10)) 3628800 >>> factorial(Integer(10)) 3628800 >>> # needs sage.symbolic >>> falling_factorial(Integer(10), RR('3.0')) 720.000000000000 >>> falling_factorial(Integer(10), RR('3.3')) 1310.11633396601 >>> a = falling_factorial(Integer(1) + I, I); a gamma(I + 2) >>> CC(a) 0.652965496420167 + 0.343065839816545*I >>> falling_factorial(Integer(1) + I, Integer(4)) 4*I + 2 >>> falling_factorial(I, Integer(4)) -10 >>> M = MatrixSpace(ZZ, Integer(4), Integer(4)) # needs sage.modules >>> A = M([Integer(1),Integer(0),Integer(1),Integer(0), Integer(1),Integer(0),Integer(1),Integer(0), Integer(1),Integer(0),Integer(10),Integer(10), Integer(1),Integer(0),Integer(1),Integer(1)]) # needs sage.modules >>> falling_factorial(A, Integer(2)) # A(A - I) # needs sage.modules [ 1 0 10 10] [ 1 0 10 10] [ 20 0 101 100] [ 2 0 11 10] >>> x = ZZ['x'].gen(0) >>> falling_factorial(x, Integer(4)) x^4 - 6*x^3 + 11*x^2 - 6*x
falling_factorial(10, 3) falling_factorial(10, 10) factorial(10) # needs sage.symbolic falling_factorial(10, RR('3.0')) falling_factorial(10, RR('3.3')) a = falling_factorial(1 + I, I); a CC(a) falling_factorial(1 + I, 4) falling_factorial(I, 4) M = MatrixSpace(ZZ, 4, 4) # needs sage.modules A = M([1,0,1,0, 1,0,1,0, 1,0,10,10, 1,0,1,1]) # needs sage.modules falling_factorial(A, 2) # A(A - I) # needs sage.modules x = ZZ['x'].0 falling_factorial(x, 4)AUTHORS:
Jaap Spies (2006-03-05)
- sage.arith.misc.four_squares(n)[source]¶
Write the integer \(n\) as a sum of four integer squares.
INPUT:
n– integer
OUTPUT: a tuple \((a,b,c,d)\) of nonnegative integers such that \(n = a^2 + b^2 + c^2 + d^2\) with \(a <= b <= c <= d\).
EXAMPLES:
sage: four_squares(3) (0, 1, 1, 1) sage: four_squares(13) (0, 0, 2, 3) sage: four_squares(130) (0, 0, 3, 11) sage: four_squares(1101011011004) (90, 102, 1220, 1049290) sage: four_squares(10^100 - 1) # needs sage.libs.pari (155024616290, 2612183768627, 14142135623730950488016887, 99999999999999999999999999999999999999999999999999) sage: for i in range(2^129, 2^129 + 10000): # long time # needs sage.libs.pari ....: S = four_squares(i) ....: assert sum(x^2 for x in S) == i
>>> from sage.all import * >>> four_squares(Integer(3)) (0, 1, 1, 1) >>> four_squares(Integer(13)) (0, 0, 2, 3) >>> four_squares(Integer(130)) (0, 0, 3, 11) >>> four_squares(Integer(1101011011004)) (90, 102, 1220, 1049290) >>> four_squares(Integer(10)**Integer(100) - Integer(1)) # needs sage.libs.pari (155024616290, 2612183768627, 14142135623730950488016887, 99999999999999999999999999999999999999999999999999) >>> for i in range(Integer(2)**Integer(129), Integer(2)**Integer(129) + Integer(10000)): # long time # needs sage.libs.pari ... S = four_squares(i) ... assert sum(x**Integer(2) for x in S) == i
four_squares(3) four_squares(13) four_squares(130) four_squares(1101011011004) four_squares(10^100 - 1) # needs sage.libs.pari for i in range(2^129, 2^129 + 10000): # long time # needs sage.libs.pari S = four_squares(i) assert sum(x^2 for x in S) == i
- sage.arith.misc.fundamental_discriminant(D)[source]¶
Return the discriminant of the quadratic extension \(K=Q(\sqrt{D})\), i.e. an integer d congruent to either 0 or 1, mod 4, and such that, at most, the only square dividing it is 4.
INPUT:
D– integer
OUTPUT: integer; the fundamental discriminant
EXAMPLES:
sage: fundamental_discriminant(102) 408 sage: fundamental_discriminant(720) 5 sage: fundamental_discriminant(2) 8
>>> from sage.all import * >>> fundamental_discriminant(Integer(102)) 408 >>> fundamental_discriminant(Integer(720)) 5 >>> fundamental_discriminant(Integer(2)) 8
fundamental_discriminant(102) fundamental_discriminant(720) fundamental_discriminant(2)
Tests with numpy and gmpy2 numbers:
sage: from numpy import int8 # needs numpy sage: fundamental_discriminant(int8(102)) # needs numpy 408 sage: from gmpy2 import mpz sage: fundamental_discriminant(mpz(102)) 408
>>> from sage.all import * >>> from numpy import int8 # needs numpy >>> fundamental_discriminant(int8(Integer(102))) # needs numpy 408 >>> from gmpy2 import mpz >>> fundamental_discriminant(mpz(Integer(102))) 408
from numpy import int8 # needs numpy fundamental_discriminant(int8(102)) # needs numpy from gmpy2 import mpz fundamental_discriminant(mpz(102))
- sage.arith.misc.gauss_sum(char_value, finite_field)[source]¶
Return the Gauss sums for a general finite field.
INPUT:
char_value– choice of multiplicative character, given by its value on thefinite_field.multiplicative_generator()finite_field– a finite field
OUTPUT:
an element of the parent ring of
char_value, that can be any field containing enough roots of unity, for example theUniversalCyclotomicField,QQbarorComplexFieldFor a finite field \(F\) of characteristic \(p\), the Gauss sum associated to a multiplicative character \(\chi\) (with values in a ring \(K\)) is defined as
\[\sum_{x \in F^{\times}} \chi(x) \zeta_p^{\operatorname{Tr} x},\]where \(\zeta_p \in K\) is a primitive root of unity of order \(p\) and Tr is the trace map from \(F\) to its prime field \(\GF{p}\).
For more info on Gauss sums, see Wikipedia article Gauss_sum.
Todo
Implement general Gauss sums for an arbitrary pair
(multiplicative_character, additive_character)EXAMPLES:
sage: # needs sage.libs.pari sage.rings.number_field sage: from sage.arith.misc import gauss_sum sage: F = GF(5); q = 5 sage: zq = UniversalCyclotomicField().zeta(q - 1) sage: L = [gauss_sum(zq**i, F) for i in range(5)]; L [-1, E(20)^4 + E(20)^13 - E(20)^16 - E(20)^17, E(5) - E(5)^2 - E(5)^3 + E(5)^4, E(20)^4 - E(20)^13 - E(20)^16 + E(20)^17, -1] sage: [g*g.conjugate() for g in L] [1, 5, 5, 5, 1] sage: # needs sage.libs.pari sage.rings.number_field sage: F = GF(11**2); q = 11**2 sage: zq = UniversalCyclotomicField().zeta(q - 1) sage: g = gauss_sum(zq**4, F) sage: g*g.conjugate() 121
>>> from sage.all import * >>> # needs sage.libs.pari sage.rings.number_field >>> from sage.arith.misc import gauss_sum >>> F = GF(Integer(5)); q = Integer(5) >>> zq = UniversalCyclotomicField().zeta(q - Integer(1)) >>> L = [gauss_sum(zq**i, F) for i in range(Integer(5))]; L [-1, E(20)^4 + E(20)^13 - E(20)^16 - E(20)^17, E(5) - E(5)^2 - E(5)^3 + E(5)^4, E(20)^4 - E(20)^13 - E(20)^16 + E(20)^17, -1] >>> [g*g.conjugate() for g in L] [1, 5, 5, 5, 1] >>> # needs sage.libs.pari sage.rings.number_field >>> F = GF(Integer(11)**Integer(2)); q = Integer(11)**Integer(2) >>> zq = UniversalCyclotomicField().zeta(q - Integer(1)) >>> g = gauss_sum(zq**Integer(4), F) >>> g*g.conjugate() 121
# needs sage.libs.pari sage.rings.number_field from sage.arith.misc import gauss_sum F = GF(5); q = 5 zq = UniversalCyclotomicField().zeta(q - 1) L = [gauss_sum(zq**i, F) for i in range(5)]; L [g*g.conjugate() for g in L] # needs sage.libs.pari sage.rings.number_field F = GF(11**2); q = 11**2 zq = UniversalCyclotomicField().zeta(q - 1) g = gauss_sum(zq**4, F) g*g.conjugate()
See also
sage.rings.padics.misc.gauss_sum()for a \(p\)-adic versionsage.modular.dirichlet.DirichletCharacter.gauss_sum()for prime finite fieldssage.modular.dirichlet.DirichletCharacter.gauss_sum_numerical()for prime finite fields
- sage.arith.misc.gcd(a, b=None, **kwargs)[source]¶
Return the greatest common divisor of
aandb.If
ais a list andbis omitted, return instead the greatest common divisor of all elements ofa.INPUT:
a,b– two elements of a ring with gcd ora– list or tuple of elements of a ring with gcd
Additional keyword arguments are passed to the respectively called methods.
OUTPUT:
The given elements are first coerced into a common parent. Then, their greatest common divisor in that common parent is returned.
EXAMPLES:
sage: GCD(97,100) 1 sage: GCD(97*10^15, 19^20*97^2) 97 sage: GCD(2/3, 4/5) 2/15 sage: GCD([2,4,6,8]) 2 sage: GCD(srange(0,10000,10)) # fast !! 10
>>> from sage.all import * >>> GCD(Integer(97),Integer(100)) 1 >>> GCD(Integer(97)*Integer(10)**Integer(15), Integer(19)**Integer(20)*Integer(97)**Integer(2)) 97 >>> GCD(Integer(2)/Integer(3), Integer(4)/Integer(5)) 2/15 >>> GCD([Integer(2),Integer(4),Integer(6),Integer(8)]) 2 >>> GCD(srange(Integer(0),Integer(10000),Integer(10))) # fast !! 10
GCD(97,100) GCD(97*10^15, 19^20*97^2) GCD(2/3, 4/5) GCD([2,4,6,8]) GCD(srange(0,10000,10)) # fast !!
Note that to take the gcd of \(n\) elements for \(n \not= 2\) you must put the elements into a list by enclosing them in
[..]. Before Issue #4988 the following wrongly returned 3 since the third parameter was just ignored:sage: gcd(3, 6, 2) Traceback (most recent call last): ... TypeError: ...gcd() takes ... sage: gcd([3, 6, 2]) 1
>>> from sage.all import * >>> gcd(Integer(3), Integer(6), Integer(2)) Traceback (most recent call last): ... TypeError: ...gcd() takes ... >>> gcd([Integer(3), Integer(6), Integer(2)]) 1
gcd(3, 6, 2) gcd([3, 6, 2])
Similarly, giving just one element (which is not a list) gives an error:
sage: gcd(3) Traceback (most recent call last): ... TypeError: 'sage.rings.integer.Integer' object is not iterable
>>> from sage.all import * >>> gcd(Integer(3)) Traceback (most recent call last): ... TypeError: 'sage.rings.integer.Integer' object is not iterable
gcd(3)
By convention, the gcd of the empty list is (the integer) 0:
sage: gcd([]) 0 sage: type(gcd([])) <class 'sage.rings.integer.Integer'>
>>> from sage.all import * >>> gcd([]) 0 >>> type(gcd([])) <class 'sage.rings.integer.Integer'>
gcd([]) type(gcd([]))
- sage.arith.misc.get_gcd(order)[source]¶
Return the fastest gcd function for integers of size no larger than order.
EXAMPLES:
sage: sage.arith.misc.get_gcd(4000) <built-in method gcd_int of sage.rings.fast_arith.arith_int object at ...> sage: sage.arith.misc.get_gcd(400000) <built-in method gcd_longlong of sage.rings.fast_arith.arith_llong object at ...> sage: sage.arith.misc.get_gcd(4000000000) <function gcd at ...>
>>> from sage.all import * >>> sage.arith.misc.get_gcd(Integer(4000)) <built-in method gcd_int of sage.rings.fast_arith.arith_int object at ...> >>> sage.arith.misc.get_gcd(Integer(400000)) <built-in method gcd_longlong of sage.rings.fast_arith.arith_llong object at ...> >>> sage.arith.misc.get_gcd(Integer(4000000000)) <function gcd at ...>
sage.arith.misc.get_gcd(4000) sage.arith.misc.get_gcd(400000) sage.arith.misc.get_gcd(4000000000)
- sage.arith.misc.get_inverse_mod(order)[source]¶
Return the fastest inverse_mod function for integers of size no larger than order.
EXAMPLES:
sage: sage.arith.misc.get_inverse_mod(6000) <built-in method inverse_mod_int of sage.rings.fast_arith.arith_int object at ...> sage: sage.arith.misc.get_inverse_mod(600000) <built-in method inverse_mod_longlong of sage.rings.fast_arith.arith_llong object at ...> sage: sage.arith.misc.get_inverse_mod(6000000000) <function inverse_mod at ...>
>>> from sage.all import * >>> sage.arith.misc.get_inverse_mod(Integer(6000)) <built-in method inverse_mod_int of sage.rings.fast_arith.arith_int object at ...> >>> sage.arith.misc.get_inverse_mod(Integer(600000)) <built-in method inverse_mod_longlong of sage.rings.fast_arith.arith_llong object at ...> >>> sage.arith.misc.get_inverse_mod(Integer(6000000000)) <function inverse_mod at ...>
sage.arith.misc.get_inverse_mod(6000) sage.arith.misc.get_inverse_mod(600000) sage.arith.misc.get_inverse_mod(6000000000)
- sage.arith.misc.hilbert_conductor(a, b)[source]¶
Return the product of all (finite) primes where the Hilbert symbol is -1.
This is the (reduced) discriminant of the quaternion algebra \((a,b)\) over \(\QQ\).
INPUT:
a,b– integers
OUTPUT: squarefree positive integer
EXAMPLES:
sage: # needs sage.libs.pari sage: hilbert_conductor(-1, -1) 2 sage: hilbert_conductor(-1, -11) 11 sage: hilbert_conductor(-2, -5) 5 sage: hilbert_conductor(-3, -17) 17
>>> from sage.all import * >>> # needs sage.libs.pari >>> hilbert_conductor(-Integer(1), -Integer(1)) 2 >>> hilbert_conductor(-Integer(1), -Integer(11)) 11 >>> hilbert_conductor(-Integer(2), -Integer(5)) 5 >>> hilbert_conductor(-Integer(3), -Integer(17)) 17
# needs sage.libs.pari hilbert_conductor(-1, -1) hilbert_conductor(-1, -11) hilbert_conductor(-2, -5) hilbert_conductor(-3, -17)
Tests with numpy and gmpy2 numbers:
sage: from numpy import int8 # needs numpy sage: hilbert_conductor(int8(-3), int8(-17)) # needs numpy sage.libs.pari 17 sage: from gmpy2 import mpz sage: hilbert_conductor(mpz(-3), mpz(-17)) # needs sage.libs.pari 17
>>> from sage.all import * >>> from numpy import int8 # needs numpy >>> hilbert_conductor(int8(-Integer(3)), int8(-Integer(17))) # needs numpy sage.libs.pari 17 >>> from gmpy2 import mpz >>> hilbert_conductor(mpz(-Integer(3)), mpz(-Integer(17))) # needs sage.libs.pari 17
from numpy import int8 # needs numpy hilbert_conductor(int8(-3), int8(-17)) # needs numpy sage.libs.pari from gmpy2 import mpz hilbert_conductor(mpz(-3), mpz(-17)) # needs sage.libs.pari
AUTHOR:
Gonzalo Tornaria (2009-03-02)
- sage.arith.misc.hilbert_conductor_inverse(d)[source]¶
Finds a pair of integers \((a,b)\) such that
hilbert_conductor(a,b) == d.The quaternion algebra \((a,b)\) over \(\QQ\) will then have (reduced) discriminant \(d\).
INPUT:
d– square-free positive integer
OUTPUT: pair of integers
EXAMPLES:
sage: # needs sage.libs.pari sage: hilbert_conductor_inverse(2) (-1, -1) sage: hilbert_conductor_inverse(3) (-1, -3) sage: hilbert_conductor_inverse(6) (-1, 3) sage: hilbert_conductor_inverse(30) (-3, -10) sage: hilbert_conductor_inverse(4) Traceback (most recent call last): ... ValueError: d needs to be squarefree sage: hilbert_conductor_inverse(-1) Traceback (most recent call last): ... ValueError: d needs to be positive
>>> from sage.all import * >>> # needs sage.libs.pari >>> hilbert_conductor_inverse(Integer(2)) (-1, -1) >>> hilbert_conductor_inverse(Integer(3)) (-1, -3) >>> hilbert_conductor_inverse(Integer(6)) (-1, 3) >>> hilbert_conductor_inverse(Integer(30)) (-3, -10) >>> hilbert_conductor_inverse(Integer(4)) Traceback (most recent call last): ... ValueError: d needs to be squarefree >>> hilbert_conductor_inverse(-Integer(1)) Traceback (most recent call last): ... ValueError: d needs to be positive
# needs sage.libs.pari hilbert_conductor_inverse(2) hilbert_conductor_inverse(3) hilbert_conductor_inverse(6) hilbert_conductor_inverse(30) hilbert_conductor_inverse(4) hilbert_conductor_inverse(-1)
AUTHOR:
Gonzalo Tornaria (2009-03-02)
- sage.arith.misc.hilbert_symbol(a, b, p, algorithm='pari')[source]¶
Return 1 if \(ax^2 + by^2\) \(p\)-adically represents a nonzero square, otherwise returns \(-1\). If either a or b is 0, returns 0.
INPUT:
a,b– integersp– integer; either prime or -1 (which represents the archimedean place)algorithm– string'pari'– (default) use the PARI C library'direct'– use a Python implementation'all'– use both PARI and direct and check thatthe results agree, then return the common answer
OUTPUT: integer (0, -1, or 1)
EXAMPLES:
sage: # needs sage.libs.pari sage: hilbert_symbol(-1, -1, -1, algorithm='all') -1 sage: hilbert_symbol(2, 3, 5, algorithm='all') 1 sage: hilbert_symbol(4, 3, 5, algorithm='all') 1 sage: hilbert_symbol(0, 3, 5, algorithm='all') 0 sage: hilbert_symbol(-1, -1, 2, algorithm='all') -1 sage: hilbert_symbol(1, -1, 2, algorithm='all') 1 sage: hilbert_symbol(3, -1, 2, algorithm='all') -1 sage: hilbert_symbol(QQ(-1)/QQ(4), -1, 2) == -1 # needs sage.libs.pari True sage: hilbert_symbol(QQ(-1)/QQ(4), -1, 3) == 1 # needs sage.libs.pari True
>>> from sage.all import * >>> # needs sage.libs.pari >>> hilbert_symbol(-Integer(1), -Integer(1), -Integer(1), algorithm='all') -1 >>> hilbert_symbol(Integer(2), Integer(3), Integer(5), algorithm='all') 1 >>> hilbert_symbol(Integer(4), Integer(3), Integer(5), algorithm='all') 1 >>> hilbert_symbol(Integer(0), Integer(3), Integer(5), algorithm='all') 0 >>> hilbert_symbol(-Integer(1), -Integer(1), Integer(2), algorithm='all') -1 >>> hilbert_symbol(Integer(1), -Integer(1), Integer(2), algorithm='all') 1 >>> hilbert_symbol(Integer(3), -Integer(1), Integer(2), algorithm='all') -1 >>> hilbert_symbol(QQ(-Integer(1))/QQ(Integer(4)), -Integer(1), Integer(2)) == -Integer(1) # needs sage.libs.pari True >>> hilbert_symbol(QQ(-Integer(1))/QQ(Integer(4)), -Integer(1), Integer(3)) == Integer(1) # needs sage.libs.pari True
# needs sage.libs.pari hilbert_symbol(-1, -1, -1, algorithm='all') hilbert_symbol(2, 3, 5, algorithm='all') hilbert_symbol(4, 3, 5, algorithm='all') hilbert_symbol(0, 3, 5, algorithm='all') hilbert_symbol(-1, -1, 2, algorithm='all') hilbert_symbol(1, -1, 2, algorithm='all') hilbert_symbol(3, -1, 2, algorithm='all') hilbert_symbol(QQ(-1)/QQ(4), -1, 2) == -1 # needs sage.libs.pari hilbert_symbol(QQ(-1)/QQ(4), -1, 3) == 1 # needs sage.libs.pari
Tests with numpy and gmpy2 numbers:
sage: from numpy import int8 # needs numpy sage: hilbert_symbol(int8(2), int8(3), int8(5), algorithm='all') # needs numpy sage.libs.pari 1 sage: from gmpy2 import mpz sage: hilbert_symbol(mpz(2), mpz(3), mpz(5), algorithm='all') # needs sage.libs.pari 1
>>> from sage.all import * >>> from numpy import int8 # needs numpy >>> hilbert_symbol(int8(Integer(2)), int8(Integer(3)), int8(Integer(5)), algorithm='all') # needs numpy sage.libs.pari 1 >>> from gmpy2 import mpz >>> hilbert_symbol(mpz(Integer(2)), mpz(Integer(3)), mpz(Integer(5)), algorithm='all') # needs sage.libs.pari 1
from numpy import int8 # needs numpy hilbert_symbol(int8(2), int8(3), int8(5), algorithm='all') # needs numpy sage.libs.pari from gmpy2 import mpz hilbert_symbol(mpz(2), mpz(3), mpz(5), algorithm='all') # needs sage.libs.pari
AUTHORS:
William Stein and David Kohel (2006-01-05)
- sage.arith.misc.integer_ceil(x)[source]¶
Return the ceiling of x.
EXAMPLES:
sage: integer_ceil(5.4) 6 sage: integer_ceil(x) # needs sage.symbolic Traceback (most recent call last): ... NotImplementedError: computation of ceil of x not implemented
>>> from sage.all import * >>> integer_ceil(RealNumber('5.4')) 6 >>> integer_ceil(x) # needs sage.symbolic Traceback (most recent call last): ... NotImplementedError: computation of ceil of x not implemented
integer_ceil(5.4) integer_ceil(x) # needs sage.symbolic
Tests with numpy and gmpy2 numbers:
sage: from numpy import float32 # needs numpy sage: integer_ceil(float32(5.4)) # needs numpy 6 sage: from gmpy2 import mpfr sage: integer_ceil(mpfr(5.4)) 6
>>> from sage.all import * >>> from numpy import float32 # needs numpy >>> integer_ceil(float32(RealNumber('5.4'))) # needs numpy 6 >>> from gmpy2 import mpfr >>> integer_ceil(mpfr(RealNumber('5.4'))) 6
from numpy import float32 # needs numpy integer_ceil(float32(5.4)) # needs numpy from gmpy2 import mpfr integer_ceil(mpfr(5.4))
- sage.arith.misc.integer_floor(x)[source]¶
Return the largest integer \(\leq x\).
INPUT:
x– an object that has a floor method or is coercible to integer
OUTPUT: integer
EXAMPLES:
sage: integer_floor(5.4) 5 sage: integer_floor(float(5.4)) 5 sage: integer_floor(-5/2) -3 sage: integer_floor(RDF(-5/2)) -3 sage: integer_floor(x) # needs sage.symbolic Traceback (most recent call last): ... NotImplementedError: computation of floor of x not implemented
>>> from sage.all import * >>> integer_floor(RealNumber('5.4')) 5 >>> integer_floor(float(RealNumber('5.4'))) 5 >>> integer_floor(-Integer(5)/Integer(2)) -3 >>> integer_floor(RDF(-Integer(5)/Integer(2))) -3 >>> integer_floor(x) # needs sage.symbolic Traceback (most recent call last): ... NotImplementedError: computation of floor of x not implemented
integer_floor(5.4) integer_floor(float(5.4)) integer_floor(-5/2) integer_floor(RDF(-5/2)) integer_floor(x) # needs sage.symbolic
Tests with numpy and gmpy2 numbers:
sage: from numpy import float32 # needs numpy sage: integer_floor(float32(5.4)) # needs numpy 5 sage: from gmpy2 import mpfr sage: integer_floor(mpfr(5.4)) 5
>>> from sage.all import * >>> from numpy import float32 # needs numpy >>> integer_floor(float32(RealNumber('5.4'))) # needs numpy 5 >>> from gmpy2 import mpfr >>> integer_floor(mpfr(RealNumber('5.4'))) 5
from numpy import float32 # needs numpy integer_floor(float32(5.4)) # needs numpy from gmpy2 import mpfr integer_floor(mpfr(5.4))
- sage.arith.misc.integer_trunc(i)[source]¶
Truncate to the integer closer to zero.
EXAMPLES:
sage: from sage.arith.misc import integer_trunc as trunc sage: trunc(-3/2), trunc(-1), trunc(-1/2), trunc(0), trunc(1/2), trunc(1), trunc(3/2) (-1, -1, 0, 0, 0, 1, 1) sage: isinstance(trunc(3/2), Integer) True
>>> from sage.all import * >>> from sage.arith.misc import integer_trunc as trunc >>> trunc(-Integer(3)/Integer(2)), trunc(-Integer(1)), trunc(-Integer(1)/Integer(2)), trunc(Integer(0)), trunc(Integer(1)/Integer(2)), trunc(Integer(1)), trunc(Integer(3)/Integer(2)) (-1, -1, 0, 0, 0, 1, 1) >>> isinstance(trunc(Integer(3)/Integer(2)), Integer) True
from sage.arith.misc import integer_trunc as trunc trunc(-3/2), trunc(-1), trunc(-1/2), trunc(0), trunc(1/2), trunc(1), trunc(3/2) isinstance(trunc(3/2), Integer)
- sage.arith.misc.inverse_mod(a, m)[source]¶
The inverse of the ring element a modulo m.
If no special inverse_mod is defined for the elements, it tries to coerce them into integers and perform the inversion there
sage: inverse_mod(7, 1) 0 sage: inverse_mod(5, 14) 3 sage: inverse_mod(3, -5) 2
>>> from sage.all import * >>> inverse_mod(Integer(7), Integer(1)) 0 >>> inverse_mod(Integer(5), Integer(14)) 3 >>> inverse_mod(Integer(3), -Integer(5)) 2
inverse_mod(7, 1) inverse_mod(5, 14) inverse_mod(3, -5)
Tests with numpy and mpz numbers:
sage: from numpy import int8 # needs numpy sage: inverse_mod(int8(5), int8(14)) # needs numpy 3 sage: from gmpy2 import mpz sage: inverse_mod(mpz(5), mpz(14)) 3
>>> from sage.all import * >>> from numpy import int8 # needs numpy >>> inverse_mod(int8(Integer(5)), int8(Integer(14))) # needs numpy 3 >>> from gmpy2 import mpz >>> inverse_mod(mpz(Integer(5)), mpz(Integer(14))) 3
from numpy import int8 # needs numpy inverse_mod(int8(5), int8(14)) # needs numpy from gmpy2 import mpz inverse_mod(mpz(5), mpz(14))
- sage.arith.misc.is_power_of_two(n)[source]¶
Return whether \(n\) is a power of 2.
INPUT:
n– integer
OUTPUT: boolean
EXAMPLES:
sage: is_power_of_two(1024) True sage: is_power_of_two(1) True sage: is_power_of_two(24) False sage: is_power_of_two(0) False sage: is_power_of_two(-4) False
>>> from sage.all import * >>> is_power_of_two(Integer(1024)) True >>> is_power_of_two(Integer(1)) True >>> is_power_of_two(Integer(24)) False >>> is_power_of_two(Integer(0)) False >>> is_power_of_two(-Integer(4)) False
is_power_of_two(1024) is_power_of_two(1) is_power_of_two(24) is_power_of_two(0) is_power_of_two(-4)
Tests with numpy and gmpy2 numbers:
sage: from numpy import int8 # needs numpy sage: is_power_of_two(int8(16)) # needs numpy True sage: is_power_of_two(int8(24)) # needs numpy False sage: from gmpy2 import mpz sage: is_power_of_two(mpz(16)) True sage: is_power_of_two(mpz(24)) False
>>> from sage.all import * >>> from numpy import int8 # needs numpy >>> is_power_of_two(int8(Integer(16))) # needs numpy True >>> is_power_of_two(int8(Integer(24))) # needs numpy False >>> from gmpy2 import mpz >>> is_power_of_two(mpz(Integer(16))) True >>> is_power_of_two(mpz(Integer(24))) False
from numpy import int8 # needs numpy is_power_of_two(int8(16)) # needs numpy is_power_of_two(int8(24)) # needs numpy from gmpy2 import mpz is_power_of_two(mpz(16)) is_power_of_two(mpz(24))
- sage.arith.misc.is_prime(n)[source]¶
Determine whether \(n\) is a prime element of its parent ring.
INPUT:
n– the object for which to determine primality
Exceptional special cases:
For integers, determine whether \(n\) is a positive prime.
For number fields except \(\QQ\), determine whether \(n\) is a prime element of the maximal order.
ALGORITHM:
For integers, this function uses a provable primality test or a strong pseudo-primality test depending on the global
arithmetic proof flag.EXAMPLES:
sage: is_prime(389) True sage: is_prime(2000) False sage: is_prime(2) True sage: is_prime(-1) False sage: is_prime(1) False sage: is_prime(-2) False
>>> from sage.all import * >>> is_prime(Integer(389)) True >>> is_prime(Integer(2000)) False >>> is_prime(Integer(2)) True >>> is_prime(-Integer(1)) False >>> is_prime(Integer(1)) False >>> is_prime(-Integer(2)) False
is_prime(389) is_prime(2000) is_prime(2) is_prime(-1) is_prime(1) is_prime(-2)
sage: a = 2**2048 + 981 sage: is_prime(a) # not tested - takes ~ 1min sage: proof.arithmetic(False) sage: is_prime(a) # instantaneous! # needs sage.libs.pari True sage: proof.arithmetic(True)
>>> from sage.all import * >>> a = Integer(2)**Integer(2048) + Integer(981) >>> is_prime(a) # not tested - takes ~ 1min >>> proof.arithmetic(False) >>> is_prime(a) # instantaneous! # needs sage.libs.pari True >>> proof.arithmetic(True)
a = 2**2048 + 981 is_prime(a) # not tested - takes ~ 1min proof.arithmetic(False) is_prime(a) # instantaneous! # needs sage.libs.pari proof.arithmetic(True)
- sage.arith.misc.is_prime_power(n, get_data=False)[source]¶
Test whether
nis a positive power of a prime number.This function simply calls the method
Integer.is_prime_power()of Integers.INPUT:
n– integerget_data– if set toTrue, return a pair(p,k)such that this integer equalsp^kinstead ofTrueor(self,0)instead ofFalse
EXAMPLES:
sage: # needs sage.libs.pari sage: is_prime_power(389) True sage: is_prime_power(2000) False sage: is_prime_power(2) True sage: is_prime_power(1024) True sage: is_prime_power(1024, get_data=True) (2, 10)
>>> from sage.all import * >>> # needs sage.libs.pari >>> is_prime_power(Integer(389)) True >>> is_prime_power(Integer(2000)) False >>> is_prime_power(Integer(2)) True >>> is_prime_power(Integer(1024)) True >>> is_prime_power(Integer(1024), get_data=True) (2, 10)
# needs sage.libs.pari is_prime_power(389) is_prime_power(2000) is_prime_power(2) is_prime_power(1024) is_prime_power(1024, get_data=True)
The same results can be obtained with:
sage: # needs sage.libs.pari sage: 389.is_prime_power() True sage: 2000.is_prime_power() False sage: 2.is_prime_power() True sage: 1024.is_prime_power() True sage: 1024.is_prime_power(get_data=True) (2, 10)
>>> from sage.all import * >>> # needs sage.libs.pari >>> Integer(389).is_prime_power() True >>> Integer(2000).is_prime_power() False >>> Integer(2).is_prime_power() True >>> Integer(1024).is_prime_power() True >>> Integer(1024).is_prime_power(get_data=True) (2, 10)
# needs sage.libs.pari 389.is_prime_power() 2000.is_prime_power() 2.is_prime_power() 1024.is_prime_power() 1024.is_prime_power(get_data=True)
- sage.arith.misc.is_pseudoprime(n)[source]¶
Test whether
nis a pseudo-prime.The result is NOT proven correct - this is a pseudo-primality test!.
INPUT:
n– integer
Note
We do not consider negatives of prime numbers as prime.
EXAMPLES:
sage: # needs sage.libs.pari sage: is_pseudoprime(389) True sage: is_pseudoprime(2000) False sage: is_pseudoprime(2) True sage: is_pseudoprime(-1) False sage: factor(-6) -1 * 2 * 3 sage: is_pseudoprime(1) False sage: is_pseudoprime(-2) False
>>> from sage.all import * >>> # needs sage.libs.pari >>> is_pseudoprime(Integer(389)) True >>> is_pseudoprime(Integer(2000)) False >>> is_pseudoprime(Integer(2)) True >>> is_pseudoprime(-Integer(1)) False >>> factor(-Integer(6)) -1 * 2 * 3 >>> is_pseudoprime(Integer(1)) False >>> is_pseudoprime(-Integer(2)) False
# needs sage.libs.pari is_pseudoprime(389) is_pseudoprime(2000) is_pseudoprime(2) is_pseudoprime(-1) factor(-6) is_pseudoprime(1) is_pseudoprime(-2)
- sage.arith.misc.is_pseudoprime_power(n, get_data=False)[source]¶
Test if
nis a power of a pseudoprime.The result is NOT proven correct - this IS a pseudo-primality test!. Note that a prime power is a positive power of a prime number so that 1 is not a prime power.
INPUT:
n– integerget_data– boolean (default:False); instead of a boolean return a pair \((p,k)\) so thatnequals \(p^k\) and \(p\) is a pseudoprime or \((n,0)\) otherwise
EXAMPLES:
sage: # needs sage.libs.pari sage: is_pseudoprime_power(389) True sage: is_pseudoprime_power(2000) False sage: is_pseudoprime_power(2) True sage: is_pseudoprime_power(1024) True sage: is_pseudoprime_power(-1) False sage: is_pseudoprime_power(1) False sage: is_pseudoprime_power(997^100) True
>>> from sage.all import * >>> # needs sage.libs.pari >>> is_pseudoprime_power(Integer(389)) True >>> is_pseudoprime_power(Integer(2000)) False >>> is_pseudoprime_power(Integer(2)) True >>> is_pseudoprime_power(Integer(1024)) True >>> is_pseudoprime_power(-Integer(1)) False >>> is_pseudoprime_power(Integer(1)) False >>> is_pseudoprime_power(Integer(997)**Integer(100)) True
# needs sage.libs.pari is_pseudoprime_power(389) is_pseudoprime_power(2000) is_pseudoprime_power(2) is_pseudoprime_power(1024) is_pseudoprime_power(-1) is_pseudoprime_power(1) is_pseudoprime_power(997^100)
Use of the get_data keyword:
sage: # needs sage.libs.pari sage: is_pseudoprime_power(3^1024, get_data=True) (3, 1024) sage: is_pseudoprime_power(2^256, get_data=True) (2, 256) sage: is_pseudoprime_power(31, get_data=True) (31, 1) sage: is_pseudoprime_power(15, get_data=True) (15, 0)
>>> from sage.all import * >>> # needs sage.libs.pari >>> is_pseudoprime_power(Integer(3)**Integer(1024), get_data=True) (3, 1024) >>> is_pseudoprime_power(Integer(2)**Integer(256), get_data=True) (2, 256) >>> is_pseudoprime_power(Integer(31), get_data=True) (31, 1) >>> is_pseudoprime_power(Integer(15), get_data=True) (15, 0)
# needs sage.libs.pari is_pseudoprime_power(3^1024, get_data=True) is_pseudoprime_power(2^256, get_data=True) is_pseudoprime_power(31, get_data=True) is_pseudoprime_power(15, get_data=True)
Tests with numpy and gmpy2 numbers:
sage: from numpy import int16 # needs numpy sage: is_pseudoprime_power(int16(1024)) # needs numpy sage.libs.pari True sage: from gmpy2 import mpz sage: is_pseudoprime_power(mpz(1024)) True
>>> from sage.all import * >>> from numpy import int16 # needs numpy >>> is_pseudoprime_power(int16(Integer(1024))) # needs numpy sage.libs.pari True >>> from gmpy2 import mpz >>> is_pseudoprime_power(mpz(Integer(1024))) True
from numpy import int16 # needs numpy is_pseudoprime_power(int16(1024)) # needs numpy sage.libs.pari from gmpy2 import mpz is_pseudoprime_power(mpz(1024))
- sage.arith.misc.is_square(n, root=False)[source]¶
Return whether or not \(n\) is square.
If \(n\) is a square also return the square root. If \(n\) is not square, also return
None.INPUT:
n– integerroot– whether or not to also return a square root (default:False)
OUTPUT:
bool– whether or not a squareobject– (optional) an actual square if found, andNoneotherwise
EXAMPLES:
sage: is_square(2) False sage: is_square(4) True sage: is_square(2.2) True sage: is_square(-2.2) False sage: is_square(CDF(-2.2)) # needs sage.rings.complex_double True sage: is_square((x-1)^2) # needs sage.symbolic Traceback (most recent call last): ... NotImplementedError: is_square() not implemented for non-constant or relational elements of Symbolic Ring
>>> from sage.all import * >>> is_square(Integer(2)) False >>> is_square(Integer(4)) True >>> is_square(RealNumber('2.2')) True >>> is_square(-RealNumber('2.2')) False >>> is_square(CDF(-RealNumber('2.2'))) # needs sage.rings.complex_double True >>> is_square((x-Integer(1))**Integer(2)) # needs sage.symbolic Traceback (most recent call last): ... NotImplementedError: is_square() not implemented for non-constant or relational elements of Symbolic Ring
is_square(2) is_square(4) is_square(2.2) is_square(-2.2) is_square(CDF(-2.2)) # needs sage.rings.complex_double is_square((x-1)^2) # needs sage.symbolic
sage: is_square(4, True) (True, 2)
>>> from sage.all import * >>> is_square(Integer(4), True) (True, 2)
is_square(4, True)
Tests with numpy and gmpy2 numbers:
sage: from numpy import int8 # needs numpy sage: is_square(int8(4)) # needs numpy True sage: from gmpy2 import mpz sage: is_square(mpz(4)) True
>>> from sage.all import * >>> from numpy import int8 # needs numpy >>> is_square(int8(Integer(4))) # needs numpy True >>> from gmpy2 import mpz >>> is_square(mpz(Integer(4))) True
from numpy import int8 # needs numpy is_square(int8(4)) # needs numpy from gmpy2 import mpz is_square(mpz(4))
Tests with Polynomial:
sage: R.<v> = LaurentPolynomialRing(QQ, 'v') sage: H = IwahoriHeckeAlgebra('A3', v**2) # needs sage.combinat sage.modules sage: R.<a,b,c,d> = QQ[] sage: p = a*b + c*d*a*d*a + 5 sage: is_square(p**2) True
>>> from sage.all import * >>> R = LaurentPolynomialRing(QQ, 'v', names=('v',)); (v,) = R._first_ngens(1) >>> H = IwahoriHeckeAlgebra('A3', v**Integer(2)) # needs sage.combinat sage.modules >>> R = QQ['a, b, c, d']; (a, b, c, d,) = R._first_ngens(4) >>> p = a*b + c*d*a*d*a + Integer(5) >>> is_square(p**Integer(2)) True
R.<v> = LaurentPolynomialRing(QQ, 'v') H = IwahoriHeckeAlgebra('A3', v**2) # needs sage.combinat sage.modules R.<a,b,c,d> = QQ[] p = a*b + c*d*a*d*a + 5 is_square(p**2)
- sage.arith.misc.is_squarefree(n)[source]¶
Test whether
nis square free.EXAMPLES:
sage: is_squarefree(100) # needs sage.libs.pari False sage: is_squarefree(101) # needs sage.libs.pari True sage: R = ZZ['x'] sage: x = R.gen() sage: is_squarefree((x^2+x+1) * (x-2)) # needs sage.libs.pari True sage: is_squarefree((x-1)**2 * (x-3)) # needs sage.libs.pari False sage: # needs sage.rings.number_field sage.symbolic sage: O = ZZ[sqrt(-1)] sage: I = O.gen(1) sage: is_squarefree(I + 1) True sage: is_squarefree(O(2)) False sage: O(2).factor() (-I) * (I + 1)^2
>>> from sage.all import * >>> is_squarefree(Integer(100)) # needs sage.libs.pari False >>> is_squarefree(Integer(101)) # needs sage.libs.pari True >>> R = ZZ['x'] >>> x = R.gen() >>> is_squarefree((x**Integer(2)+x+Integer(1)) * (x-Integer(2))) # needs sage.libs.pari True >>> is_squarefree((x-Integer(1))**Integer(2) * (x-Integer(3))) # needs sage.libs.pari False >>> # needs sage.rings.number_field sage.symbolic >>> O = ZZ[sqrt(-Integer(1))] >>> I = O.gen(Integer(1)) >>> is_squarefree(I + Integer(1)) True >>> is_squarefree(O(Integer(2))) False >>> O(Integer(2)).factor() (-I) * (I + 1)^2
is_squarefree(100) # needs sage.libs.pari is_squarefree(101) # needs sage.libs.pari R = ZZ['x'] x = R.gen() is_squarefree((x^2+x+1) * (x-2)) # needs sage.libs.pari is_squarefree((x-1)**2 * (x-3)) # needs sage.libs.pari # needs sage.rings.number_field sage.symbolic O = ZZ[sqrt(-1)] I = O.gen(1) is_squarefree(I + 1) is_squarefree(O(2)) O(2).factor()
This method fails on domains which are not Unique Factorization Domains:
sage: O = ZZ[sqrt(-5)] # needs sage.rings.number_field sage.symbolic sage: a = O.gen(1) # needs sage.rings.number_field sage.symbolic sage: is_squarefree(a - 3) # needs sage.rings.number_field sage.symbolic Traceback (most recent call last): ... ArithmeticError: non-principal ideal in factorization
>>> from sage.all import * >>> O = ZZ[sqrt(-Integer(5))] # needs sage.rings.number_field sage.symbolic >>> a = O.gen(Integer(1)) # needs sage.rings.number_field sage.symbolic >>> is_squarefree(a - Integer(3)) # needs sage.rings.number_field sage.symbolic Traceback (most recent call last): ... ArithmeticError: non-principal ideal in factorization
O = ZZ[sqrt(-5)] # needs sage.rings.number_field sage.symbolic a = O.gen(1) # needs sage.rings.number_field sage.symbolic is_squarefree(a - 3) # needs sage.rings.number_field sage.symbolic
Tests with numpy and gmpy2 numbers:
sage: # needs sage.libs.pari sage: from numpy import int8 # needs numpy sage: is_squarefree(int8(100)) # needs numpy False sage: is_squarefree(int8(101)) # needs numpy True sage: from gmpy2 import mpz sage: is_squarefree(mpz(100)) False sage: is_squarefree(mpz(101)) True
>>> from sage.all import * >>> # needs sage.libs.pari >>> from numpy import int8 # needs numpy >>> is_squarefree(int8(Integer(100))) # needs numpy False >>> is_squarefree(int8(Integer(101))) # needs numpy True >>> from gmpy2 import mpz >>> is_squarefree(mpz(Integer(100))) False >>> is_squarefree(mpz(Integer(101))) True
# needs sage.libs.pari from numpy import int8 # needs numpy is_squarefree(int8(100)) # needs numpy is_squarefree(int8(101)) # needs numpy from gmpy2 import mpz is_squarefree(mpz(100)) is_squarefree(mpz(101))
- sage.arith.misc.jacobi_symbol(a, b)[source]¶
The Jacobi symbol of integers \(a\) and \(b\), where \(b\) is odd.
Note
The
kronecker_symbol()command extends the Jacobi symbol to all integers \(b\).If
\(b = p_1^{e_1} * ... * p_r^{e_r}\)
then
\((a|b) = (a|p_1)^{e_1} ... (a|p_r)^{e_r}\)
where \((a|p_j)\) are Legendre Symbols.
INPUT:
a– integerb– odd integer
EXAMPLES:
sage: jacobi_symbol(10,777) -1 sage: jacobi_symbol(10,5) 0 sage: jacobi_symbol(10,2) Traceback (most recent call last): ... ValueError: second input must be odd, 2 is not odd
>>> from sage.all import * >>> jacobi_symbol(Integer(10),Integer(777)) -1 >>> jacobi_symbol(Integer(10),Integer(5)) 0 >>> jacobi_symbol(Integer(10),Integer(2)) Traceback (most recent call last): ... ValueError: second input must be odd, 2 is not odd
jacobi_symbol(10,777) jacobi_symbol(10,5) jacobi_symbol(10,2)
Tests with numpy and gmpy2 numbers:
sage: from numpy import int16 # needs numpy sage: jacobi_symbol(int16(10), int16(777)) # needs numpy -1 sage: from gmpy2 import mpz sage: jacobi_symbol(mpz(10),mpz(777)) -1
>>> from sage.all import * >>> from numpy import int16 # needs numpy >>> jacobi_symbol(int16(Integer(10)), int16(Integer(777))) # needs numpy -1 >>> from gmpy2 import mpz >>> jacobi_symbol(mpz(Integer(10)),mpz(Integer(777))) -1
from numpy import int16 # needs numpy jacobi_symbol(int16(10), int16(777)) # needs numpy from gmpy2 import mpz jacobi_symbol(mpz(10),mpz(777))
- sage.arith.misc.kronecker(x, y)[source]¶
The Kronecker symbol \((x|y)\).
INPUT:
x– integery– integer
OUTPUT: integer
EXAMPLES:
sage: kronecker_symbol(13,21) -1 sage: kronecker_symbol(101,4) 1
>>> from sage.all import * >>> kronecker_symbol(Integer(13),Integer(21)) -1 >>> kronecker_symbol(Integer(101),Integer(4)) 1
kronecker_symbol(13,21) kronecker_symbol(101,4)
This is also available as
kronecker():sage: kronecker(3,5) -1 sage: kronecker(3,15) 0 sage: kronecker(2,15) 1 sage: kronecker(-2,15) -1 sage: kronecker(2/3,5) 1
>>> from sage.all import * >>> kronecker(Integer(3),Integer(5)) -1 >>> kronecker(Integer(3),Integer(15)) 0 >>> kronecker(Integer(2),Integer(15)) 1 >>> kronecker(-Integer(2),Integer(15)) -1 >>> kronecker(Integer(2)/Integer(3),Integer(5)) 1
kronecker(3,5) kronecker(3,15) kronecker(2,15) kronecker(-2,15) kronecker(2/3,5)
Tests with numpy and gmpy2 numbers:
sage: from numpy import int8 # needs numpy sage: kronecker_symbol(int8(13),int8(21)) # needs numpy -1 sage: from gmpy2 import mpz sage: kronecker_symbol(mpz(13),mpz(21)) -1
>>> from sage.all import * >>> from numpy import int8 # needs numpy >>> kronecker_symbol(int8(Integer(13)),int8(Integer(21))) # needs numpy -1 >>> from gmpy2 import mpz >>> kronecker_symbol(mpz(Integer(13)),mpz(Integer(21))) -1
from numpy import int8 # needs numpy kronecker_symbol(int8(13),int8(21)) # needs numpy from gmpy2 import mpz kronecker_symbol(mpz(13),mpz(21))
- sage.arith.misc.kronecker_symbol(x, y)[source]¶
The Kronecker symbol \((x|y)\).
INPUT:
x– integery– integer
OUTPUT: integer
EXAMPLES:
sage: kronecker_symbol(13,21) -1 sage: kronecker_symbol(101,4) 1
>>> from sage.all import * >>> kronecker_symbol(Integer(13),Integer(21)) -1 >>> kronecker_symbol(Integer(101),Integer(4)) 1
kronecker_symbol(13,21) kronecker_symbol(101,4)
This is also available as
kronecker():sage: kronecker(3,5) -1 sage: kronecker(3,15) 0 sage: kronecker(2,15) 1 sage: kronecker(-2,15) -1 sage: kronecker(2/3,5) 1
>>> from sage.all import * >>> kronecker(Integer(3),Integer(5)) -1 >>> kronecker(Integer(3),Integer(15)) 0 >>> kronecker(Integer(2),Integer(15)) 1 >>> kronecker(-Integer(2),Integer(15)) -1 >>> kronecker(Integer(2)/Integer(3),Integer(5)) 1
kronecker(3,5) kronecker(3,15) kronecker(2,15) kronecker(-2,15) kronecker(2/3,5)
Tests with numpy and gmpy2 numbers:
sage: from numpy import int8 # needs numpy sage: kronecker_symbol(int8(13),int8(21)) # needs numpy -1 sage: from gmpy2 import mpz sage: kronecker_symbol(mpz(13),mpz(21)) -1
>>> from sage.all import * >>> from numpy import int8 # needs numpy >>> kronecker_symbol(int8(Integer(13)),int8(Integer(21))) # needs numpy -1 >>> from gmpy2 import mpz >>> kronecker_symbol(mpz(Integer(13)),mpz(Integer(21))) -1
from numpy import int8 # needs numpy kronecker_symbol(int8(13),int8(21)) # needs numpy from gmpy2 import mpz kronecker_symbol(mpz(13),mpz(21))
- sage.arith.misc.legendre_symbol(x, p)[source]¶
The Legendre symbol \((x|p)\), for \(p\) prime.
Note
The
kronecker_symbol()command extends the Legendre symbol to composite moduli and \(p=2\).INPUT:
x– integerp– odd prime number
EXAMPLES:
sage: legendre_symbol(2,3) -1 sage: legendre_symbol(1,3) 1 sage: legendre_symbol(1,2) Traceback (most recent call last): ... ValueError: p must be odd sage: legendre_symbol(2,15) Traceback (most recent call last): ... ValueError: p must be a prime sage: kronecker_symbol(2,15) 1 sage: legendre_symbol(2/3,7) -1
>>> from sage.all import * >>> legendre_symbol(Integer(2),Integer(3)) -1 >>> legendre_symbol(Integer(1),Integer(3)) 1 >>> legendre_symbol(Integer(1),Integer(2)) Traceback (most recent call last): ... ValueError: p must be odd >>> legendre_symbol(Integer(2),Integer(15)) Traceback (most recent call last): ... ValueError: p must be a prime >>> kronecker_symbol(Integer(2),Integer(15)) 1 >>> legendre_symbol(Integer(2)/Integer(3),Integer(7)) -1
legendre_symbol(2,3) legendre_symbol(1,3) legendre_symbol(1,2) legendre_symbol(2,15) kronecker_symbol(2,15) legendre_symbol(2/3,7)
Tests with numpy and gmpy2 numbers:
sage: from numpy import int8 # needs numpy sage: legendre_symbol(int8(2), int8(3)) # needs numpy -1 sage: from gmpy2 import mpz sage: legendre_symbol(mpz(2),mpz(3)) -1
>>> from sage.all import * >>> from numpy import int8 # needs numpy >>> legendre_symbol(int8(Integer(2)), int8(Integer(3))) # needs numpy -1 >>> from gmpy2 import mpz >>> legendre_symbol(mpz(Integer(2)),mpz(Integer(3))) -1
from numpy import int8 # needs numpy legendre_symbol(int8(2), int8(3)) # needs numpy from gmpy2 import mpz legendre_symbol(mpz(2),mpz(3))
- sage.arith.misc.mqrr_rational_reconstruction(u, m, T)[source]¶
Maximal Quotient Rational Reconstruction.
For research purposes only - this is pure Python, so slow.
INPUT:
u,m,T– integers such that \(m > u \ge 0\), \(T > 0\)
OUTPUT:
Either integers \(n,d\) such that \(d>0\), \(\mathop{\mathrm{gcd}}(n,d)=1\), \(n/d=u \bmod m\), and \(T \cdot d \cdot |n| < m\), or
None.Reference: Monagan, Maximal Quotient Rational Reconstruction: An Almost Optimal Algorithm for Rational Reconstruction (page 11)
This algorithm is probabilistic.
EXAMPLES:
sage: mqrr_rational_reconstruction(21, 3100, 13) (21, 1)
>>> from sage.all import * >>> mqrr_rational_reconstruction(Integer(21), Integer(3100), Integer(13)) (21, 1)
mqrr_rational_reconstruction(21, 3100, 13)
Tests with numpy and gmpy2 numbers:
sage: from numpy import int16 # needs numpy sage: mqrr_rational_reconstruction(int16(21), int16(3100), int16(13)) # needs numpy (21, 1) sage: from gmpy2 import mpz sage: mqrr_rational_reconstruction(mpz(21), mpz(3100), mpz(13)) (21, 1)
>>> from sage.all import * >>> from numpy import int16 # needs numpy >>> mqrr_rational_reconstruction(int16(Integer(21)), int16(Integer(3100)), int16(Integer(13))) # needs numpy (21, 1) >>> from gmpy2 import mpz >>> mqrr_rational_reconstruction(mpz(Integer(21)), mpz(Integer(3100)), mpz(Integer(13))) (21, 1)
from numpy import int16 # needs numpy mqrr_rational_reconstruction(int16(21), int16(3100), int16(13)) # needs numpy from gmpy2 import mpz mqrr_rational_reconstruction(mpz(21), mpz(3100), mpz(13))
- sage.arith.misc.multinomial(*ks)[source]¶
Return the multinomial coefficient.
INPUT:
either an arbitrary number of integer arguments \(k_1,\dots,k_n\)
or an iterable (e.g. a list) of integers \([k_1,\dots,k_n]\)
OUTPUT: the integer:
\[\binom{k_1 + \cdots + k_n}{k_1, \cdots, k_n} =\frac{\left(\sum_{i=1}^n k_i\right)!}{\prod_{i=1}^n k_i!} = \prod_{i=1}^n \binom{\sum_{j=1}^i k_j}{k_i}\]EXAMPLES:
sage: multinomial(0, 0, 2, 1, 0, 0) 3 sage: multinomial([0, 0, 2, 1, 0, 0]) 3 sage: multinomial(3, 2) 10 sage: multinomial(2^30, 2, 1) 618970023101454657175683075 sage: multinomial([2^30, 2, 1]) 618970023101454657175683075 sage: multinomial(Composition([1, 3])) 4 sage: multinomial(Partition([4, 2])) # needs sage.combinat 15
>>> from sage.all import * >>> multinomial(Integer(0), Integer(0), Integer(2), Integer(1), Integer(0), Integer(0)) 3 >>> multinomial([Integer(0), Integer(0), Integer(2), Integer(1), Integer(0), Integer(0)]) 3 >>> multinomial(Integer(3), Integer(2)) 10 >>> multinomial(Integer(2)**Integer(30), Integer(2), Integer(1)) 618970023101454657175683075 >>> multinomial([Integer(2)**Integer(30), Integer(2), Integer(1)]) 618970023101454657175683075 >>> multinomial(Composition([Integer(1), Integer(3)])) 4 >>> multinomial(Partition([Integer(4), Integer(2)])) # needs sage.combinat 15
multinomial(0, 0, 2, 1, 0, 0) multinomial([0, 0, 2, 1, 0, 0]) multinomial(3, 2) multinomial(2^30, 2, 1) multinomial([2^30, 2, 1]) multinomial(Composition([1, 3])) multinomial(Partition([4, 2])) # needs sage.combinat
AUTHORS:
Gabriel Ebner
- sage.arith.misc.multinomial_coefficients(m, n)[source]¶
Return a dictionary containing pairs \(\{(k_1, k_2, ..., k_m) : C_{k, n}\}\) where \(C_{k, n}\) are multinomial coefficients such that \(n = k_1 + k_2 + ...+ k_m\).
INPUT:
m– integern– integer
OUTPUT: dictionary
EXAMPLES:
sage: sorted(multinomial_coefficients(2, 5).items()) [((0, 5), 1), ((1, 4), 5), ((2, 3), 10), ((3, 2), 10), ((4, 1), 5), ((5, 0), 1)]
>>> from sage.all import * >>> sorted(multinomial_coefficients(Integer(2), Integer(5)).items()) [((0, 5), 1), ((1, 4), 5), ((2, 3), 10), ((3, 2), 10), ((4, 1), 5), ((5, 0), 1)]
sorted(multinomial_coefficients(2, 5).items())
Notice that these are the coefficients of \((x+y)^5\):
sage: R.<x,y> = QQ[] sage: (x+y)^5 x^5 + 5*x^4*y + 10*x^3*y^2 + 10*x^2*y^3 + 5*x*y^4 + y^5
>>> from sage.all import * >>> R = QQ['x, y']; (x, y,) = R._first_ngens(2) >>> (x+y)**Integer(5) x^5 + 5*x^4*y + 10*x^3*y^2 + 10*x^2*y^3 + 5*x*y^4 + y^5
R.<x,y> = QQ[] (x+y)^5
sage: sorted(multinomial_coefficients(3, 2).items()) [((0, 0, 2), 1), ((0, 1, 1), 2), ((0, 2, 0), 1), ((1, 0, 1), 2), ((1, 1, 0), 2), ((2, 0, 0), 1)]
>>> from sage.all import * >>> sorted(multinomial_coefficients(Integer(3), Integer(2)).items()) [((0, 0, 2), 1), ((0, 1, 1), 2), ((0, 2, 0), 1), ((1, 0, 1), 2), ((1, 1, 0), 2), ((2, 0, 0), 1)]
sorted(multinomial_coefficients(3, 2).items())
ALGORITHM: The algorithm we implement for computing the multinomial coefficients is based on the following result:
\[\binom{n}{k_1, \cdots, k_m} = \frac{k_1+1}{n-k_1}\sum_{i=2}^m \binom{n}{k_1+1, \cdots, k_i-1, \cdots}\]e.g.:
sage: k = (2, 4, 1, 0, 2, 6, 0, 0, 3, 5, 7, 1) # random value sage: n = sum(k) sage: s = 0 sage: for i in range(1, len(k)): ....: ki = list(k) ....: ki[0] += 1 ....: ki[i] -= 1 ....: s += multinomial(n, *ki) sage: multinomial(n, *k) == (k[0] + 1) / (n - k[0]) * s True
>>> from sage.all import * >>> k = (Integer(2), Integer(4), Integer(1), Integer(0), Integer(2), Integer(6), Integer(0), Integer(0), Integer(3), Integer(5), Integer(7), Integer(1)) # random value >>> n = sum(k) >>> s = Integer(0) >>> for i in range(Integer(1), len(k)): ... ki = list(k) ... ki[Integer(0)] += Integer(1) ... ki[i] -= Integer(1) ... s += multinomial(n, *ki) >>> multinomial(n, *k) == (k[Integer(0)] + Integer(1)) / (n - k[Integer(0)]) * s True
k = (2, 4, 1, 0, 2, 6, 0, 0, 3, 5, 7, 1) # random value n = sum(k) s = 0 for i in range(1, len(k)): ki = list(k) ki[0] += 1 ki[i] -= 1 s += multinomial(n, *ki) multinomial(n, *k) == (k[0] + 1) / (n - k[0]) * s
- sage.arith.misc.next_prime(n, proof=None)[source]¶
The next prime greater than the integer \(n\). If \(n\) is prime, then this function does not return \(n\), but the next prime after \(n\). If the optional argument proof is
False, this function only returns a pseudo-prime, as defined by the PARI nextprime function. If it isNone, uses the global default (seesage.structure.proof.proof)INPUT:
n– integerproof– boolean orNone(default:None)
EXAMPLES:
sage: # needs sage.libs.pari sage: next_prime(-100) 2 sage: next_prime(1) 2 sage: next_prime(2) 3 sage: next_prime(3) 5 sage: next_prime(4) 5
>>> from sage.all import * >>> # needs sage.libs.pari >>> next_prime(-Integer(100)) 2 >>> next_prime(Integer(1)) 2 >>> next_prime(Integer(2)) 3 >>> next_prime(Integer(3)) 5 >>> next_prime(Integer(4)) 5
# needs sage.libs.pari next_prime(-100) next_prime(1) next_prime(2) next_prime(3) next_prime(4)
Notice that the next_prime(5) is not 5 but 7.
sage: next_prime(5) # needs sage.libs.pari 7 sage: next_prime(2004) # needs sage.libs.pari 2011
>>> from sage.all import * >>> next_prime(Integer(5)) # needs sage.libs.pari 7 >>> next_prime(Integer(2004)) # needs sage.libs.pari 2011
next_prime(5) # needs sage.libs.pari next_prime(2004) # needs sage.libs.pari
- sage.arith.misc.next_prime_power(n)[source]¶
Return the smallest prime power greater than
n.Note that if
nis a prime power, then this function does not returnn, but the next prime power aftern.This function just calls the method
Integer.next_prime_power()of Integers.See also
EXAMPLES:
sage: # needs sage.libs.pari sage: next_prime_power(1) 2 sage: next_prime_power(2) 3 sage: next_prime_power(10) 11 sage: next_prime_power(7) 8 sage: next_prime_power(99) 101
>>> from sage.all import * >>> # needs sage.libs.pari >>> next_prime_power(Integer(1)) 2 >>> next_prime_power(Integer(2)) 3 >>> next_prime_power(Integer(10)) 11 >>> next_prime_power(Integer(7)) 8 >>> next_prime_power(Integer(99)) 101
# needs sage.libs.pari next_prime_power(1) next_prime_power(2) next_prime_power(10) next_prime_power(7) next_prime_power(99)
The same results can be obtained with:
sage: 1.next_prime_power() 2 sage: 2.next_prime_power() 3 sage: 10.next_prime_power() 11
>>> from sage.all import * >>> Integer(1).next_prime_power() 2 >>> Integer(2).next_prime_power() 3 >>> Integer(10).next_prime_power() 11
1.next_prime_power() 2.next_prime_power() 10.next_prime_power()
Note that \(2\) is the smallest prime power:
sage: next_prime_power(-10) 2 sage: next_prime_power(0) 2
>>> from sage.all import * >>> next_prime_power(-Integer(10)) 2 >>> next_prime_power(Integer(0)) 2
next_prime_power(-10) next_prime_power(0)
- sage.arith.misc.next_probable_prime(n)[source]¶
Return the next probable prime after self, as determined by PARI.
INPUT:
n– integer
EXAMPLES:
sage: # needs sage.libs.pari sage: next_probable_prime(-100) 2 sage: next_probable_prime(19) 23 sage: next_probable_prime(int(999999999)) 1000000007 sage: next_probable_prime(2^768) 1552518092300708935148979488462502555256886017116696611139052038026050952686376886330878408828646477950487730697131073206171580044114814391444287275041181139204454976020849905550265285631598444825262999193716468750892846853816058039
>>> from sage.all import * >>> # needs sage.libs.pari >>> next_probable_prime(-Integer(100)) 2 >>> next_probable_prime(Integer(19)) 23 >>> next_probable_prime(int(Integer(999999999))) 1000000007 >>> next_probable_prime(Integer(2)**Integer(768)) 1552518092300708935148979488462502555256886017116696611139052038026050952686376886330878408828646477950487730697131073206171580044114814391444287275041181139204454976020849905550265285631598444825262999193716468750892846853816058039
# needs sage.libs.pari next_probable_prime(-100) next_probable_prime(19) next_probable_prime(int(999999999)) next_probable_prime(2^768)
- sage.arith.misc.nth_prime(n)[source]¶
Return the \(n\)-th prime number (1-indexed, so that 2 is the 1st prime).
INPUT:
n– positive integer
OUTPUT: the \(n\)-th prime number
EXAMPLES:
sage: nth_prime(3) # needs sage.libs.pari 5 sage: nth_prime(10) # needs sage.libs.pari 29 sage: nth_prime(10^7) # needs sage.libs.pari 179424673
>>> from sage.all import * >>> nth_prime(Integer(3)) # needs sage.libs.pari 5 >>> nth_prime(Integer(10)) # needs sage.libs.pari 29 >>> nth_prime(Integer(10)**Integer(7)) # needs sage.libs.pari 179424673
nth_prime(3) # needs sage.libs.pari nth_prime(10) # needs sage.libs.pari nth_prime(10^7) # needs sage.libs.pari
sage: nth_prime(0) Traceback (most recent call last): ... ValueError: nth prime meaningless for nonpositive n (=0)
>>> from sage.all import * >>> nth_prime(Integer(0)) Traceback (most recent call last): ... ValueError: nth prime meaningless for nonpositive n (=0)
nth_prime(0)
- sage.arith.misc.number_of_divisors(n)[source]¶
Return the number of divisors of the integer \(n\).
INPUT:
n– nonzero integer
OUTPUT: integer; the number of divisors of \(n\)
EXAMPLES:
sage: number_of_divisors(100) # needs sage.libs.pari 9 sage: number_of_divisors(-720) # needs sage.libs.pari 30
>>> from sage.all import * >>> number_of_divisors(Integer(100)) # needs sage.libs.pari 9 >>> number_of_divisors(-Integer(720)) # needs sage.libs.pari 30
number_of_divisors(100) # needs sage.libs.pari number_of_divisors(-720) # needs sage.libs.pari
Tests with numpy and gmpy2 numbers:
sage: from numpy import int8 # needs numpy sage: number_of_divisors(int8(100)) # needs numpy sage.libs.pari 9 sage: from gmpy2 import mpz sage: number_of_divisors(mpz(100)) # needs sage.libs.pari 9
>>> from sage.all import * >>> from numpy import int8 # needs numpy >>> number_of_divisors(int8(Integer(100))) # needs numpy sage.libs.pari 9 >>> from gmpy2 import mpz >>> number_of_divisors(mpz(Integer(100))) # needs sage.libs.pari 9
from numpy import int8 # needs numpy number_of_divisors(int8(100)) # needs numpy sage.libs.pari from gmpy2 import mpz number_of_divisors(mpz(100)) # needs sage.libs.pari
- sage.arith.misc.odd_part(n)[source]¶
The odd part of the integer \(n\). This is \(n / 2^v\), where \(v = \mathrm{valuation}(n,2)\).
EXAMPLES:
sage: odd_part(5) 5 sage: odd_part(4) 1 sage: odd_part(factorial(31)) 122529844256906551386796875
>>> from sage.all import * >>> odd_part(Integer(5)) 5 >>> odd_part(Integer(4)) 1 >>> odd_part(factorial(Integer(31))) 122529844256906551386796875
odd_part(5) odd_part(4) odd_part(factorial(31))
Tests with numpy and gmpy2 numbers:
sage: from numpy import int8 # needs numpy sage: odd_part(int8(5)) # needs numpy 5 sage: from gmpy2 import mpz sage: odd_part(mpz(5)) 5
>>> from sage.all import * >>> from numpy import int8 # needs numpy >>> odd_part(int8(Integer(5))) # needs numpy 5 >>> from gmpy2 import mpz >>> odd_part(mpz(Integer(5))) 5
from numpy import int8 # needs numpy odd_part(int8(5)) # needs numpy from gmpy2 import mpz odd_part(mpz(5))
- sage.arith.misc.power_mod(a, n, m)[source]¶
Return the \(n\)-th power of \(a\) modulo \(m\), where \(a\) and \(m\) are elements of a ring that implements the modulo operator
%.ALGORITHM: square-and-multiply
EXAMPLES:
sage: power_mod(2, 388, 389) 1 sage: power_mod(2, 390, 391) 285 sage: power_mod(2, -1, 7) 4 sage: power_mod(11, 1, 7) 4
>>> from sage.all import * >>> power_mod(Integer(2), Integer(388), Integer(389)) 1 >>> power_mod(Integer(2), Integer(390), Integer(391)) 285 >>> power_mod(Integer(2), -Integer(1), Integer(7)) 4 >>> power_mod(Integer(11), Integer(1), Integer(7)) 4
power_mod(2, 388, 389) power_mod(2, 390, 391) power_mod(2, -1, 7) power_mod(11, 1, 7)
This function works for fairly general rings:
sage: R.<x> = ZZ[] sage: power_mod(3*x, 10, 7) 4*x^10 sage: power_mod(-3*x^2 + 4, 7, 2*x^3 - 5) x^14 + x^8 + x^6 + x^3 + 962509*x^2 - 791910*x - 698281
>>> from sage.all import * >>> R = ZZ['x']; (x,) = R._first_ngens(1) >>> power_mod(Integer(3)*x, Integer(10), Integer(7)) 4*x^10 >>> power_mod(-Integer(3)*x**Integer(2) + Integer(4), Integer(7), Integer(2)*x**Integer(3) - Integer(5)) x^14 + x^8 + x^6 + x^3 + 962509*x^2 - 791910*x - 698281
R.<x> = ZZ[] power_mod(3*x, 10, 7) power_mod(-3*x^2 + 4, 7, 2*x^3 - 5)
- sage.arith.misc.previous_prime(n)[source]¶
The largest prime < n. The result is provably correct. If n <= 1, this function raises a
ValueError.EXAMPLES:
sage: # needs sage.libs.pari sage: previous_prime(10) 7 sage: previous_prime(7) 5 sage: previous_prime(8) 7 sage: previous_prime(7) 5 sage: previous_prime(5) 3 sage: previous_prime(3) 2 sage: previous_prime(2) Traceback (most recent call last): ... ValueError: no previous prime sage: previous_prime(1) Traceback (most recent call last): ... ValueError: no previous prime sage: previous_prime(-20) Traceback (most recent call last): ... ValueError: no previous prime
>>> from sage.all import * >>> # needs sage.libs.pari >>> previous_prime(Integer(10)) 7 >>> previous_prime(Integer(7)) 5 >>> previous_prime(Integer(8)) 7 >>> previous_prime(Integer(7)) 5 >>> previous_prime(Integer(5)) 3 >>> previous_prime(Integer(3)) 2 >>> previous_prime(Integer(2)) Traceback (most recent call last): ... ValueError: no previous prime >>> previous_prime(Integer(1)) Traceback (most recent call last): ... ValueError: no previous prime >>> previous_prime(-Integer(20)) Traceback (most recent call last): ... ValueError: no previous prime
# needs sage.libs.pari previous_prime(10) previous_prime(7) previous_prime(8) previous_prime(7) previous_prime(5) previous_prime(3) previous_prime(2) previous_prime(1) previous_prime(-20)
- sage.arith.misc.previous_prime_power(n)[source]¶
Return the largest prime power smaller than
n.The result is provably correct. If
nis smaller or equal than2this function raises an error.This function simply call the method
Integer.previous_prime_power()of Integers.See also
EXAMPLES:
sage: # needs sage.libs.pari sage: previous_prime_power(3) 2 sage: previous_prime_power(10) 9 sage: previous_prime_power(7) 5 sage: previous_prime_power(127) 125
>>> from sage.all import * >>> # needs sage.libs.pari >>> previous_prime_power(Integer(3)) 2 >>> previous_prime_power(Integer(10)) 9 >>> previous_prime_power(Integer(7)) 5 >>> previous_prime_power(Integer(127)) 125
# needs sage.libs.pari previous_prime_power(3) previous_prime_power(10) previous_prime_power(7) previous_prime_power(127)
The same results can be obtained with:
sage: # needs sage.libs.pari sage: 3.previous_prime_power() 2 sage: 10.previous_prime_power() 9 sage: 7.previous_prime_power() 5 sage: 127.previous_prime_power() 125
>>> from sage.all import * >>> # needs sage.libs.pari >>> Integer(3).previous_prime_power() 2 >>> Integer(10).previous_prime_power() 9 >>> Integer(7).previous_prime_power() 5 >>> Integer(127).previous_prime_power() 125
# needs sage.libs.pari 3.previous_prime_power() 10.previous_prime_power() 7.previous_prime_power() 127.previous_prime_power()
Input less than or equal to \(2\) raises errors:
sage: previous_prime_power(2) Traceback (most recent call last): ... ValueError: no prime power less than 2 sage: previous_prime_power(-10) Traceback (most recent call last): ... ValueError: no prime power less than 2
>>> from sage.all import * >>> previous_prime_power(Integer(2)) Traceback (most recent call last): ... ValueError: no prime power less than 2 >>> previous_prime_power(-Integer(10)) Traceback (most recent call last): ... ValueError: no prime power less than 2
previous_prime_power(2) previous_prime_power(-10)
sage: n = previous_prime_power(2^16 - 1) # needs sage.libs.pari sage: while is_prime(n): # needs sage.libs.pari ....: n = previous_prime_power(n) sage: factor(n) # needs sage.libs.pari 251^2
>>> from sage.all import * >>> n = previous_prime_power(Integer(2)**Integer(16) - Integer(1)) # needs sage.libs.pari >>> while is_prime(n): # needs sage.libs.pari ... n = previous_prime_power(n) >>> factor(n) # needs sage.libs.pari 251^2
n = previous_prime_power(2^16 - 1) # needs sage.libs.pari while is_prime(n): # needs sage.libs.pari n = previous_prime_power(n) factor(n) # needs sage.libs.pari
- sage.arith.misc.prime_divisors(n)[source]¶
Return the list of prime divisors (up to units) of this element of a unique factorization domain.
INPUT:
n– any object which can be decomposed into prime factors
OUTPUT:
A list of prime factors of
n. For integers, this list is sorted in increasing order.EXAMPLES:
Prime divisors of positive integers:
sage: prime_divisors(1) [] sage: prime_divisors(100) [2, 5] sage: prime_divisors(2004) [2, 3, 167]
>>> from sage.all import * >>> prime_divisors(Integer(1)) [] >>> prime_divisors(Integer(100)) [2, 5] >>> prime_divisors(Integer(2004)) [2, 3, 167]
prime_divisors(1) prime_divisors(100) prime_divisors(2004)
If
nis negative, we do not include -1 among the prime divisors, since -1 is not a prime number:sage: prime_divisors(-100) [2, 5]
>>> from sage.all import * >>> prime_divisors(-Integer(100)) [2, 5]
prime_divisors(-100)
For polynomials we get all irreducible factors:
sage: R.<x> = PolynomialRing(QQ) sage: prime_divisors(x^12 - 1) # needs sage.libs.pari [x - 1, x + 1, x^2 - x + 1, x^2 + 1, x^2 + x + 1, x^4 - x^2 + 1]
>>> from sage.all import * >>> R = PolynomialRing(QQ, names=('x',)); (x,) = R._first_ngens(1) >>> prime_divisors(x**Integer(12) - Integer(1)) # needs sage.libs.pari [x - 1, x + 1, x^2 - x + 1, x^2 + 1, x^2 + x + 1, x^4 - x^2 + 1]
R.<x> = PolynomialRing(QQ) prime_divisors(x^12 - 1) # needs sage.libs.pari
Tests with numpy and gmpy2 numbers:
sage: from numpy import int8 # needs numpy sage: prime_divisors(int8(-100)) # needs numpy [2, 5] sage: from gmpy2 import mpz sage: prime_divisors(mpz(-100)) [2, 5]
>>> from sage.all import * >>> from numpy import int8 # needs numpy >>> prime_divisors(int8(-Integer(100))) # needs numpy [2, 5] >>> from gmpy2 import mpz >>> prime_divisors(mpz(-Integer(100))) [2, 5]
from numpy import int8 # needs numpy prime_divisors(int8(-100)) # needs numpy from gmpy2 import mpz prime_divisors(mpz(-100))
- sage.arith.misc.prime_factors(n)[source]¶
Return the list of prime divisors (up to units) of this element of a unique factorization domain.
INPUT:
n– any object which can be decomposed into prime factors
OUTPUT:
A list of prime factors of
n. For integers, this list is sorted in increasing order.EXAMPLES:
Prime divisors of positive integers:
sage: prime_divisors(1) [] sage: prime_divisors(100) [2, 5] sage: prime_divisors(2004) [2, 3, 167]
>>> from sage.all import * >>> prime_divisors(Integer(1)) [] >>> prime_divisors(Integer(100)) [2, 5] >>> prime_divisors(Integer(2004)) [2, 3, 167]
prime_divisors(1) prime_divisors(100) prime_divisors(2004)
If
nis negative, we do not include -1 among the prime divisors, since -1 is not a prime number:sage: prime_divisors(-100) [2, 5]
>>> from sage.all import * >>> prime_divisors(-Integer(100)) [2, 5]
prime_divisors(-100)
For polynomials we get all irreducible factors:
sage: R.<x> = PolynomialRing(QQ) sage: prime_divisors(x^12 - 1) # needs sage.libs.pari [x - 1, x + 1, x^2 - x + 1, x^2 + 1, x^2 + x + 1, x^4 - x^2 + 1]
>>> from sage.all import * >>> R = PolynomialRing(QQ, names=('x',)); (x,) = R._first_ngens(1) >>> prime_divisors(x**Integer(12) - Integer(1)) # needs sage.libs.pari [x - 1, x + 1, x^2 - x + 1, x^2 + 1, x^2 + x + 1, x^4 - x^2 + 1]
R.<x> = PolynomialRing(QQ) prime_divisors(x^12 - 1) # needs sage.libs.pari
Tests with numpy and gmpy2 numbers:
sage: from numpy import int8 # needs numpy sage: prime_divisors(int8(-100)) # needs numpy [2, 5] sage: from gmpy2 import mpz sage: prime_divisors(mpz(-100)) [2, 5]
>>> from sage.all import * >>> from numpy import int8 # needs numpy >>> prime_divisors(int8(-Integer(100))) # needs numpy [2, 5] >>> from gmpy2 import mpz >>> prime_divisors(mpz(-Integer(100))) [2, 5]
from numpy import int8 # needs numpy prime_divisors(int8(-100)) # needs numpy from gmpy2 import mpz prime_divisors(mpz(-100))
- sage.arith.misc.prime_powers(start, stop=None)[source]¶
List of all positive primes powers between
startandstop-1, inclusive. If the second argument is omitted, returns the prime powers up to the first argument.INPUT:
start– integer; if two inputs are given, a lower bound for the returned set of prime powers. If this is the only input, then it is an upper bound.stop– integer (default:None); an upper bound for the returned set of prime powers
OUTPUT:
The set of all prime powers between
startandstopor, if only one argument is passed, the set of all prime powers between 1 andstart. The number \(n\) is a prime power if \(n=p^k\), where \(p\) is a prime number and \(k\) is a positive integer. Thus, \(1\) is not a prime power.EXAMPLES:
sage: prime_powers(20) # needs sage.libs.pari [2, 3, 4, 5, 7, 8, 9, 11, 13, 16, 17, 19] sage: len(prime_powers(1000)) # needs sage.libs.pari 193 sage: len(prime_range(1000)) # needs sage.libs.pari 168 sage: # needs sage.libs.pari sage: a = [z for z in range(95, 1234) if is_prime_power(z)] sage: b = prime_powers(95, 1234) sage: len(b) 194 sage: len(a) 194 sage: a[:10] [97, 101, 103, 107, 109, 113, 121, 125, 127, 128] sage: b[:10] [97, 101, 103, 107, 109, 113, 121, 125, 127, 128] sage: a == b True sage: prime_powers(100) == [i for i in range(100) if is_prime_power(i)] # needs sage.libs.pari True sage: prime_powers(10, 7) [] sage: prime_powers(-5) [] sage: prime_powers(-1, 3) # needs sage.libs.pari [2]
>>> from sage.all import * >>> prime_powers(Integer(20)) # needs sage.libs.pari [2, 3, 4, 5, 7, 8, 9, 11, 13, 16, 17, 19] >>> len(prime_powers(Integer(1000))) # needs sage.libs.pari 193 >>> len(prime_range(Integer(1000))) # needs sage.libs.pari 168 >>> # needs sage.libs.pari >>> a = [z for z in range(Integer(95), Integer(1234)) if is_prime_power(z)] >>> b = prime_powers(Integer(95), Integer(1234)) >>> len(b) 194 >>> len(a) 194 >>> a[:Integer(10)] [97, 101, 103, 107, 109, 113, 121, 125, 127, 128] >>> b[:Integer(10)] [97, 101, 103, 107, 109, 113, 121, 125, 127, 128] >>> a == b True >>> prime_powers(Integer(100)) == [i for i in range(Integer(100)) if is_prime_power(i)] # needs sage.libs.pari True >>> prime_powers(Integer(10), Integer(7)) [] >>> prime_powers(-Integer(5)) [] >>> prime_powers(-Integer(1), Integer(3)) # needs sage.libs.pari [2]
prime_powers(20) # needs sage.libs.pari len(prime_powers(1000)) # needs sage.libs.pari len(prime_range(1000)) # needs sage.libs.pari # needs sage.libs.pari a = [z for z in range(95, 1234) if is_prime_power(z)] b = prime_powers(95, 1234) len(b) len(a) a[:10] b[:10] a == b prime_powers(100) == [i for i in range(100) if is_prime_power(i)] # needs sage.libs.pari prime_powers(10, 7) prime_powers(-5) prime_powers(-1, 3) # needs sage.libs.pari
- sage.arith.misc.prime_to_m_part(n, m)[source]¶
Return the prime-to-\(m\) part of \(n\).
This is the largest divisor of \(n\) that is coprime to \(m\).
INPUT:
n– integer (nonzero)m– integer
OUTPUT: integer
EXAMPLES:
sage: prime_to_m_part(240,2) 15 sage: prime_to_m_part(240,3) 80 sage: prime_to_m_part(240,5) 48 sage: prime_to_m_part(43434,20) 21717
>>> from sage.all import * >>> prime_to_m_part(Integer(240),Integer(2)) 15 >>> prime_to_m_part(Integer(240),Integer(3)) 80 >>> prime_to_m_part(Integer(240),Integer(5)) 48 >>> prime_to_m_part(Integer(43434),Integer(20)) 21717
prime_to_m_part(240,2) prime_to_m_part(240,3) prime_to_m_part(240,5) prime_to_m_part(43434,20)
Note that integers also have a method with the same name:
sage: 240.prime_to_m_part(2) 15
>>> from sage.all import * >>> Integer(240).prime_to_m_part(Integer(2)) 15
240.prime_to_m_part(2)
Tests with numpy and gmpy2 numbers:
sage: from numpy import int16 # needs numpy sage: prime_to_m_part(int16(240), int16(2)) # needs numpy 15 sage: from gmpy2 import mpz sage: prime_to_m_part(mpz(240), mpz(2)) 15
>>> from sage.all import * >>> from numpy import int16 # needs numpy >>> prime_to_m_part(int16(Integer(240)), int16(Integer(2))) # needs numpy 15 >>> from gmpy2 import mpz >>> prime_to_m_part(mpz(Integer(240)), mpz(Integer(2))) 15
from numpy import int16 # needs numpy prime_to_m_part(int16(240), int16(2)) # needs numpy from gmpy2 import mpz prime_to_m_part(mpz(240), mpz(2))
- sage.arith.misc.primes(start=2, stop=None, proof=None)[source]¶
Return an iterator over all primes between
startandstop-1, inclusive. This is much slower thanprime_range(), but potentially uses less memory. As withnext_prime(), the optional argumentproofcontrols whether the numbers returned are guaranteed to be prime or not.This command is like the Python 3
range()command, except it only iterates over primes. In some cases it is better to useprimes()thanprime_range(), becauseprimes()does not build a list of all primes in the range in memory all at once. However, it is potentially much slower since it simply calls thenext_prime()function repeatedly, andnext_prime()is slow.INPUT:
start– integer (default: 2); lower bound for the primesstop– integer (or infinity); upper (open) bound for the primesproof– boolean orNone(default:None); ifTrue, the function yields only proven primes. IfFalse, the function uses a pseudo-primality test, which is much faster for really big numbers but does not provide a proof of primality. IfNone, uses the global default (seesage.structure.proof.proof)
OUTPUT: an iterator over primes from
starttostop-1, inclusiveEXAMPLES:
sage: # needs sage.libs.pari sage: for p in primes(5, 10): ....: print(p) 5 7 sage: list(primes(13)) [2, 3, 5, 7, 11] sage: list(primes(10000000000, 10000000100)) [10000000019, 10000000033, 10000000061, 10000000069, 10000000097] sage: max(primes(10^100, 10^100+10^4, proof=False)) 10000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000009631 sage: next(p for p in primes(10^20, infinity) if is_prime(2*p+1)) 100000000000000001243
>>> from sage.all import * >>> # needs sage.libs.pari >>> for p in primes(Integer(5), Integer(10)): ... print(p) 5 7 >>> list(primes(Integer(13))) [2, 3, 5, 7, 11] >>> list(primes(Integer(10000000000), Integer(10000000100))) [10000000019, 10000000033, 10000000061, 10000000069, 10000000097] >>> max(primes(Integer(10)**Integer(100), Integer(10)**Integer(100)+Integer(10)**Integer(4), proof=False)) 10000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000009631 >>> next(p for p in primes(Integer(10)**Integer(20), infinity) if is_prime(Integer(2)*p+Integer(1))) 100000000000000001243
# needs sage.libs.pari for p in primes(5, 10): print(p) list(primes(13)) list(primes(10000000000, 10000000100)) max(primes(10^100, 10^100+10^4, proof=False)) next(p for p in primes(10^20, infinity) if is_prime(2*p+1))
- sage.arith.misc.primes_first_n(n, leave_pari=False)[source]¶
Return the first \(n\) primes.
INPUT:
n– nonnegative integer
OUTPUT: list of the first \(n\) prime numbers
EXAMPLES:
sage: primes_first_n(10) # needs sage.libs.pari [2, 3, 5, 7, 11, 13, 17, 19, 23, 29] sage: len(primes_first_n(1000)) # needs sage.libs.pari 1000 sage: primes_first_n(0) []
>>> from sage.all import * >>> primes_first_n(Integer(10)) # needs sage.libs.pari [2, 3, 5, 7, 11, 13, 17, 19, 23, 29] >>> len(primes_first_n(Integer(1000))) # needs sage.libs.pari 1000 >>> primes_first_n(Integer(0)) []
primes_first_n(10) # needs sage.libs.pari len(primes_first_n(1000)) # needs sage.libs.pari primes_first_n(0)
- sage.arith.misc.primitive_root(n, check=True)[source]¶
Return a positive integer that generates the multiplicative group of integers modulo \(n\), if one exists; otherwise, raise a
ValueError.A primitive root exists if \(n=4\) or \(n=p^k\) or \(n=2p^k\), where \(p\) is an odd prime and \(k\) is a nonnegative number.
INPUT:
n– nonzero integercheck– boolean (default:True); ifFalse, then \(n\) is assumed to be a positive integer possessing a primitive root, and behavior is undefined otherwise.
OUTPUT:
A primitive root of \(n\). If \(n\) is prime, this is the smallest primitive root.
EXAMPLES:
sage: # needs sage.libs.pari sage: primitive_root(23) 5 sage: primitive_root(-46) 5 sage: primitive_root(25) 2 sage: print([primitive_root(p) for p in primes(100)]) [1, 2, 2, 3, 2, 2, 3, 2, 5, 2, 3, 2, 6, 3, 5, 2, 2, 2, 2, 7, 5, 3, 2, 3, 5] sage: primitive_root(8) Traceback (most recent call last): ... ValueError: no primitive root
>>> from sage.all import * >>> # needs sage.libs.pari >>> primitive_root(Integer(23)) 5 >>> primitive_root(-Integer(46)) 5 >>> primitive_root(Integer(25)) 2 >>> print([primitive_root(p) for p in primes(Integer(100))]) [1, 2, 2, 3, 2, 2, 3, 2, 5, 2, 3, 2, 6, 3, 5, 2, 2, 2, 2, 7, 5, 3, 2, 3, 5] >>> primitive_root(Integer(8)) Traceback (most recent call last): ... ValueError: no primitive root
# needs sage.libs.pari primitive_root(23) primitive_root(-46) primitive_root(25) print([primitive_root(p) for p in primes(100)]) primitive_root(8)
Note
It takes extra work to check if \(n\) has a primitive root; to avoid this, use
check=False, which may slightly speed things up (but could also result in undefined behavior). For example, the second call below is an order of magnitude faster than the first:sage: n = 10^50 + 151 # a prime sage: primitive_root(n) # needs sage.libs.pari 11 sage: primitive_root(n, check=False) # needs sage.libs.pari 11
>>> from sage.all import * >>> n = Integer(10)**Integer(50) + Integer(151) # a prime >>> primitive_root(n) # needs sage.libs.pari 11 >>> primitive_root(n, check=False) # needs sage.libs.pari 11
n = 10^50 + 151 # a prime primitive_root(n) # needs sage.libs.pari primitive_root(n, check=False) # needs sage.libs.pari
- sage.arith.misc.quadratic_residues(n)[source]¶
Return a sorted list of all squares modulo the integer \(n\) in the range \(0\leq x < |n|\).
EXAMPLES:
sage: quadratic_residues(11) [0, 1, 3, 4, 5, 9] sage: quadratic_residues(1) [0] sage: quadratic_residues(2) [0, 1] sage: quadratic_residues(8) [0, 1, 4] sage: quadratic_residues(-10) [0, 1, 4, 5, 6, 9] sage: v = quadratic_residues(1000); len(v) 159
>>> from sage.all import * >>> quadratic_residues(Integer(11)) [0, 1, 3, 4, 5, 9] >>> quadratic_residues(Integer(1)) [0] >>> quadratic_residues(Integer(2)) [0, 1] >>> quadratic_residues(Integer(8)) [0, 1, 4] >>> quadratic_residues(-Integer(10)) [0, 1, 4, 5, 6, 9] >>> v = quadratic_residues(Integer(1000)); len(v) 159
quadratic_residues(11) quadratic_residues(1) quadratic_residues(2) quadratic_residues(8) quadratic_residues(-10) v = quadratic_residues(1000); len(v)
Tests with numpy and gmpy2 numbers:
sage: from numpy import int8 # needs numpy sage: quadratic_residues(int8(11)) # needs numpy [0, 1, 3, 4, 5, 9] sage: from gmpy2 import mpz sage: quadratic_residues(mpz(11)) [0, 1, 3, 4, 5, 9]
>>> from sage.all import * >>> from numpy import int8 # needs numpy >>> quadratic_residues(int8(Integer(11))) # needs numpy [0, 1, 3, 4, 5, 9] >>> from gmpy2 import mpz >>> quadratic_residues(mpz(Integer(11))) [0, 1, 3, 4, 5, 9]
from numpy import int8 # needs numpy quadratic_residues(int8(11)) # needs numpy from gmpy2 import mpz quadratic_residues(mpz(11))
- sage.arith.misc.radical(n, *args, **kwds)[source]¶
Return the product of the prime divisors of n.
This calls
n.radical(*args, **kwds).EXAMPLES:
sage: radical(2 * 3^2 * 5^5) 30 sage: radical(0) Traceback (most recent call last): ... ArithmeticError: radical of 0 is not defined sage: K.<i> = QuadraticField(-1) # needs sage.rings.number_field sage: radical(K(2)) # needs sage.rings.number_field i + 1
>>> from sage.all import * >>> radical(Integer(2) * Integer(3)**Integer(2) * Integer(5)**Integer(5)) 30 >>> radical(Integer(0)) Traceback (most recent call last): ... ArithmeticError: radical of 0 is not defined >>> K = QuadraticField(-Integer(1), names=('i',)); (i,) = K._first_ngens(1)# needs sage.rings.number_field >>> radical(K(Integer(2))) # needs sage.rings.number_field i + 1
radical(2 * 3^2 * 5^5) radical(0) K.<i> = QuadraticField(-1) # needs sage.rings.number_field radical(K(2)) # needs sage.rings.number_field
Tests with numpy and gmpy2 numbers:
sage: from numpy import int8 # needs numpy sage: radical(int8(50)) # needs numpy 10 sage: from gmpy2 import mpz sage: radical(mpz(50)) 10
>>> from sage.all import * >>> from numpy import int8 # needs numpy >>> radical(int8(Integer(50))) # needs numpy 10 >>> from gmpy2 import mpz >>> radical(mpz(Integer(50))) 10
from numpy import int8 # needs numpy radical(int8(50)) # needs numpy from gmpy2 import mpz radical(mpz(50))
- sage.arith.misc.random_prime(n, proof=None, lbound=2)[source]¶
Return a random prime \(p\) between
lboundand \(n\).The returned prime \(p\) satisfies
lbound\(\leq p \leq n\).The returned prime \(p\) is chosen uniformly at random from the set of prime numbers less than or equal to \(n\).
INPUT:
n– integer \(\geq 2\)proof– boolean orNone(default:None); ifFalse, the function uses a pseudo-primality test, which is much faster for really big numbers but does not provide a proof of primality. IfNone, uses the global default (seesage.structure.proof.proof)lbound– integer; \(\geq 2\), lower bound for the chosen primes
EXAMPLES:
sage: # needs sage.libs.pari sage: p = random_prime(100000) sage: p.is_prime() True sage: p <= 100000 True sage: random_prime(2) 2
>>> from sage.all import * >>> # needs sage.libs.pari >>> p = random_prime(Integer(100000)) >>> p.is_prime() True >>> p <= Integer(100000) True >>> random_prime(Integer(2)) 2
# needs sage.libs.pari p = random_prime(100000) p.is_prime() p <= 100000 random_prime(2)
Here we generate a random prime between 100 and 200:
sage: p = random_prime(200, lbound=100) sage: p.is_prime() True sage: 100 <= p <= 200 True
>>> from sage.all import * >>> p = random_prime(Integer(200), lbound=Integer(100)) >>> p.is_prime() True >>> Integer(100) <= p <= Integer(200) True
p = random_prime(200, lbound=100) p.is_prime() 100 <= p <= 200
If all we care about is finding a pseudo prime, then we can pass in
proof=Falsesage: p = random_prime(200, proof=False, lbound=100) # needs sage.libs.pari sage: p.is_pseudoprime() # needs sage.libs.pari True sage: 100 <= p <= 200 True
>>> from sage.all import * >>> p = random_prime(Integer(200), proof=False, lbound=Integer(100)) # needs sage.libs.pari >>> p.is_pseudoprime() # needs sage.libs.pari True >>> Integer(100) <= p <= Integer(200) True
p = random_prime(200, proof=False, lbound=100) # needs sage.libs.pari p.is_pseudoprime() # needs sage.libs.pari 100 <= p <= 200
AUTHORS:
Jon Hanke (2006-08-08): with standard Stein cleanup
Jonathan Bober (2007-03-17)
- sage.arith.misc.rational_reconstruction(a, m, algorithm='fast')[source]¶
This function tries to compute \(x/y\), where \(x/y\) is a rational number in lowest terms such that the reduction of \(x/y\) modulo \(m\) is equal to \(a\) and the absolute values of \(x\) and \(y\) are both \(\le \sqrt{m/2}\). If such \(x/y\) exists, that pair is unique and this function returns it. If no such pair exists, this function raises
ZeroDivisionError.An efficient algorithm for computing rational reconstruction is very similar to the extended Euclidean algorithm. For more details, see Knuth, Vol 2, 3rd ed, pages 656-657.
INPUT:
a– integerm– a modulusalgorithm– string (default:'fast')'fast'– a fast implementation using direct GMP library callsin Cython
OUTPUT:
Numerator and denominator \(n\), \(d\) of the unique rational number \(r=n/d\), if it exists, with \(n\) and \(|d| \le \sqrt{N/2}\). Return \((0,0)\) if no such number exists.
The algorithm for rational reconstruction is described (with a complete nontrivial proof) on pages 656-657 of Knuth, Vol 2, 3rd ed. as the solution to exercise 51 on page 379. See in particular the conclusion paragraph right in the middle of page 657, which describes the algorithm thus:
This discussion proves that the problem can be solved efficiently by applying Algorithm 4.5.2X with \(u=m\) and \(v=a\), but with the following replacement for step X2: If \(v3 \le \sqrt{m/2}\), the algorithm terminates. The pair \((x,y)=(|v2|,v3*\mathrm{sign}(v2))\) is then the unique solution, provided that \(x\) and \(y\) are coprime and \(x \le \sqrt{m/2}\); otherwise there is no solution. (Alg 4.5.2X is the extended Euclidean algorithm.)
Knuth remarks that this algorithm is due to Wang, Kornerup, and Gregory from around 1983.
EXAMPLES:
sage: m = 100000 sage: (119*inverse_mod(53,m))%m 11323 sage: rational_reconstruction(11323,m) 119/53
>>> from sage.all import * >>> m = Integer(100000) >>> (Integer(119)*inverse_mod(Integer(53),m))%m 11323 >>> rational_reconstruction(Integer(11323),m) 119/53
m = 100000 (119*inverse_mod(53,m))%m rational_reconstruction(11323,m)
sage: rational_reconstruction(400,1000) Traceback (most recent call last): ... ArithmeticError: rational reconstruction of 400 (mod 1000) does not exist
>>> from sage.all import * >>> rational_reconstruction(Integer(400),Integer(1000)) Traceback (most recent call last): ... ArithmeticError: rational reconstruction of 400 (mod 1000) does not exist
rational_reconstruction(400,1000)
sage: rational_reconstruction(3, 292393) 3 sage: a = Integers(292393)(45/97); a 204977 sage: rational_reconstruction(a, 292393, algorithm='fast') 45/97 sage: rational_reconstruction(293048, 292393) Traceback (most recent call last): ... ArithmeticError: rational reconstruction of 655 (mod 292393) does not exist sage: rational_reconstruction(0, 0) Traceback (most recent call last): ... ZeroDivisionError: rational reconstruction with zero modulus sage: rational_reconstruction(0, 1, algorithm='foobar') Traceback (most recent call last): ... ValueError: unknown algorithm 'foobar'
>>> from sage.all import * >>> rational_reconstruction(Integer(3), Integer(292393)) 3 >>> a = Integers(Integer(292393))(Integer(45)/Integer(97)); a 204977 >>> rational_reconstruction(a, Integer(292393), algorithm='fast') 45/97 >>> rational_reconstruction(Integer(293048), Integer(292393)) Traceback (most recent call last): ... ArithmeticError: rational reconstruction of 655 (mod 292393) does not exist >>> rational_reconstruction(Integer(0), Integer(0)) Traceback (most recent call last): ... ZeroDivisionError: rational reconstruction with zero modulus >>> rational_reconstruction(Integer(0), Integer(1), algorithm='foobar') Traceback (most recent call last): ... ValueError: unknown algorithm 'foobar'
rational_reconstruction(3, 292393) a = Integers(292393)(45/97); a rational_reconstruction(a, 292393, algorithm='fast') rational_reconstruction(293048, 292393) rational_reconstruction(0, 0) rational_reconstruction(0, 1, algorithm='foobar')
Tests with numpy and gmpy2 numbers:
sage: from numpy import int32 # needs numpy sage: rational_reconstruction(int32(3), int32(292393)) # needs numpy 3 sage: from gmpy2 import mpz sage: rational_reconstruction(mpz(3), mpz(292393)) 3
>>> from sage.all import * >>> from numpy import int32 # needs numpy >>> rational_reconstruction(int32(Integer(3)), int32(Integer(292393))) # needs numpy 3 >>> from gmpy2 import mpz >>> rational_reconstruction(mpz(Integer(3)), mpz(Integer(292393))) 3
from numpy import int32 # needs numpy rational_reconstruction(int32(3), int32(292393)) # needs numpy from gmpy2 import mpz rational_reconstruction(mpz(3), mpz(292393))
- sage.arith.misc.rising_factorial(x, a)[source]¶
Return the rising factorial \((x)^a\).
The notation in the literature is a mess: often \((x)^a\), but there are many other notations: GKP: Concrete Mathematics uses \(x^{\overline{a}}\).
The rising factorial is also known as the Pochhammer symbol, see Maple and Mathematica.
Definition: for integer \(a \ge 0\) we have \(x(x+1) \cdots (x+a-1)\). In all other cases we use the GAMMA-function: \(\frac {\Gamma(x+a)} {\Gamma(x)}\).
INPUT:
x– element of a ringa– nonnegative integer orx,a– any numbers
OUTPUT: the rising factorial
See also
EXAMPLES:
sage: rising_factorial(10,3) 1320 sage: # needs sage.symbolic sage: rising_factorial(10, RR('3.0')) 1320.00000000000 sage: rising_factorial(10, RR('3.3')) 2826.38895824964 sage: a = rising_factorial(1+I, I); a gamma(2*I + 1)/gamma(I + 1) sage: CC(a) 0.266816390637832 + 0.122783354006372*I sage: a = rising_factorial(I, 4); a -10 sage: x = polygen(ZZ) sage: rising_factorial(x, 4) x^4 + 6*x^3 + 11*x^2 + 6*x
>>> from sage.all import * >>> rising_factorial(Integer(10),Integer(3)) 1320 >>> # needs sage.symbolic >>> rising_factorial(Integer(10), RR('3.0')) 1320.00000000000 >>> rising_factorial(Integer(10), RR('3.3')) 2826.38895824964 >>> a = rising_factorial(Integer(1)+I, I); a gamma(2*I + 1)/gamma(I + 1) >>> CC(a) 0.266816390637832 + 0.122783354006372*I >>> a = rising_factorial(I, Integer(4)); a -10 >>> x = polygen(ZZ) >>> rising_factorial(x, Integer(4)) x^4 + 6*x^3 + 11*x^2 + 6*x
rising_factorial(10,3) # needs sage.symbolic rising_factorial(10, RR('3.0')) rising_factorial(10, RR('3.3')) a = rising_factorial(1+I, I); a CC(a) a = rising_factorial(I, 4); a x = polygen(ZZ) rising_factorial(x, 4)AUTHORS:
Jaap Spies (2006-03-05)
- sage.arith.misc.smooth_part(x, base)[source]¶
Given an element
xof a Euclidean domain and a factor basebase, return aFactorizationobject corresponding to the largest divisor ofxthat splits completely overbase.The factor base can be specified in the following ways:
A sequence of elements.
A
ProductTreebuilt from such a sequence. (Caching the tree in the caller will speed things up if this function is called multiple times with the same factor base.)
EXAMPLES:
sage: from sage.arith.misc import smooth_part sage: from sage.rings.generic import ProductTree sage: smooth_part(10^77+1, primes(1000)) 11^2 * 23 * 463 sage: tree = ProductTree(primes(1000)) sage: smooth_part(10^77+1, tree) 11^2 * 23 * 463 sage: smooth_part(10^99+1, tree) 7 * 11^2 * 13 * 19 * 23
>>> from sage.all import * >>> from sage.arith.misc import smooth_part >>> from sage.rings.generic import ProductTree >>> smooth_part(Integer(10)**Integer(77)+Integer(1), primes(Integer(1000))) 11^2 * 23 * 463 >>> tree = ProductTree(primes(Integer(1000))) >>> smooth_part(Integer(10)**Integer(77)+Integer(1), tree) 11^2 * 23 * 463 >>> smooth_part(Integer(10)**Integer(99)+Integer(1), tree) 7 * 11^2 * 13 * 19 * 23
from sage.arith.misc import smooth_part from sage.rings.generic import ProductTree smooth_part(10^77+1, primes(1000)) tree = ProductTree(primes(1000)) smooth_part(10^77+1, tree) smooth_part(10^99+1, tree)
- sage.arith.misc.sort_complex_numbers_for_display(nums)[source]¶
Given a list of complex numbers (or a list of tuples, where the first element of each tuple is a complex number), we sort the list in a “pretty” order.
Real numbers (with a zero imaginary part) come before complex numbers, and are sorted. Complex numbers are sorted by their real part unless their real parts are quite close, in which case they are sorted by their imaginary part.
This is not a useful function mathematically (not least because there is no principled way to determine whether the real components should be treated as equal or not). It is called by various polynomial root-finders; its purpose is to make doctest printing more reproducible.
We deliberately choose a cumbersome name for this function to discourage use, since it is mathematically meaningless.
EXAMPLES:
sage: # needs sage.rings.complex_double sage: import sage.arith.misc sage: sort_c = sort_complex_numbers_for_display sage: nums = [CDF(i) for i in range(3)] sage: for i in range(3): ....: nums.append(CDF(i + RDF.random_element(-3e-11, 3e-11), ....: RDF.random_element())) ....: nums.append(CDF(i + RDF.random_element(-3e-11, 3e-11), ....: RDF.random_element())) sage: shuffle(nums) sage: nums = sort_c(nums) sage: for i in range(len(nums)): ....: if nums[i].imag(): ....: first_non_real = i ....: break ....: else: ....: first_non_real = len(nums) sage: assert first_non_real >= 3 sage: for i in range(first_non_real - 1): ....: assert nums[i].real() <= nums[i + 1].real() sage: def truncate(n): ....: if n.real() < 1e-10: ....: return 0 ....: else: ....: return n.real().n(digits=9) sage: for i in range(first_non_real, len(nums)-1): ....: assert truncate(nums[i]) <= truncate(nums[i + 1]) ....: if truncate(nums[i]) == truncate(nums[i + 1]): ....: assert nums[i].imag() <= nums[i+1].imag()
>>> from sage.all import * >>> # needs sage.rings.complex_double >>> import sage.arith.misc >>> sort_c = sort_complex_numbers_for_display >>> nums = [CDF(i) for i in range(Integer(3))] >>> for i in range(Integer(3)): ... nums.append(CDF(i + RDF.random_element(-RealNumber('3e-11'), RealNumber('3e-11')), ... RDF.random_element())) ... nums.append(CDF(i + RDF.random_element(-RealNumber('3e-11'), RealNumber('3e-11')), ... RDF.random_element())) >>> shuffle(nums) >>> nums = sort_c(nums) >>> for i in range(len(nums)): ... if nums[i].imag(): ... first_non_real = i ... break ... else: ... first_non_real = len(nums) >>> assert first_non_real >= Integer(3) >>> for i in range(first_non_real - Integer(1)): ... assert nums[i].real() <= nums[i + Integer(1)].real() >>> def truncate(n): ... if n.real() < RealNumber('1e-10'): ... return Integer(0) ... else: ... return n.real().n(digits=Integer(9)) >>> for i in range(first_non_real, len(nums)-Integer(1)): ... assert truncate(nums[i]) <= truncate(nums[i + Integer(1)]) ... if truncate(nums[i]) == truncate(nums[i + Integer(1)]): ... assert nums[i].imag() <= nums[i+Integer(1)].imag()
# needs sage.rings.complex_double import sage.arith.misc sort_c = sort_complex_numbers_for_display nums = [CDF(i) for i in range(3)] for i in range(3): nums.append(CDF(i + RDF.random_element(-3e-11, 3e-11), RDF.random_element())) nums.append(CDF(i + RDF.random_element(-3e-11, 3e-11), RDF.random_element())) shuffle(nums) nums = sort_c(nums) for i in range(len(nums)): if nums[i].imag(): first_non_real = i break else: first_non_real = len(nums) assert first_non_real >= 3 for i in range(first_non_real - 1): assert nums[i].real() <= nums[i + 1].real() def truncate(n): if n.real() < 1e-10: return 0 else: return n.real().n(digits=9) for i in range(first_non_real, len(nums)-1): assert truncate(nums[i]) <= truncate(nums[i + 1]) if truncate(nums[i]) == truncate(nums[i + 1]): assert nums[i].imag() <= nums[i+1].imag()
- sage.arith.misc.squarefree_divisors(x)[source]¶
Return an iterator over the squarefree divisors (up to units) of this ring element.
Depends on the output of the prime_divisors function.
Squarefree divisors of an integer are not necessarily yielded in increasing order.
INPUT:
x– an element of any ring for which the prime_divisors function works
EXAMPLES:
Integers with few prime divisors:
sage: list(squarefree_divisors(7)) [1, 7] sage: list(squarefree_divisors(6)) [1, 2, 3, 6] sage: list(squarefree_divisors(12)) [1, 2, 3, 6]
>>> from sage.all import * >>> list(squarefree_divisors(Integer(7))) [1, 7] >>> list(squarefree_divisors(Integer(6))) [1, 2, 3, 6] >>> list(squarefree_divisors(Integer(12))) [1, 2, 3, 6]
list(squarefree_divisors(7)) list(squarefree_divisors(6)) list(squarefree_divisors(12))
Squarefree divisors are not yielded in increasing order:
sage: list(squarefree_divisors(30)) [1, 2, 3, 6, 5, 10, 15, 30]
>>> from sage.all import * >>> list(squarefree_divisors(Integer(30))) [1, 2, 3, 6, 5, 10, 15, 30]
list(squarefree_divisors(30))
- sage.arith.misc.subfactorial(n)[source]¶
Subfactorial or rencontres numbers, or derangements: number of permutations of \(n\) elements with no fixed points.
INPUT:
n– nonnegative integer
OUTPUT: integer
EXAMPLES:
sage: subfactorial(0) 1 sage: subfactorial(1) 0 sage: subfactorial(8) 14833
>>> from sage.all import * >>> subfactorial(Integer(0)) 1 >>> subfactorial(Integer(1)) 0 >>> subfactorial(Integer(8)) 14833
subfactorial(0) subfactorial(1) subfactorial(8)
Tests with numpy and gmpy2 numbers:
sage: from numpy import int8 # needs numpy sage: subfactorial(int8(8)) # needs numpy 14833 sage: from gmpy2 import mpz sage: subfactorial(mpz(8)) 14833
>>> from sage.all import * >>> from numpy import int8 # needs numpy >>> subfactorial(int8(Integer(8))) # needs numpy 14833 >>> from gmpy2 import mpz >>> subfactorial(mpz(Integer(8))) 14833
from numpy import int8 # needs numpy subfactorial(int8(8)) # needs numpy from gmpy2 import mpz subfactorial(mpz(8))
AUTHORS:
Jaap Spies (2007-01-23)
- sage.arith.misc.sum_of_k_squares(k, n)[source]¶
Write the integer \(n\) as a sum of \(k\) integer squares if possible; otherwise raise a
ValueError.INPUT:
k– nonnegative integern– integer
OUTPUT: a tuple \((x_1, ..., x_k)\) of nonnegative integers such that their squares sum to \(n\).
EXAMPLES:
sage: sum_of_k_squares(2, 9634) (15, 97) sage: sum_of_k_squares(3, 9634) (0, 15, 97) sage: sum_of_k_squares(4, 9634) (1, 2, 5, 98) sage: sum_of_k_squares(5, 9634) (0, 1, 2, 5, 98) sage: sum_of_k_squares(6, 11^1111 - 1) # needs sage.libs.pari (19215400822645944253860920437586326284, 37204645194585992174252915693267578306, 3473654819477394665857484221256136567800161086815834297092488779216863122, 5860191799617673633547572610351797996721850737768032876360978911074629287841061578270832330322236796556721252602860754789786937515870682024273948, 20457423294558182494001919812379023992538802203730791019728543439765347851316366537094696896669915675685581905102118246887673397020172285247862426612188418787649371716686651256443143210952163970564228423098202682066311189439731080552623884051737264415984619097656479060977602722566383385989, 311628095411678159849237738619458396497534696043580912225334269371611836910345930320700816649653412141574887113710604828156159177769285115652741014638785285820578943010943846225597311231847997461959204894255074229895666356909071243390280307709880906261008237873840245959883405303580405277298513108957483306488193844321589356441983980532251051786704380984788999660195252373574924026139168936921591652831237741973242604363696352878914129671292072201700073286987126265965322808664802662993006926302359371379531571194266134916767573373504566621665949840469229781956838744551367172353) sage: sum_of_k_squares(7, 0) (0, 0, 0, 0, 0, 0, 0) sage: sum_of_k_squares(30,999999) (0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 2, 3, 7, 44, 999) sage: sum_of_k_squares(1, 9) (3,) sage: sum_of_k_squares(1, 10) Traceback (most recent call last): ... ValueError: 10 is not a sum of 1 square sage: sum_of_k_squares(1, -10) Traceback (most recent call last): ... ValueError: -10 is not a sum of 1 square sage: sum_of_k_squares(0, 9) Traceback (most recent call last): ... ValueError: 9 is not a sum of 0 squares sage: sum_of_k_squares(0, 0) () sage: sum_of_k_squares(7, -1) Traceback (most recent call last): ... ValueError: -1 is not a sum of 7 squares sage: sum_of_k_squares(-1, 0) Traceback (most recent call last): ... ValueError: k = -1 must be nonnegative
>>> from sage.all import * >>> sum_of_k_squares(Integer(2), Integer(9634)) (15, 97) >>> sum_of_k_squares(Integer(3), Integer(9634)) (0, 15, 97) >>> sum_of_k_squares(Integer(4), Integer(9634)) (1, 2, 5, 98) >>> sum_of_k_squares(Integer(5), Integer(9634)) (0, 1, 2, 5, 98) >>> sum_of_k_squares(Integer(6), Integer(11)**Integer(1111) - Integer(1)) # needs sage.libs.pari (19215400822645944253860920437586326284, 37204645194585992174252915693267578306, 3473654819477394665857484221256136567800161086815834297092488779216863122, 5860191799617673633547572610351797996721850737768032876360978911074629287841061578270832330322236796556721252602860754789786937515870682024273948, 20457423294558182494001919812379023992538802203730791019728543439765347851316366537094696896669915675685581905102118246887673397020172285247862426612188418787649371716686651256443143210952163970564228423098202682066311189439731080552623884051737264415984619097656479060977602722566383385989, 311628095411678159849237738619458396497534696043580912225334269371611836910345930320700816649653412141574887113710604828156159177769285115652741014638785285820578943010943846225597311231847997461959204894255074229895666356909071243390280307709880906261008237873840245959883405303580405277298513108957483306488193844321589356441983980532251051786704380984788999660195252373574924026139168936921591652831237741973242604363696352878914129671292072201700073286987126265965322808664802662993006926302359371379531571194266134916767573373504566621665949840469229781956838744551367172353) >>> sum_of_k_squares(Integer(7), Integer(0)) (0, 0, 0, 0, 0, 0, 0) >>> sum_of_k_squares(Integer(30),Integer(999999)) (0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 2, 3, 7, 44, 999) >>> sum_of_k_squares(Integer(1), Integer(9)) (3,) >>> sum_of_k_squares(Integer(1), Integer(10)) Traceback (most recent call last): ... ValueError: 10 is not a sum of 1 square >>> sum_of_k_squares(Integer(1), -Integer(10)) Traceback (most recent call last): ... ValueError: -10 is not a sum of 1 square >>> sum_of_k_squares(Integer(0), Integer(9)) Traceback (most recent call last): ... ValueError: 9 is not a sum of 0 squares >>> sum_of_k_squares(Integer(0), Integer(0)) () >>> sum_of_k_squares(Integer(7), -Integer(1)) Traceback (most recent call last): ... ValueError: -1 is not a sum of 7 squares >>> sum_of_k_squares(-Integer(1), Integer(0)) Traceback (most recent call last): ... ValueError: k = -1 must be nonnegative
sum_of_k_squares(2, 9634) sum_of_k_squares(3, 9634) sum_of_k_squares(4, 9634) sum_of_k_squares(5, 9634) sum_of_k_squares(6, 11^1111 - 1) # needs sage.libs.pari sum_of_k_squares(7, 0) sum_of_k_squares(30,999999) sum_of_k_squares(1, 9) sum_of_k_squares(1, 10) sum_of_k_squares(1, -10) sum_of_k_squares(0, 9) sum_of_k_squares(0, 0) sum_of_k_squares(7, -1) sum_of_k_squares(-1, 0)
Tests with numpy and gmpy2 numbers:
sage: from numpy import int16 # needs numpy sage: sum_of_k_squares(int16(2), int16(9634)) # needs numpy (15, 97) sage: from gmpy2 import mpz sage: sum_of_k_squares(mpz(2), mpz(9634)) (15, 97)
>>> from sage.all import * >>> from numpy import int16 # needs numpy >>> sum_of_k_squares(int16(Integer(2)), int16(Integer(9634))) # needs numpy (15, 97) >>> from gmpy2 import mpz >>> sum_of_k_squares(mpz(Integer(2)), mpz(Integer(9634))) (15, 97)
from numpy import int16 # needs numpy sum_of_k_squares(int16(2), int16(9634)) # needs numpy from gmpy2 import mpz sum_of_k_squares(mpz(2), mpz(9634))
- sage.arith.misc.three_squares(n)[source]¶
Write the integer \(n\) as a sum of three integer squares if possible; otherwise raise a
ValueError.INPUT:
n– integer
OUTPUT: a tuple \((a,b,c)\) of nonnegative integers such that \(n = a^2 + b^2 + c^2\) with \(a <= b <= c\).
EXAMPLES:
sage: three_squares(389) (1, 8, 18) sage: three_squares(946) (9, 9, 28) sage: three_squares(2986) (3, 24, 49) sage: three_squares(7^100) (0, 0, 1798465042647412146620280340569649349251249) sage: three_squares(11^111 - 1) # needs sage.libs.pari (616274160655975340150706442680, 901582938385735143295060746161, 6270382387635744140394001363065311967964099981788593947233) sage: three_squares(7 * 2^41) # needs sage.libs.pari (1048576, 2097152, 3145728) sage: three_squares(7 * 2^42) Traceback (most recent call last): ... ValueError: 30786325577728 is not a sum of 3 squares sage: three_squares(0) (0, 0, 0) sage: three_squares(-1) Traceback (most recent call last): ... ValueError: -1 is not a sum of 3 squares
>>> from sage.all import * >>> three_squares(Integer(389)) (1, 8, 18) >>> three_squares(Integer(946)) (9, 9, 28) >>> three_squares(Integer(2986)) (3, 24, 49) >>> three_squares(Integer(7)**Integer(100)) (0, 0, 1798465042647412146620280340569649349251249) >>> three_squares(Integer(11)**Integer(111) - Integer(1)) # needs sage.libs.pari (616274160655975340150706442680, 901582938385735143295060746161, 6270382387635744140394001363065311967964099981788593947233) >>> three_squares(Integer(7) * Integer(2)**Integer(41)) # needs sage.libs.pari (1048576, 2097152, 3145728) >>> three_squares(Integer(7) * Integer(2)**Integer(42)) Traceback (most recent call last): ... ValueError: 30786325577728 is not a sum of 3 squares >>> three_squares(Integer(0)) (0, 0, 0) >>> three_squares(-Integer(1)) Traceback (most recent call last): ... ValueError: -1 is not a sum of 3 squares
three_squares(389) three_squares(946) three_squares(2986) three_squares(7^100) three_squares(11^111 - 1) # needs sage.libs.pari three_squares(7 * 2^41) # needs sage.libs.pari three_squares(7 * 2^42) three_squares(0) three_squares(-1)
ALGORITHM:
- sage.arith.misc.trial_division(n, bound=None)[source]¶
Return the smallest prime divisor less than or equal to
boundof the positive integer \(n\), or \(n\) if there is no such prime. If the optional argument bound is omitted, then bound \(\leq n\).INPUT:
n– positive integerbound– (optional) positive integer
OUTPUT: a prime
p=boundthat divides \(n\), or \(n\) if there is no such primeEXAMPLES:
sage: trial_division(15) 3 sage: trial_division(91) 7 sage: trial_division(11) 11 sage: trial_division(387833, 300) 387833 sage: # 300 is not big enough to split off a sage: # factor, but 400 is. sage: trial_division(387833, 400) 389
>>> from sage.all import * >>> trial_division(Integer(15)) 3 >>> trial_division(Integer(91)) 7 >>> trial_division(Integer(11)) 11 >>> trial_division(Integer(387833), Integer(300)) 387833 >>> # 300 is not big enough to split off a >>> # factor, but 400 is. >>> trial_division(Integer(387833), Integer(400)) 389
trial_division(15) trial_division(91) trial_division(11) trial_division(387833, 300) # 300 is not big enough to split off a # factor, but 400 is. trial_division(387833, 400)
Tests with numpy and gmpy2 numbers:
sage: from numpy import int8 # needs numpy sage: trial_division(int8(91)) # needs numpy 7 sage: from gmpy2 import mpz sage: trial_division(mpz(91)) 7
>>> from sage.all import * >>> from numpy import int8 # needs numpy >>> trial_division(int8(Integer(91))) # needs numpy 7 >>> from gmpy2 import mpz >>> trial_division(mpz(Integer(91))) 7
from numpy import int8 # needs numpy trial_division(int8(91)) # needs numpy from gmpy2 import mpz trial_division(mpz(91))
- sage.arith.misc.two_squares(n)[source]¶
Write the integer \(n\) as a sum of two integer squares if possible; otherwise raise a
ValueError.INPUT:
n– integer
OUTPUT: a tuple \((a,b)\) of nonnegative integers such that \(n = a^2 + b^2\) with \(a <= b\).
EXAMPLES:
sage: two_squares(389) (10, 17) sage: two_squares(21) Traceback (most recent call last): ... ValueError: 21 is not a sum of 2 squares sage: two_squares(21^2) (0, 21) sage: a, b = two_squares(100000000000000000129); a, b # needs sage.libs.pari (4418521500, 8970878873) sage: a^2 + b^2 # needs sage.libs.pari 100000000000000000129 sage: two_squares(2^222 + 1) # needs sage.libs.pari (253801659504708621991421712450521, 2583712713213354898490304645018692) sage: two_squares(0) (0, 0) sage: two_squares(-1) Traceback (most recent call last): ... ValueError: -1 is not a sum of 2 squares
>>> from sage.all import * >>> two_squares(Integer(389)) (10, 17) >>> two_squares(Integer(21)) Traceback (most recent call last): ... ValueError: 21 is not a sum of 2 squares >>> two_squares(Integer(21)**Integer(2)) (0, 21) >>> a, b = two_squares(Integer(100000000000000000129)); a, b # needs sage.libs.pari (4418521500, 8970878873) >>> a**Integer(2) + b**Integer(2) # needs sage.libs.pari 100000000000000000129 >>> two_squares(Integer(2)**Integer(222) + Integer(1)) # needs sage.libs.pari (253801659504708621991421712450521, 2583712713213354898490304645018692) >>> two_squares(Integer(0)) (0, 0) >>> two_squares(-Integer(1)) Traceback (most recent call last): ... ValueError: -1 is not a sum of 2 squares
two_squares(389) two_squares(21) two_squares(21^2) a, b = two_squares(100000000000000000129); a, b # needs sage.libs.pari a^2 + b^2 # needs sage.libs.pari two_squares(2^222 + 1) # needs sage.libs.pari two_squares(0) two_squares(-1)
ALGORITHM:
- sage.arith.misc.valuation(m, *args, **kwds)[source]¶
Return the valuation of
m.This function simply calls the m.valuation() method. See the documentation of m.valuation() for a more precise description.
Note that the use of this functions is discouraged as it is better to use m.valuation() directly.
Note
This is not always a valuation in the mathematical sense. For more information see: sage.rings.finite_rings.integer_mod.IntegerMod_int.valuation
EXAMPLES:
sage: valuation(512,2) 9 sage: valuation(1,2) 0 sage: valuation(5/9, 3) -2
>>> from sage.all import * >>> valuation(Integer(512),Integer(2)) 9 >>> valuation(Integer(1),Integer(2)) 0 >>> valuation(Integer(5)/Integer(9), Integer(3)) -2
valuation(512,2) valuation(1,2) valuation(5/9, 3)
Valuation of 0 is defined, but valuation with respect to 0 is not:
sage: valuation(0,7) +Infinity sage: valuation(3,0) Traceback (most recent call last): ... ValueError: You can only compute the valuation with respect to a integer larger than 1.
>>> from sage.all import * >>> valuation(Integer(0),Integer(7)) +Infinity >>> valuation(Integer(3),Integer(0)) Traceback (most recent call last): ... ValueError: You can only compute the valuation with respect to a integer larger than 1.
valuation(0,7) valuation(3,0)
Here are some other examples:
sage: valuation(100,10) 2 sage: valuation(200,10) 2 sage: valuation(243,3) 5 sage: valuation(243*10007,3) 5 sage: valuation(243*10007,10007) 1 sage: y = QQ['y'].gen() sage: valuation(y^3, y) 3 sage: x = QQ[['x']].gen() sage: valuation((x^3-x^2)/(x-4)) 2 sage: valuation(4r,2r) 2 sage: valuation(1r,1r) Traceback (most recent call last): ... ValueError: You can only compute the valuation with respect to a integer larger than 1. sage: from numpy import int16 # needs numpy sage: valuation(int16(512), int16(2)) # needs numpy 9 sage: from gmpy2 import mpz sage: valuation(mpz(512), mpz(2)) 9
>>> from sage.all import * >>> valuation(Integer(100),Integer(10)) 2 >>> valuation(Integer(200),Integer(10)) 2 >>> valuation(Integer(243),Integer(3)) 5 >>> valuation(Integer(243)*Integer(10007),Integer(3)) 5 >>> valuation(Integer(243)*Integer(10007),Integer(10007)) 1 >>> y = QQ['y'].gen() >>> valuation(y**Integer(3), y) 3 >>> x = QQ[['x']].gen() >>> valuation((x**Integer(3)-x**Integer(2))/(x-Integer(4))) 2 >>> valuation(4,2) 2 >>> valuation(1,1) Traceback (most recent call last): ... ValueError: You can only compute the valuation with respect to a integer larger than 1. >>> from numpy import int16 # needs numpy >>> valuation(int16(Integer(512)), int16(Integer(2))) # needs numpy 9 >>> from gmpy2 import mpz >>> valuation(mpz(Integer(512)), mpz(Integer(2))) 9
valuation(100,10) valuation(200,10) valuation(243,3) valuation(243*10007,3) valuation(243*10007,10007) y = QQ['y'].gen() valuation(y^3, y) x = QQ[['x']].gen() valuation((x^3-x^2)/(x-4)) valuation(4r,2r) valuation(1r,1r) from numpy import int16 # needs numpy valuation(int16(512), int16(2)) # needs numpy from gmpy2 import mpz valuation(mpz(512), mpz(2))
- sage.arith.misc.xgcd(a, b=None)[source]¶
Return the greatest common divisor and the Bézout coefficients of the input arguments.
When both
aandbare given, then return a triple(g,s,t)such that \(g = s\cdot a+t\cdot b = \gcd(a,b)\). When onlyais given, then return a tuplerof lengthlen(a) + 1such that \(r_0 = \sum_{i = 0}^{len(a) - 1} r_{i + 1}a_i = gcd(a_0, \dots, a_{len(a) - 1})\)Note
One exception is if the elements are not in a principal ideal domain (see Wikipedia article Principal_ideal_domain), e.g., they are both polynomials over the integers. Then this function can’t in general return
(g,s,t)orras above, since they need not exist. Instead, over the integers, whenaandbare given, we first multiply \(g\) by a divisor of the resultant of \(a/g\) and \(b/g\), up to sign.INPUT:
One of the following:
a, b– integers or more generally, element of a ring for which the xgcd make sense (e.g. a field or univariate polynomials).a– a list or tuple of at least two integers or more generally, elements of a ring which the xgcd make sense.
OUTPUT:
One of the following:
g, s, t– when two inputsa, bare given. They satisfy \(g = s\cdot a + t\cdot b\).r– a tuple, when onlyais given (andb = None). Its first entryr[0]is the gcd of the inputs, and has length one longer than the length ofa. Its entries satisfy \(r_0 = \sum_{i = 0}^{len(a) - 1} r_{i + 1}a_i\).
Note
There is no guarantee that the returned cofactors (s and t) are minimal.
EXAMPLES:
sage: xgcd(56, 44) (4, 4, -5) sage: 4*56 + (-5)*44 4 sage: xgcd([56, 44]) (4, 4, -5) sage: r = xgcd([30, 105, 70, 42]); r (1, -255, 85, -17, -2) sage: (-255)*30 + 85*105 + (-17)*70 + (-2)*42 1 sage: xgcd([]) (0,) sage: xgcd([42]) (42, 1) sage: g, a, b = xgcd(5/1, 7/1); g, a, b (1, 3, -2) sage: a*(5/1) + b*(7/1) == g True sage: x = polygen(QQ) sage: xgcd(x^3 - 1, x^2 - 1) (x - 1, 1, -x) sage: g, a, b, c = xgcd([x^4 - x, x^6 - 1, x^4 - 1]); g, a, b, c (x - 1, x^3, -x, 1) sage: a*(x^4 - x) + b*(x^6 - 1) + c*(x^4 - 1) == g True sage: K.<g> = NumberField(x^2 - 3) # needs sage.rings.number_field sage: g.xgcd(g + 2) # needs sage.rings.number_field (1, 1/3*g, 0) sage: # needs sage.rings.number_field sage: R.<a,b> = K[] sage: S.<y> = R.fraction_field()[] sage: xgcd(y^2, a*y + b) (1, a^2/b^2, ((-a)/b^2)*y + 1/b) sage: xgcd((b+g)*y^2, (a+g)*y + b) (1, (a^2 + (2*g)*a + 3)/(b^3 + g*b^2), ((-a + (-g))/b^2)*y + 1/b)
>>> from sage.all import * >>> xgcd(Integer(56), Integer(44)) (4, 4, -5) >>> Integer(4)*Integer(56) + (-Integer(5))*Integer(44) 4 >>> xgcd([Integer(56), Integer(44)]) (4, 4, -5) >>> r = xgcd([Integer(30), Integer(105), Integer(70), Integer(42)]); r (1, -255, 85, -17, -2) >>> (-Integer(255))*Integer(30) + Integer(85)*Integer(105) + (-Integer(17))*Integer(70) + (-Integer(2))*Integer(42) 1 >>> xgcd([]) (0,) >>> xgcd([Integer(42)]) (42, 1) >>> g, a, b = xgcd(Integer(5)/Integer(1), Integer(7)/Integer(1)); g, a, b (1, 3, -2) >>> a*(Integer(5)/Integer(1)) + b*(Integer(7)/Integer(1)) == g True >>> x = polygen(QQ) >>> xgcd(x**Integer(3) - Integer(1), x**Integer(2) - Integer(1)) (x - 1, 1, -x) >>> g, a, b, c = xgcd([x**Integer(4) - x, x**Integer(6) - Integer(1), x**Integer(4) - Integer(1)]); g, a, b, c (x - 1, x^3, -x, 1) >>> a*(x**Integer(4) - x) + b*(x**Integer(6) - Integer(1)) + c*(x**Integer(4) - Integer(1)) == g True >>> K = NumberField(x**Integer(2) - Integer(3), names=('g',)); (g,) = K._first_ngens(1)# needs sage.rings.number_field >>> g.xgcd(g + Integer(2)) # needs sage.rings.number_field (1, 1/3*g, 0) >>> # needs sage.rings.number_field >>> R = K['a, b']; (a, b,) = R._first_ngens(2) >>> S = R.fraction_field()['y']; (y,) = S._first_ngens(1) >>> xgcd(y**Integer(2), a*y + b) (1, a^2/b^2, ((-a)/b^2)*y + 1/b) >>> xgcd((b+g)*y**Integer(2), (a+g)*y + b) (1, (a^2 + (2*g)*a + 3)/(b^3 + g*b^2), ((-a + (-g))/b^2)*y + 1/b)
xgcd(56, 44) 4*56 + (-5)*44 xgcd([56, 44]) r = xgcd([30, 105, 70, 42]); r (-255)*30 + 85*105 + (-17)*70 + (-2)*42 xgcd([]) xgcd([42]) g, a, b = xgcd(5/1, 7/1); g, a, b a*(5/1) + b*(7/1) == g x = polygen(QQ) xgcd(x^3 - 1, x^2 - 1) g, a, b, c = xgcd([x^4 - x, x^6 - 1, x^4 - 1]); g, a, b, c a*(x^4 - x) + b*(x^6 - 1) + c*(x^4 - 1) == g K.<g> = NumberField(x^2 - 3) # needs sage.rings.number_field g.xgcd(g + 2) # needs sage.rings.number_field # needs sage.rings.number_field R.<a,b> = K[] S.<y> = R.fraction_field()[] xgcd(y^2, a*y + b) xgcd((b+g)*y^2, (a+g)*y + b)
Here is an example of a xgcd for two polynomials over the integers, where the linear combination is not the gcd but the gcd multiplied by the resultant:
sage: R.<x> = ZZ[] sage: gcd(2*x*(x-1), x^2) x sage: xgcd(2*x*(x-1), x^2) (2*x, -1, 2) sage: (2*(x-1)).resultant(x) # needs sage.libs.pari 2
>>> from sage.all import * >>> R = ZZ['x']; (x,) = R._first_ngens(1) >>> gcd(Integer(2)*x*(x-Integer(1)), x**Integer(2)) x >>> xgcd(Integer(2)*x*(x-Integer(1)), x**Integer(2)) (2*x, -1, 2) >>> (Integer(2)*(x-Integer(1))).resultant(x) # needs sage.libs.pari 2
R.<x> = ZZ[] gcd(2*x*(x-1), x^2) xgcd(2*x*(x-1), x^2) (2*(x-1)).resultant(x) # needs sage.libs.pari
Tests with numpy and gmpy2 types:
sage: from numpy import int8 # needs numpy sage: xgcd(4, int8(8)) # needs numpy (4, 1, 0) sage: xgcd(int8(4), int8(8)) # needs numpy (4, 1, 0) sage: xgcd([int8(4), int8(8), int(10)]) # needs numpy (2, -2, 0, 1) sage: from gmpy2 import mpz sage: xgcd(mpz(4), mpz(8)) (4, 1, 0) sage: xgcd(4, mpz(8)) (4, 1, 0) sage: xgcd([4, mpz(8), mpz(10)]) (2, -2, 0, 1)
>>> from sage.all import * >>> from numpy import int8 # needs numpy >>> xgcd(Integer(4), int8(Integer(8))) # needs numpy (4, 1, 0) >>> xgcd(int8(Integer(4)), int8(Integer(8))) # needs numpy (4, 1, 0) >>> xgcd([int8(Integer(4)), int8(Integer(8)), int(Integer(10))]) # needs numpy (2, -2, 0, 1) >>> from gmpy2 import mpz >>> xgcd(mpz(Integer(4)), mpz(Integer(8))) (4, 1, 0) >>> xgcd(Integer(4), mpz(Integer(8))) (4, 1, 0) >>> xgcd([Integer(4), mpz(Integer(8)), mpz(Integer(10))]) (2, -2, 0, 1)
from numpy import int8 # needs numpy xgcd(4, int8(8)) # needs numpy xgcd(int8(4), int8(8)) # needs numpy xgcd([int8(4), int8(8), int(10)]) # needs numpy from gmpy2 import mpz xgcd(mpz(4), mpz(8)) xgcd(4, mpz(8)) xgcd([4, mpz(8), mpz(10)])
- sage.arith.misc.xkcd(n='')[source]¶
This function is similar to the xgcd function, but behaves in a completely different way.
See https://xkcd.com/json.html for more details.
INPUT:
n– integer (optional)
OUTPUT: a fragment of HTML
EXAMPLES:
sage: xkcd(353) # optional - internet <h1>Python</h1><img src="https://imgs.xkcd.com/comics/python.png" title="I wrote 20 short programs in Python yesterday. It was wonderful. Perl, I'm leaving you."><div>Source: <a href="http://xkcd.com/353" target="_blank">http://xkcd.com/353</a></div>
>>> from sage.all import * >>> xkcd(Integer(353)) # optional - internet <h1>Python</h1><img src="https://imgs.xkcd.com/comics/python.png" title="I wrote 20 short programs in Python yesterday. It was wonderful. Perl, I'm leaving you."><div>Source: <a href="http://xkcd.com/353" target="_blank">http://xkcd.com/353</a></div>
xkcd(353) # optional - internet
- sage.arith.misc.xlcm(m, n)[source]¶
Extended lcm function: given two positive integers \(m,n\), returns a triple \((l,m_1,n_1)\) such that \(l=\mathop{\mathrm{lcm}}(m,n)=m_1 \cdot n_1\) where \(m_1|m\), \(n_1|n\) and \(\gcd(m_1,n_1)=1\), all with no factorization.
Used to construct an element of order \(l\) from elements of orders \(m,n\) in any group: see sage/groups/generic.py for examples.
EXAMPLES:
sage: xlcm(120,36) (360, 40, 9)
>>> from sage.all import * >>> xlcm(Integer(120),Integer(36)) (360, 40, 9)
xlcm(120,36)