cone/faces.py: add preliminary dual_face() operation.

diff --git a/mjo/cone/faces.py b/mjo/cone/faces.py
index 8e5f2436c47802f4e8c61fbb2166f44d7370537d..a6820164999e55d6396d5d1b081419e57f5c791b 100644 (file)
--- a/mjo/cone/faces.py
+++ b/mjo/cone/faces.py
@@ -1,5 +1,6 @@
 from sage.all import *
 
+
 def face_generated_by(K,S):
     r"""
     Return the intersection of all faces of ``K`` that contain ``S``.
@@ -130,3 +131,68 @@ def face_generated_by(K,S):
         raise ValueError('S is not a subset of the cone')
     else:
         return face_lattice.sorted(candidates)[0]
+
+
+def dual_face(K,F):
+    r"""
+    Return the dual face of ``F`` with respect to the cone ``K``.
+
+    OUTPUT:
+
+    A face of ``K.dual()``.
+
+    REFERENCES:
+
+    .. [Tam] Bit-Shun Tam. On the duality operator of a convex cone. Linear
+       Algebra and its Applications, 64:33-56, 1985, doi:10.1016/0024-3795(85)
+       90265-4.
+
+    SETUP::
+
+        sage: from mjo.cone.faces import dual_face
+
+    EXAMPLES:
+
+    The dual face of the first standard basis vector in three dimensions
+    is the face generated by the other two standard basis vectors::
+
+        sage: K = Cone([(1,0,0),(0,1,0),(0,0,1)])
+        sage: F = Cone([(1,0,0)])
+        sage: dual_face(K,F).rays()
+        M(0, 0, 1),
+        M(0, 1, 0)
+        in 3-d lattice M
+
+    TESTS:
+
+    The dual face of ``K`` with respect to itself should be the
+    lineality space of its dual [Tam]_::
+
+        sage: set_random_seed()
+        sage: K = random_cone(max_ambient_dim=8, max_rays=10)
+        sage: K_dual = K.dual()
+        sage: lKd_gens = [ dir*l for dir in [1,-1] for l in K_dual.lines() ]
+        sage: linspace_K_dual = Cone(lKd_gens, K_dual.lattice())
+        sage: dual_face(K,K).is_equivalent(linspace_K_dual)
+        True
+
+    If ``K`` is proper, then the dual face of its trivial face is the
+    dual of ``K`` [Tam]_::
+
+        sage: set_random_seed()
+        sage: K = random_cone(max_ambient_dim=8,
+        ....:                 max_rays=10,
+        ....:                 strictly_convex=True,
+        ....:                 solid=True)
+        sage: L = K.lattice()
+        sage: trivial_face = Cone([L.zero()], L)
+        sage: dual_face(K,trivial_face).is_equivalent(K.dual())
+        True
+
+    """
+    # Ensure that F is actually a face of K before continuing.
+    if not F.is_face_of(K):
+        raise ValueError("%s is not a face of %s" % (F,K))
+
+    span_F = Cone([c*g for c in [1,-1] for g in F], F.lattice())
+    return K.dual().intersection(span_F.dual())