Bases: sage.structure.sage_object.SageObject
The projection of a Polyhedron.
This class keeps track of the necessary data to plot the input polyhedron.
Convert a coordinate vector to its internal index.
EXAMPLES:
sage: p = polytopes.n_cube(3)
sage: proj = p.projection()
sage: proj.coord_index_of(vector((1,1,1)))
7
Convert list of coordinate vectors to the corresponding list of internal indices.
EXAMPLES:
sage: p = polytopes.n_cube(3)
sage: proj = p.projection()
sage: proj.coord_indices_of([vector((1,1,1)),vector((1,-1,1))])
[7, 5]
Given a list of indices, return the projected coordinates.
EXAMPLES:
sage: p = polytopes.n_simplex(4).projection()
sage: p.coordinates_of([1])
[[0, -81649/100000, 7217/25000, 22361/100000]]
Return the identity projection of the polyhedron.
EXAMPLES:
sage: p = polytopes.icosahedron()
sage: from sage.geometry.polyhedron.plot import Projection
sage: pproj = Projection(p)
sage: ppid = pproj.identity()
sage: ppid.dimension
3
Return the filled interior (a polygon) of a polyhedron in 2d.
EXAMPLES:
sage: cps = [i^3 for i in srange(-2,2,1/5)]
sage: p = Polyhedron(vertices = [[(t^2-1)/(t^2+1),2*t/(t^2+1)] for t in cps])
sage: proj = p.projection()
sage: filled_poly = proj.render_fill_2d()
sage: filled_poly.axes_width()
0.8
Return the outline (edges) of a polyhedron in 2d.
EXAMPLES:
sage: penta = polytopes.regular_polygon(5)
sage: outline = penta.projection().render_outline_2d()
sage: outline._objects[0]
Line defined by 2 points
Return the points of a polyhedron in 2d.
EXAMPLES:
sage: hex = polytopes.regular_polygon(6)
sage: proj = hex.projection()
sage: hex_points = proj.render_points_2d()
sage: hex_points._objects
[Point set defined by 6 point(s)]
Return solid 3d rendering of a 3d polytope.
EXAMPLES:
sage: p = polytopes.n_cube(3).projection()
sage: p_solid = p.render_solid_3d(opacity = .7)
sage: type(p_solid)
<class 'sage.plot.plot3d.base.Graphics3dGroup'>
Return the 3d rendering of the vertices.
EXAMPLES:
sage: p = polytopes.cross_polytope(3)
sage: proj = p.projection()
sage: verts = proj.render_vertices_3d()
sage: verts.bounding_box()
((-1.0, -1.0, -1.0), (1.0, 1.0, 1.0))
Return the 3d wireframe rendering.
EXAMPLES:
sage: cube = polytopes.n_cube(3)
sage: cube_proj = cube.projection()
sage: wire = cube_proj.render_wireframe_3d()
sage: print wire.tachyon().split('\n')[77] # for testing
FCylinder base -1.0 1.0 -1.0 apex -1.0 -1.0 -1.0 rad 0.005 texture...
Return the Schlegel projection.
The vertices are normalized to the unit sphere, and stereographically projected from a point slightly outside of the sphere.
INPUT:
- projection_direction - The direction of the Schlegel projection. By default, the vector consisting of the first n primes is chosen.
EXAMPLES:
sage: cube4 = polytopes.n_cube(4)
sage: from sage.geometry.polyhedron.plot import Projection
sage: Projection(cube4).schlegel([1,0,0,0])
The projection of a polyhedron into 3 dimensions
sage: _.show()
TESTS:
sage: Projection(cube4).schlegel()
The projection of a polyhedron into 3 dimensions
Return the stereographic projection.
INPUT:
EXAMPLES:
sage: from sage.geometry.polyhedron.plot import Projection
sage: proj = Projection(polytopes.buckyball()) #long time
sage: proj #long time
The projection of a polyhedron into 3 dimensions
sage: proj.stereographic([5,2,3]).show() #long time
sage: Projection( polytopes.twenty_four_cell() ).stereographic([2,0,0,0])
The projection of a polyhedron into 3 dimensions
The Schlegel projection from the given input point.
EXAMPLES:
sage: from sage.geometry.polyhedron.plot import ProjectionFuncSchlegel
sage: proj = ProjectionFuncSchlegel([2,2,2])
sage: proj(vector([1.1,1.1,1.11]))[0]
0.0302...
The stereographic (or perspective) projection.
EXAMPLES:
sage: from sage.geometry.polyhedron.plot import ProjectionFuncStereographic
sage: cube = polytopes.n_cube(3).vertices()
sage: proj = ProjectionFuncStereographic([1.2, 3.4, 5.6])
sage: ppoints = [proj(vector(x)) for x in cube]
sage: ppoints[0]
(-0.0486511..., 0.0859565...)
Return the vertices/rays in cyclic order if possible.
NOTES:
This works if and only if each vertex/ray is adjacent to exactly two others. For example, any 2-dimensional polyhedron satisfies this.
See vertex_adjacency_matrix() for a discussion of “adjacent”.
EXAMPLES:
sage: from sage.geometry.polyhedron.plot import cyclic_sort_vertices_2d
sage: square = Polyhedron([[1,0],[-1,0],[0,1],[0,-1]])
sage: vertices = [v for v in square.vertex_generator()]
sage: vertices
[A vertex at (-1, 0),
A vertex at (0, -1),
A vertex at (0, 1),
A vertex at (1, 0)]
sage: cyclic_sort_vertices_2d(vertices)
[A vertex at (1, 0),
A vertex at (0, -1),
A vertex at (-1, 0),
A vertex at (0, 1)]
Rays are allowed, too:
sage: P = Polyhedron(vertices=[(0, 1), (1, 0), (2, 0), (3, 0), (4, 1)], rays=[(0,1)])
sage: P.adjacency_matrix()
[0 1 0 1 0]
[1 0 1 0 0]
[0 1 0 0 1]
[1 0 0 0 1]
[0 0 1 1 0]
sage: cyclic_sort_vertices_2d(P.Vrepresentation())
[A vertex at (3, 0),
A vertex at (1, 0),
A vertex at (0, 1),
A ray in the direction (0, 1),
A vertex at (4, 1)]
sage: P = Polyhedron(vertices=[(0, 1), (1, 0), (2, 0), (3, 0), (4, 1)], rays=[(0,1), (1,1)])
sage: P.adjacency_matrix()
[0 1 0 0 0]
[1 0 1 0 0]
[0 1 0 0 1]
[0 0 0 0 1]
[0 0 1 1 0]
sage: cyclic_sort_vertices_2d(P.Vrepresentation())
[A ray in the direction (1, 1),
A vertex at (3, 0),
A vertex at (1, 0),
A vertex at (0, 1),
A ray in the direction (0, 1)]
sage: P = Polyhedron(vertices=[(1,2)], rays=[(0,1)], lines=[(1,0)])
sage: P.adjacency_matrix()
[0 0 1]
[0 0 0]
[1 0 0]
sage: cyclic_sort_vertices_2d(P.Vrepresentation())
[A vertex at (0, 2),
A line in the direction (1, 0),
A ray in the direction (0, 1)]
The identity projection.
EXAMPLES:
sage: from sage.geometry.polyhedron.plot import projection_func_identity
sage: projection_func_identity((1,2,3))
[1, 2, 3]
Return 2d rendering of the projection of a polyhedron into 2-dimensional ambient space.
EXAMPLES:
sage: p1 = Polyhedron(vertices=[[1,1]], rays=[[1,1]])
sage: q1 = p1.projection()
sage: p2 = Polyhedron(vertices=[[1,0], [0,1], [0,0]])
sage: q2 = p2.projection()
sage: p3 = Polyhedron(vertices=[[1,2]])
sage: q3 = p3.projection()
sage: p4 = Polyhedron(vertices=[[2,0]], rays=[[1,-1]], lines=[[1,1]])
sage: q4 = p4.projection()
sage: q1.show() + q2.show() + q3.show() + q4.show()
sage: from sage.geometry.polyhedron.plot import render_2d
sage: q = render_2d(p1.projection())
sage: q._objects
[Point set defined by 1 point(s),
Arrow from (1.0,1.0) to (2.0,2.0),
Polygon defined by 3 points]
Return 3d rendering of a polyhedron projected into 3-dimensional ambient space.
NOTE:
This method, render_3d, is used in the show() method of a polyhedron if it is in 3 dimensions.
EXAMPLES:
sage: p1 = Polyhedron(vertices=[[1,1,1]], rays=[[1,1,1]])
sage: p2 = Polyhedron(vertices=[[2,0,0], [0,2,0], [0,0,2]])
sage: p3 = Polyhedron(vertices=[[1,0,0], [0,1,0], [0,0,1]], rays=[[-1,-1,-1]])
sage: p1.projection().show() + p2.projection().show() + p3.projection().show()
It correctly handles various degenerate cases:
sage: Polyhedron(lines=[[1,0,0],[0,1,0],[0,0,1]]).show() # whole space
sage: Polyhedron(vertices=[[1,1,1]], rays=[[1,0,0]], lines=[[0,1,0],[0,0,1]]).show() # half space
sage: Polyhedron(vertices=[[1,1,1]], lines=[[0,1,0],[0,0,1]]).show() # R^2 in R^3
sage: Polyhedron(rays=[[0,1,0],[0,0,1]], lines=[[1,0,0]]).show() # quadrant wedge in R^2
sage: Polyhedron(rays=[[0,1,0]], lines=[[1,0,0]]).show() # upper half plane in R^3
sage: Polyhedron(lines=[[1,0,0]]).show() # R^1 in R^2
sage: Polyhedron(rays=[[0,1,0]]).show() # Half-line in R^3
sage: Polyhedron(vertices=[[1,1,1]]).show() # point in R^3
Return a 3d rendering of the Schlegel projection of a 4d polyhedron projected into 3-dimensional space.
NOTES:
The show() method of Polyhedron() uses this to draw itself if the ambient dimension is 4.
INPUT:
EXAMPLES:
sage: poly = polytopes.twenty_four_cell()
sage: poly
A 4-dimensional polyhedron in QQ^4 defined as the convex hull of 24 vertices
sage: poly.show()
sage: poly.show(projection_direction=[2,5,11,17])
sage: type( poly.show() )
<class 'sage.plot.plot3d.base.Graphics3dGroup'>
TESTS:
sage: from sage.geometry.polyhedron.plot import render_4d
sage: p = polytopes.n_cube(4)
sage: q = render_4d(p)
sage: tach_str = q.tachyon()
sage: tach_str.count('FCylinder')
32