Genus¶
-
sage.quadratic_forms.genera.genus.
Genus
(A)¶ Given a nonsingular symmetric matrix \(A\), return the genus of \(A\).
INPUT:
- \(A\) – a symmetric matrix with coefficients in \(\ZZ\)
OUTPUT:
A
GenusSymbol_global_ring
object, encoding the Conway-Sloane genus symbol of the quadratic form whose Gram matrix is \(A\).EXAMPLES:
sage: from sage.quadratic_forms.genera.genus import Genus sage: A = Matrix(ZZ, 2, 2, [1,1,1,2]) sage: Genus(A) Genus of [1 1] [1 2]
-
class
sage.quadratic_forms.genera.genus.
GenusSymbol_global_ring
(A, max_elem_divisors=None)¶ Bases:
object
This represents a collection of local genus symbols (at primes) and signature information which represent the genus of a non-degenerate integral lattice.
-
determinant
()¶ Returns the determinant of this genus, where the determinant is the Hessian determinant of the quadratic form whose Gram matrix is the Gram matrix giving rise to this global genus symbol.
INPUT:
None
OUTPUT:
an integer
EXAMPLES:
sage: from sage.quadratic_forms.genera.genus import GenusSymbol_global_ring sage: A = DiagonalQuadraticForm(ZZ, [1,-2,3,4]).Hessian_matrix() sage: GS = GenusSymbol_global_ring(A) sage: GS.determinant() -384
-
signature_pair_of_matrix
()¶ Returns the signature pair (p, n) of the (non-degenerate) global genus symbol, where p is the number of positive eigenvalues and n is the number of negative eigenvalues.
INPUT:
None
OUTPUT:
a pair of integers (p, n) each >= 0
EXAMPLES:
sage: from sage.quadratic_forms.genera.genus import GenusSymbol_global_ring sage: A = DiagonalQuadraticForm(ZZ, [1,-2,3,4,8,-11]).Hessian_matrix() sage: GS = GenusSymbol_global_ring(A) sage: GS.signature_pair_of_matrix() (4, 2)
-
-
class
sage.quadratic_forms.genera.genus.
Genus_Symbol_p_adic_ring
(prime, symbol, check=True)¶ Bases:
object
Local genus symbol over a p-adic ring.
-
canonical_symbol
()¶ Return (and cache) the canonical p-adic genus symbol. This is only really affects the 2-adic symbol, since when p > 2 the symbol is already canonical.
INPUT:
None
OUTPUT:
a list of lists of integers
EXAMPLES:
sage: from sage.quadratic_forms.genera.genus import p_adic_symbol sage: from sage.quadratic_forms.genera.genus import Genus_Symbol_p_adic_ring sage: A = Matrix(ZZ, 2, 2, [1,1,1,2]) sage: p = 2 sage: G2 = Genus_Symbol_p_adic_ring(p, p_adic_symbol(A, p, 2)); G2 Genus symbol at 2 : [[0, 2, 1, 1, 2]] sage: G2.canonical_symbol() [[0, 2, 1, 1, 2]] sage: A = Matrix(ZZ, 2, 2, [1,0,0,2]) sage: p = 2 sage: G2 = Genus_Symbol_p_adic_ring(p, p_adic_symbol(A, p, 2)); G2 Genus symbol at 2 : [[0, 1, 1, 1, 1], [1, 1, 1, 1, 1]] sage: G2.canonical_symbol() ## Oddity fusion occurred here! [[0, 1, 1, 1, 2], [1, 1, 1, 1, 0]] sage: A = DiagonalQuadraticForm(ZZ, [1,2,3,4]).Hessian_matrix() sage: p = 2 sage: G2 = Genus_Symbol_p_adic_ring(p, p_adic_symbol(A, p, 2)); G2 Genus symbol at 2 : [[1, 2, 3, 1, 4], [2, 1, 1, 1, 1], [3, 1, 1, 1, 1]] sage: G2.canonical_symbol() ## Oddity fusion occurred here! [[1, 2, -1, 1, 6], [2, 1, 1, 1, 0], [3, 1, 1, 1, 0]] sage: A = Matrix(ZZ, 2, 2, [2,1,1,2]) sage: p = 2 sage: G2 = Genus_Symbol_p_adic_ring(p, p_adic_symbol(A, p, 2)); G2 Genus symbol at 2 : [[0, 2, 3, 0, 0]] sage: G2.canonical_symbol() [[0, 2, -1, 0, 0]] sage: A = DiagonalQuadraticForm(ZZ, [1,2,3,4]).Hessian_matrix() sage: p = 3 sage: G3 = Genus_Symbol_p_adic_ring(p, p_adic_symbol(A, p, 2)); G3 Genus symbol at 3 : [[0, 3, 1], [1, 1, -1]] sage: G3.canonical_symbol() [[0, 3, 1], [1, 1, -1]]
Note
See Conway-Sloane 3rd edition, pp. 381-382 for definitions and examples.
Todo
Add an example where sign walking occurs!
-
compartments
()¶ Compute the indices for each of the compartments in this local genus symbol if it is associated to the prime p=2 (and raise an error for all other primes).
INPUT:
None
OUTPUT:
a list of integers >= 0
EXAMPLES:
sage: from sage.quadratic_forms.genera.genus import p_adic_symbol sage: from sage.quadratic_forms.genera.genus import Genus_Symbol_p_adic_ring sage: A = DiagonalQuadraticForm(ZZ, [1,2,3,4]).Hessian_matrix() sage: p = 2 sage: G2 = Genus_Symbol_p_adic_ring(p, p_adic_symbol(A, p, 2)); G2 Genus symbol at 2 : [[1, 2, 3, 1, 4], [2, 1, 1, 1, 1], [3, 1, 1, 1, 1]] sage: G2.compartments() [[0, 1, 2]]
-
determinant
()¶ Returns the (p-part of the) determinant (square-class) of the Hessian matrix of the quadratic form (given by regarding the integral symmetric matrix which generated this genus symbol as the Gram matrix of Q) associated to this local genus symbol.
INPUT:
None
OUTPUT:
an integer
EXAMPLES:
sage: from sage.quadratic_forms.genera.genus import p_adic_symbol sage: from sage.quadratic_forms.genera.genus import Genus_Symbol_p_adic_ring sage: A = DiagonalQuadraticForm(ZZ, [1,2,3,4]).Hessian_matrix() sage: p = 2 sage: G2 = Genus_Symbol_p_adic_ring(p, p_adic_symbol(A, p, 2)); G2 Genus symbol at 2 : [[1, 2, 3, 1, 4], [2, 1, 1, 1, 1], [3, 1, 1, 1, 1]] sage: G2.determinant() 128 sage: A = DiagonalQuadraticForm(ZZ, [1,2,3,4]).Hessian_matrix() sage: p = 3 sage: G3 = Genus_Symbol_p_adic_ring(p, p_adic_symbol(A, p, 2)); G3 Genus symbol at 3 : [[0, 3, 1], [1, 1, -1]] sage: G3.determinant() 3
-
dimension
()¶ Returns the dimension of a quadratic form associated to this genus symbol.
INPUT:
None
OUTPUT:
an integer >= 0
EXAMPLES:
sage: from sage.quadratic_forms.genera.genus import p_adic_symbol sage: from sage.quadratic_forms.genera.genus import Genus_Symbol_p_adic_ring sage: A = DiagonalQuadraticForm(ZZ, [1,2,3,4]).Hessian_matrix() sage: p = 2 sage: G2 = Genus_Symbol_p_adic_ring(p, p_adic_symbol(A, p, 2)); G2 Genus symbol at 2 : [[1, 2, 3, 1, 4], [2, 1, 1, 1, 1], [3, 1, 1, 1, 1]] sage: G2.dimension() 4 sage: A = DiagonalQuadraticForm(ZZ, [1,2,3,4]).Hessian_matrix() sage: p = 3 sage: G3 = Genus_Symbol_p_adic_ring(p, p_adic_symbol(A, p, 2)); G3 Genus symbol at 3 : [[0, 3, 1], [1, 1, -1]] sage: G3.dimension() 4
-
excess
()¶ Returns the p-excess of the quadratic form whose Hessian matrix is the symmetric matrix A. When p = 2 the p-excess is called the oddity.
WARNING/NOTE: This normalization seems non-standard, and we should review this entire class to make sure that we have our doubling conventions straight throughout!
REFERENCE:
Conway and Sloane Book, 3rd edition, pp 370-371.
INPUT:
None
OUTPUT:
an integer
EXAMPLES:
sage: from sage.quadratic_forms.genera.genus import p_adic_symbol sage: from sage.quadratic_forms.genera.genus import Genus_Symbol_p_adic_ring sage: AC = diagonal_matrix(ZZ, [1,3,-3]) sage: p=2; Genus_Symbol_p_adic_ring(p, p_adic_symbol(AC, p, 2)).excess() 1 sage: p=3; Genus_Symbol_p_adic_ring(p, p_adic_symbol(AC, p, 2)).excess() 0 sage: p=5; Genus_Symbol_p_adic_ring(p, p_adic_symbol(AC, p, 2)).excess() 0 sage: p=7; Genus_Symbol_p_adic_ring(p, p_adic_symbol(AC, p, 2)).excess() 0 sage: p=11; Genus_Symbol_p_adic_ring(p, p_adic_symbol(AC, p, 2)).excess() 0 sage: AC = 2 * diagonal_matrix(ZZ, [1,3,-3]) sage: p=2; Genus_Symbol_p_adic_ring(p, p_adic_symbol(AC, p, 2)).excess() 1 sage: p=3; Genus_Symbol_p_adic_ring(p, p_adic_symbol(AC, p, 2)).excess() 0 sage: p=5; Genus_Symbol_p_adic_ring(p, p_adic_symbol(AC, p, 2)).excess() 0 sage: p=7; Genus_Symbol_p_adic_ring(p, p_adic_symbol(AC, p, 2)).excess() 0 sage: p=11; Genus_Symbol_p_adic_ring(p, p_adic_symbol(AC, p, 2)).excess() 0 sage: A = 2*diagonal_matrix(ZZ, [1,2,3,4]) sage: p=2; Genus_Symbol_p_adic_ring(p, p_adic_symbol(A, p, 2)).excess() 2 sage: p=3; Genus_Symbol_p_adic_ring(p, p_adic_symbol(A, p, 2)).excess() 6 sage: p=5; Genus_Symbol_p_adic_ring(p, p_adic_symbol(A, p, 2)).excess() 0 sage: p=7; Genus_Symbol_p_adic_ring(p, p_adic_symbol(A, p, 2)).excess() 0 sage: p=11; Genus_Symbol_p_adic_ring(p, p_adic_symbol(A, p, 2)).excess() 0
-
number_of_blocks
()¶ Returns the number of positive dimensional symbols/Jordan blocks
INPUT:
None
OUTPUT:
integer >= 0
EXAMPLES:
sage: from sage.quadratic_forms.genera.genus import p_adic_symbol sage: from sage.quadratic_forms.genera.genus import Genus_Symbol_p_adic_ring sage: A = DiagonalQuadraticForm(ZZ, [1,2,3,4]).Hessian_matrix() sage: p = 2 sage: G2 = Genus_Symbol_p_adic_ring(p, p_adic_symbol(A, p, 2)); G2 Genus symbol at 2 : [[1, 2, 3, 1, 4], [2, 1, 1, 1, 1], [3, 1, 1, 1, 1]] sage: G2.number_of_blocks() 3 sage: A = DiagonalQuadraticForm(ZZ, [1,2,3,4]).Hessian_matrix() sage: p = 3 sage: G3 = Genus_Symbol_p_adic_ring(p, p_adic_symbol(A, p, 2)); G3 Genus symbol at 3 : [[0, 3, 1], [1, 1, -1]] sage: G3.number_of_blocks() 2
-
rank
()¶ Returns the dimension of a quadratic form associated to this genus symbol.
Todo
DELETE THIS METHOD IN FAVOR OF THE dimension() METHOD BELOW!
INPUT:
None
OUTPUT:
an integer >= 0
EXAMPLES:
sage: from sage.quadratic_forms.genera.genus import p_adic_symbol sage: from sage.quadratic_forms.genera.genus import Genus_Symbol_p_adic_ring sage: A = DiagonalQuadraticForm(ZZ, [1,2,3,4]).Hessian_matrix() sage: p = 2 sage: G2 = Genus_Symbol_p_adic_ring(p, p_adic_symbol(A, p, 2)); G2 Genus symbol at 2 : [[1, 2, 3, 1, 4], [2, 1, 1, 1, 1], [3, 1, 1, 1, 1]] sage: G2.rank() 4 sage: A = DiagonalQuadraticForm(ZZ, [1,2,3,4]).Hessian_matrix() sage: p = 3 sage: G3 = Genus_Symbol_p_adic_ring(p, p_adic_symbol(A, p, 2)); G3 Genus symbol at 3 : [[0, 3, 1], [1, 1, -1]] sage: G3.rank() 4
-
symbol_tuple_list
()¶ Returns the underlying list of lists of integers defining the genus symbol.
INPUT:
None
OUTPUT:
list of lists of integers
EXAMPLES:
sage: from sage.quadratic_forms.genera.genus import p_adic_symbol sage: from sage.quadratic_forms.genera.genus import Genus_Symbol_p_adic_ring sage: A = DiagonalQuadraticForm(ZZ, [1,2,3,4]).Hessian_matrix() sage: p = 3 sage: G3 = Genus_Symbol_p_adic_ring(p, p_adic_symbol(A, p, 2)); G3 Genus symbol at 3 : [[0, 3, 1], [1, 1, -1]] sage: G3.symbol_tuple_list() [[0, 3, 1], [1, 1, -1]] sage: type(G3.symbol_tuple_list()) <type 'list'> sage: A = DiagonalQuadraticForm(ZZ, [1,2,3,4]).Hessian_matrix() sage: p = 2 sage: G2 = Genus_Symbol_p_adic_ring(p, p_adic_symbol(A, p, 2)); G2 Genus symbol at 2 : [[1, 2, 3, 1, 4], [2, 1, 1, 1, 1], [3, 1, 1, 1, 1]] sage: G2.symbol_tuple_list() [[1, 2, 3, 1, 4], [2, 1, 1, 1, 1], [3, 1, 1, 1, 1]] sage: type(G2.symbol_tuple_list()) <type 'list'>
-
trains
()¶ Compute the indices for each of the trains in this local genus symbol if it is associated to the prime p=2 (and raise an error for all other primes).
INPUT:
None
OUTPUT:
a list of integers >= 0
EXAMPLES:
sage: from sage.quadratic_forms.genera.genus import p_adic_symbol sage: from sage.quadratic_forms.genera.genus import Genus_Symbol_p_adic_ring sage: A = DiagonalQuadraticForm(ZZ, [1,2,3,4]).Hessian_matrix() sage: p = 2 sage: G2 = Genus_Symbol_p_adic_ring(p, p_adic_symbol(A, p, 2)); G2 Genus symbol at 2 : [[1, 2, 3, 1, 4], [2, 1, 1, 1, 1], [3, 1, 1, 1, 1]] sage: G2.trains() [[0, 1, 2]]
-
-
sage.quadratic_forms.genera.genus.
LocalGenusSymbol
(A, p)¶ Given a nonsingular symmetric matrix A, return the local symbol of A at the prime p.
INPUT:
- A – a symmetric matrix with coefficients in ZZ
- p – an integer prime p > 0
OUTPUT:
A Genus_Symbol_p_adic_ring object, encoding the Conway-Sloane genus symbol at p of the quadratic form whose Gram matrix is A.
EXAMPLES:
sage: from sage.quadratic_forms.genera.genus import LocalGenusSymbol sage: A = Matrix(ZZ, 2, 2, [1,1,1,2]) sage: LocalGenusSymbol(A, 2) Genus symbol at 2 : [[0, 2, 1, 1, 2]] sage: LocalGenusSymbol(A, 3) Genus symbol at 3 : [[0, 2, 1]] sage: A = Matrix(ZZ, 2, 2, [1,0,0,2]) sage: LocalGenusSymbol(A, 2) Genus symbol at 2 : [[0, 1, 1, 1, 1], [1, 1, 1, 1, 1]] sage: LocalGenusSymbol(A, 3) Genus symbol at 3 : [[0, 2, -1]]
-
sage.quadratic_forms.genera.genus.
basis_complement
(B)¶ Given an echelonized basis matrix (over a field), calculate a matrix whose rows form a basis complement (to the rows of B).
INPUT:
- B – matrix over a field in row echelon form
OUTPUT:
a rectangular matrix over a field
EXAMPLES:
sage: from sage.quadratic_forms.genera.genus import basis_complement sage: A = Matrix(ZZ, 2, 2, [1,1,1,1]) sage: B = A.kernel().echelonized_basis_matrix(); B [ 1 -1] sage: basis_complement(B) [0 1]
-
sage.quadratic_forms.genera.genus.
canonical_2_adic_compartments
(genus_symbol_quintuple_list)¶ Given a 2-adic local symbol (as the underlying list of quintuples) this returns a list of lists of indices of the genus_symbol_quintuple_list which are in the same compartment. A compartment is defined to be a maximal interval of Jordan components all (scaled) of type I (i.e. odd).
INPUT:
- genus_symbol_quintuple_list – a quintuple of integers (with certain restrictions).
OUTPUT:
a list of lists of integers.
EXAMPLES:
sage: from sage.quadratic_forms.genera.genus import LocalGenusSymbol sage: from sage.quadratic_forms.genera.genus import canonical_2_adic_compartments sage: A = Matrix(ZZ, 2, 2, [1,1,1,2]) sage: G2 = LocalGenusSymbol(A, 2); G2 Genus symbol at 2 : [[0, 2, 1, 1, 2]] sage: canonical_2_adic_compartments(G2.symbol_tuple_list()) [[0]] sage: A = Matrix(ZZ, 2, 2, [1,0,0,2]) sage: G2 = LocalGenusSymbol(A, 2); G2 Genus symbol at 2 : [[0, 1, 1, 1, 1], [1, 1, 1, 1, 1]] sage: canonical_2_adic_compartments(G2.symbol_tuple_list()) [[0, 1]] sage: A = DiagonalQuadraticForm(ZZ, [1,2,3,4]).Hessian_matrix() sage: G2 = LocalGenusSymbol(A, 2); G2 Genus symbol at 2 : [[1, 2, 3, 1, 4], [2, 1, 1, 1, 1], [3, 1, 1, 1, 1]] sage: canonical_2_adic_compartments(G2.symbol_tuple_list()) [[0, 1, 2]] sage: A = Matrix(ZZ, 2, 2, [2,1,1,2]) sage: G2 = LocalGenusSymbol(A, 2); G2 Genus symbol at 2 : [[0, 2, 3, 0, 0]] sage: canonical_2_adic_compartments(G2.symbol_tuple_list()) ## No compartments here! []
NOTES:
See Conway-Sloane 3rd edition, pp. 381-382 for definitions and examples.
-
sage.quadratic_forms.genera.genus.
canonical_2_adic_reduction
(genus_symbol_quintuple_list)¶ Given a 2-adic local symbol (as the underlying list of quintuples) this returns a canonical 2-adic symbol (again as a raw list of quintuples of integers) which has at most one minus sign per train and this sign appears on the smallest dimensional Jordan component in each train. This results from applying the “sign-walking” and “oddity fusion” equivalences.
INPUT:
- genus_symbol_quintuple_list – a quintuple of integers (with certain restrictions).
- compartments – a list of lists of distinct integers (optional)
OUTPUT:
a list of lists of distinct integers.
EXAMPLES:
sage: from sage.quadratic_forms.genera.genus import LocalGenusSymbol sage: from sage.quadratic_forms.genera.genus import canonical_2_adic_reduction sage: A = Matrix(ZZ, 2, 2, [1,1,1,2]) sage: G2 = LocalGenusSymbol(A, 2); G2 Genus symbol at 2 : [[0, 2, 1, 1, 2]] sage: canonical_2_adic_reduction(G2.symbol_tuple_list()) [[0, 2, 1, 1, 2]] sage: A = Matrix(ZZ, 2, 2, [1,0,0,2]) sage: G2 = LocalGenusSymbol(A, 2); G2 Genus symbol at 2 : [[0, 1, 1, 1, 1], [1, 1, 1, 1, 1]] sage: canonical_2_adic_reduction(G2.symbol_tuple_list()) ## Oddity fusion occurred here! [[0, 1, 1, 1, 2], [1, 1, 1, 1, 0]] sage: A = DiagonalQuadraticForm(ZZ, [1,2,3,4]).Hessian_matrix() sage: G2 = LocalGenusSymbol(A, 2); G2 Genus symbol at 2 : [[1, 2, 3, 1, 4], [2, 1, 1, 1, 1], [3, 1, 1, 1, 1]] sage: canonical_2_adic_reduction(G2.symbol_tuple_list()) ## Oddity fusion occurred here! [[1, 2, -1, 1, 6], [2, 1, 1, 1, 0], [3, 1, 1, 1, 0]] sage: A = Matrix(ZZ, 2, 2, [2,1,1,2]) sage: G2 = LocalGenusSymbol(A, 2); G2 Genus symbol at 2 : [[0, 2, 3, 0, 0]] sage: canonical_2_adic_reduction(G2.symbol_tuple_list()) [[0, 2, -1, 0, 0]]
Note
See Conway-Sloane 3rd edition, pp. 381-382 for definitions and examples.
Todo
Add an example where sign walking occurs!
-
sage.quadratic_forms.genera.genus.
canonical_2_adic_trains
(genus_symbol_quintuple_list, compartments=None)¶ Given a 2-adic local symbol (as the underlying list of quintuples) this returns a list of lists of indices of the genus_symbol_quintuple_list which are in the same train. A train is defined to be a maximal interval of Jordan components so that at least one of each adjacent pair (allowing zero-dimensional Jordan components) is (scaled) of type I (i.e. odd).
INPUT:
- genus_symbol_quintuple_list – a quintuple of integers (with certain restrictions).
- compartments – a list of lists of distinct integers (optional)
OUTPUT:
a list of lists of distinct integers.
EXAMPLES:
sage: from sage.quadratic_forms.genera.genus import LocalGenusSymbol sage: from sage.quadratic_forms.genera.genus import canonical_2_adic_compartments sage: from sage.quadratic_forms.genera.genus import canonical_2_adic_trains sage: A = Matrix(ZZ, 2, 2, [1,1,1,2]) sage: G2 = LocalGenusSymbol(A, 2); G2 Genus symbol at 2 : [[0, 2, 1, 1, 2]] sage: c = canonical_2_adic_compartments(G2.symbol_tuple_list()); c [[0]] sage: canonical_2_adic_trains(G2.symbol_tuple_list(), c) [[0]] sage: A = Matrix(ZZ, 2, 2, [1,0,0,2]) sage: G2 = LocalGenusSymbol(A, 2); G2 Genus symbol at 2 : [[0, 1, 1, 1, 1], [1, 1, 1, 1, 1]] sage: c = canonical_2_adic_compartments(G2.symbol_tuple_list()); c [[0, 1]] sage: canonical_2_adic_trains(G2.symbol_tuple_list(), c) [[0, 1]] sage: canonical_2_adic_trains(G2.symbol_tuple_list()) [[0, 1]] sage: A = DiagonalQuadraticForm(ZZ, [1,2,3,4]).Hessian_matrix() sage: G2 = LocalGenusSymbol(A, 2); G2 Genus symbol at 2 : [[1, 2, 3, 1, 4], [2, 1, 1, 1, 1], [3, 1, 1, 1, 1]] sage: c = canonical_2_adic_compartments(G2.symbol_tuple_list()); c [[0, 1, 2]] sage: canonical_2_adic_trains(G2.symbol_tuple_list()) [[0, 1, 2]] sage: A = Matrix(ZZ, 2, 2, [2,1,1,2]) sage: G2 = LocalGenusSymbol(A, 2); G2 Genus symbol at 2 : [[0, 2, 3, 0, 0]] sage: c = canonical_2_adic_compartments(G2.symbol_tuple_list()); c ## No compartments here! [] sage: canonical_2_adic_trains(G2.symbol_tuple_list()) []
Note
See Conway-Sloane 3rd edition, pp. 381-382 for definitions and examples.
Todo
Add a non-trivial example in the doctest here!
-
sage.quadratic_forms.genera.genus.
is_2_adic_genus
(genus_symbol_quintuple_list)¶ Given a 2-adic local symbol (as the underlying list of quintuples) check whether it is the 2-adic symbol of a 2-adic form.
INPUT:
- genus_symbol_quintuple_list – a quintuple of integers (with certain restrictions).
OUTPUT:
boolean
EXAMPLES:
sage: from sage.quadratic_forms.genera.genus import LocalGenusSymbol sage: A = Matrix(ZZ, 2, 2, [1,1,1,2]) sage: G2 = LocalGenusSymbol(A, 2) sage: is_2_adic_genus(G2.symbol_tuple_list()) True sage: A = Matrix(ZZ, 2, 2, [1,1,1,2]) sage: G3 = LocalGenusSymbol(A, 3) sage: is_2_adic_genus(G3.symbol_tuple_list()) ## This raises an error Traceback (most recent call last): ... TypeError: The genus symbols are not quintuples, so it's not a genus symbol for the prime p=2. sage: A = Matrix(ZZ, 2, 2, [1,0,0,2]) sage: G2 = LocalGenusSymbol(A, 2) sage: is_2_adic_genus(G2.symbol_tuple_list()) True
-
sage.quadratic_forms.genera.genus.
is_GlobalGenus
(G)¶ Given a genus symbol G (specified by a collection of local symbols), return True in G represents the genus of a global quadratic form or lattice.
INPUT:
- G – GenusSymbol_global_ring object
OUTPUT:
boolean
EXAMPLES:
sage: from sage.quadratic_forms.genera.genus import GenusSymbol_global_ring sage: from sage.quadratic_forms.genera.genus import Genus, is_GlobalGenus sage: A = Matrix(ZZ, 2, 2, [1,1,1,2]) sage: G = Genus(A) sage: is_GlobalGenus(G) True
-
sage.quadratic_forms.genera.genus.
is_even_matrix
(A)¶ Determines if the integral symmetric matrix A is even (i.e. represents only even numbers). If not, then it returns the index of an odd diagonal entry. If it is even, then we return the index -1.
INPUT:
- A – symmetric integer matrix
OUTPUT:
a pair of the form (boolean, integer)
EXAMPLES:
sage: from sage.quadratic_forms.genera.genus import is_even_matrix sage: A = Matrix(ZZ, 2, 2, [1,1,1,1]) sage: is_even_matrix(A) (False, 0) sage: A = Matrix(ZZ, 2, 2, [2,1,1,2]) sage: is_even_matrix(A) (True, -1)
-
sage.quadratic_forms.genera.genus.
p_adic_symbol
(A, p, val)¶ Given a symmetric matrix A and prime p, return the genus symbol at p.
val = valuation of the maximal elementary divisor of A needed to obtain enough precision calculation is modulo p to the val+3
Todo
Some description of the definition of the genus symbol.
INPUT:
- A – symmetric matrix with integer coefficients
- p – prime number > 0
- val – integer >= 0
OUTPUT:
a list of lists of integers
EXAMPLES:
sage: from sage.quadratic_forms.genera.genus import p_adic_symbol sage: A = DiagonalQuadraticForm(ZZ, [1,2,3,4]).Hessian_matrix() sage: p_adic_symbol(A, 2, 2) [[1, 2, 3, 1, 4], [2, 1, 1, 1, 1], [3, 1, 1, 1, 1]] sage: p_adic_symbol(A, 3, 1) [[0, 3, 1], [1, 1, -1]]
-
sage.quadratic_forms.genera.genus.
signature_pair_of_matrix
(A)¶ Computes the signature pair (p, n) of a non-degenerate symmetric matrix, where
- p = number of positive eigenvalues of A
- n = number of negative eigenvalues of A
INPUT:
- A – symmetric matrix (assumed to be non-degenerate)
OUTPUT:
a pair (tuple) of integers.
EXAMPLES:
sage: from sage.quadratic_forms.genera.genus import signature_pair_of_matrix sage: A = Matrix(ZZ, 2, 2, [-1,0,0,3]) sage: signature_pair_of_matrix(A) (1, 1) sage: A = Matrix(ZZ, 2, 2, [-1,1,1,7]) sage: signature_pair_of_matrix(A) (1, 1) sage: A = Matrix(ZZ, 2, 2, [3,1,1,7]) sage: signature_pair_of_matrix(A) (2, 0) sage: A = Matrix(ZZ, 2, 2, [-3,1,1,-11]) sage: signature_pair_of_matrix(A) (0, 2) sage: A = Matrix(ZZ, 2, 2, [1,1,1,1]) sage: signature_pair_of_matrix(A) Traceback (most recent call last): ... ArithmeticError: given matrix is not invertible
-
sage.quadratic_forms.genera.genus.
split_odd
(A)¶ Given a non-degenerate Gram matrix A (mod 8), return a splitting [u] + B such that u is odd and B is not even.
INPUT:
- A – an odd symmetric matrix with integer coefficients (which admits a splitting as above).
OUTPUT:
a pair (u, B) consisting of an odd integer u and an odd integral symmetric matrix B.
EXAMPLES:
sage: from sage.quadratic_forms.genera.genus import is_even_matrix sage: from sage.quadratic_forms.genera.genus import split_odd sage: A = Matrix(ZZ, 2, 2, [1,2,2,3]) sage: is_even_matrix(A) (False, 0) sage: split_odd(A) (1, [-1]) sage: A = Matrix(ZZ, 2, 2, [1,2,2,5]) sage: split_odd(A) (1, [1]) sage: A = Matrix(ZZ, 2, 2, [1,1,1,1]) sage: is_even_matrix(A) (False, 0) sage: split_odd(A) ## This fails because no such splitting exists. =( Traceback (most recent call last): ... RuntimeError: The matrix A does not admit a non-even splitting. sage: A = Matrix(ZZ, 2, 2, [1,2,2,6]) sage: split_odd(A) ## This fails because no such splitting exists. =( Traceback (most recent call last): ... RuntimeError: The matrix A does not admit a non-even splitting.
-
sage.quadratic_forms.genera.genus.
trace_diag_mod_8
(A)¶ Return the trace of the diagonalised form of A of an integral symmetric matrix which is diagonalizable mod 8. (Note that since the Jordan decomposition into blocks of size <= 2 is not unique here, this is not the same as saying that A is always diagonal in any 2-adic Jordan decomposition!)
INPUT:
- A – symmetric matrix with coefficients in Z which is odd in Z/2Z and has determinant not divisible by 8.
OUTPUT:
an integer
EXAMPLES:
sage: from sage.quadratic_forms.genera.genus import is_even_matrix sage: from sage.quadratic_forms.genera.genus import split_odd sage: from sage.quadratic_forms.genera.genus import trace_diag_mod_8 sage: A = Matrix(ZZ, 2, 2, [1,2,2,3]) sage: is_even_matrix(A) (False, 0) sage: split_odd(A) (1, [-1]) sage: trace_diag_mod_8(A) 0 sage: A = Matrix(ZZ, 2, 2, [1,2,2,5]) sage: split_odd(A) (1, [1]) sage: trace_diag_mod_8(A) 2
-
sage.quadratic_forms.genera.genus.
two_adic_symbol
(A, val)¶ Given a symmetric matrix A and prime p, return the genus symbol at p.
val = valuation of maximal 2-elementary divisor
The genus symbol of a component 2^m*f is of the form (m,n,s,d[,o]), where
- m = valuation of the component
- n = dimension of f
- d = det(f) in {1,3,5,7}
- s = 0 (or 1) if even (or odd)
- o = oddity of f (= 0 if s = 0) in Z/8Z
INPUT:
- A – symmetric matrix with integer coefficients
- val – integer >=0
OUTPUT:
a list of lists of integers (representing a Conway-Sloane 2-adic symbol)
EXAMPLES:
sage: from sage.quadratic_forms.genera.genus import two_adic_symbol sage: A = diagonal_matrix(ZZ, [1,2,3,4]) sage: two_adic_symbol(A, 2) [[0, 2, 3, 1, 4], [1, 1, 1, 1, 1], [2, 1, 1, 1, 1]]