Merge #22081

Author: Matthias Koeppe

src/doc/en/reference/references/index.rst

@@ -1384,6 +1384,10 @@ REFERENCES:
 .. [CFI1992] Cai, JY., Fürer, M. & Immerman, N. Combinatorica (1992) 12: 389.
              :doi:`10.1007/BF01305232`
 
+.. [CFKP1997] James W. Cannon, William J. Floyd, Richard Kenyon and Walter R. Parry.
+              *Hyperbolic Geometry*. Flavors of Geometry, MSRI Publications,
+              Volume 31, 1997.
+
 .. [CFKPR2010] Ioannis Caragiannis, Afonso Ferreira, Christos Kaklamanis,
                Stéphane Pérennes, Hervé Rivano.
                *Fractional Path Coloring in Bounded Degree Trees with
@@ -5308,6 +5312,9 @@ REFERENCES:
              Variable I. Lecture Notes in Mathematics,
              vol 320. (1973), Springer, Berlin, Heidelberg
 
+.. [Sta1993] Stahl, Saul: *The Poincaré Half-plane: A Gateway to Modern Geometry*
+             Jones & Bartlett Learning, 1993. ISBN 086720298X, 9780867202984
+
 .. [Sta1995] \R. P. Stanley, *A symmetric function generalization
              of the chromatic polynomial of a graph*, Adv. Math., ***111***
              no.1 (1995), 166-194. :doi:`10.1006/aima.1995.1020`.

src/sage/geometry/hyperbolic_space/hyperbolic_geodesic.py

@@ -91,6 +91,8 @@
 from sage.misc.lazy_import import lazy_import
 lazy_import("sage.plot.arc", "arc")
 lazy_import("sage.plot.line", "line")
+lazy_import("sage.plot.arc", "arc")
+lazy_import("sage.plot.bezier_path", "bezier_path")
 from sage.symbolic.constants import pi
 from sage.modules.free_module_element import vector
 from sage.matrix.constructor import matrix
@@ -696,15 +698,27 @@ def complete(self):
         ::
 
             sage: KM = H.KM()
-            sage: KM.get_geodesic((0,0), (0, 1/2)).complete()
+            sage: KM.get_geodesic((0,0), (0,1/2)).complete()
             Geodesic in KM from (0, -1) to (0, 1)
 
         .. PLOT::
 
-            g = HyperbolicPlane().KM().get_geodesic((0.0,0.0), (0.0, 0.5))
+            g = HyperbolicPlane().KM().get_geodesic(CC(0,0), CC(0, 0.5))
+            h = g.complete()
+            sphinx_plot(g.plot()+h.plot(linestyle='dashed'))
+
+        ::
+
+            sage: KM.get_geodesic(-I, 1).complete()
+            Geodesic in KM from -I to 1
+
+        .. PLOT::
+
+            g = HyperbolicPlane().KM().get_geodesic(CC(0,-1), CC(1, 0))
             h = g.complete()
             sphinx_plot(g.plot()+h.plot(linestyle='dashed'))
 
+
         ::
 
             sage: HM = H.HM()
@@ -1186,12 +1200,15 @@ def plot(self, boundary=True, **options):
             elif end_2 == CC(infinity):
                 end_2 = (real(end_1), (imag(end_1) + 10))
                 end_1 = (real(end_1), imag(end_1))
-            pic = line((end_1, end_2), **opts)
+            else:
+                end_1 = (real(end_1), imag(end_1))
+                end_2 = (real(end_2), imag(end_2))
+            pic = bezier_path([[end_1, end_2]], **opts)
             if boundary:
                 cent = min(bd_1, bd_2)
                 bd_dict = {'bd_min': cent - 3, 'bd_max': cent + 3}
                 bd_pic = self._model.get_background_graphic(**bd_dict)
-                pic = bd_pic + pic
+                pic += bd_pic
             return pic
         else:
             center = (bd_1 + bd_2) / 2  # Circle center
@@ -1213,7 +1230,7 @@ def plot(self, boundary=True, **options):
                 bd_dict = {'bd_min': midpoint - length, 'bd_max': midpoint +
                            length}
                 bd_pic = self._model.get_background_graphic(**bd_dict)
-                pic = bd_pic + pic
+                pic += bd_pic
             return pic
 
     def ideal_endpoints(self):
@@ -1244,10 +1261,10 @@ def ideal_endpoints(self):
         M = self._model
         # infinity is the first endpoint, so the other ideal endpoint
         # is just the real part of the second coordinate
-        if start == infinity:
+        if CC(start).is_infinity():
             return [M.get_point(start), M.get_point(x2)]
         # Same idea as above
-        if end == infinity:
+        if CC(end).is_infinity():
             return [M.get_point(x1), M.get_point(end)]
         # We could also have a vertical line with two interior points
         if x1 == x2:
@@ -2174,11 +2191,9 @@ class HyperbolicGeodesicPD(HyperbolicGeodesic):
          g = PD.get_geodesic(I,-I/2)
          h = PD.get_geodesic(-0.5+I*0.5,0.5+I*0.5)
          sphinx_plot(g.plot()+h.plot(color='green'))
-
     """
 
     def plot(self, boundary=True, **options):
-
         r"""
         Plot ``self``.
 
@@ -2235,7 +2250,7 @@ def plot(self, boundary=True, **options):
         bd_1, bd_2 = [CC(k.coordinates()) for k in self.ideal_endpoints()]
         # Check to see if it's a line
         if abs(bd_1 + bd_2) < EPSILON:
-            pic = line([end_1, end_2], **opts)
+            pic = bezier_path([[(real(end_1), imag(end_1)), (real(end_2), imag(end_2))]], **opts)
         else:
             # If we are here, we know it's not a line
             # So we compute the center and radius of the circle
@@ -2276,20 +2291,18 @@ class HyperbolicGeodesicKM(HyperbolicGeodesic):
     EXAMPLES::
 
         sage: KM = HyperbolicPlane().KM()
-        sage: g = KM.get_geodesic(KM.get_point((0.1,0.9)), KM.get_point((-0.1,-0.9)))
         sage: g = KM.get_geodesic((0.1,0.9),(-0.1,-0.9))
         sage: h = KM.get_geodesic((-0.707106781,-0.707106781),(0.707106781,-0.707106781))
         sage: P = g.plot(color='orange')+h.plot(); P  # optional - sage.plot
         Graphics object consisting of 4 graphics primitives
 
-
     .. PLOT::
 
         KM = HyperbolicPlane().KM()
-        g = KM.get_geodesic((0.1,0.9),
-                            (-0.1,-0.9))
-        h = KM.get_geodesic((-0.707106781,-0.707106781),
-                            (0.707106781,-0.707106781))
+        g = KM.get_geodesic(CC(0.1,0.9),
+                            CC(-0.1,-0.9))
+        h = KM.get_geodesic(CC(-0.707106781,-0.707106781),
+                            CC(0.707106781,-0.707106781))
         sphinx_plot(g.plot(color='orange')+h.plot())
 
     """
@@ -2300,19 +2313,24 @@ def plot(self, boundary=True, **options):
 
         EXAMPLES::
 
-            sage: HyperbolicPlane().KM().get_geodesic((0,0), (1,0)).plot()  # optional - sage.plot
+            sage: HyperbolicPlane().KM().get_geodesic(0, 1).plot()  # optional - sage.plot
             Graphics object consisting of 2 graphics primitives
 
         .. PLOT::
 
             KM = HyperbolicPlane().KM()
-            sphinx_plot(KM.get_geodesic((0,0), (1,0)).plot())
+            sphinx_plot(KM.get_geodesic(CC(0,0), CC(1,0)).plot())
 
         """
         opts = {'axes': False, 'aspect_ratio': 1}
         opts.update(self.graphics_options())
         opts.update(options)
-        pic = line([k.coordinates() for k in self.endpoints()], **opts)
+        def map_pt(pt):
+            if pt in CC:
+                return CC(pt)
+            return CC(*pt)
+        end_1, end_2 = [map_pt(k.coordinates()) for k in self.endpoints()]
+        pic = bezier_path([[(real(end_1), imag(end_1)), (real(end_2), imag(end_2))]], **opts)
         if boundary:
             pic += self._model.get_background_graphic()
         return pic
@@ -2334,7 +2352,6 @@ class HyperbolicGeodesicHM(HyperbolicGeodesic):
 
     EXAMPLES::
 
-        sage: from sage.geometry.hyperbolic_space.hyperbolic_geodesic import *
         sage: HM = HyperbolicPlane().HM()
         sage: p1 = HM.get_point((4, -4, sqrt(33)))
         sage: p2 = HM.get_point((-3,-3,sqrt(19)))
@@ -2350,6 +2367,51 @@ class HyperbolicGeodesicHM(HyperbolicGeodesic):
         sphinx_plot(g.plot(color='blue'))
 
     """
+    def _plot_vertices(self, points=75):
+        r"""
+        Return ``self`` plotting vertices in `\RR^3`.
+
+        Auxiliary function needed to plot polygons.
+
+        EXAMPLES::
+
+            sage: HM = HyperbolicPlane().HM()
+            sage: p1 = HM.get_point((4, -4, sqrt(33)))
+            sage: p2 = HM.get_point((-3,-3,sqrt(19)))
+            sage: g = HM.get_geodesic(p1, p2)
+            sage: g._plot_vertices(5)
+            [(4.0, -4.0, 5.744562646538029),
+             (1.3632131724692114, -1.6370738298435326, 2.353372235315133),
+             (0.13856858387448234, -0.9699800871213693, 1.4000223647674197),
+             (-0.9425338542843988, -1.3076813974501533, 1.8969450977056184),
+             (-3.0, -3.0, 4.358898943540652)]
+        """
+        from sage.plot.misc import setup_for_eval_on_grid
+        from sage.arith.srange import xsrange
+
+        x = SR.var('x')
+        v1, u2 = [vector(k.coordinates()) for k in self.endpoints()]
+        # Lorentzian Gram Shmidt.  The original vectors will be
+        # u1, u2 and the orthogonal ones will be v1, v2.  Except
+        # v1 = u1, and I don't want to declare another variable,
+        # hence the odd naming convention above.
+        # We need the Lorentz dot product of v1 and u2.
+        v1_ldot_u2 = u2[0]*v1[0] + u2[1]*v1[1] - u2[2]*v1[2]
+        v2 = u2 + v1_ldot_u2 * v1
+        v2_norm = sqrt(v2[0]**2 + v2[1]**2 - v2[2]**2)
+        v2 = v2 / v2_norm
+        v2_ldot_u2 = u2[0]*v2[0] + u2[1]*v2[1] - u2[2]*v2[2]
+        # Now v1 and v2 are Lorentz orthogonal, and |v1| = -1, |v2|=1
+        # That is, v1 is unit timelike and v2 is unit spacelike.
+        # This means that cosh(x)*v1 + sinh(x)*v2 is unit timelike.
+        hyperbola = tuple(cosh(x)*v1 + sinh(x)*v2)
+        endtime = arcsinh(v2_ldot_u2)
+        # mimic the function _parametric_plot3d_curve using a bezier3d instead of a line3d
+        # this is required in order to be able to plot hyperbolic polygons whithin the plot library
+        g, ranges = setup_for_eval_on_grid(hyperbola, [(x, 0, endtime)], points)
+        f_x, f_y, f_z = g
+        points = [(f_x(u), f_y(u), f_z(u)) for u in xsrange(*ranges[0], include_endpoint=True)]
+        return points
 
     def plot(self, show_hyperboloid=True, **graphics_options):
         r"""