return t
-def nonlineal_part(K):
+def nonlineal_part(K, orthogonal=False):
r"""
- Return the pointed component of this cone.
-
- Every closed convex cone can be decomposed into the direct sum of
- a vector space (its :meth:`linear_subspace`) and a pointed closed
- convex cone that lives in the orthogonal complement of the
- lineality space. The straightforward way to obtain the generators
- of the pointed component is through orthogonal projection, but we
- can't (directly) orthogonally project over the rationals.
-
- This function takes a different approach, positively scaling the
- generators that live outside of the lineality space, and then
- removing their lineal parts. The resulting generators (along with
- those for the lineality space) still generate the original cone
- because the original generating set can be recovered by adding
- back the positive multiples of the orthogonal basis for the
- lineality space.
+ Return a pointed cone that is complementary to this cone's
+ :meth:`linear_subspace`.
+
+ The largest :meth:`linear_subspace` contained in a given convex
+ cone is commonly known as its *lineality space*. Keeping in mind
+ that vector subspaces are also convex cones, Wets and Witzgall
+ showed that a conic generating set for the lineality space can be
+ obtained from that of the original cone by selecting those that
+ live in the lineality space [WetsWitz1967]_. Correspondingly, Wets
+ and Witzgall call that subset the *lineal part* of the generating
+ set. They call the remaining generators the *conical part*, but
+ this is too confusing; we shall refer to it as the *nonlineal
+ part*.
+
+ It is also known that every closed convex cone can be decomposed
+ as the orthogonal direct sum of its lineality space and a pointed
+ closed convex cone [StoerWitz1970]_. (In the decomposition, the
+ pointed component lives in the orthogonal complement of the
+ lineality space.) A (minimal) generating set for the pointed
+ component can be obtained by orthogonally projecting the nonlineal
+ part of the given (minimal) generating set [MayerWarth2024]_.
+
+ From that, it can be shown that the nonlineal part (before the
+ orthogonal projection) already generates a pointed convex cone. So
+ in either case (projection or no) we express the given cone as a
+ sum of its lineality space and a pointed cone. Abusing the
+ notation, we call that pointed cone the *nonlineal part* of the
+ cone. This method returns it, with the ``orthgonal`` argument
+ determining whether or not we "project" the nonlineal part of the
+ generating set onto the orthogonal complement of the lineality
+ space.
+
+ .. NOTE:
+
+ If we forego the orthogonal projection, the sum of cones
+ need not be direct.
+
+ The straightforward way to obtain the generators of the pointed
+ component (when ``orthogonal`` is ``True`` would be to just do the
+ orthogonal projection. That may however not be possible over the
+ rationals, and our convex cones have rational generators. So we
+ take a different approach: positively scaling the generators that
+ live outside of the lineality space, and then removing their
+ lineal parts. The resulting generators (along with those for the
+ lineality space) still generate the original cone because the
+ original generating set can be recovered by adding back the
+ positive multiples of the orthogonal basis for the lineality
+ space. They are orthogonal to the lineality space by construction
+ (think Gram-Schmidt).
+
+ INPUT:
+
+ - ``orthogonal`` -- boolean (default: False); whether or not
+ the resulting cone should be orthogonal to the lineality
+ space.
+
+ OUTPUT:
+
+ A pointed convex subcone that, when added to the lineality space
+ of the given cone, recovers it. Mayer and von der Warth show that
+ the sum is the same in either case, but to show that the result is
+ pointed when ``orthogonal`` is ``False`` requires another step
+ [MayerWarth2024]_.
+
+ REFERENCES:
+
+ .. [WetsWitz1967] Roger Wets and Christoph Witzgall.
+ Algorithms for frames and lineality spaces of cones.
+ Journal of Research of the National Bureau of Standards,
+ Section B: Mathematics and Mathematical Physics,
+ 71B(1):1-7, 1967. :doi:`10.6028/jres.071B.001`.
+ .. [StoerWitz1970] Josef Stoer and Christoph Witzgall.
+ Convexity and Optimization in Finite Dimensions I, vol.
+ 163 of Grundlehren der mathematischen Wissenschaften.
+ Springer-Verlag, New York, 1970. ISBN 9783642462184,
+ :doi:`10.1007/978-3-642-46216-0`.
+ - [MayerWarth2024] Matthias Georg Mayer and Fabian von der Warth.
+ An algorithm for minimum cardinality generators of cones.
+ :arxiv:`2412.01451`.
SETUP::
EXAMPLES:
- It's difficult to see this function working because Sage (via PPL)
- minimizes the generating sets of its cones by default. In many
- cases, this minimization process will do the same thing that this
- method does, subtracting out the lineal portions of the generators
- that do not live in the lineality space. Here is one example of
- the ell-one cone living above a line::
+ It's difficult to see ``orthogonal=True`` working because Sage
+ (via PPL) minimizes the generating sets of its cones by
+ default. In many cases, this minimization process will do the same
+ thing that our orthogonalization does, subtracting out the lineal
+ portions of the generators that do not live in the lineality
+ space. Here is one example of the ell-one cone living above a
+ line::
sage: N = ToricLattice(3)
sage: G = [ (1,0,0), (-1,0,0),
....: (1,0,1), (-1,0,1), (0,1,1), (0,-1,1) ]
sage: K = Cone(list(map(N,G)), check=False, normalize=False)
- sage: P = nonlineal_part(K)
+ sage: P = nonlineal_part(K, orthogonal=True)
sage: P.rays()
N(0, -1, 1),
N(0, 1, 1)
sage: K.ray(0)*M(K.lines()[0].dense_coefficient_list())
-146
- But we get it right::
+ But we do it with ``orthogonal=True``::
- sage: P = nonlineal_part(K)
+ sage: P = nonlineal_part(K, orthogonal=True)
sage: P.rays()
N(7, 223, -146, -123)
in 4-d lattice N
sage: P.ray(0)*M(K.lines()[0].dense_coefficient_list())
0
- TESTS::
+ TESTS:
+
+ Verify that the sum of the lineality space and the nonlineal part
+ is the original cone for both values of ``orthogonal``. When
+ ``orthogonal`` is ``True``, we verify that orthogonality::
sage: K = random_cone(strictly_convex=False)
sage: lat = K.lattice()
- sage: P = nonlineal_part(K)
sage: V = lat.vector_space()
+ sage: linspace_rays = [ c*l for l in K.lines() for c in [-1,1] ]
+ sage: P = nonlineal_part(K)
+ sage: J = Cone(list(P.rays()) + linspace_rays, lattice=lat)
+ sage: J.is_equivalent(K)
+ True
+ sage: P = nonlineal_part(K, orthogonal=True)
sage: all( V(p).inner_product(V(l)).is_zero()
....: for p in P
....: for l in K.lines())
True
- sage: linspace_rays = [ c*l for l in K.lines() for c in [-1,1] ]
sage: J = Cone(list(P.rays()) + linspace_rays, lattice=lat)
sage: J.is_equivalent(K)
True
The nonlineal part of a pointed cone is itself::
sage: K = random_cone(strictly_convex=True)
- sage: nonlineal_part(K) is K
+ sage: nonlineal_part(K, orthogonal=True) is K
+ True
+ sage: nonlineal_part(K, orthogonal=False) is K
True
Verify the assumption that :meth:`lines` is always orthogonal,
sage: K = random_cone(strictly_convex=False)
sage: nonlineal_part(K).is_strictly_convex()
+ sage: nonlineal_part(K, orthogonal=True).is_strictly_convex()
True
sage: V = K.lattice().vector_space()
sage: L = [V(l) for l in K.lines()]
....: if i != j )
True
+ An example where the sum is not direct. The idea is that
+ ``[l1,l2]`` and ``[h1,h2,h3,h4]`` begin as orthogonal generating
+ sets for the lineality space and pointed component of a cone
+ ``K``. By adding ``l1`` and ``l2`` to the ``hN``, we obtain a new
+ set ``[g1,g2,g3,g4]``. The set ``G = [l1,l2,g1,g2,g3,g4]`` should
+ still generate the same ``K``, and its lineality space has not
+ changed, but the nonlineal part of ``G``, namely the ``gN``, now
+ have a span that intersects the lineality space nontrivially; in
+ other words, the sum `K == Cone([l1,l2]) + Cone([g1,g2,g3,g4])``
+ is not direct.
+
+ By default PPL will minimize away the lineal part of the ``gN``,
+ so we have to check ourselves that the set ``G`` is a minimal set
+ of generators for ``K``, and then construct a new cone using
+ ``check=False`` and ``normalize=False``. We begin by declaring
+ the various generators we'll use::
+
+ sage: V = QQ^4
+ sage: l1 = V([0,0,0,1])
+ sage: l2 = V([0,0,0,-1])
+ sage: h1 = V([1,0,1,0])
+ sage: h2 = V([-1,0,1,0])
+ sage: h3 = V([0,1,1,0])
+ sage: h4 = V([0,-1,1,0])
+ sage: g1 = h1 + l1
+ sage: g2 = h2 + l1
+ sage: g3 = h3 + l2
+ sage: g4 = h4 + l2
+
+ Build the cone ``K``::
+
+ sage: H = [l1,l2,h1,h2,h3,h4]
+ sage: K = Cone(H)
+
+ Next we confirm that ``G`` is also minimal generating set for
+ ``K``::
+
+ sage: G = [l1,l2,g1,g2,g3,g4]
+ sage: any( g in Cone(G_bar)
+ ....: for g in G
+ ....: if (G_bar := [x for x in G if not x == g]) )
+ False
+
+ Build a new (equivalent) cone from ``G``::
+
+ sage: J = Cone(G, check=False, normalize=False, lattice=K.lattice())
+ sage: list(J.rays()) == G
+ True
+ sage: K.is_equivalent(J)
+ True
+
+ Finally, the sum is not direct::
+
+ sage: L = J.linear_subspace(); L
+ Vector space of degree 4 and dimension 1 over Rational Field
+ Basis matrix:
+ [0 0 0 1]
+ sage: P = nonlinear_part(J)
+ sage: list(P.rays()) == [g1,g2,g3,g4]
+ True
+ sage: L.intersection(span(J.rays())).dimension()
+ 1
+
"""
+ lat = K.lattice()
+
+ if not orthogonal:
+ # Open question: is this generating set guaranteed to be
+ # minimal?
+ return Cone([r for r in K if not r in K.lines()], lattice=lat)
+
# We need to compute inner products and doing it in the ambient
# vector space (rather than in a more picky lattice) is the
# easiest way to do it.
- lat = K.lattice()
V = lat.vector_space()
linspace_basis = [V(l) for l in K.lines()]
if not linspace_basis: