Integer Lattices

An integer lattice $L$ is a finitely generated $\mathbb{Z}$-submodule of a quadratic vector space $V = \mathbb{Q}^n$ over the rational numbers. Integer lattices are also known as quadratic forms over the integers. We will refer to them as $\mathbb{Z}$-lattices.

A $\mathbb{Z}$-lattice $L$ has the type ZLat. It is given in terms of its ambient quadratic space $V$ together with a basis matrix $B$ whose rows span $L$, i.e. $L = \mathbb{Z}^r B$ where $r$ is the ($\mathbb{Z}$-module) rank of $L$.

To access $V$ and $B$ see ambient_space(L::ZLat) and basis_matrix(L::ZLat).

Creation of integer lattices

From a gram matrix

ZlatticeMethod
Zlattice([B::MatElem]; gram) -> ZLat

Return the Z-lattice with basis matrix $B$ inside the quadratic space with Gram matrix gram.

If the keyword gram is not specified, the Gram matrix is the identity matrix. If $B$ is not specified, the basis matrix is the identity matrix.

Examples

julia> L = Zlattice(matrix(QQ, 2, 2, [1//2, 0, 0, 2]));

julia> gram_matrix(L) == matrix(QQ, 2, 2, [1//4, 0, 0, 4])
true

julia> L = Zlattice(gram = matrix(ZZ, [2 -1; -1 2]));

julia> gram_matrix(L) == matrix(ZZ, [2 -1; -1 2])
true

In a quadratic space

latticeMethod
lattice(V::AbsSpace, basis::MatElem ; check::Bool = true) -> AbsLat

Given an ambient space V and a matrix basis, return the lattice spanned by the rows of basis inside V. If V is hermitian (resp. quadratic) then the output is a hermitian (resp. quadratic) lattice.

By default, basis is checked to be of full rank. This test can be disabled by setting check to false.

Special lattices

root_latticeMethod
root_lattice(R::Symbol, n::Int)

Return the root lattice of type R given by :A, :D or :E with parameter n.

hyperbolic_plane_latticeMethod
hyperbolic_plane_lattice(n::RationalUnion = 1)

Return the hyperbolic plane with intersection form of scale n, that is, the unique (up to isometry) even unimodular hyperbolic $\mathbb Z$-lattice of rank 2, rescaled by n.

Examples

julia> L = hyperbolic_plane_lattice(6);

julia> gram_matrix(L)
[0   6]
[6   0]

julia> L = hyperbolic_plane_lattice(ZZ(-13));

julia> gram_matrix(L)
[  0   -13]
[-13     0]
ZlatticeMethod
Zlattice(S::Symbol, n::RationalUnion = 1) -> Zlat

Given S = :H or S = :U, return a $\mathbb Z$-lattice admitting $n*J_2$ as Gram matrix in some basis, where $J_2$ is the 2-by-2 matrix with 0's on the main diagonal and 1's elsewhere.

leech_latticeFunction
leech_lattice()

Return the Leech lattice.

leech_lattice(niemeier_lattice::ZLat) -> Leech, neighbor vector, index

Return a triple L, v, h where L is the Leech lattice.

L is an h-neighbor of the Niemeier lattice N with respect to v. This means that L / L ∩ N ≅ ℤ / h ℤ. Here h is the Coxeter number of the Niemeier lattice.

This implements the 23 holy constructions of the Leech lattice in J. H. Conway, N. J. A. Sloane (1999).

Examples

julia> R = Zlattice(gram=2 * identity_matrix(ZZ, 24));

julia> N = maximal_even_lattice(R) # Some Niemeier lattice
Quadratic lattice of rank 24 and degree 24 over the rationals

julia> minimum(N)
2

julia> det(N)
1

julia> L, v, h = leech_lattice(N);

julia> minimum(L)
4

julia> det(L)
1

julia> h == index(L, intersect(L, N))
true

We illustrate how the Leech lattice is constructed from N, h and v.

julia> Zmodh = ResidueRing(ZZ, h);

julia> V = ambient_space(N);

julia> vG = map_entries(x->Zmodh(ZZ(x)), inner_product(V, v, basis_matrix(N)));

julia> LN = transpose(lift(kernel(vG)[2]))*basis_matrix(N); # vectors whose inner product with `v` is divisible by `h`.

julia> lattice(V, LN) == intersect(L, N)
true

julia> gensL = vcat(LN, 1//h * v);

julia> lattice(V, gensL, isbasis=false) == L
true

From a genus

Integer lattices can be created as representatives of a genus. See (representative(L::ZGenus))

Rescaling the Quadratic Form

rescaleMethod
rescale(L::ZLat, r::RationalUnion) -> ZLat

Return the lattice L in the quadratic space with form r \Phi.

Examples

This can be useful to apply methods intended for positive definite lattices.

julia> L = Zlattice(gram=ZZ[-1 0; 0 -1])
Quadratic lattice of rank 2 and degree 2 over the rationals

julia> shortest_vectors(rescale(L, -1))
2-element Vector{Vector{fmpz}}:
 [0, 1]
 [1, 0]

Attributes

ambient_spaceMethod
ambient_space(L::AbsLat) -> AbsSpace

Return the ambient space of the lattice L. If the ambient space is not known, an error is raised.

basis_matrixMethod
basis_matrix(L::ZLat)

Return the basis matrix $B$ of the integer lattice $L$.

The lattice is given by the row span of $B$ seen inside of the ambient quadratic space of $L$.

gram_matrixMethod
gram_matrix(L::ZLat) -> fmpq_mat

Return the gram matrix of $L$.

Examples

julia> L = Zlattice(matrix(ZZ, [2 0; -1 2]));

julia> gram_matrix(L)
[ 4   -2]
[-2    5]
rational_spanMethod
rational_span(L::ZLat) -> QuadSpace

Return the rational span of $L$, which is the quadratic space with Gram matrix equal to gram_matrix(L).

Examples

julia> L = Zlattice(matrix(ZZ, [2 0; -1 2]));

julia> rational_span(L)
Quadratic space over
Rational Field
with Gram matrix
[4 -2; -2 5]
base_ringMethod
base_ring(M::PMat)

The PMat $M$ defines an $R$-module for some maximal order $R$. This function returns the $R$ that was used to defined $M$.

base_ring(L::AbsLat) -> Ring

Return the order over which the lattice L is defined.

base_ring(I::MPolyIdeal)

Return the ambient ring of I.

Examples

julia> R, (x, y) = PolynomialRing(QQ, ["x", "y"])
(Multivariate Polynomial Ring in x, y over Rational Field, fmpq_mpoly[x, y])

julia> I = ideal(R, [x, y])^2
ideal(x^2, x*y, y^2)

julia> base_ring(I)
Multivariate Polynomial Ring in x, y over Rational Field
source
base_ring(X::AbsSpec)

On an affine scheme $X/𝕜$ over $𝕜$ this returns the ring $𝕜$.

Examples

julia> X = affine_space(QQ,3)
Spec of Multivariate Polynomial Ring in x1, x2, x3 over Rational Field

julia> base_ring(X)
Rational Field
source
base_fieldMethod
base_field(E::EllCrv) -> Field

Return the base field over which E is defined.

base_field(C::HypellCrv) -> Field

Return the base field over which C is defined.

base_field(L::AbsLat) -> Field

Return the algebra over which the rational span of the lattice L is defined.

Invariants

rankMethod
rank(L::AbsLat) -> Int

Return the rank of the underlying module of the lattice L.

detMethod
det(L::ZLat) -> fmpq

Return the determinant of the gram matrix of L.

scaleMethod
scale(L::ZLat) -> fmpq

Return the scale of L.

The scale of L is defined as the positive generator of the $\mathbb Z$-ideal generated by $\{\Phi(x, y) : x, y \in L\}$.

normMethod
norm(L::ZLat) -> fmpq

Return the norm of L.

The norm of L is defined as the positive generator of the $\mathbb Z$- ideal generated by $\{\Phi(x,x) : x \in L\}$.

isevenMethod
iseven(L::ZLat) -> Bool

Return whether L is even.

An integer lattice L in the rational quadratic space $(V,\Phi)$ is called even if $\Phi(x,x) \in 2\mathbb{Z}$ for all $x in L$.

is_integralMethod
is_integral(L::AbsLat) -> Bool

Return whether the lattice L is integral.

is_primary_with_primeMethod
is_primary_with_prime(L::ZLat) -> Bool, fmpz

Given a $\mathbb Z$-lattice L, return whether L is primary, that is whether L is integral and its discriminant group (see discriminant_group) is a p-group for some prime number p. In case it is, p is also returned as second output.

Note that for unimodular lattices, this function returns (true, 1). If the lattice is not primary, the second return value is -1 by default.

is_primaryMethod
is_primary(L::ZLat, p::Union{Integer, fmpz}) -> Bool

Given an integral $\mathbb Z$-lattice L and a prime number p, return whether L is p-primary, that is whether its discriminant group (see discriminant_group) is a p-group.

is_elementary_with_primeMethod
is_elementary_with_prime(L::ZLat) -> Bool, fmpz

Given a $\mathbb Z$-lattice L, return whether L is elementary, that is whether L is integral and its discriminant group (see discriminant_group) is an elemenentary p-group for some prime number p. In case it is, p is also returned as second output.

Note that for unimodular lattices, this function returns (true, 1). If the lattice is not elementary, the second return value is -1 by default.

is_elementaryMethod
is_elementary(L::ZLat, p::Union{Integer, fmpz}) -> Bool

Given an integral $\mathbb Z$-lattice L and a prime number p, return whether L is p-elementary, that is whether its discriminant group (see discriminant_group) is an elementary p-group.

The Genus

For an integral lattice The genus of an integer lattice collects its local invariants. genus(::ZLat)

massMethod
mass(L::ZLat) -> fmpq

Return the mass of the genus of L.

genus_representativesMethod
genus_representatives(L::ZLat) -> Vector{ZLat}

Return representatives for the isometry classes in the genus of L.

Real invariants

signature_tupleMethod
signature_tuple(L::ZLat) -> Tuple{Int,Int,Int}

Return the number of (positive, zero, negative) inertia of L.

is_positive_definiteMethod
is_positive_definite(L::AbsLat) -> Bool

Return whether the rational span of the lattice L is positive definite.

is_negative_definiteMethod
is_negative_definite(L::AbsLat) -> Bool

Return whether the rational span of the lattice L is negative definite.

is_definiteMethod
is_definite(L::AbsLat) -> Bool

Return whether the rational span of the lattice L is definite.

Isometries

automorphism_group_generatorsMethod
automorphism_group_generators(E::EllCrv) -> Vector{EllCrvIso}

Return generators of the automorphism group of $E$.

automorphism_group_generators(L::AbsLat; ambient_representation::Bool = true)
                                                      -> Vector{MatElem}

Given a definite lattice L, return generators for the automorphism group of L. If ambient_representation == true (the default), the transformations are represented with respect to the ambient space of L. Otherwise, the transformations are represented with respect to the (pseudo-)basis of L.

automorphism_group_orderMethod
automorphism_group_order(L::AbsLat) -> Int

Given a definite lattice L, return the order of the automorphism group of L.

is_isometricMethod
is_isometric(L::AbsLat, M::AbsLat) -> Bool

Return whether the lattices L and M are isometric.

is_locally_isometricMethod
is_locally_isometric(L::ZLat, M::ZLat, p::Int) -> Bool

Return whether L and M are isometric over the p-adic integers.

i.e. whether $L \otimes \Z_p \cong M\otimes \Z_p$.

Root lattices

root_lattice_recognitionMethod
root_lattice_recognition(L::ZLat)

Return the ADE type of the root sublattice of L.

Input:

L – a definite and integral $\mathbb{Z}$-lattice.

Output:

Two lists, the first one conaining the ADE types and the second one the irreducible root sublattices.

For more recognizable gram matrices use root_lattice_recognition_fundamental.

Examples

julia> L = Zlattice(gram=ZZ[4  0 0  0 3  0 3  0;
                            0 16 8 12 2 12 6 10;
                            0  8 8  6 2  8 4  5;
                            0 12 6 10 2  9 5  8;
                            3  2 2  2 4  2 4  2;
                            0 12 8  9 2 12 6  9;
                            3  6 4  5 4  6 6  5;
                            0 10 5  8 2  9 5  8])
Quadratic lattice of rank 8 and degree 8 over the rationals

julia> R = root_lattice_recognition(L)
([(:A, 1), (:D, 6)], ZLat[Quadratic lattice of rank 1 and degree 8 over the rationals, Quadratic lattice of rank 6 and degree 8 over the rationals])
root_lattice_recognition_fundamentalMethod
root_lattice_recognition_fundamental(L::ZLat)

Return the ADE type of the root sublattice of L as well as the corresponding irreducible root sublattices with basis given by a fundamental root system.

Input:

L – a definite and integral $\mathbb Z$-lattice.

Output:

  • the root sublattice, with basis given by a fundamental root system
  • the ADE types
  • a Vector consisting of the irreducible root sublattices.

Examples

julia> L = Zlattice(gram=ZZ[4  0 0  0 3  0 3  0;
                            0 16 8 12 2 12 6 10;
                            0  8 8  6 2  8 4  5;
                            0 12 6 10 2  9 5  8;
                            3  2 2  2 4  2 4  2;
                            0 12 8  9 2 12 6  9;
                            3  6 4  5 4  6 6  5;
                            0 10 5  8 2  9 5  8])
Quadratic lattice of rank 8 and degree 8 over the rationals

julia> R = root_lattice_recognition_fundamental(L);

julia> gram_matrix(R[1])
[2    0    0    0    0    0    0]
[0    2    0   -1    0    0    0]
[0    0    2   -1    0    0    0]
[0   -1   -1    2   -1    0    0]
[0    0    0   -1    2   -1    0]
[0    0    0    0   -1    2   -1]
[0    0    0    0    0   -1    2]
ADE_typeMethod
ADE_type(G::MatrixElem) -> Tuple{Symbol,Int64}

Return the type of the irreducible root lattice with gram matrix G.

See also root_lattice_recognition.

Examples

julia> Hecke.ADE_type(gram_matrix(root_lattice(:A,3)))
(:A, 3)
coxeter_numberMethod
coxeter_number(ADE::Symbol, n)

Return the Coxeter number of the corresponding ADE root lattice.

If $L$ is a root lattice and $R$ its set of roots, then the Coxeter number $h$ is $|R|/n$ where n is the rank of $L$.

julia> coxeter_number(:D, 4)
6
highest_rootMethod
highest_root(ADE::Symbol, n) -> fmpz_mat

Return coordinates of the highest root of root_lattice(ADE, n).

julia> highest_root(:E, 6)
[1   2   3   2   1   2]

Module operations

Most module operations assume that the lattices live in the same ambient space. For instance only lattices in the same ambient space compare.

==Method

Return true if both lattices have the same ambient quadratic space and the same underlying module.

is_sublatticeMethod
is_sublattice(L::AbsLat, M::AbsLat) -> Bool

Return whether M is a sublattice of the lattice L.

is_sublattice_with_relationsMethod
is_sublattice_with_relations(M::ZLat, N::ZLat) -> Bool, fmpq_mat

Returns whether $N$ is a sublattice of $M$. In this case, the second return value is a matrix $B$ such that $B B_M = B_N$, where $B_M$ and $B_N$ are the basis matrices of $M$ and $N$ respectively.

+Method
+(L::AbsLat, M::AbsLat) -> AbsLat

Return the sum of the lattices L and M.

The lattices L and M must have the same ambient space.

*Method
*(a::RationalUnion, L::ZLat) -> ZLat

Returns the lattice $aM$ inside the ambient space of $M$.

intersectMethod
intersect(L::AbsLat, M::AbsLat) -> AbsLat

Return the intersection of the lattices L and M.

The lattices L and M must have the same ambient space.

inMethod
Base.in(v::Vector, L::ZLat) -> Bool

Return whether the vector v lies in the lattice L.

inMethod
Base.in(v::fmpq_mat, L::ZLat) -> Bool

Return whether the row span of v lies in the lattice L.

primitive_closureMethod
primitive_closure(M::ZLat, N::ZLat) -> ZLat

Given two $\mathbb Z$-lattices M and N with $N \subseteq \mathbb{Q} M$, return the primitive closure $M \cap \mathbb{Q} N$ of N in M.

Examples

julia> M = root_lattice(:D, 6);

julia> N = lattice_in_same_ambient_space(M, 3*basis_matrix(M)[1,:]);

julia> basis_matrix(N)
[3   0   0   0   0   0]

julia> N2 = primitive_closure(M, N)
Quadratic lattice of rank 1 and degree 6 over the rationals

julia> basis_matrix(N2)
[1   0   0   0   0   0]

julia> M2 = primitive_closure(dual(M), M);

julia> is_integral(M2)
false
is_primitiveMethod
is_primitive(M::ZLat, N::ZLat) -> Bool

Given two $\mathbb Z$-lattices $N \subseteq M$, return whether N is a primitive sublattice of M.

Examples

julia> U = hyperbolic_plane_lattice(3);

julia> bU = basis_matrix(U);

julia> e1, e2 = bU[1,:], bU[2,:]
([1 0], [0 1])

julia> N = lattice_in_same_ambient_space(U, e1 + e2)
Quadratic lattice of rank 1 and degree 2 over the rationals

julia> is_primitive(U, N)
true

julia> M = root_lattice(:A, 3);

julia> f = matrix(QQ, 3, 3, [0 1 1; -1 -1 -1; 1 1 0]);

julia> N = kernel_lattice(M, f+1)
Quadratic lattice of rank 1 and degree 3 over the rationals

julia> is_primitive(M, N)
true
is_primitiveMethod
is_primitive(L::ZLat, v::Union{Vector, fmpq_mat}) -> Bool

Return whether the vector v is primitive in L.

A vector v in a $\mathbb Z$-lattice L is called primitive if for all w in L such that $v = dw$ for some integer d, then $d = \pm 1$.

divisibilityMethod
divisibility(L::ZLat, v::Union{Vector, fmpq_mat}) -> fmpq

Return the divisibility of v with respect to L.

For a vector v in the ambient quadratic space $(V, \Phi)$ of L, we call the divisibility of v with the respect to L the non-negative generator of the fractional $\mathbb Z$-ideal $\Phi(v, L)$.

Embeddings

Orthogonal sublattices

orthogonal_sumMethod
orthogonal_sum(L1::ZLat, L2::ZLat)

Return the orthogonal direct sum of the lattices L1 and L2.

It lives in the orthogonal direct sum of their ambient spaces.

orthogonal_submoduleMethod
orthogonal_submodule(L::ZLat, S::ZLat) -> ZLat

Return the largest submodule of L orthogonal to S.

direct_sumMethod
direct_sum(x::Vararg{ZLat}) -> ZLat, Vector{AbsSpaceMor}, Vector{AbsSpaceMor}
direct_sum(x::Vector{ZLat}) -> ZLat, Vector{AbsSpaceMor}, Vector{AbsSpaceMor}

Given a collection of $\mathbb Z$-lattices $L_1, \ldots, L_n$, return their complete direct sum $L := L_1 \oplus \ldots \oplus L_n$, together with the injections $L_i \to L$ and the projections $L \to L_i$ (seen as maps between the corresponding ambient spaces).

irreducible_componentsMethod
irreducible_components(L::ZLat)

Return the irreducible components $L_i$ of the positive definite lattice $L$.

This yields a maximal orthogonal splitting of L as

\[L = \bigoplus_i L_i.\]

Dual lattice

dualMethod
dual(L::AbsLat) -> AbsLat

Return the dual lattice of the lattice L.

Discriminant group

See discriminant_group(L::ZLat).

Overlattices

glue_mapMethod
glue_map(L::ZLat, S::ZLat, R::ZLat; check=true)
                       -> Tuple{TorQuadModMor, TorQuadModMor, TorQuadModMor}

Given three integral $\mathbb Z$-lattices L, S and R, with S and R primitive sublattices of L and such that the sum of the ranks of S and R is equal to the rank of L, return the glue map $\gamma$ of the primitive extension $S+R \subseteq L$, as well as the inclusion maps of the domain and codomain of $\gamma$ into the respective discriminant groups of S and R.

Example

julia> M = root_lattice(:E,8);

julia> f = matrix(QQ, 8, 8, [-1 -1  0  0  0  0  0  0;
                              1  0  0  0  0  0  0  0;
                              0  1  1  0  0  0  0  0;
                              0  0  0  1  0  0  0  0;
                              0  0  0  0  1  0  0  0;
                              0  0  0  0  0  1  1  0;
                             -2 -4 -6 -5 -4 -3 -2 -3;
                              0  0  0  0  0  0  0  1]);

julia> S = kernel_lattice(M ,f-1)
Quadratic lattice of rank 4 and degree 8 over the rationals

julia> R = kernel_lattice(M , f^2+f+1)
Quadratic lattice of rank 4 and degree 8 over the rationals

julia> glue, iS, iR = glue_map(M, S, R)
(Map with following data
Domain:
=======
TorQuadMod [4//3 0; 0 4//3]
Codomain:
=========
TorQuadMod [2//3 0; 0 2//3], Map with following data
Domain:
=======
TorQuadMod [4//3 0; 0 4//3]
Codomain:
=========
TorQuadMod [4//3 2//3; 2//3 2//3], Map with following data
Domain:
=======
TorQuadMod [2//3 0; 0 2//3]
Codomain:
=========
TorQuadMod [2//3 1//3; 1//3 4//3])

julia> is_bijective(glue)
true
overlatticeMethod
overlattice(glue_map::TorQuadModMor) -> ZLat

Given the glue map of a primitive extension of $\mathbb Z$-lattices $S+R \subseteq L$, return L.

Example

julia> M = root_lattice(:E,8);

julia> f = matrix(QQ, 8, 8, [ 1  0  0  0  0  0  0  0;
                              0  1  0  0  0  0  0  0;
                              1  2  4  4  3  2  1  2;
                             -2 -4 -6 -5 -4 -3 -2 -3;
                              2  4  6  4  3  2  1  3;
                             -1 -2 -3 -2 -1  0  0 -2;
                              0  0  0  0  0 -1  0  0;
                             -1 -2 -3 -3 -2 -1  0 -1]);

julia> S = kernel_lattice(M ,f-1)
Quadratic lattice of rank 4 and degree 8 over the rationals

julia> R = kernel_lattice(M , f^4+f^3+f^2+f+1)
Quadratic lattice of rank 4 and degree 8 over the rationals

julia> glue, iS, iR = glue_map(M, S, R);

julia> overlattice(glue) == M
true
local_modificationMethod
local_modification(M::ZLat, L::ZLat, p)

Return a local modification of M that matches L at p.

INPUT:

  • $M$ – a \mathbb{Z}_p-maximal lattice
  • $L$ – the a lattice isomorphic to M over \QQ_p
  • $p$ – a prime number

OUTPUT:

an integral lattice M' in the ambient space of M such that M and M' are locally equal at all completions except at p where M' is locally isometric to the lattice L.

maximal_integral_latticeMethod
maximal_integral_lattice(L::AbsLat) -> AbsLat

Given a lattice L, return a lattice M in the ambient space of L which is maximal integral and which contains L.

Sublattices defined by endomorphisms

kernel_latticeMethod
kernel_lattice(L::ZLat, f::MatElem;
               ambient_representation::Bool = true) -> ZLat

Given a $\mathbf{Z}$-lattice $L$ and a matrix $f$ inducing an endomorphism of $L$, return $\ker(f)$ is a sublattice of $L$.

If ambient_representation is true (the default), the endomorphism is represented with respect to the ambient space of $L$. Otherwise, the endomorphism is represented with respect to the basis of $L$.

invariant_latticeMethod
invariant_lattice(L::ZLat, G::Vector{MatElem};
                  ambient_representation::Bool = true) -> ZLat

Given a $\mathbf{Z}$-lattice $L$ and a list of matrices $G$ inducing endomorphisms of $L$, return the lattice $L^G$, consisting of elements fixed by $G$.

If ambient_representation is true (the default), the endomorphism is represented with respect to the ambient space of $L$. Otherwise, the endomorphism is represented with respect to the basis of $L$.

Computing embeddings

embedFunction

embed(S::ZLat, G::Genus, primitive=true) -> Bool, embedding

Return a (primitive) embedding of the integral lattice S into some lattice in the genus of G.

julia> G = genera((8,0), 1, even=true)[1];

julia> L, S, i = embed(root_lattice(:A,5), G);
embed_in_unimodularMethod
embed_in_unimodular(S::ZLat, pos, neg, primitive=true, even=true) -> Bool, L, S', iS, iR

Return a (primitive) embedding of the integral lattice S into some (even) unimodular lattice of signature (pos, neg).

For now this works only for even lattices.

julia> NS = orthogonal_sum(Zlattice(:U), rescale(root_lattice(:A, 16), -1))[1];

julia> LK3, iNS, i = embed_in_unimodular(NS, 3, 19);

julia> genus(LK3)
ZGenus
Signature: (3, 19)
Genus symbol at 2:   1^22

julia> iNS
Quadratic lattice of rank 18 and degree 22 over the rationals

julia> is_primitive(LK3, iNS)
true

LLL, Short and Close Vectors

LLL and indefinite LLL

lllMethod
lll(L::ZLat, same_ambient::Bool = true) -> ZLat

Given an integral $\mathbb Z$-lattice L with basis matrix B, compute a basis C of L such that the gram matrix $G_C$ of L with respect to C is LLL-reduced.

By default, it creates the lattice in the same ambient space as L. This can be disabled by setting same_ambient = false. Works with both definite and indefinite lattices.

Short Vectors

short_vectorsFunction
short_vectors(L::ZLat, [lb = 0], ub, [elem_type = fmpz]; check::Bool = true)
                                   -> Vector{Tuple{Vector{elem_type}, fmpq}}

Returns all tuples (v, n) such that n = v G v^t satisfies lb <= n <= ub, where G is the Gram matrix of L and v is non-zero.

Note that the vectors are computed up to sign (so only one of v and -v appears).

It is assumed and checked that L is positive definite.

See also short_vectors_iterator for an iterator version.

shortest_vectorsFunction
shortest_vectors(L::ZLat, [elem_type = fmpz]; check::Bool = true)
                                           -> fmpq, Vector{elem_type}, fmpq}

Returns the list of shortest non-zero vectors. Note that the vectors are computed up to sign (so only one of v and -v appears).

It is assumed and checked that L is positive definite.

See also minimum.

short_vectors_iteratorFunction
short_vectors_iterator(L::ZLat, [lb = 0], ub,
                       [elem_type = fmpz]; check::Bool = true)
                                -> Tuple{Vector{elem_type}, fmpq} (iterator)

Returns an iterator for all tuples (v, n) such that n = v G v^t satisfies lb <= n <= ub, where G is the Gram matrix of L and v is non-zero.

Note that the vectors are computed up to sign (so only one of v and -v appears).

It is assumed and checked that L is positive definite.

See also short_vectors.

minimumMethod
minimum(L::ZLat)

Return the minimum squared length among the non-zero vectors in L.

kissing_numberMethod
kissing_number(L::ZLat)

Return the Kissing number of the sphere packing defined by L.

This is the number of non-overlapping spheres touching any other given sphere.

Close Vectors

close_vectorsMethod
close_vectors(L:ZLat, v:Vector, [lb,], ub; check::Bool = false)
                                        -> Vector{Tuple{Vector{Int}}, fmpq}

Return all tuples (x, d) where x is an element of L such that d = b(v - x, v - x) <= ub. If lb is provided, then also lb <= d.

If filter is not nothing, then only those x with filter(x) evaluating to true are returned.

By default, it will be checked whether L is positive definite. This can be disabled setting check = false.

Both input and output are with respect to the basis matrix of L.

Examples

julia> L = Zlattice(matrix(QQ, 2, 2, [1, 0, 0, 2]));

julia> close_vectors(L, [1, 1], 1)
3-element Vector{Tuple{Vector{fmpz}, fmpq}}:
 ([2, 1], 1)
 ([0, 1], 1)
 ([1, 1], 0)

julia> close_vectors(L, [1, 1], 1, 1)
2-element Vector{Tuple{Vector{fmpz}, fmpq}}:
 ([2, 1], 1)
 ([0, 1], 1)