# HG changeset patch
# User Andrey Novoseltsev <novoselt@gmail.com>
# Date 1276029901 21600
# Node ID 125efc24ad6b1db4e2e95c1ae54458f8c9307dfc
# Parent 09dbbfb5c3f65c62405778b8d5c4b15e5da91771
Trac 9188: Bug in facet_normal of lattice polytopes  reviewer changes.
Added abs(...) to dual scaling and prettified doctests.
diff r 09dbbfb5c3f6 r 125efc24ad6b sage/geometry/lattice_polytope.py
a

b


521  521  self._shift_vector = p0 
522  522  if compute_vertices: 
523  523  self._vertices = self._embed(H_polytope.vertices()) 
524   
 524  # In order to use facet normals obtained from subpolytopes, we 
 525  # need the following (see Trac #9188). 
525  526  M = self._embedding_matrix 
526   basis = matrix( M.columns() + M.integer_kernel().basis() ).transpose() 
527   self._dual_embedding_scale = basis.det() 
 527  # Basis for the ambient space with spanned subspace in front 
 528  basis = M.columns() + M.integer_kernel().basis() 
 529  # Let's represent it as columns of a matrix 
 530  basis = matrix(basis).transpose() 
 531  # Absolute value helps to keep normals "inner" 
 532  self._dual_embedding_scale = abs(basis.det()) 
528  533  dualbasis = matrix(ZZ, self._dual_embedding_scale * basis.inverse()) 
529  534  self._dual_embedding_matrix = dualbasis.submatrix(0,0,M.ncols()) 
530  535  
… 
… 

1513  1518  sage: p.facet_constant(0) 
1514  1519  1 
1515  1520  
1516   A lattice polytope of strictly smaller dimension (=2) than the 
1517   ambient dimension(=4):: 
1518   
1519   sage: LatticePolytope(matrix([(1, 1, 0), (1, 1, 0), (1, 1, 0), (3, 3, 0)])) 
 1521  This is a 2dimensional lattice polytope in a 4dimensional space:: 
 1522  
 1523  sage: m = matrix([(1,1,1,3), (1,1,1,3), (0,0,0,0)]) 
 1524  sage: p = LatticePolytope(m.transpose()) 
 1525  sage: p 
1520  1526  A lattice polytope: 2dimensional, 3 vertices. 
1521   sage: lp = LatticePolytope(matrix([(1, 1, 0), (1, 1, 0), (1, 1, 0), (3, 3, 0)])) 
1522   sage: lp.vertices() 
 1527  sage: p.vertices() 
1523  1528  [ 1 1 0] 
1524  1529  [1 1 0] 
1525  1530  [ 1 1 0] 
1526  1531  [ 3 3 0] 
1527   sage: matrix([[ lp.facet_normal(i)*lp.vertex(j) + lp.facet_constant(i) for i in range(0,3)] for j in range(0,3)]) 
1528   [22 0 0] 
1529   [ 0 22 0] 
1530   [ 0 0 11] 
 1532  sage: fns = [p.facet_normal(i) for i in range(p.nfacets())] 
 1533  sage: fns 
 1534  [(11, 1, 1, 3), (11, 1, 1, 3), (0, 1, 1, 3)] 
 1535  sage: fcs = [p.facet_constant(i) for i in range(p.nfacets())] 
 1536  sage: fcs 
 1537  [0, 0, 11] 
 1538  
 1539  Now we manually compute the distance matrix of this polytope. Since it 
 1540  is a triangle, each line (corresponding to a facet) should have two 
 1541  zeros (vertices of the corresponding facet) and one positive number 
 1542  (since our normals are inner):: 
 1543  
 1544  sage: matrix([[fns[i] * p.vertex(j) + fcs[i] 
 1545  ... for j in range(p.nvertices())] 
 1546  ... for i in range(p.nfacets())]) 
 1547  [22 0 0] 
 1548  [ 0 22 0] 
 1549  [ 0 0 11] 
1531  1550  """ 
1532  1551  if self.is_reflexive(): 
1533  1552  return 1 
1534  1553  elif self.ambient_dim() == self.dim(): 
1535  1554  return self._facet_constants[i] 
1536  1555  else: 
1537   return (self._sublattice_polytope.facet_constant(i) * self._dual_embedding_scale 
1538    self.facet_normal(i) * self._shift_vector ) 
 1556  return (self._sublattice_polytope.facet_constant(i) 
 1557  * self._dual_embedding_scale 
 1558   self.facet_normal(i) * self._shift_vector) 
1539  1559  
1540  1560  def facet_normal(self, i): 
1541  1561  r""" 
1542  1562  Return the inner normal to the ``i``th facet of this polytope. 
1543  1563  
1544  1564  If this polytope is not fulldimensional, facet normals will be 
1545   parallel to the affine subspace spanned by this polytope. 
 1565  orthogonal to the integer kernel of the affine subspace spanned by 
 1566  this polytope. 
1546  1567  
1547  1568  INPUT: 
1548  1569  
… 
… 

1588  1609  sage: p.facet_constant(0) 
1589  1610  1 
1590  1611  
1591   Here is an example of a 3dimensional polytope in 4d:: 
1592   
1593   sage: lp = LatticePolytope(matrix([[1,1,1,0],[1,1,1,0],[1,1,1,0],[3,3,3,0]])) 
1594   sage: lp.vertices() 
1595   [ 1 1 1 0] 
1596   [ 1 1 1 0] 
1597   [ 1 1 1 0] 
1598   [ 3 3 3 0] 
1599   sage: ker = lp.vertices().integer_kernel().matrix() 
 1612  Here is an example of a 3dimensional polytope in a 4dimensional 
 1613  space:: 
 1614  
 1615  sage: m = matrix([(0,0,0,0), (1,1,1,3), (1,1,1,3), (1,1,1,3)]) 
 1616  sage: p = LatticePolytope(m.transpose()) 
 1617  sage: p.vertices() 
 1618  [ 0 1 1 1] 
 1619  [ 0 1 1 1] 
 1620  [ 0 1 1 1] 
 1621  [ 0 3 3 3] 
 1622  sage: ker = p.vertices().integer_kernel().matrix() 
1600  1623  sage: ker 
1601  1624  [ 0 0 3 1] 
1602   sage: [ ker * lp.facet_normal(i) for i in range(0,4) ] 
 1625  sage: [ker * p.facet_normal(i) for i in range(p.nfacets())] 
1603  1626  [(0), (0), (0), (0)] 
1604   sage: matrix([lp.facet_normal(i) for i in range(0,4)]) * lp.vertices() 
1605   [ 0 0 20 0] 
1606   [ 0 20 0 0] 
1607   [ 10 10 10 0] 
1608   [20 0 0 0] 
1609   sage: matrix([[ lp.facet_normal(i)*lp.vertex(j) + lp.facet_constant(i) for i in range(0,4)] for j in range(0,4)]) 
1610   [ 0 0 0 20] 
1611   [ 0 20 0 0] 
1612   [20 0 0 0] 
1613   [ 0 0 10 0] 
 1627  
 1628  Now we manually compute the distance matrix of this polytope. Since it 
 1629  is a simplex, each line (corresponding to a facet) should consist of 
 1630  zeros (indicating generating vertices of the corresponding facet) and 
 1631  a single positive number (since our normals are inner):: 
 1632  
 1633  sage: matrix([[p.facet_normal(i) * p.vertex(j) 
 1634  ... + p.facet_constant(i) 
 1635  ... for j in range(p.nvertices())] 
 1636  ... for i in range(p.nfacets())]) 
 1637  [ 0 0 0 20] 
 1638  [ 0 0 20 0] 
 1639  [ 0 20 0 0] 
 1640  [10 0 0 0] 
1614  1641  """ 
1615  1642  if self.is_reflexive(): 
1616  1643  return self.polar().vertex(i) 