Commits (1)
  • Jan Pöschko's avatar
    Author: Jan Pöschko <jan@poeschko.com> · 7af2600b
    Jan Pöschko authored
     Date:   Thu Aug 23 23:01:44 2012 +0200
    
         doc
    
     Author: Jan Pöschko <jan@poeschko.com>
     Date:   Thu Aug 23 22:31:53 2012 +0200
    
         cleanup, copyright
    
     Author: Jan Pöschko <jan@poeschko.com>
     Date:   Thu Aug 23 21:46:23 2012 +0200
    
         plotting
    
     Author: Jan Pöschko <jan@poeschko.com>
     Date:   Thu Aug 23 21:24:57 2012 +0200
    
         add convex_volume to Polyhedron_base
    
     Author: Jan Pöschko <jan@poeschko.com>
     Date:   Thu Aug 23 13:22:11 2012 +0200
    
         random and special lattices, separate constructor file
    
     Author: Jan Pöschko <jan@poeschko.com>
     Date:   Thu Aug 23 05:41:26 2012 +0200
    
         tests (including volume of polyhedron)
    
     Author: Jan Pöschko <jan@poeschko.com>
     Date:   Thu Aug 23 04:53:15 2012 +0200
    
         fix dependencies of lrs_volume
    
     Author: Jan Pöschko <jan@poeschko.com>
     Date:   Thu Aug 23 04:50:30 2012 +0200
    
         original version
    
     Author: Jan Pöschko <jan@poeschko.com>
     Date:   Thu Aug 23 04:14:55 2012 +0200
    
         cleanup & doc
    
     Author: Jan Pöschko <jan@poeschko.com>
     Date:   Thu Aug 23 03:34:42 2012 +0200
    
         new class structure
    
     Author: Jan Pöschko <jan@poeschko.com>
     Date:   Thu Aug 23 02:03:59 2012 +0200
    
         allow non-square result of LLL_gram
    
     Author: Jan Pöschko <jan@poeschko.com>
     Date:   Thu Aug 23 02:02:55 2012 +0200
    
         original version
    
     Author: Jan Pöschko <jan@poeschko.com>
     Date:   Thu Aug 23 00:46:48 2012 +0200
    
         LLL_gram for rational matrix
    
     Author: Jan Pöschko <jan@poeschko.com>
     Date:   Thu Aug 23 00:27:26 2012 +0200
    
         ** instead of ^
    
     Author: Jan Pöschko <jan@poeschko.com>
     Date:   Fri Aug 17 12:45:05 2012 +0200
    
         inner_product_matrix patch
    
     Author: Jan Pöschko <jan@poeschko.com>
     Date:   Fri Aug 17 12:41:15 2012 +0200
    
         original version
    
     Merge: 79d4703 6e36031
     Author: Jan Pöschko <jan@poeschko.com>
     Date:   Fri Aug 17 00:01:12 2012 +0200
    
         Merge branch 'master' of github.com:poeschko/sage
    
         Conflicts:
         	devel/sage-lattices/sage/lattices/lattice.py
    
     Author: Jan Pöschko <jan@poeschko.com>
     Date:   Thu Aug 16 23:39:10 2012 +0200
    
         docs
    
     Author: Jan Pöschko <jan@poeschko.com>
     Date:   Thu Aug 16 23:39:10 2012 +0200
    
         docs
    
     Author: Jan Pöschko <jan@poeschko.com>
     Date:   Thu Aug 16 21:45:32 2012 +0200
    
         update readme
    
     Author: Jan Pöschko <jan@poeschko.com>
     Date:   Thu Aug 16 21:41:38 2012 +0200
    
         cholesky instead of deprecated cholesky_decomposition
    
     Author: Jan Pöschko <jan@poeschko.com>
     Date:   Thu Aug 16 21:38:35 2012 +0200
    
         integrate patch to #12051 (rational LLL)
    
     Author: Jan Pöschko <jan@poeschko.com>
     Date:   Thu Aug 16 21:36:24 2012 +0200
    
         original version of matrix_rational_dense.pyx
    
     Author: Jan Pöschko <jan@poeschko.com>
     Date:   Thu Aug 16 21:08:10 2012 +0200
    
         added .gitignore
    
     Author: Jan Pöschko <jan@poeschko.com>
     Date:   Thu Aug 16 21:06:09 2012 +0200
    
         rebase to Sage 5.2
    
     Author: Jan Pöschko <jan@poeschko.com>
     Date:   Wed Aug 8 01:18:41 2012 +0200
    
         documentation
    
     Author: Jan Pöschko <jan@poeschko.com>
     Date:   Wed Aug 8 00:42:57 2012 +0200
    
         cleanup
    
     Author: Jan Pöschko <jan@poeschko.com>
     Date:   Wed Aug 8 00:37:18 2012 +0200
    
         closest vector
    
     Author: Jan Pöschko <jan@poeschko.com>
     Date:   Tue Jul 10 16:04:59 2012 +0200
    
         improved documentation
    
     Author: Jan Pöschko <jan@poeschko.com>
     Date:   Tue Jul 10 15:52:02 2012 +0200
    
         allow application of matrix() to a matrix
    
     Author: Jan Pöschko <jan@poeschko.com>
     Date:   Tue Jul 10 15:50:18 2012 +0200
    
         constructing lattices via quadratic form
    
     Author: Jan Pöschko <jan@poeschko.com>
     Date:   Sun Jul 8 18:19:06 2012 +0200
    
         timing as6 lattice
    
     Author: Jan Pöschko <jan@poeschko.com>
     Date:   Sun Jul 8 17:19:54 2012 +0200
    
         allow non-square lattice bases in Voronoi cell computation
    
     Author: Jan Pöschko <jan@poeschko.com>
     Date:   Wed Jul 4 21:24:26 2012 +0200
    
         perform actual cuts accumulated in the end
    
     Author: Jan Pöschko <jan@poeschko.com>
     Date:   Wed Jul 4 20:59:35 2012 +0200
    
         cleanup
    
     Author: Jan Pöschko <jan@poeschko.com>
     Date:   Wed Jul 4 20:50:49 2012 +0200
    
         fixed errors in implementation of algorithm
    
     Author: Jan Pöschko <jan@poeschko.com>
     Date:   Wed Jul 4 18:25:23 2012 +0200
    
         first implementation of diamond cutting using Polyhedron
    
     Author: Jan Pöschko <jan@poeschko.com>
     Date:   Wed Jun 20 04:57:03 2012 +0200
    
         removed basis_inner_product_matrix again (= gram_matrix)
    
     Author: Jan Pöschko <jan@poeschko.com>
     Date:   Wed Jun 20 04:55:29 2012 +0200
    
         discriminant
    
     Author: Jan Pöschko <jan@poeschko.com>
     Date:   Wed Jun 20 04:43:30 2012 +0200
    
         first implementation of shortest_vectors
    
     Author: Jan Pöschko <jan@poeschko.com>
     Date:   Wed Jun 20 03:40:03 2012 +0200
    
         LLL basis reduction
    
     Author: Jan Pöschko <jan@poeschko.com>
     Date:   Wed Jun 20 02:31:40 2012 +0200
    
         reverted prevention of coercion from reals to rationals
    
     Author: Jan Pöschko <jan@poeschko.com>
     Date:   Mon Jun 18 22:21:22 2012 +0200
    
         use TestSuite
    
     Author: Jan Pöschko <jan@poeschko.com>
     Date:   Mon Jun 4 02:25:06 2012 +0200
    
         added documentation
    
     Author: Jan Pöschko <jan@poeschko.com>
     Date:   Mon Jun 4 01:19:49 2012 +0200
    
         doc
    
     Author: Jan Pöschko <jan@poeschko.com>
     Date:   Mon Jun 4 00:29:01 2012 +0200
    
         more tests
    
     Author: Jan Pöschko <jan@poeschko.com>
     Date:   Mon Jun 4 00:24:00 2012 +0200
    
         cleanup and doc
    
     Author: Jan Pöschko <jan@poeschko.com>
     Date:   Sun Jun 3 23:44:29 2012 +0200
    
         Z-lattices over R^n and C^n
    
     Author: Jan Pöschko <jan@poeschko.com>
     Date:   Sun Jun 3 00:29:15 2012 +0200
    
         refactored class structure
    
     Author: Jan Pöschko <jan@poeschko.com>
     Date:   Fri Jun 1 23:27:36 2012 +0200
    
         removed unnecessary code
    
     Author: Jan Pöschko <jan@poeschko.com>
     Date:   Fri Jun 1 23:26:23 2012 +0200
    
         basic class structure
    
     Author: Jan Pöschko <jan@poeschko.com>
     Date:   Sun May 27 12:14:51 2012 +0200
    
         base class
    
     Author: Jan Pöschko <jan@poeschko.com>
     Date:   Sat May 26 04:58:16 2012 +0200
    
         removed README.txt
    
     Author: Jan Pöschko <jan@poeschko.com>
     Date:   Sat May 26 04:55:35 2012 +0200
    
         updated readme
    
     Author: Jan Pöschko <jan@poeschko.com>
     Date:   Sat May 26 04:51:19 2012 +0200
    
         lattices files
    
     Author: Jan Pöschko <jan@poeschko.com>
     Date:   Sat May 26 04:48:55 2012 +0200
    
         readme and copying
    
     Author: Jan Pöschko <jan@poeschko.com>
     Date:   Fri May 25 19:26:16 2012 +0200
    
         readme
    
     Author: Jan Pöschko <jan@poeschko.com>
     Date:   Fri May 25 14:40:16 2012 +0200
    
         first commit
    7af2600b
......@@ -84,7 +84,8 @@ Groups, Monoids, Matrices, Modules
* :doc:`Monoids <monoids/index>`
* :doc:`Matrices and Spaces of Matrices <matrices/index>`
* :doc:`Modules <modules/index>`
* :doc:`Lattices <lattices/index>`
Geometry and Topology
---------------------
......
../conf_sub.py
\ No newline at end of file
.. _ch:lattices:
Lattices
========
.. toctree::
:maxdepth: 2
sage/lattices/constructor
sage/lattices/lattice
sage/lattices/special_lattices
sage/lattices/plot
......@@ -168,6 +168,8 @@ sage.misc.lazy_import.lazy_import('sage.sandpiles.all', '*', globals())
from sage.tensor.all import *
from sage.lattices.all import *
from sage.matroids.all import *
# Lazily import notebook functions and interacts (#15335)
......
from constructor import Lattice, random_lattice
from special_lattices import special_lattice, special_lattice_names
from plot import plot_lattice
\ No newline at end of file
"""
Lattice constructor
Lattices are created using the `Lattice` factory function, given
either a basis, an inner product matrix, a quadratic form,
or an order in an algebraic number field.
AUTHORS:
- Jan Poeschko (2012-08-10): initial version
"""
#*****************************************************************************
# Copyright (C) 2012 Jan Poeschko <jan@poeschko.com>
#
# Distributed under the terms of the GNU General Public License (GPL)
# as published by the Free Software Foundation; either version 2 of
# the License, or (at your option) any later version.
# http://www.gnu.org/licenses/
#*****************************************************************************
from sage.rings.integer_ring import ZZ
from sage.matrix.constructor import matrix, random_matrix
from lattice import Lattice_with_basis, Lattice_ZZ
def Lattice(basis=None, inner_product_matrix=None, quadratic_form=None, algebraic_order=None,
base_ring=ZZ, **kwargs):
"""
The :func:`Lattice` function creates lattices using a given base,
an inner product matrix, an underlying quadratic form,
or an order in an algebraic number field.
INPUT:
- ``basis``
- ``inner_product_matrix``
- ``quadratic_form``
- ``algebraic_order``
- ``base_ring``
OUTPUT:
A lattice.
EXAMPLES::
sage: Lattice([[2, 0, 0], [0, 1, 0]])
ZZ-lattice of degree 3 and rank 2
Inner product matrix:
[1 0]
[0 4]
Basis matrix:
[ 0 1 0]
[-2 0 0]
A lattice can be specified by a quadratic form::
sage: Lattice(quadratic_form=QuadraticForm(ZZ, 3, [1,2,3,4,5,6]))
Lattice of degree 3 and rank 3 over Integer Ring
Inner product matrix:
[ 1 1 3/2]
[ 1 4 5/2]
[3/2 5/2 6]
Basis matrix:
[ 1 0 0]
[ 1 1.732050807568878? 0]
[ 3/2 0.5773502691896258? 1.848422751068237?]
It is also possible to specify an order in an algebraic number field::
sage: K.<a> = NumberField(x ^ 2 + 1)
sage: O = K.order(a)
sage: Lattice(algebraic_order=O)
ZZ-lattice of degree 2 and rank 2
Inner product matrix:
[1 0]
[0 1]
Basis matrix:
[1 0]
[0 1]
It is an error to specify an empty basis::
sage: Lattice([])
Traceback (most recent call last):
...
ValueError: basis must not be empty
"""
if quadratic_form is not None:
inner_product_matrix = quadratic_form.Gram_matrix_rational()
if algebraic_order is not None:
basis = algebraic_order.free_module().basis()
if basis is not None:
if not basis:
raise ValueError("basis must not be empty")
basis = matrix(basis)
rank = basis.dimensions()[0]
if inner_product_matrix is None:
inner_product_matrix = basis * basis.transpose()
elif inner_product_matrix is not None:
inner_product_matrix = matrix(inner_product_matrix)
rank = inner_product_matrix.dimensions()[0]
basis = inner_product_matrix.cholesky()
else:
raise TypeError("basis or inner_product_matrix must be given")
if base_ring == ZZ and inner_product_matrix.base_ring() == ZZ:
return Lattice_ZZ(rank, inner_product_matrix, basis, **kwargs)
else:
return Lattice_with_basis(base_ring, rank, inner_product_matrix, basis)
def random_lattice(dimension):
"""
Construct a random ZZ-lattice of a given dimension.
INPUT:
- ``dimension`` -- dimension of the constructed lattice.
OUTPUT:
A lattice with integer basis vectors.
EXAMPLES::
sage: random_lattice(3)
ZZ-lattice of degree 3 and rank 3
Inner product matrix:
[ 2 0 0]
[ 0 66 -22]
[ 0 -22 4246]
Basis matrix:
[ 0 1 -1]
[-8 1 1]
[14 45 45]
sage: random_lattice(100).dimension()
100
"""
basis = random_matrix(ZZ, dimension, dimension)
return Lattice(basis)
"""
Diamond cutting implementation
AUTHORS:
- Jan Poeschko (2012-07-02): initial version
"""
#*****************************************************************************
# Copyright (C) 2012 Jan Poeschko <jan@poeschko.com>
#
# Distributed under the terms of the GNU General Public License (GPL)
# as published by the Free Software Foundation; either version 2 of
# the License, or (at your option) any later version.
# http://www.gnu.org/licenses/
#*****************************************************************************
from sage.geometry.polyhedron.constructor import Polyhedron
from sage.matrix.constructor import matrix, identity_matrix
from sage.modules.free_module_element import vector
from sage.rings.rational_field import QQ
from math import sqrt, floor, ceil
def plane_inequality(v):
"""
Return the inequality for points on the same side as the origin
with respect to the plane through v normal to v.
sage: from sage.lattices.diamond_cutting import plane_inequality
sage: ieq = plane_inequality([1, -1]); ieq
[2, -1, 1]
sage: ieq[0] + vector(ieq[1:]) * vector([1, -1])
0
"""
v = vector(v)
c = -v * v
if c < 0:
c, v = -c, -v
return [c] + list(v)
def jacobi(M):
"""
Compute the Cholesky/Jacobi decomposition of M.
EXAMPLES::
sage: from sage.lattices.diamond_cutting import jacobi
sage: jacobi(identity_matrix(3) * 4)
[4 0 0]
[0 4 0]
[0 0 4]
"""
dim = M.dimensions()
assert dim[0] == dim[1]
dim = dim[0]
q = [list(row) for row in M]
for i in range(dim - 1):
for j in range(i + 1, dim):
q[j][i] = q[i][j]
q[i][j] = q[i][j] / q[i][i]
for k in range(i + 1, dim):
for l in range(k, dim):
q[k][l] = q[k][l] - q[k][i] * q[i][l]
for i in range(1, dim):
for j in range(i):
q[i][j] = 0
return matrix(q)
def diamond_cut(V, GM, C, debug=False):
"""
Perform diamond cutting on polyhedron V with basis matrix GM and radius C.
INPUT:
- ``V``: polyhedron to cut from;
- ``GM``: half of the basis matrix of the lattice;
- ``C``: radius to use in cutting algorithm;
- ``debug``: whether to print debug information.
OUTPUT:
A :class:``Polyhedron`` instance.
EXAMPLES::
sage: from sage.lattices.diamond_cutting import diamond_cut
sage: V = Polyhedron([[0], [2]])
sage: GM = matrix([2])
sage: V = diamond_cut(V, GM, 4)
sage: V.vertices()
(A vertex at (2), A vertex at (0))
"""
# coerce to floats
GM = GM.N()
C = float(C)
if debug:
print "Cut\n%s\nwith radius %s" % (GM, C)
dim = GM.dimensions()
assert dim[0] == dim[1]
dim = dim[0]
T = [0] * dim
U = [0] * dim
x = [0] * dim
L = [0] * dim
# calculate the Gram matrix
q = matrix([[sum(GM[i][k] * GM[j][k] for k in range(dim)) for j in range(dim)] for i in range(dim)])
if debug:
print "q:\n%s" % q.N()
# apply Cholesky/Jacobi decomposition
q = jacobi(q)
if debug:
print "q:\n%s" % q.N()
i = dim - 1
T[i] = C
U[i] = 0
new_dimension = True
cut_count = 0
inequalities = []
while True:
if debug:
print "Dimension: %d" % i
if new_dimension:
Z = sqrt(T[i] / q[i][i])
if debug:
print "Z: %s" % Z
L[i] = int(floor(Z - U[i]))
if debug:
print "L: %s" % L
x[i] = int(ceil(-Z - U[i]) - 1)
new_dimension = False
x[i] += 1
if debug:
print "x: %s" % x
if x[i] > L[i]:
i += 1
elif i > 0:
T[i - 1] = T[i] - q[i][i] * (x[i] + U[i]) ** 2
i -= 1
U[i] = 0
for j in range(i + 1, dim):
U[i] += q[i][j] * x[j]
new_dimension = True
else:
if all(elmt == 0 for elmt in x):
break
hv = [0] * dim
for k in range(dim):
for j in range(dim):
hv[k] += x[j] * GM[j][k]
hv = vector(hv)
for hv in [hv, -hv]:
cut_count += 1
if debug:
print "\n%d) Cut using normal vector %s" % (cut_count, hv)
hv = [QQ(round(elmt, 6)) for elmt in hv]
inequalities.append(plane_inequality(hv))
#cut = Polyhedron(ieqs=[plane_inequality(hv)])
#V = V.intersection(cut)
if debug:
print "Final cut"
cut = Polyhedron(ieqs=inequalities)
V = V.intersection(cut)
if debug:
print "End"
return V
def calculate_voronoi_cell(basis, radius=None, debug=False):
"""
Calculate the Voronoi cell of the lattice defined by basis
INPUT:
- ``basis``: embedded basis matrix of the lattice;
- ``radius``: radius of basis vectors to consider;
- ``debug``: whether to print debug information.
OUTPUT:
A :class:``Polyhedron`` instance.
EXAMPLES::
sage: from sage.lattices.diamond_cutting import calculate_voronoi_cell
sage: V = calculate_voronoi_cell(matrix([[1, 0], [0, 1]]))
sage: V.volume()
1
"""
dim = basis.dimensions()
artificial_length = None
if dim[0] < dim[1]:
# introduce "artificial" basis points (representing infinity)
artificial_length = ceil(max(abs(v) for v in basis)) * 2
additional_vectors = identity_matrix(dim[1]) * artificial_length
basis = basis.stack(additional_vectors)
# LLL-reduce to get quadratic matrix
basis = basis.LLL()
basis = matrix([v for v in basis if v])
dim = basis.dimensions()
assert dim[0] == dim[1]
basis = basis / 2
ieqs = []
for v in basis:
ieqs.append(plane_inequality(v))
ieqs.append(plane_inequality(-v))
Q = Polyhedron(ieqs=ieqs)
# twice the length of longest vertex in Q is a safe choice
if radius is None:
radius = 2 * max(abs(v) ** 2 for v in basis)
V = diamond_cut(Q, basis, radius, debug=debug)
if artificial_length is not None:
# remove inequalities introduced by artificial basis points
H = V.Hrepresentation()
H = [v for v in H if all(not V._is_zero(v.A() * w / 2 - v.b() and
not V._is_zero(v.A() * (-w) / 2 - v.b())) for w in additional_vectors)]
V = Polyhedron(ieqs=H)
return V
This diff is collapsed.
"""
Lattice plot functionality
AUTHORS:
- Jan Poeschko (2012-08-15): initial version
EXAMPLES::
sage: L = special_lattice('SimpleCubic')
sage: g = plot_lattice(L)
You can also add the Voronoi cell to the plot::
sage: V = L.voronoi_cell()
sage: g += V.plot()
sage: g.show(viewer='tachyon')
"""
#*****************************************************************************
# Copyright (C) 2012 Jan Poeschko <jan@poeschko.com>
#
# Distributed under the terms of the GNU General Public License (GPL)
# as published by the Free Software Foundation; either version 2 of
# the License, or (at your option) any later version.
# http://www.gnu.org/licenses/
#*****************************************************************************
from sage import plot
def subsets(s):
"""
Generate all subsets of a list of elements.
INPUT:
- ``s`` -- a list of elements.
OUTPUT:
A generator iterating through all subsets of `s`.
EXAMPLES::
sage: from sage.lattices.plot import subsets
sage: list(subsets([1, 2, 3]))
[[], [1], [2], [1, 2], [3], [1, 3], [2, 3], [1, 2, 3]]
"""
n_elems = len(s)
n_subsets = 2**len(s)
for i in range(0,n_subsets):
sub = []
for j in range(0,n_elems):
if (i >> j & 1):
sub.append(s[j])
yield sub
def plot_lattice(L, **options):
"""
Plot a lattice.
The sum of each subset of basis vectors (including their negated counterparts)
is plotted.
INPUT:
- ``L`` -- a `Lattice` instance;
- ``options`` - options to be passed to `list_plot`.
OUTPUT:
Graphics containing the plot of the lattice.
EXAMPLES::
sage: L = Lattice([[2, 0], [0, 3]])
sage: plot_lattice(L)
Three-dimensional lattices are supported, too::
sage: L = Lattice([[1, 0, 0], [0, 2, 0], [0, 0, 3]])
sage: plot_lattice(L)
Higher dimensions are not supported, though::
sage: L = random_lattice(4)
sage: plot_lattice(L)
Traceback (most recent call last):
...
ValueError: only 2-dimensional and 3-dimensional lattices can be plotted
"""
dim = L.dimension()
basis = L.embedded_basis()
basis = sum(([-v, v] for v in basis), [])
points = [sum(v, L(0)) for v in subsets(basis)]
points = [tuple(v.list()) for v in points]
points = set(points)
if dim not in (2, 3):
raise ValueError("only 2-dimensional and 3-dimensional lattices can be plotted")
return plot.point.points(points)
"""
Special named lattices
This module provides a constructor for some special named lattices.
AUTHORS:
- Jan Poeschko (2012-08-10): initial version
"""
#*****************************************************************************
# Copyright (C) 2012 Jan Poeschko <jan@poeschko.com>
#
# Distributed under the terms of the GNU General Public License (GPL)
# as published by the Free Software Foundation; either version 2 of
# the License, or (at your option) any later version.
# http://www.gnu.org/licenses/
#*****************************************************************************
def special_lattice(name):
"""
Construct a special named lattice.
INPUT:
- ``name`` -- the name of the lattice to construct (as a string).
The list of supported names is returned by :func:`special_lattice_names`.
OUTPUT:
A lattice.
EXAMPLES::
sage: special_lattice('SquareLattice').embedded_basis()
[(1, 0), (0, 1)]
Consider the Leech lattice::
sage: leech = special_lattice('Leech')
sage: leech
ZZ-lattice of degree 24 and rank 24
Inner product matrix:
24 x 24 dense matrix over Integer Ring
Basis matrix:
24 x 24 dense matrix over Algebraic Field
The Leech lattice is unimodular::
sage: leech.is_unimodular()
True
"""
from constructor import Lattice
basis = _basis.get(name)
if basis is not None:
return Lattice(basis=basis, reduce=False)
ipm = _inner_product_matrices.get(name)
if ipm is not None:
return Lattice(inner_product_matrix=ipm)
raise ValueError("unkown special lattice: %s" % name)
def special_lattice_names():
"""
Return the list of supported special named lattices.
OUTPUT:
A list of strings containing the supported names.
EXAMPLES::
sage: special_lattice_names()
['BodyCenteredCubic', 'FaceCenteredCubic', 'HexagonalLattice', 'Leech', 'SimpleCubic', 'SimpleHexagonal', 'SimpleOrthorhombic', 'SquareLattice']
"""
return sorted(set(_inner_product_matrices.keys() + _basis.keys()))
# named lattices given by basis vectors
_basis = {
'BodyCenteredCubic': [
[2, 0, 0],
[0, 2, 0],
[1, 1, 1]
],
'FaceCenteredCubic': [
[-1, -1, 0],
[1, -1, 0],
[0, 1, -1]
],
'HexagonalLattice': [
[-1, 1],
[0, -1]
],
'SimpleCubic': [
[1, 0, 0],
[0, 1, 0],
[0, 0, 1]
],
'SquareLattice': [
[1, 0],
[0, 1]
],
}
# named lattices given by their inner product matrix
_inner_product_matrices = {
'Leech': [
[8, 4, 4, 4, 4, 4, 4, 2, 4, 4, 4, 2, 4, 2, 2, 2, 4, 2, 2, 2, 0, 0, 0, -3],
[4, 4, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 1, 1, 2, 1, 1, 2, 1, 0, 0, -1],
[4, 2, 4, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 1, 2, 1, 2, 2, 1, 1, 1, 0, 0, -1],
[4, 2, 2, 4, 2, 2, 2, 2, 2, 2, 2, 2, 2, 1, 1, 2, 2, 1, 2, 1, 1, 0, 0, -1],
[4, 2, 2, 2, 4, 2, 2, 2, 2, 2, 2, 1, 2, 2, 2, 2, 2, 2, 2, 2, 1, 0, 0, -1],
[4, 2, 2, 2, 2, 4, 2, 2, 2, 2, 2, 1, 2, 2, 1, 1, 2, 1, 2, 1, 0, 0, 0, -1],
[4, 2, 2, 2, 2, 2, 4, 2, 2, 2, 2, 1, 2, 1, 2, 1, 2, 1, 1, 2, 0, 0, 0, -1],
[2, 2, 2, 2, 2, 2, 2, 4, 1, 1, 1, 2, 1, 2, 2, 2, 1, 2, 2, 2, 2, 0, 0, 1],
[4, 2, 2, 2, 2, 2, 2, 1, 4, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 1, 1, 1, -1],
[4, 2, 2, 2, 2, 2, 2, 1, 2, 4, 2, 2, 2, 2, 1, 1, 2, 2, 1, 1, 0, 1, 0, -1],
[4, 2, 2, 2, 2, 2, 2, 1, 2, 2, 4, 2, 2, 1, 2, 1, 2, 1, 2, 1, 0, 0, 1, -1],
[2, 2, 2, 2, 1, 1, 1, 2, 2, 2, 2, 4, 1, 2, 2, 2, 1, 2, 2, 2, 2, 1, 1, 1],
[4, 2, 2, 2, 2, 2, 2, 1, 2, 2, 2, 1, 4, 2, 2, 2, 2, 1, 1, 1, 1, 1, 1, -1],
[2, 2, 1, 1, 2, 2, 1, 2, 2, 2, 1, 2, 2, 4, 2, 2, 1, 2, 2, 2, 2, 2, 1, 1],
[2, 1, 2, 1, 2, 1, 2, 2, 2, 1, 2, 2, 2, 2, 4, 2, 1, 2, 2, 2, 2, 1, 2, 1],
[2, 1, 1, 2, 2, 1, 1, 2, 2, 1, 1, 2, 2, 2, 2, 4, 1, 2, 2, 2, 2, 1, 1, 1],
[4, 2, 2, 2, 2, 2, 2, 1, 2, 2, 2, 1, 2, 1, 1, 1, 4, 2, 2, 2, 1, 1, 1, -1],
[2, 1, 2, 1, 2, 1, 1, 2, 2, 2, 1, 2, 1, 2, 2, 2, 2, 4, 2, 2, 2, 2, 1, 1],
[2, 1, 1, 2, 2, 2, 1, 2, 2, 1, 2, 2, 1, 2, 2, 2, 2, 2, 4, 2, 2, 1, 2, 1],
[2, 2, 1, 1, 2, 1, 2, 2, 2, 1, 1, 2, 1, 2, 2, 2, 2, 2, 2, 4, 2, 1, 1, 1],
[0, 1, 1, 1, 1, 0, 0, 2, 1, 0, 0, 2, 1, 2, 2, 2, 1, 2, 2, 2, 4, 2, 2, 2],
[0, 0, 0, 0, 0, 0, 0, 0, 1, 1, 0, 1, 1, 2, 1, 1, 1, 2, 1, 1, 2, 4, 2, 2],
[0, 0, 0, 0, 0, 0, 0, 0, 1, 0, 1, 1, 1, 1, 2, 1, 1, 1, 2, 1, 2, 2, 4, 2],
[-3, -1, -1, -1, -1, -1, -1, 1, -1, -1, -1, 1, -1, 1, 1, 1, -1, 1, 1, 1, 2, 2, 2, 4]
],
'SimpleHexagonal': [
[1, 0, 0],
[0, 2, 1],
[0, 1, 2]
],
'SimpleOrthorhombic': [
[1, 0, 0],
[0, 2, 0],
[0, 0, 3]
],
}
"""
Internal test script
This is not meant to appear in the final Sage release.
"""
from sage.all import *
from diamond_cutting import *
import time
"""L = Lattice([[2, 0], [0, 3]])
g = plot_lattice(L)
g.show()
L = special_lattice('BodyCenteredCubic')
print L.embedded_basis()
plot_lattice(L).show(viewer='tachyon')
"""
L = special_lattice('SimpleCubic')
g = plot_lattice(L)
V = L.voronoi_cell()
g += V.plot()
g.show(viewer='tachyon')
g.save('lattice.png')
#L = Lattice([[GF(3)(1), 0], [0, 1]])
#print L
#print L.zero() + L.an_element()
#L_complex = Lattice([[i, 0], [0, 1]])
#L = Lattice([[2, 0], [0, 1]])
#print L([2, 1])
#L = Lattice([[6, 1], [9, 0]])
#print L
#L.voronoi_cell()
"""L = Lattice([[1, 2, 0], [4, 0, 6], [0, 8, 9]])
#V = L.voronoi_cell()
V = calculate_voronoi_cell(L.embedded_basis_matrix(), radius=None, debug=True)
print V
print V.Vrepresentation()
print L.embedded_basis_matrix()
print L.gram_matrix()
print L.discriminant()
print V.volume()
print V.lrs_volume()
L = Lattice(matrix([[1, 2, 3, 4], [5, 6, 7, 8], [9, 10, 11, 12], [13, 14, 15, 16]]))
#V = L.voronoi_cell(radius=100)
V = calculate_voronoi_cell(L.embedded_basis_matrix(), radius=847, debug=True)
print V.volume()
print V.lrs_volume()
print L.discriminant()
"""
cr63 = """
4.47 0.00 0.00 0.00 0.00 0.00
-3.35 2.96 0.00 0.00 0.00 0.00
1.34 -3.55 2.37 0.00 0.00 0.00
-0.22 1.77 -3.55 2.05 0.00 0.00
-0.22 -0.59 1.77 -3.76 1.53 0.00
1.34 1.18 0.59 3.07 -2.29 1.32
"""
cr6 = """
2.4494897427832E+00 0.0000000000000E+00 0.0000000000000E+00 0.0000000000000E+00 0.0000000000000E+00 0.0000000000000E+00
-1.6329931618555E+00 1.8257418583506E+00 0.0000000000000E+00 0.0000000000000E+00 0.0000000000000E+00 0.0000000000000E+00
4.0824829046386E-01 -1.8257418583506E+00 1.5811388300842E+00 0.0000000000000E+00 0.0000000000000E+00 0.0000000000000E+00
0.0000000000000E+00 5.4772255750517E-01 -1.8973665961010E+00 1.4491376746189E+00 0.0000000000000E+00 0.0000000000000E+00
0.0000000000000E+00 0.0000000000000E+00 6.3245553203368E-01 -1.9321835661586E+00 1.3662601021279E+00 0.0000000000000E+00
4.0824829046386E-01 3.6514837167011E-01 3.1622776601684E-01 9.6609178307930E-01 -1.7078251276599E+00 1.3228756555323E+00
"""
m2a = """
2 0
0 1
"""
m2b = """
-2 -2
2 -2
"""
es6 = """
1.41 0.00 0.00 0.00 0.00 0.00
-0.71 1.22 0.00 0.00 0.00 0.00
0.00 -0.82 1.15 0.00 0.00 0.00
0.00 0.00 -0.87 1.12 0.00 0.00
0.00 -0.82 -0.58 -0.45 0.37 0.00
0.00 0.00 -0.87 -0.67 0.55 0.71
"""
lv6 = """
1.4135770e+000 0.0000000e+000 0.0000000e+000 0.0000000e+000 0.0000000e+000 0.0000000e+000
3.7281308e-001 1.3635287e+000 0.0000000e+000 0.0000000e+000 0.0000000e+000 0.0000000e+000
-2.9499631e-001 -2.2516690e-001 1.4037368e+000 0.0000000e+000 0.0000000e+000 0.0000000e+000
-3.7281308e-001 -2.8456344e-001 -8.7491668e-001 1.0063572e+000 0.0000000e+000 0.0000000e+000
3.7281308e-001 2.8456344e-001 -1.7307180e-001 -4.5556135e-001 1.2412671e+000 0.0000000e+000
-7.4569690e-001 -5.6918088e-001 3.4619130e-001 -5.5058489e-001 -6.4855953e-001 6.2011907e-001
"""
matas6 = """
2.4494897427832E+00 0.0000000000000E+00 0.0000000000000E+00 0.0000000000000E+00 0.0000000000000E+00 0.0000000000000E+00
-4.0824829046386E-01 2.4152294576982E+00 0.0000000000000E+00 0.0000000000000E+00 0.0000000000000E+00 0.0000000000000E+00
-4.0824829046386E-01 -4.8304589153965E-01 2.3664319132398E+00 0.0000000000000E+00 0.0000000000000E+00 0.0000000000000E+00
-4.0824829046386E-01 -4.8304589153965E-01 -5.9160797830996E-01 2.2912878474779E+00 0.0000000000000E+00 0.0000000000000E+00
-4.0824829046386E-01 -4.8304589153965E-01 -5.9160797830996E-01 -7.6376261582597E-01 2.1602468994693E+00 0.0000000000000E+00
-4.0824829046386E-01 -4.8304589153965E-01 -5.9160797830996E-01 -7.6376261582597E-01 -1.0801234497346E+00 1.8708286933870E+00
"""
def volume_simplex(S):
n = len(S)
V = matrix([vector(list(S[k] - S[0])) for k in range(1, n)]).transpose()
return abs(V.determinant()) / factorial(n - 1)
def volume(P):
from scipy.spatial import Delaunay
if hasattr(P, 'vertices'):
V = P.vertices()
else:
V = P
D = Delaunay(V)
vol = 0
for simplex in D.vertices:
vol += volume_simplex([D.points[index] for index in simplex])
return vol
def test_lattice(data):
M = []
for line in data.splitlines():
row = []
for value in line.split(' '):
value = value.strip()
if value:
value = QQ(float(value))
row.append(value)
if row:
M.append(row)
M = matrix(M)
start = time.clock()
V = calculate_voronoi_cell(M, radius=10, debug=True)
stop = time.clock()
print V.Vrepresentation()
print V.Hrepresentation()
print "Computed Voronoi cell in %f seconds" % (stop - start)
print "Volume: %f" % volume(V)
L = Lattice(M)
print "Det: %f" % L.determinant()
#test_lattice(cr63)
#test_lattice("""
#2 543 2 123
#432 425 123 5
#43 21 8 0
#45 7 0 1
#""")
m3 = """
1 2 0
0 3 4
3 1 0
"""
#test_lattice(m3)
#test_lattice(m2b)
#test_lattice(cr6)
#test_lattice(lv6)
#test_lattice(matas6)
#test_lattice("1 0 0 0 0\n0 2 1 0 0")
#L = Lattice([[1, 0], [2, 0], [0, 2]])
#print L.voronoi_cell().Vrepresentation()
#L = Lattice(quadratic_form=[[2,0], [0,2]])
#L = Lattice([[-1, -1, 0], [1, -1, 0], [0, 1, -1]])
#V = L.voronoi_cell()
#print V.lrs_volume()
"""L = Lattice(quadratic_form=[[1,0], [0,1]])
print L.closest_vector([0.4, 0.2])
print L.closest_vector([1.5, 2.5])
print L.closest_vector([1, 1])
print L.closest_vector([6, 0.6])
L = Lattice([[1,0,0], [0,1,0]])
print L.closest_vector([4, 2.6, 3.6])"""
"""print Lattice(basis=[[1, 0], [0, 1]])
print Lattice(basis=[[3, 2, 0], [1, 0, 1]])
print Lattice(basis=[[ZZ(1) / ZZ(2), 0], [0, 1]])
print Lattice(basis=[[0.5, 0], [0, 1]])
print Lattice(basis=[[1, 0], [0, 1]], base_ring=Zp(2))
L = Lattice([[6, 1], [9, 0]])
print L
print L.gram_matrix()
print L.determinant()
print L.discriminant()
L = Lattice([[2, 0], [0, 3]])
print L.shortest_vectors()
leech_lattice = [
[8, 4, 4, 4, 4, 4, 4, 2, 4, 4, 4, 2, 4, 2, 2, 2, 4, 2, 2, 2, 0, 0, 0, -3],
[4, 4, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 1, 1, 2, 1, 1, 2, 1, 0, 0, -1],
[4, 2, 4, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 1, 2, 1, 2, 2, 1, 1, 1, 0, 0, -1],
[4, 2, 2, 4, 2, 2, 2, 2, 2, 2, 2, 2, 2, 1, 1, 2, 2, 1, 2, 1, 1, 0, 0, -1],
[4, 2, 2, 2, 4, 2, 2, 2, 2, 2, 2, 1, 2, 2, 2, 2, 2, 2, 2, 2, 1, 0, 0, -1],
[4, 2, 2, 2, 2, 4, 2, 2, 2, 2, 2, 1, 2, 2, 1, 1, 2, 1, 2, 1, 0, 0, 0, -1],
[4, 2, 2, 2, 2, 2, 4, 2, 2, 2, 2, 1, 2, 1, 2, 1, 2, 1, 1, 2, 0, 0, 0, -1],
[2, 2, 2, 2, 2, 2, 2, 4, 1, 1, 1, 2, 1, 2, 2, 2, 1, 2, 2, 2, 2, 0, 0, 1],
[4, 2, 2, 2, 2, 2, 2, 1, 4, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 1, 1, 1, -1],
[4, 2, 2, 2, 2, 2, 2, 1, 2, 4, 2, 2, 2, 2, 1, 1, 2, 2, 1, 1, 0, 1, 0, -1],
[4, 2, 2, 2, 2, 2, 2, 1, 2, 2, 4, 2, 2, 1, 2, 1, 2, 1, 2, 1, 0, 0, 1, -1],
[2, 2, 2, 2, 1, 1, 1, 2, 2, 2, 2, 4, 1, 2, 2, 2, 1, 2, 2, 2, 2, 1, 1, 1],
[4, 2, 2, 2, 2, 2, 2, 1, 2, 2, 2, 1, 4, 2, 2, 2, 2, 1, 1, 1, 1, 1, 1, -1],
[2, 2, 1, 1, 2, 2, 1, 2, 2, 2, 1, 2, 2, 4, 2, 2, 1, 2, 2, 2, 2, 2, 1, 1],
[2, 1, 2, 1, 2, 1, 2, 2, 2, 1, 2, 2, 2, 2, 4, 2, 1, 2, 2, 2, 2, 1, 2, 1],
[2, 1, 1, 2, 2, 1, 1, 2, 2, 1, 1, 2, 2, 2, 2, 4, 1, 2, 2, 2, 2, 1, 1, 1],
[4, 2, 2, 2, 2, 2, 2, 1, 2, 2, 2, 1, 2, 1, 1, 1, 4, 2, 2, 2, 1, 1, 1, -1],
[2, 1, 2, 1, 2, 1, 1, 2, 2, 2, 1, 2, 1, 2, 2, 2, 2, 4, 2, 2, 2, 2, 1, 1],
[2, 1, 1, 2, 2, 2, 1, 2, 2, 1, 2, 2, 1, 2, 2, 2, 2, 2, 4, 2, 2, 1, 2, 1],
[2, 2, 1, 1, 2, 1, 2, 2, 2, 1, 1, 2, 1, 2, 2, 2, 2, 2, 2, 4, 2, 1, 1, 1],
[0, 1, 1, 1, 1, 0, 0, 2, 1, 0, 0, 2, 1, 2, 2, 2, 1, 2, 2, 2, 4, 2, 2, 2],
[0, 0, 0, 0, 0, 0, 0, 0, 1, 1, 0, 1, 1, 2, 1, 1, 1, 2, 1, 1, 2, 4, 2, 2],
[0, 0, 0, 0, 0, 0, 0, 0, 1, 0, 1, 1, 1, 1, 2, 1, 1, 1, 2, 1, 2, 2, 4, 2],
[-3, -1, -1, -1, -1, -1, -1, 1, -1, -1, -1, 1, -1, 1, 1, 1, -1, 1, 1, 1, 2, 2, 2, 4]
]
L = Lattice(inner_product_matrix=leech_lattice)
print L.shortest_vectors_count()
L = Lattice([[1, 0, 0], [0, 1, 0], [0, 0, 1]])
print L.shortest_vectors()
"""
"""
#GM = matrix([[0, 3], [3, -1]])
L = Lattice([[1, 0], [0, 1]])
GM = L.basis()
V = calculate_voronoi_cell(matrix(GM))
print V.Hrepresentation()
print V.Vrepresentation()
"""
......@@ -2418,20 +2418,17 @@ cdef class Matrix_integer_dense(matrix_dense.Matrix_dense): # dense or sparse
# LLL
####################################################################################
def LLL_gram(self):
"""
LLL reduction of the lattice whose gram matrix is self.
"""LLL reduction of the lattice whose gram matrix is self.
INPUT:
- ``M`` - gram matrix of a definite quadratic form
- ``M`` - gram matrix of a definite quadratic form
OUTPUT:
- ``U`` - unimodular transformation matrix such that
U.transpose() \* M \* U is LLL-reduced.
- ``U`` - unimodular transformation matrix such that U.transpose() * M *
U is LLL-reduced.
ALGORITHM: Use PARI
......
......@@ -2696,6 +2696,46 @@ cdef class Matrix_rational_dense(matrix_dense.Matrix_dense):
mpq_init(v._entries[j]); mpq_set(v._entries[j], self._matrix[i][j])
return v
def LLL(self, *args, **kwargs):
"""Return an LLL reduced or approximated LLL reduced lattice for
``self`` interpreted as a lattice.
For details, see
:meth:`sage.matrix.matrix_integer_dense.Matrix_integer_dense.LLL`.
EXAMPLE::
sage: A = Matrix(QQ, 3, 3, [1/n for n in range(1, 10)])
sage: A.LLL()
[ 1/28 -1/40 -1/18]
[ 1/28 -1/40 1/18]
[ 0 -3/40 0]
"""
A, d = self._clear_denom()
return A.LLL(*args, **kwargs) / d
def LLL_gram(self, *args, **kwargs):
"""LLL reduction of the lattice whose gram matrix is ``self``.
For details, see
:meth:`sage.matrix.matrix_integer_dense.Matrix_integer_dense.LLL_gram`.
EXAMPLES::
sage: M = Matrix(ZZ, 2, 2, [5,3,3,2]) / 7 ; M
[5/7 3/7]
[3/7 2/7]
sage: U = M.LLL_gram(); U
[-1 1]
[ 1 -2]
sage: U.transpose() * M * U
[1/7 0]
[ 0 1/7]
"""
A, d = self._clear_denom()
return A.LLL_gram(*args, **kwargs)
def column(self, Py_ssize_t i, from_list=False):
"""
Return the i-th column of this matrix as a dense vector.
......
......@@ -357,7 +357,7 @@ class FreeQuadraticModule_generic(free_module.FreeModule_generic):
r = n//2
else:
r = (n-1)//2
return (-1)^r*self.gram_matrix().determinant()
return (-1)**r * self.gram_matrix().determinant()
def gram_matrix(self):
"""
......