doubly_nonnegative.py: Begin implementing the extreme vector algorithms.

diff --git a/mjo/cone/doubly_nonnegative.py b/mjo/cone/doubly_nonnegative.py
index f705e5e98b083d6cf85fac39148d171cdf452d99..b43c974b2076a5dbda14ba43f679fd8caee30327 100644 (file)
--- a/mjo/cone/doubly_nonnegative.py
+++ b/mjo/cone/doubly_nonnegative.py
@@ -29,7 +29,7 @@ def is_doubly_nonnegative(A):
 
     INPUT:
 
-      - ``A`` - The matrix in question
+    - ``A`` - The matrix in question
 
     OUTPUT:
 
@@ -67,9 +67,81 @@ def is_doubly_nonnegative(A):
 
 
 
+def has_admissible_extreme_rank(A):
+    """
+    The extreme matrices of the doubly-nonnegative cone have some
+    restrictions on their ranks. This function checks to see whether or
+    not ``A`` could be extreme based on its rank.
+
+    INPUT:
+
+    - ``A`` - The matrix in question
+
+    OUTPUT:
+
+    ``False`` if the rank of ``A`` precludes it from being an extreme
+    matrix of the doubly-nonnegative cone, ``True`` otherwise.
+
+    REFERENCE:
+
+    Hamilton-Jester, Crista Lee; Li, Chi-Kwong. Extreme Vectors of
+    Doubly Nonnegative Matrices. Rocky Mountain Journal of Mathematics
+    26 (1996), no. 4, 1371--1383. doi:10.1216/rmjm/1181071993.
+    http://projecteuclid.org/euclid.rmjm/1181071993.
+
+    EXAMPLES:
+
+    The zero matrix has rank zero, which is admissible::
+
+    sage: A = zero_matrix(QQ, 5, 5)
+    sage: has_admissible_extreme_rank(A)
+    True
+
+    """
+    if not A.is_symmetric():
+        raise ValueError('The matrix ``A`` must be symmetric.')
+
+    r = rank(A)
+    n = A.nrows() # Columns would work, too, since ``A`` is symmetric.
+
+    if r == 0:
+        # Zero is in the doubly-nonnegative cone.
+        return True
+
+    # See Theorem 3.1 in the cited reference.
+    if r == 2:
+        return False
+
+    if n.mod(2) == 0:
+        # n is even
+        return r <= max(1, n-3)
+    else:
+        # n is odd
+        return r <= max(1, n-2)
+
+
 def is_extreme_doubly_nonnegative(A):
     """
     Returns ``True`` if the given matrix is an extreme matrix of the
     doubly-nonnegative cone, and ``False`` otherwise.
+
+    EXAMPLES:
+
+    The zero matrix is an extreme matrix::
+
+        sage: A = zero_matrix(QQ, 5, 5)
+        sage: is_extreme_doubly_nonnegative(A)
+        True
+
     """
+
+    r = A.rank()
+
+    if r == 0:
+        # Short circuit, we know the zero matrix is extreme.
+        return True
+
+    if not is_admissible_extreme_rank(r):
+        return False
+
     raise NotImplementedError()