Lists of Manin symbols over number fields, elements of \(\mathbb{P}^1(R/N)\)¶
Lists of elements of \(\mathbb{P}^1(R/N)\) where \(R\) is the ring of integers of a number field \(K\) and \(N\) is an integral ideal.
AUTHORS:
- Maite Aranes (2009): Initial version 
EXAMPLES:
We define a P1NFList:
sage: x = polygen(QQ, 'x')
sage: k.<a> = NumberField(x^3 + 11)
sage: N = k.ideal(5, a^2 - a + 1)
sage: P = P1NFList(N); P
The projective line over
 the ring of integers modulo the Fractional ideal (5, a^2 - a + 1)
>>> from sage.all import *
>>> x = polygen(QQ, 'x')
>>> k = NumberField(x**Integer(3) + Integer(11), names=('a',)); (a,) = k._first_ngens(1)
>>> N = k.ideal(Integer(5), a**Integer(2) - a + Integer(1))
>>> P = P1NFList(N); P
The projective line over
 the ring of integers modulo the Fractional ideal (5, a^2 - a + 1)
List operations with the P1NFList:
sage: len(P)
26
sage: [p for p in P]
[M-symbol (0: 1) of level Fractional ideal (5, a^2 - a + 1),
...
M-symbol (1: 2*a^2 + 2*a) of level Fractional ideal (5, a^2 - a + 1)]
>>> from sage.all import *
>>> len(P)
26
>>> [p for p in P]
[M-symbol (0: 1) of level Fractional ideal (5, a^2 - a + 1),
...
M-symbol (1: 2*a^2 + 2*a) of level Fractional ideal (5, a^2 - a + 1)]
The elements of the P1NFList are M-symbols:
sage: type(P[2])
<class 'sage.modular.modsym.p1list_nf.MSymbol'>
>>> from sage.all import *
>>> type(P[Integer(2)])
<class 'sage.modular.modsym.p1list_nf.MSymbol'>
Definition of MSymbols:
sage: alpha = MSymbol(N, 3, a^2); alpha
M-symbol (3: a^2) of level Fractional ideal (5, a^2 - a + 1)
>>> from sage.all import *
>>> alpha = MSymbol(N, Integer(3), a**Integer(2)); alpha
M-symbol (3: a^2) of level Fractional ideal (5, a^2 - a + 1)
Find the index of the class of an M-Symbol \((c: d)\) in the list:
sage: i = P.index(alpha)
sage: P[i].c*alpha.d - P[i].d*alpha.c in N
True
>>> from sage.all import *
>>> i = P.index(alpha)
>>> P[i].c*alpha.d - P[i].d*alpha.c in N
True
Lift an MSymbol to a matrix in \(SL(2, R)\):
sage: alpha = MSymbol(N, a + 2, 3*a^2)
sage: alpha.lift_to_sl2_Ok()
[-a - 1, 15*a^2 - 38*a + 86, a + 2, -a^2 + 9*a - 19]
sage: Ok = k.ring_of_integers()
sage: M = Matrix(Ok, 2, alpha.lift_to_sl2_Ok())
sage: det(M)
1
sage: M[1][1] - alpha.d in N
True
>>> from sage.all import *
>>> alpha = MSymbol(N, a + Integer(2), Integer(3)*a**Integer(2))
>>> alpha.lift_to_sl2_Ok()
[-a - 1, 15*a^2 - 38*a + 86, a + 2, -a^2 + 9*a - 19]
>>> Ok = k.ring_of_integers()
>>> M = Matrix(Ok, Integer(2), alpha.lift_to_sl2_Ok())
>>> det(M)
1
>>> M[Integer(1)][Integer(1)] - alpha.d in N
True
Lift an MSymbol from P1NFList to a matrix in \(SL(2, R)\)
sage: P[3]
M-symbol (1: -2*a) of level Fractional ideal (5, a^2 - a + 1)
sage: P.lift_to_sl2_Ok(3)
[0, -1, 1, -2*a]
>>> from sage.all import *
>>> P[Integer(3)]
M-symbol (1: -2*a) of level Fractional ideal (5, a^2 - a + 1)
>>> P.lift_to_sl2_Ok(Integer(3))
[0, -1, 1, -2*a]
- class sage.modular.modsym.p1list_nf.MSymbol(N, c, d=None, check=True)[source]¶
- Bases: - SageObject- The constructor for an M-symbol over a number field. - INPUT: - N– integral ideal (the modulus or level)
- c– integral element of the underlying number field or an MSymbol of level N
- d– (optional) when present, it must be an integral element such that \(\langle c\rangle + \langle d\rangle + N = R\), where \(R\) is the corresponding ring of integers
- check– boolean (default:- True); if- check=Falsethe constructor does not check the condition \(\langle c\rangle + \langle d\rangle + N = R\)
 - OUTPUT: - An M-symbol modulo the given ideal \(N\), i.e. an element of the projective line \(\mathbb{P}^1(R/N)\), where \(R\) is the ring of integers of the underlying number field. - EXAMPLES: - sage: x = polygen(QQ, 'x') sage: k.<a> = NumberField(x^3 + 11) sage: N = k.ideal(a + 1, 2) sage: MSymbol(N, 3, a^2 + 1) M-symbol (3: a^2 + 1) of level Fractional ideal (2, a + 1) - >>> from sage.all import * >>> x = polygen(QQ, 'x') >>> k = NumberField(x**Integer(3) + Integer(11), names=('a',)); (a,) = k._first_ngens(1) >>> N = k.ideal(a + Integer(1), Integer(2)) >>> MSymbol(N, Integer(3), a**Integer(2) + Integer(1)) M-symbol (3: a^2 + 1) of level Fractional ideal (2, a + 1) - We can give a tuple as input: - sage: MSymbol(N, (1, 0)) M-symbol (1: 0) of level Fractional ideal (2, a + 1) - >>> from sage.all import * >>> MSymbol(N, (Integer(1), Integer(0))) M-symbol (1: 0) of level Fractional ideal (2, a + 1) - We get an error if \(\langle c\rangle\), \(\langle d\rangle\) and \(N\) are not coprime: - sage: MSymbol(N, 2*a, a - 1) Traceback (most recent call last): ... ValueError: (2*a, a - 1) is not an element of P1(R/N). sage: MSymbol(N, (0, 0)) Traceback (most recent call last): ... ValueError: (0, 0) is not an element of P1(R/N). - >>> from sage.all import * >>> MSymbol(N, Integer(2)*a, a - Integer(1)) Traceback (most recent call last): ... ValueError: (2*a, a - 1) is not an element of P1(R/N). >>> MSymbol(N, (Integer(0), Integer(0))) Traceback (most recent call last): ... ValueError: (0, 0) is not an element of P1(R/N). - Saving and loading works: - sage: alpha = MSymbol(N, 3, a^2 + 1) sage: loads(dumps(alpha))==alpha True - >>> from sage.all import * >>> alpha = MSymbol(N, Integer(3), a**Integer(2) + Integer(1)) >>> loads(dumps(alpha))==alpha True - N()[source]¶
- Return the level or modulus of this MSymbol. - EXAMPLES: - sage: x = polygen(QQ, 'x') sage: k.<a> = NumberField(x^2 + 23) sage: N = k.ideal(3, a - 1) sage: alpha = MSymbol(N, 3, a) sage: alpha.N() Fractional ideal (3, 1/2*a - 1/2) - >>> from sage.all import * >>> x = polygen(QQ, 'x') >>> k = NumberField(x**Integer(2) + Integer(23), names=('a',)); (a,) = k._first_ngens(1) >>> N = k.ideal(Integer(3), a - Integer(1)) >>> alpha = MSymbol(N, Integer(3), a) >>> alpha.N() Fractional ideal (3, 1/2*a - 1/2) 
 - property c¶
- Return the first coefficient of the M-symbol. - EXAMPLES: - sage: x = polygen(QQ, 'x') sage: k.<a> = NumberField(x^3 + 11) sage: N = k.ideal(a + 1, 2) sage: alpha = MSymbol(N, 3, a^2 + 1) sage: alpha.c # indirect doctest 3 - >>> from sage.all import * >>> x = polygen(QQ, 'x') >>> k = NumberField(x**Integer(3) + Integer(11), names=('a',)); (a,) = k._first_ngens(1) >>> N = k.ideal(a + Integer(1), Integer(2)) >>> alpha = MSymbol(N, Integer(3), a**Integer(2) + Integer(1)) >>> alpha.c # indirect doctest 3 
 - property d¶
- Return the second coefficient of the M-symbol. - EXAMPLES: - sage: x = polygen(QQ, 'x') sage: k.<a> = NumberField(x^3 + 11) sage: N = k.ideal(a + 1, 2) sage: alpha = MSymbol(N, 3, a^2 + 1) sage: alpha.d # indirect doctest a^2 + 1 - >>> from sage.all import * >>> x = polygen(QQ, 'x') >>> k = NumberField(x**Integer(3) + Integer(11), names=('a',)); (a,) = k._first_ngens(1) >>> N = k.ideal(a + Integer(1), Integer(2)) >>> alpha = MSymbol(N, Integer(3), a**Integer(2) + Integer(1)) >>> alpha.d # indirect doctest a^2 + 1 
 - lift_to_sl2_Ok()[source]¶
- Lift the - MSymbolto an element of \(SL(2, O_k)\), where \(O_k\) is the ring of integers of the corresponding number field.- OUTPUT: - A list of integral elements \([a, b, c', d']\) that are the entries of a \(2\times 2\) matrix with determinant 1. The lower two entries are congruent (modulo the level) to the coefficients \(c\), \(d\) of the - MSymbol- self.- EXAMPLES: - sage: x = polygen(QQ, 'x') sage: k.<a> = NumberField(x^2 + 23) sage: N = k.ideal(3, a - 1) sage: alpha = MSymbol(N, 3*a + 1, a) sage: alpha.lift_to_sl2_Ok() [0, -1, 1, a] - >>> from sage.all import * >>> x = polygen(QQ, 'x') >>> k = NumberField(x**Integer(2) + Integer(23), names=('a',)); (a,) = k._first_ngens(1) >>> N = k.ideal(Integer(3), a - Integer(1)) >>> alpha = MSymbol(N, Integer(3)*a + Integer(1), a) >>> alpha.lift_to_sl2_Ok() [0, -1, 1, a] 
 - normalize(with_scalar=False)[source]¶
- Return a normalized - MSymbol(a canonical representative of an element of \(\mathbb{P}^1(R/N)\)) equivalent to- self.- INPUT: - with_scalar– boolean (default:- False)
 - OUTPUT: - (only if - with_scalar=True) a transforming scalar \(u\), such that \((u*c', u*d')\) is congruent to \((c: d)\) (mod \(N\)), where \((c: d)\) are the coefficients of- selfand \(N\) is the level.
- a normalized - MSymbol\((c': d')\) equivalent to- self.
 - EXAMPLES: - sage: x = polygen(QQ, 'x') sage: k.<a> = NumberField(x^2 + 23) sage: N = k.ideal(3, a - 1) sage: alpha1 = MSymbol(N, 3, a); alpha1 M-symbol (3: a) of level Fractional ideal (3, 1/2*a - 1/2) sage: alpha1.normalize() M-symbol (0: 1) of level Fractional ideal (3, 1/2*a - 1/2) sage: alpha2 = MSymbol(N, 4, a + 1) sage: alpha2.normalize() M-symbol (1: -a) of level Fractional ideal (3, 1/2*a - 1/2) - >>> from sage.all import * >>> x = polygen(QQ, 'x') >>> k = NumberField(x**Integer(2) + Integer(23), names=('a',)); (a,) = k._first_ngens(1) >>> N = k.ideal(Integer(3), a - Integer(1)) >>> alpha1 = MSymbol(N, Integer(3), a); alpha1 M-symbol (3: a) of level Fractional ideal (3, 1/2*a - 1/2) >>> alpha1.normalize() M-symbol (0: 1) of level Fractional ideal (3, 1/2*a - 1/2) >>> alpha2 = MSymbol(N, Integer(4), a + Integer(1)) >>> alpha2.normalize() M-symbol (1: -a) of level Fractional ideal (3, 1/2*a - 1/2) - We get the scaling factor by setting - with_scalar=True:- sage: alpha1.normalize(with_scalar=True) (a, M-symbol (0: 1) of level Fractional ideal (3, 1/2*a - 1/2)) sage: r, beta1 = alpha1.normalize(with_scalar=True) sage: r*beta1.c - alpha1.c in N True sage: r*beta1.d - alpha1.d in N True sage: r, beta2 = alpha2.normalize(with_scalar=True) sage: r*beta2.c - alpha2.c in N True sage: r*beta2.d - alpha2.d in N True - >>> from sage.all import * >>> alpha1.normalize(with_scalar=True) (a, M-symbol (0: 1) of level Fractional ideal (3, 1/2*a - 1/2)) >>> r, beta1 = alpha1.normalize(with_scalar=True) >>> r*beta1.c - alpha1.c in N True >>> r*beta1.d - alpha1.d in N True >>> r, beta2 = alpha2.normalize(with_scalar=True) >>> r*beta2.c - alpha2.c in N True >>> r*beta2.d - alpha2.d in N True 
 - tuple()[source]¶
- Return the - MSymbolas a list \((c, d)\).- EXAMPLES: - sage: x = polygen(QQ, 'x') sage: k.<a> = NumberField(x^2 + 23) sage: N = k.ideal(3, a - 1) sage: alpha = MSymbol(N, 3, a); alpha M-symbol (3: a) of level Fractional ideal (3, 1/2*a - 1/2) sage: alpha.tuple() (3, a) - >>> from sage.all import * >>> x = polygen(QQ, 'x') >>> k = NumberField(x**Integer(2) + Integer(23), names=('a',)); (a,) = k._first_ngens(1) >>> N = k.ideal(Integer(3), a - Integer(1)) >>> alpha = MSymbol(N, Integer(3), a); alpha M-symbol (3: a) of level Fractional ideal (3, 1/2*a - 1/2) >>> alpha.tuple() (3, a) 
 
- class sage.modular.modsym.p1list_nf.P1NFList(N)[source]¶
- Bases: - SageObject- The class for \(\mathbb{P}^1(R/N)\), the projective line modulo \(N\), where \(R\) is the ring of integers of a number field \(K\) and \(N\) is an integral ideal. - INPUT: - N– integral ideal (the modulus or level)
 - OUTPUT: - A - P1NFListobject representing \(\mathbb{P}^1(R/N)\).- EXAMPLES: - sage: x = polygen(QQ, 'x') sage: k.<a> = NumberField(x^3 + 11) sage: N = k.ideal(5, a + 1) sage: P = P1NFList(N); P The projective line over the ring of integers modulo the Fractional ideal (5, a + 1) - >>> from sage.all import * >>> x = polygen(QQ, 'x') >>> k = NumberField(x**Integer(3) + Integer(11), names=('a',)); (a,) = k._first_ngens(1) >>> N = k.ideal(Integer(5), a + Integer(1)) >>> P = P1NFList(N); P The projective line over the ring of integers modulo the Fractional ideal (5, a + 1) - Saving and loading works. - sage: loads(dumps(P)) == P True - >>> from sage.all import * >>> loads(dumps(P)) == P True - N()[source]¶
- Return the level or modulus of this - P1NFList.- EXAMPLES: - sage: x = polygen(QQ, 'x') sage: k.<a> = NumberField(x^2 + 31) sage: N = k.ideal(5, a + 3) sage: P = P1NFList(N) sage: P.N() Fractional ideal (5, 1/2*a + 3/2) - >>> from sage.all import * >>> x = polygen(QQ, 'x') >>> k = NumberField(x**Integer(2) + Integer(31), names=('a',)); (a,) = k._first_ngens(1) >>> N = k.ideal(Integer(5), a + Integer(3)) >>> P = P1NFList(N) >>> P.N() Fractional ideal (5, 1/2*a + 3/2) 
 - apply_J_epsilon(i, e1, e2=1)[source]¶
- Apply the matrix \(J_{\epsilon}\) = [e1, 0, 0, e2] to the \(i\)-th M-Symbol of the list. - e1, e2 are units of the underlying number field. - INPUT: - i– integer
- e1– unit
- e2– unit (default: 1)
 - OUTPUT: - integer – the index of the M-Symbol obtained by the right action of the matrix \(J_{\epsilon}\) = [e1, 0, 0, e2] on the \(i\)-th M-Symbol. - EXAMPLES: - sage: x = polygen(QQ, 'x') sage: k.<a> = NumberField(x^3 + 11) sage: N = k.ideal(5, a + 1) sage: P = P1NFList(N) sage: u = k.unit_group().gens_values(); u [-1, 2*a^2 + 4*a - 1] sage: P.apply_J_epsilon(4, -1) 2 sage: P.apply_J_epsilon(4, u[0], u[1]) 1 - >>> from sage.all import * >>> x = polygen(QQ, 'x') >>> k = NumberField(x**Integer(3) + Integer(11), names=('a',)); (a,) = k._first_ngens(1) >>> N = k.ideal(Integer(5), a + Integer(1)) >>> P = P1NFList(N) >>> u = k.unit_group().gens_values(); u [-1, 2*a^2 + 4*a - 1] >>> P.apply_J_epsilon(Integer(4), -Integer(1)) 2 >>> P.apply_J_epsilon(Integer(4), u[Integer(0)], u[Integer(1)]) 1 - sage: k.<a> = NumberField(x^4 + 13*x - 7) sage: N = k.ideal(a + 1) sage: P = P1NFList(N) sage: u = k.unit_group().gens_values(); u [-1, -a^3 - a^2 - a - 12, -a^3 - 3*a^2 + 1] sage: P.apply_J_epsilon(3, u[2]^2)==P.apply_J_epsilon(P.apply_J_epsilon(3, u[2]),u[2]) True - >>> from sage.all import * >>> k = NumberField(x**Integer(4) + Integer(13)*x - Integer(7), names=('a',)); (a,) = k._first_ngens(1) >>> N = k.ideal(a + Integer(1)) >>> P = P1NFList(N) >>> u = k.unit_group().gens_values(); u [-1, -a^3 - a^2 - a - 12, -a^3 - 3*a^2 + 1] >>> P.apply_J_epsilon(Integer(3), u[Integer(2)]**Integer(2))==P.apply_J_epsilon(P.apply_J_epsilon(Integer(3), u[Integer(2)]),u[Integer(2)]) True 
 - apply_S(i)[source]¶
- Applies the matrix \(S\) = [0, -1, 1, 0] to the \(i\)-th M-Symbol of the list. - INPUT: - i– integer
 - OUTPUT: - integer – the index of the M-Symbol obtained by the right action of the matrix \(S\) = [0, -1, 1, 0] on the \(i\)-th M-Symbol. - EXAMPLES: - sage: x = polygen(QQ, 'x') sage: k.<a> = NumberField(x^3 + 11) sage: N = k.ideal(5, a + 1) sage: P = P1NFList(N) sage: j = P.apply_S(P.index_of_normalized_pair(1, 0)) sage: P[j] M-symbol (0: 1) of level Fractional ideal (5, a + 1) - >>> from sage.all import * >>> x = polygen(QQ, 'x') >>> k = NumberField(x**Integer(3) + Integer(11), names=('a',)); (a,) = k._first_ngens(1) >>> N = k.ideal(Integer(5), a + Integer(1)) >>> P = P1NFList(N) >>> j = P.apply_S(P.index_of_normalized_pair(Integer(1), Integer(0))) >>> P[j] M-symbol (0: 1) of level Fractional ideal (5, a + 1) - We test that S has order 2: - sage: j = randint(0,len(P)-1) sage: P.apply_S(P.apply_S(j))==j True - >>> from sage.all import * >>> j = randint(Integer(0),len(P)-Integer(1)) >>> P.apply_S(P.apply_S(j))==j True 
 - apply_TS(i)[source]¶
- Applies the matrix \(TS\) = [1, -1, 0, 1] to the \(i\)-th M-Symbol of the list. - INPUT: - i– integer
 - OUTPUT: - integer – the index of the M-Symbol obtained by the right action of the matrix \(TS\) = [1, -1, 0, 1] on the \(i\)-th M-Symbol. - EXAMPLES: - sage: x = polygen(QQ, 'x') sage: k.<a> = NumberField(x^3 + 11) sage: N = k.ideal(5, a + 1) sage: P = P1NFList(N) sage: P.apply_TS(3) 2 - >>> from sage.all import * >>> x = polygen(QQ, 'x') >>> k = NumberField(x**Integer(3) + Integer(11), names=('a',)); (a,) = k._first_ngens(1) >>> N = k.ideal(Integer(5), a + Integer(1)) >>> P = P1NFList(N) >>> P.apply_TS(Integer(3)) 2 - We test that \(TS\) has order 3: - sage: j = randint(0,len(P)-1) sage: P.apply_TS(P.apply_TS(P.apply_TS(j)))==j True - >>> from sage.all import * >>> j = randint(Integer(0),len(P)-Integer(1)) >>> P.apply_TS(P.apply_TS(P.apply_TS(j)))==j True 
 - apply_T_alpha(i, alpha=1)[source]¶
- Applies the matrix \(T_{alpha}\) = [1, \(alpha\), 0, 1] to the \(i\)-th M-Symbol of the list. - INPUT: - i– integer
- alpha– (default: 1) element of the corresponding ring of integers
 - OUTPUT: - integer – the index of the M-Symbol obtained by the right action of the matrix \(T_{alpha}\) = [1, \(alpha\), 0, 1] on the \(i\)-th M-Symbol. - EXAMPLES: - sage: x = polygen(QQ, 'x') sage: k.<a> = NumberField(x^3 + 11) sage: N = k.ideal(5, a + 1) sage: P = P1NFList(N) sage: P.apply_T_alpha(4, a^ 2 - 2) 3 - >>> from sage.all import * >>> x = polygen(QQ, 'x') >>> k = NumberField(x**Integer(3) + Integer(11), names=('a',)); (a,) = k._first_ngens(1) >>> N = k.ideal(Integer(5), a + Integer(1)) >>> P = P1NFList(N) >>> P.apply_T_alpha(Integer(4), a** Integer(2) - Integer(2)) 3 - We test that \(T_a*T_b = T_{(a+b)}\): - sage: P.apply_T_alpha(3, a^2 - 2)==P.apply_T_alpha(P.apply_T_alpha(3,a^2),-2) True - >>> from sage.all import * >>> P.apply_T_alpha(Integer(3), a**Integer(2) - Integer(2))==P.apply_T_alpha(P.apply_T_alpha(Integer(3),a**Integer(2)),-Integer(2)) True 
 - index(c, d=None, with_scalar=False)[source]¶
- Return the index of the class of the pair \((c, d)\) in the fixed list of representatives of \(\mathbb{P}^1(R/N)\). - INPUT: - c– integral element of the corresponding number field, or an- MSymbol
- d– (optional) when present, it must be an integral element of the number field such that \((c, d)\) defines an M-symbol of level \(N\)
- with_scalar– boolean (default:- False)
 - OUTPUT: - u– the normalizing scalar (only if- with_scalar=True)
- i– the index of \((c, d)\) in the list
 - EXAMPLES: - sage: x = polygen(QQ, 'x') sage: k.<a> = NumberField(x^2 + 31) sage: N = k.ideal(5, a + 3) sage: P = P1NFList(N) sage: P.index(3,a) 5 sage: P[5]==MSymbol(N, 3, a).normalize() True - >>> from sage.all import * >>> x = polygen(QQ, 'x') >>> k = NumberField(x**Integer(2) + Integer(31), names=('a',)); (a,) = k._first_ngens(1) >>> N = k.ideal(Integer(5), a + Integer(3)) >>> P = P1NFList(N) >>> P.index(Integer(3),a) 5 >>> P[Integer(5)]==MSymbol(N, Integer(3), a).normalize() True - We can give an - MSymbolas input:- sage: alpha = MSymbol(N, 3, a) sage: P.index(alpha) 5 - >>> from sage.all import * >>> alpha = MSymbol(N, Integer(3), a) >>> P.index(alpha) 5 - We cannot look for the class of an - MSymbolof a different level:- sage: M = k.ideal(a + 1) sage: beta = MSymbol(M, 0, 1) sage: P.index(beta) Traceback (most recent call last): ... ValueError: The MSymbol is of a different level - >>> from sage.all import * >>> M = k.ideal(a + Integer(1)) >>> beta = MSymbol(M, Integer(0), Integer(1)) >>> P.index(beta) Traceback (most recent call last): ... ValueError: The MSymbol is of a different level - If we are interested in the transforming scalar: - sage: alpha = MSymbol(N, 3, a) sage: P.index(alpha, with_scalar=True) (-a, 5) sage: u, i = P.index(alpha, with_scalar=True) sage: (u*P[i].c - alpha.c in N) and (u*P[i].d - alpha.d in N) True - >>> from sage.all import * >>> alpha = MSymbol(N, Integer(3), a) >>> P.index(alpha, with_scalar=True) (-a, 5) >>> u, i = P.index(alpha, with_scalar=True) >>> (u*P[i].c - alpha.c in N) and (u*P[i].d - alpha.d in N) True 
 - index_of_normalized_pair(c, d=None)[source]¶
- Return the index of the class \((c, d)\) in the fixed list of representatives of \(\mathbb(P)^1(R/N)\). - INPUT: - c– integral element of the corresponding number field, or a normalized- MSymbol
- d– (optional) when present, it must be an integral element of the number field such that \((c, d)\) defines a normalized M-symbol of level \(N\)
 - OUTPUT: - i– the index of \((c, d)\) in the list- EXAMPLES: - sage: x = polygen(QQ, 'x') sage: k.<a> = NumberField(x^2 + 31) sage: N = k.ideal(5, a + 3) sage: P = P1NFList(N) sage: P.index_of_normalized_pair(1, 0) 3 sage: j = randint(0,len(P)-1) sage: P.index_of_normalized_pair(P[j])==j True - >>> from sage.all import * >>> x = polygen(QQ, 'x') >>> k = NumberField(x**Integer(2) + Integer(31), names=('a',)); (a,) = k._first_ngens(1) >>> N = k.ideal(Integer(5), a + Integer(3)) >>> P = P1NFList(N) >>> P.index_of_normalized_pair(Integer(1), Integer(0)) 3 >>> j = randint(Integer(0),len(P)-Integer(1)) >>> P.index_of_normalized_pair(P[j])==j True 
 - lift_to_sl2_Ok(i)[source]¶
- Lift the \(i\)-th element of this - P1NFListto an element of \(SL(2, R)\), where \(R\) is the ring of integers of the corresponding number field.- INPUT: - i– integer (index of the element to lift)
 - OUTPUT: - If the \(i\)-th element is \((c : d)\), the function returns a list of integral elements \([a, b, c', d']\) that defines a \(2\times 2\) matrix with determinant 1 and such that \(c=c'\) (mod \(N\)) and \(d=d'\) (mod \(N\)). - EXAMPLES: - sage: x = polygen(QQ, 'x') sage: k.<a> = NumberField(x^2 + 23) sage: N = k.ideal(3) sage: P = P1NFList(N) sage: len(P) 16 sage: P[5] M-symbol (1/2*a + 1/2: -a) of level Fractional ideal (3) sage: P.lift_to_sl2_Ok(5) [-a, 2*a - 2, 1/2*a + 1/2, -a] - >>> from sage.all import * >>> x = polygen(QQ, 'x') >>> k = NumberField(x**Integer(2) + Integer(23), names=('a',)); (a,) = k._first_ngens(1) >>> N = k.ideal(Integer(3)) >>> P = P1NFList(N) >>> len(P) 16 >>> P[Integer(5)] M-symbol (1/2*a + 1/2: -a) of level Fractional ideal (3) >>> P.lift_to_sl2_Ok(Integer(5)) [-a, 2*a - 2, 1/2*a + 1/2, -a] - sage: Ok = k.ring_of_integers() sage: L = [Matrix(Ok, 2, P.lift_to_sl2_Ok(i)) for i in range(len(P))] sage: all(det(L[i]) == 1 for i in range(len(L))) True - >>> from sage.all import * >>> Ok = k.ring_of_integers() >>> L = [Matrix(Ok, Integer(2), P.lift_to_sl2_Ok(i)) for i in range(len(P))] >>> all(det(L[i]) == Integer(1) for i in range(len(L))) True 
 - list()[source]¶
- Return the underlying list of this - P1NFListobject.- EXAMPLES: - sage: x = polygen(QQ, 'x') sage: k.<a> = NumberField(x^3 + 11) sage: N = k.ideal(5, a+1) sage: P = P1NFList(N) sage: type(P) <class 'sage.modular.modsym.p1list_nf.P1NFList'> sage: type(P.list()) <... 'list'> - >>> from sage.all import * >>> x = polygen(QQ, 'x') >>> k = NumberField(x**Integer(3) + Integer(11), names=('a',)); (a,) = k._first_ngens(1) >>> N = k.ideal(Integer(5), a+Integer(1)) >>> P = P1NFList(N) >>> type(P) <class 'sage.modular.modsym.p1list_nf.P1NFList'> >>> type(P.list()) <... 'list'> 
 - normalize(c, d=None, with_scalar=False)[source]¶
- Return a normalised element of \(\mathbb{P}^1(R/N)\). - INPUT: - c– integral element of the underlying number field, or an MSymbol
- d– (optional) when present, it must be an integral element of the number field such that \((c, d)\) defines an M-symbol of level \(N\)
- with_scalar– boolean (default:- False)
 - OUTPUT: - (only if - with_scalar=True) a transforming scalar \(u\), such that \((u*c', u*d')\) is congruent to \((c: d)\) (mod \(N\)).
- a normalized - MSymbol\((c': d')\) equivalent to \((c: d)\).
 - EXAMPLES: - sage: x = polygen(QQ, 'x') sage: k.<a> = NumberField(x^2 + 31) sage: N = k.ideal(5, a + 3) sage: P = P1NFList(N) sage: P.normalize(3, a) M-symbol (1: 2*a) of level Fractional ideal (5, 1/2*a + 3/2) - >>> from sage.all import * >>> x = polygen(QQ, 'x') >>> k = NumberField(x**Integer(2) + Integer(31), names=('a',)); (a,) = k._first_ngens(1) >>> N = k.ideal(Integer(5), a + Integer(3)) >>> P = P1NFList(N) >>> P.normalize(Integer(3), a) M-symbol (1: 2*a) of level Fractional ideal (5, 1/2*a + 3/2) - We can use an - MSymbolas input:- sage: alpha = MSymbol(N, 3, a) sage: P.normalize(alpha) M-symbol (1: 2*a) of level Fractional ideal (5, 1/2*a + 3/2) - >>> from sage.all import * >>> alpha = MSymbol(N, Integer(3), a) >>> P.normalize(alpha) M-symbol (1: 2*a) of level Fractional ideal (5, 1/2*a + 3/2) - If we are interested in the normalizing scalar: - sage: P.normalize(alpha, with_scalar=True) (-a, M-symbol (1: 2*a) of level Fractional ideal (5, 1/2*a + 3/2)) sage: r, beta = P.normalize(alpha, with_scalar=True) sage: (r*beta.c - alpha.c in N) and (r*beta.d - alpha.d in N) True - >>> from sage.all import * >>> P.normalize(alpha, with_scalar=True) (-a, M-symbol (1: 2*a) of level Fractional ideal (5, 1/2*a + 3/2)) >>> r, beta = P.normalize(alpha, with_scalar=True) >>> (r*beta.c - alpha.c in N) and (r*beta.d - alpha.d in N) True 
 
- sage.modular.modsym.p1list_nf.P1NFList_clear_level_cache()[source]¶
- Clear the global cache of data for the level ideals. - EXAMPLES: - sage: x = polygen(QQ, 'x') sage: k.<a> = NumberField(x^3 + 11) sage: N = k.ideal(a+1) sage: alpha = MSymbol(N, 2*a^2, 5) sage: alpha.normalize() M-symbol (-4*a^2: 5*a^2) of level Fractional ideal (a + 1) sage: sage.modular.modsym.p1list_nf._level_cache {Fractional ideal (a + 1): (...)} sage: sage.modular.modsym.p1list_nf.P1NFList_clear_level_cache() sage: sage.modular.modsym.p1list_nf._level_cache {} - >>> from sage.all import * >>> x = polygen(QQ, 'x') >>> k = NumberField(x**Integer(3) + Integer(11), names=('a',)); (a,) = k._first_ngens(1) >>> N = k.ideal(a+Integer(1)) >>> alpha = MSymbol(N, Integer(2)*a**Integer(2), Integer(5)) >>> alpha.normalize() M-symbol (-4*a^2: 5*a^2) of level Fractional ideal (a + 1) >>> sage.modular.modsym.p1list_nf._level_cache {Fractional ideal (a + 1): (...)} >>> sage.modular.modsym.p1list_nf.P1NFList_clear_level_cache() >>> sage.modular.modsym.p1list_nf._level_cache {} 
- sage.modular.modsym.p1list_nf.lift_to_sl2_Ok(N, c, d)[source]¶
- Lift a pair (c, d) to an element of \(SL(2, O_k)\), where \(O_k\) is the ring of integers of the corresponding number field. - INPUT: - N– number field ideal
- c– integral element of the number field
- d– integral element of the number field
 - OUTPUT: - A list \([a, b, c', d']\) of integral elements that are the entries of a \(2\times 2\) matrix with determinant 1. The lower two entries are congruent to \(c\), \(d\) modulo the ideal \(N\). - EXAMPLES: - sage: from sage.modular.modsym.p1list_nf import lift_to_sl2_Ok sage: x = polygen(QQ, 'x') sage: k.<a> = NumberField(x^2 + 23) sage: Ok = k.ring_of_integers() sage: N = k.ideal(3) sage: M = Matrix(Ok, 2, lift_to_sl2_Ok(N, 1, a)) sage: det(M) 1 sage: M = Matrix(Ok, 2, lift_to_sl2_Ok(N, 0, a)) sage: det(M) 1 sage: (M[1][0] in N) and (M[1][1] - a in N) True sage: M = Matrix(Ok, 2, lift_to_sl2_Ok(N, 0, 0)) Traceback (most recent call last): ... ValueError: Cannot lift (0, 0) to an element of Sl2(Ok). - >>> from sage.all import * >>> from sage.modular.modsym.p1list_nf import lift_to_sl2_Ok >>> x = polygen(QQ, 'x') >>> k = NumberField(x**Integer(2) + Integer(23), names=('a',)); (a,) = k._first_ngens(1) >>> Ok = k.ring_of_integers() >>> N = k.ideal(Integer(3)) >>> M = Matrix(Ok, Integer(2), lift_to_sl2_Ok(N, Integer(1), a)) >>> det(M) 1 >>> M = Matrix(Ok, Integer(2), lift_to_sl2_Ok(N, Integer(0), a)) >>> det(M) 1 >>> (M[Integer(1)][Integer(0)] in N) and (M[Integer(1)][Integer(1)] - a in N) True >>> M = Matrix(Ok, Integer(2), lift_to_sl2_Ok(N, Integer(0), Integer(0))) Traceback (most recent call last): ... ValueError: Cannot lift (0, 0) to an element of Sl2(Ok). - sage: k.<a> = NumberField(x^3 + 11) sage: Ok = k.ring_of_integers() sage: N = k.ideal(3, a - 1) sage: M = Matrix(Ok, 2, lift_to_sl2_Ok(N, 2*a, 0)) sage: det(M) 1 sage: (M[1][0] - 2*a in N) and (M[1][1] in N) True sage: M = Matrix(Ok, 2, lift_to_sl2_Ok(N, 4*a^2, a + 1)) sage: det(M) 1 sage: (M[1][0] - 4*a^2 in N) and (M[1][1] - (a+1) in N) True - >>> from sage.all import * >>> k = NumberField(x**Integer(3) + Integer(11), names=('a',)); (a,) = k._first_ngens(1) >>> Ok = k.ring_of_integers() >>> N = k.ideal(Integer(3), a - Integer(1)) >>> M = Matrix(Ok, Integer(2), lift_to_sl2_Ok(N, Integer(2)*a, Integer(0))) >>> det(M) 1 >>> (M[Integer(1)][Integer(0)] - Integer(2)*a in N) and (M[Integer(1)][Integer(1)] in N) True >>> M = Matrix(Ok, Integer(2), lift_to_sl2_Ok(N, Integer(4)*a**Integer(2), a + Integer(1))) >>> det(M) 1 >>> (M[Integer(1)][Integer(0)] - Integer(4)*a**Integer(2) in N) and (M[Integer(1)][Integer(1)] - (a+Integer(1)) in N) True - sage: k.<a> = NumberField(x^4 - x^3 -21*x^2 + 17*x + 133) sage: Ok = k.ring_of_integers() sage: N = k.ideal(7, a) sage: M = Matrix(Ok, 2, lift_to_sl2_Ok(N, 0, a^2 - 1)) sage: det(M) 1 sage: (M[1][0] in N) and (M[1][1] - (a^2-1) in N) True sage: M = Matrix(Ok, 2, lift_to_sl2_Ok(N, 0, 7)) Traceback (most recent call last): ... ValueError: <0> + <7> and the Fractional ideal (7, -4/7*a^3 + 13/7*a^2 + 39/7*a - 19) are not coprime. - >>> from sage.all import * >>> k = NumberField(x**Integer(4) - x**Integer(3) -Integer(21)*x**Integer(2) + Integer(17)*x + Integer(133), names=('a',)); (a,) = k._first_ngens(1) >>> Ok = k.ring_of_integers() >>> N = k.ideal(Integer(7), a) >>> M = Matrix(Ok, Integer(2), lift_to_sl2_Ok(N, Integer(0), a**Integer(2) - Integer(1))) >>> det(M) 1 >>> (M[Integer(1)][Integer(0)] in N) and (M[Integer(1)][Integer(1)] - (a**Integer(2)-Integer(1)) in N) True >>> M = Matrix(Ok, Integer(2), lift_to_sl2_Ok(N, Integer(0), Integer(7))) Traceback (most recent call last): ... ValueError: <0> + <7> and the Fractional ideal (7, -4/7*a^3 + 13/7*a^2 + 39/7*a - 19) are not coprime. 
- sage.modular.modsym.p1list_nf.make_coprime(N, c, d)[source]¶
- Return (c, d’) so d’ is congruent to d modulo N, and such that c and d’ are coprime (\(\langle c\rangle + \langle d'\rangle = R\)). - INPUT: - N– number field ideal
- c– integral element of the number field
- d– integral element of the number field
 - OUTPUT: - A pair \((c, d')\) where \(c\), \(d'\) are integral elements of the corresponding number field, with \(d'\) congruent to \(d\) mod \(N\), and such that \(\langle c\rangle + \langle d'\rangle = R\) (\(R\) being the corresponding ring of integers). - EXAMPLES: - sage: from sage.modular.modsym.p1list_nf import make_coprime sage: x = polygen(QQ, 'x') sage: k.<a> = NumberField(x^2 + 23) sage: N = k.ideal(3, a - 1) sage: c = 2*a; d = a + 1 sage: N.is_coprime(k.ideal(c, d)) True sage: k.ideal(c).is_coprime(d) False sage: c, dp = make_coprime(N, c, d) sage: k.ideal(c).is_coprime(dp) True - >>> from sage.all import * >>> from sage.modular.modsym.p1list_nf import make_coprime >>> x = polygen(QQ, 'x') >>> k = NumberField(x**Integer(2) + Integer(23), names=('a',)); (a,) = k._first_ngens(1) >>> N = k.ideal(Integer(3), a - Integer(1)) >>> c = Integer(2)*a; d = a + Integer(1) >>> N.is_coprime(k.ideal(c, d)) True >>> k.ideal(c).is_coprime(d) False >>> c, dp = make_coprime(N, c, d) >>> k.ideal(c).is_coprime(dp) True 
- sage.modular.modsym.p1list_nf.p1NFlist(N)[source]¶
- Return a list of the normalized elements of \(\mathbb{P}^1(R/N)\), where \(N\) is an integral ideal. - INPUT: - N– integral ideal (the level or modulus)
 - EXAMPLES: - sage: x = polygen(QQ, 'x') sage: k.<a> = NumberField(x^2 + 23) sage: N = k.ideal(3) sage: from sage.modular.modsym.p1list_nf import p1NFlist, psi sage: len(p1NFlist(N))==psi(N) True - >>> from sage.all import * >>> x = polygen(QQ, 'x') >>> k = NumberField(x**Integer(2) + Integer(23), names=('a',)); (a,) = k._first_ngens(1) >>> N = k.ideal(Integer(3)) >>> from sage.modular.modsym.p1list_nf import p1NFlist, psi >>> len(p1NFlist(N))==psi(N) True 
- sage.modular.modsym.p1list_nf.psi(N)[source]¶
- The index \([\Gamma : \Gamma_0(N)]\), where \(\Gamma = GL(2, R)\) for \(R\) the corresponding ring of integers, and \(\Gamma_0(N)\) standard congruence subgroup. - EXAMPLES: - sage: from sage.modular.modsym.p1list_nf import psi sage: x = polygen(QQ, 'x') sage: k.<a> = NumberField(x^2 + 23) sage: N = k.ideal(3, a - 1) sage: psi(N) 4 - >>> from sage.all import * >>> from sage.modular.modsym.p1list_nf import psi >>> x = polygen(QQ, 'x') >>> k = NumberField(x**Integer(2) + Integer(23), names=('a',)); (a,) = k._first_ngens(1) >>> N = k.ideal(Integer(3), a - Integer(1)) >>> psi(N) 4 - sage: k.<a> = NumberField(x^2 + 23) sage: N = k.ideal(5) sage: psi(N) 26 - >>> from sage.all import * >>> k = NumberField(x**Integer(2) + Integer(23), names=('a',)); (a,) = k._first_ngens(1) >>> N = k.ideal(Integer(5)) >>> psi(N) 26