X-Git-Url: http://gitweb.michael.orlitzky.com/?a=blobdiff_plain;f=mjo%2Fcone%2Fcone.py;h=daeeb8936edb5416b269a131afe95c286da75598;hb=3da1b3bbfcc95758bcce3596a86b0697765bd177;hp=e9d70e42ef61e0b44d492e6d46687f383ab8b1c6;hpb=14571032f67138970cf702181970ecd1290f6ef2;p=sage.d.git diff --git a/mjo/cone/cone.py b/mjo/cone/cone.py index e9d70e4..daeeb89 100644 --- a/mjo/cone/cone.py +++ b/mjo/cone/cone.py @@ -40,14 +40,14 @@ def is_lyapunov_like(L,K): The identity is always Lyapunov-like in a nontrivial space:: sage: set_random_seed() - sage: K = random_cone(min_ambient_dim = 1, max_ambient_dim = 8) + sage: K = random_cone(min_ambient_dim=1, max_ambient_dim=8) sage: L = identity_matrix(K.lattice_dim()) sage: is_lyapunov_like(L,K) True As is the "zero" transformation:: - sage: K = random_cone(min_ambient_dim = 1, max_ambient_dim = 8) + sage: K = random_cone(min_ambient_dim=1, max_ambient_dim=8) sage: R = K.lattice().vector_space().base_ring() sage: L = zero_matrix(R, K.lattice_dim()) sage: is_lyapunov_like(L,K) @@ -56,7 +56,7 @@ def is_lyapunov_like(L,K): Everything in ``K.lyapunov_like_basis()`` should be Lyapunov-like on ``K``:: - sage: K = random_cone(min_ambient_dim = 1, max_ambient_dim = 6) + sage: K = random_cone(min_ambient_dim=1, max_ambient_dim=6) sage: all([ is_lyapunov_like(L,K) for L in K.lyapunov_like_basis() ]) True @@ -99,28 +99,87 @@ def random_element(K): sage: random_element(K) (0, 0, 0) + A random element of the nonnegative orthant should have all + components nonnegative:: + + sage: set_random_seed() + sage: K = Cone([(1,0,0),(0,1,0),(0,0,1)]) + sage: all([ x >= 0 for x in random_element(K) ]) + True + TESTS: - Any cone should contain an element of itself:: + Any cone should contain a random element of itself:: sage: set_random_seed() - sage: K = random_cone(max_rays = 8) + sage: K = random_cone(max_ambient_dim=8) sage: K.contains(random_element(K)) True + A strictly convex cone contains no lines, and thus no negative + multiples of any of its elements besides zero:: + + sage: set_random_seed() + sage: K = random_cone(max_ambient_dim=8, strictly_convex=True) + sage: x = random_element(K) + sage: x.is_zero() or not K.contains(-x) + True + + The sum of random elements of a cone lies in the cone:: + + sage: set_random_seed() + sage: K = random_cone(max_ambient_dim=8) + sage: K.contains(sum([random_element(K) for i in range(10)])) + True + """ V = K.lattice().vector_space() - F = V.base_ring() - coefficients = [ F.random_element().abs() for i in range(K.nrays()) ] - vector_gens = map(V, K.rays()) - scaled_gens = [ coefficients[i]*vector_gens[i] - for i in range(len(vector_gens)) ] + scaled_gens = [ V.base_field().random_element().abs()*V(r) for r in K ] # Make sure we return a vector. Without the coercion, we might # return ``0`` when ``K`` has no rays. - v = V(sum(scaled_gens)) - return v + return V(sum(scaled_gens)) + + +def pointed_decomposition(K): + """ + Every convex cone is the direct sum of a pointed cone and a linear + subspace. Return a pair ``(P,S)`` of cones such that ``P`` is + pointed, ``S`` is a subspace, and ``K`` is the direct sum of ``P`` + and ``S``. + + OUTPUT: + + An ordered pair ``(P,S)`` of closed convex polyhedral cones where + ``P`` is pointed, ``S`` is a subspace, and ``K`` is the direct sum + of ``P`` and ``S``. + + TESTS: + + A random point in the cone should belong to either the pointed + subcone ``P`` or the subspace ``S``. If the point is nonzero, it + should lie in one but not both of them:: + + sage: set_random_seed() + sage: K = random_cone(max_ambient_dim=8) + sage: (P,S) = pointed_decomposition(K) + sage: x = random_element(K) + sage: P.contains(x) or S.contains(x) + True + sage: x.is_zero() or (P.contains(x) != S.contains(x)) + True + """ + linspace_gens = [ copy(b) for b in K.linear_subspace().basis() ] + linspace_gens += [ -b for b in linspace_gens ] + + S = Cone(linspace_gens, K.lattice()) + # Since ``S`` is a subspace, its dual is its orthogonal complement + # (albeit in the wrong lattice). + S_perp = Cone(S.dual(), K.lattice()) + P = K.intersection(S_perp) + + return (P,S) def positive_operator_gens(K): r""" @@ -177,7 +236,7 @@ def positive_operator_gens(K): A positive operator on a cone should send its generators into the cone:: sage: set_random_seed() - sage: K = random_cone(max_ambient_dim = 5) + sage: K = random_cone(max_ambient_dim=5) sage: pi_of_K = positive_operator_gens(K) sage: all([K.contains(p*x) for p in pi_of_K for x in K.rays()]) True @@ -186,7 +245,7 @@ def positive_operator_gens(K): corollary in my paper:: sage: set_random_seed() - sage: K = random_cone(max_ambient_dim = 5) + sage: K = random_cone(max_ambient_dim=5) sage: n = K.lattice_dim() sage: m = K.dim() sage: l = K.lineality() @@ -197,6 +256,29 @@ def positive_operator_gens(K): sage: actual == expected True + The lineality of the cone of positive operators is given by the + corollary in my paper:: + + sage: set_random_seed() + sage: K = random_cone(max_ambient_dim=5) + sage: n = K.lattice_dim() + sage: pi_of_K = positive_operator_gens(K) + sage: L = ToricLattice(n**2) + sage: actual = Cone([p.list() for p in pi_of_K], lattice=L).lineality() + sage: expected = n**2 - K.dim()*K.dual().dim() + sage: actual == expected + True + + The cone ``K`` is proper if and only if the cone of positive + operators on ``K`` is proper:: + + sage: set_random_seed() + sage: K = random_cone(max_ambient_dim=5) + sage: pi_of_K = positive_operator_gens(K) + sage: L = ToricLattice(K.lattice_dim()**2) + sage: pi_cone = Cone([p.list() for p in pi_of_K], lattice=L) + sage: K.is_proper() == pi_cone.is_proper() + True """ # Matrices are not vectors in Sage, so we have to convert them # to vectors explicitly before we can find a basis. We need these @@ -273,7 +355,7 @@ def Z_transformation_gens(K): The Z-property is possessed by every Z-transformation:: sage: set_random_seed() - sage: K = random_cone(max_ambient_dim = 6) + sage: K = random_cone(max_ambient_dim=6) sage: Z_of_K = Z_transformation_gens(K) sage: dcs = K.discrete_complementarity_set() sage: all([(z*x).inner_product(s) <= 0 for z in Z_of_K @@ -283,12 +365,33 @@ def Z_transformation_gens(K): The lineality space of Z is LL:: sage: set_random_seed() - sage: K = random_cone(min_ambient_dim = 1, max_ambient_dim = 6) + sage: K = random_cone(min_ambient_dim=1, max_ambient_dim=6) sage: lls = span([ vector(l.list()) for l in K.lyapunov_like_basis() ]) sage: z_cone = Cone([ z.list() for z in Z_transformation_gens(K) ]) sage: z_cone.linear_subspace() == lls True + And thus, the lineality of Z is the Lyapunov rank:: + + sage: set_random_seed() + sage: K = random_cone(max_ambient_dim=6) + sage: Z_of_K = Z_transformation_gens(K) + sage: L = ToricLattice(K.lattice_dim()**2) + sage: z_cone = Cone([ z.list() for z in Z_of_K ], lattice=L) + sage: z_cone.lineality() == K.lyapunov_rank() + True + + The lineality spaces of pi-star and Z-star are equal: + + sage: set_random_seed() + sage: K = random_cone(max_ambient_dim=5) + sage: pi_of_K = positive_operator_gens(K) + sage: Z_of_K = Z_transformation_gens(K) + sage: L = ToricLattice(K.lattice_dim()**2) + sage: pi_star = Cone([p.list() for p in pi_of_K], lattice=L).dual() + sage: z_star = Cone([ z.list() for z in Z_of_K], lattice=L).dual() + sage: pi_star.linear_subspace() == z_star.linear_subspace() + True """ # Matrices are not vectors in Sage, so we have to convert them # to vectors explicitly before we can find a basis. We need these @@ -317,3 +420,18 @@ def Z_transformation_gens(K): # not cross-positive ones. M = MatrixSpace(F, n) return [ -M(v.list()) for v in Sigma_cone.rays() ] + + +def Z_cone(K): + gens = Z_transformation_gens(K) + L = None + if len(gens) == 0: + L = ToricLattice(0) + return Cone([ g.list() for g in gens ], lattice=L) + +def pi_cone(K): + gens = positive_operator_gens(K) + L = None + if len(gens) == 0: + L = ToricLattice(0) + return Cone([ g.list() for g in gens ], lattice=L)