Add the project_span() function.

diff --git a/mjo/cone/cone.py b/mjo/cone/cone.py
index 7597cd66ab5275ec914df0c4929152b2937bacf0..60f9c34ec8bc271d65812859f51ca77636c8cbbc 100644 (file)
--- a/mjo/cone/cone.py
+++ b/mjo/cone/cone.py
@@ -7,6 +7,56 @@ addsitedir(abspath('../../'))
 
 from sage.all import *
 
+def project_span(K):
+    r"""
+    Project ``K`` into its own span.
+
+    EXAMPLES::
+
+        sage: K = Cone([(1,)])
+        sage: project_span(K) == K
+        True
+
+        sage: K2 = Cone([(1,0)])
+        sage: project_span(K2).rays()
+        N(1)
+        in 1-d lattice N
+        sage: K3 = Cone([(1,0,0)])
+        sage: project_span(K3).rays()
+        N(1)
+        in 1-d lattice N
+        sage: project_span(K2) == project_span(K3)
+        True
+
+    TESTS:
+
+    The projected cone should always be solid::
+
+        sage: K = random_cone()
+        sage: K_S = project_span(K)
+        sage: K_S.is_solid()
+        True
+
+    If we do this according to our paper, then the result is proper::
+
+        sage: K = random_cone()
+        sage: K_S = project_span(K)
+        sage: P = project_span(K_S.dual()).dual()
+        sage: P.is_proper()
+        True
+
+    """
+    F = K.lattice().base_field()
+    Q = K.lattice().quotient(K.sublattice_complement())
+    vecs = [ vector(F, reversed(list(Q(r)))) for r in K.rays() ]
+
+    L = None
+    if len(vecs) == 0:
+        L = ToricLattice(0)
+
+    return Cone(vecs, lattice=L)
+
+
 def rename_lattice(L,s):
     r"""
     Change all names of the given lattice to ``s``.
@@ -458,5 +508,17 @@ def lyapunov_rank(K):
         sage: actual == expected
         True
 
+    And test with the project_span function::
+
+        sage: K = random_cone(max_dim=15)
+        sage: actual = lyapunov_rank(K)
+        sage: K_S = project_span(K)
+        sage: P = project_span(K_S.dual()).dual()
+        sage: l = K.linear_subspace().dimension()
+        sage: codim = K.lattice_dim() - K.dim()
+        sage: expected = lyapunov_rank(P) + K.dim()*(l + codim) + codim**2
+        sage: actual == expected
+        True
+
     """
     return len(LL(K))