Ticket #3813: trac_3813-final.patch

a/sage/plot/plot.py

@@ -1170 +1170 @@
-        z-coordinate. 
-
-        EXAMPLES: 
-
+ #long
-        """
-        from sage.plot.plot3d.base import Graphics3dGroup
-        g = Graphics3dGroup([g.plot3d(**kwds) for g in self.__objects])
@@ -1285 +1285 @@
-        Add grid lines at specific positions (using iterators).
-            sage: def maple_leaf(t):
-            ...     return (100/(100+(t-pi/2)^8))*(2-sin(7*t)-cos(30*t)/2)
-
-
+ #long
+ #long
-
-        Add grid lines at specific positions (using functions).
-            sage: y = x^5 + 4*x^4 - 10*x^3 - 40*x^2 + 9*x + 36
@@ -3372 +3372 @@
-
-    PLOT OPTIONS:
-    The plot options are
-
-
-
-
-
+
+
+
+
+
+
+
+
+
+
+
+
+
-        xmin -- starting x value
-        xmax -- ending x value
-        color -- an rgb-tuple (r,g,b) with each of r,g,b between 0 and 1, or
@@ -3419 +3427 @@
-        Graphics object consisting of 1 graphics primitive
-        sage: len(P)     # number of graphics primitives
-        1
-
-200
+ (random)
+6369
-        sage: P          # render
-        
-        sage: P = plot(sin, (0,10), plot_points=10); print P
-        Graphics object consisting of 1 graphics primitive
-        sage: len(P[0])  # random output
-80
+32
-        sage: P          # render
-
-    We plot with randomize=False, which makes the initial sample 
-    points evenly spaced (hence always the same).  Adaptive plotting 
-plot_divi
-plot_divi
+adaptive_recur
+adaptive_recur
-        sage: list(p[0])
-        [(0.0, 1.0), (1.0, 1.0), (2.0, 1.0), (3.0, 1.0)]
-
@@ -3446 +3454 @@
-        sage: plot([sin(n*x) for n in [1..4]], (0, pi))
-
-
-
-
-
+
+
-        sage: plot(sin(1/x), (x, -1, 1))
-
-    With the \code{plot_points} option you can increase the number
-    of sample points, to obtain a more accurate plot. 
-        sage: plot(sin(1/x), (x, -1, 1), plot_points=1000)
-
-    Note that the independent variable may be omitted if there is no
-    ambiguity:
-
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
-
-    The actual sample points are slightly randomized, so the above
-    plots may look slightly different each time you draw them.
@@ -3487 +3508 @@
-        sage: set_verbose(0)
-            
-    To plot the negative real cube root, use something like the following.
- 
+
-
-    TESTS:
-    We do not randomize the endpoints:
@@ -3533 +3554 @@
-def _plot(funcs, xrange, parametric=False,
-              polar=False, label='', randomize=True, **kwds):
-    options = {'alpha':1,'thickness':1,'rgbcolor':(0,0,1),
-plot_division':1000
-max_bend': 0.1
+adaptive_tolerance':0.01
+adaptive_recursion':5
-
-    if kwds.has_key('color') and not kwds.has_key('rgbcolor'):
-        kwds['rgbcolor'] = kwds['color']
@@ -3555 +3576 @@
-
-    plot_points = int(options.pop('plot_points'))
-    x, data = var_and_list_of_values(xrange, plot_points)
-    data = list(data)
-    xmin = data[0]
-    xmax = data[-1]
-
@@ -3564 +3584 @@
-    if isinstance(funcs, (list, tuple)) and not parametric:
-        return reduce(operator.add, (plot(f, (xmin, xmax), polar=polar, **kwds) for f in funcs))
-
-
-
-
-
+
-
-    random = current_randstate().python_random().random
-    exceptions = 0; msg=''
@@ -3576 +3593 @@
-        xi = data[i]
-        # Slightly randomize the interior sample points if
-        # randomize is true
-
-
-
-
-
-
-
-
+
+
+
-        try:
-            data[i] = (float(xi), float(f(xi)))
-            if str(data[i][1]) in ['nan', 'NaN']:
@@ -3592 +3604 @@
-                exception_indices.append(i)
-        except (ZeroDivisionError, TypeError, ValueError, OverflowError), msg:
-            sage.misc.misc.verbose("%s\nUnable to compute f(%s)"%(msg, x),1)
+
+            if i == 0:
+                for j in range(1, 99):
+                    xj = xi + delta*j/100.0
+                    try:
+                        data[i] = (float(xj), float(f(xj)))
+                        # nan != nan
+                        if data[i][1] != data[i][1]:
+                            continue
+                        break
+                    except (ZeroDivisionError, TypeError, ValueError, OverflowError), msg:
+                        pass
+                else:
+                    exceptions += 1
+                    exception_indices.append(i)
+            elif i == plot_points-1:
+                for j in range(1, 99):
+                    xj = xi - delta*j/100.0
+                    try:
+                        data[i] = (float(xj), float(f(xj)))
+                        # nan != nan
+                        if data[i][1] != data[i][1]:
+                            continue
+                        break
+                    except (ZeroDivisionError, TypeError, ValueError, OverflowError), msg:
+                        pass
+                else:
+                    exceptions += 1
+                    exception_indices.append(i)
+            else:
+                exceptions += 1
+                exception_indices.append(i)
+    
-            exceptions += 1
-            exception_indices.append(i)
-
@@ -3600 +3645 @@
-
-    # adaptive refinement
-    i, j = 0, 0
-
-
-
-
+
+
+
-    while i < len(data) - 1:
-
-
-
-
-
-
-
-
-
-
-
-
+
+
+
+
-            i += 1
+       i += 1
-
-    if (len(data) == 0 and exceptions > 0) or exceptions > 10:
-        sage.misc.misc.verbose("WARNING: When plotting, failed to evaluate function at %s points."%exceptions, level=0)
@@ -4438 +4475 @@
-                'ymin':ymin, 'ymax':ymax}
-    else:
-        return xmin, xmax, ymin, ymax
+
+def adaptive_refinement(f, p1, p2, adaptive_tolerance=0.01, adaptive_recursion=5, level=0):
+    r"""
+    The adaptive refinement algorithm for plotting a function f. See the
+    docstring for plot or PlotFactory for a description of the algorithm.
+
+    INPUT:
+        f -- a function of one variable
+        p1, p2 -- two points to refine between
+        adaptive_recursion -- (default: 5) how many levels of recursion to go 
+                              before giving up when doing adaptive refinement.
+                              Setting this to 0 disables adaptive refinement.
+        adaptive_tolerance -- (default: 0.01) how large a difference should be
+                              before the adaptive refinement code considers 
+                              it significant.  See the documentation for
+                              plot() for more information.
+    
+    OUTPUT:
+        list -- a list of points to insert between p1 and p2 to get
+                a better linear approximation between them
+
+    TESTS:
+        sage: from sage.plot.plot import adaptive_refinement
+        sage: adaptive_refinement(sin, (0,0), (pi,0), adaptive_tolerance=0.01, adaptive_recursion=0)
+        []
+        sage: adaptive_refinement(sin, (0,0), (pi,0), adaptive_tolerance=0.01)
+        [(0.125000000000000*pi, 0.38268343236508978), (0.187500000000000*pi, 0.55557023301960218), (0.250000000000000*pi, 0.707106781186547...), (0.312500000000000*pi, 0.83146961230254524), (0.375000000000000*pi, 0.92387953251128674), (0.437500000000000*pi, 0.98078528040323043), (0.500000000000000*pi, 1.0), (0.562500000000000*pi, 0.98078528040323043), (0.625000000000000*pi, 0.92387953251128674), (0.687500000000000*pi, 0.83146961230254546), (0.750000000000000*pi, 0.70710678118654757), (0.812500000000000*pi, 0.55557023301960218), (0.875000000000000*pi, 0.38268343236508989)]
+
+    This shows that lowering adaptive_tolerance and raising
+    adaptive_recursion both increase the number of subdivision points:
+
+        sage: x = var('x')
+        sage: f = sin(1/x)
+        sage: n1 = len(adaptive_refinement(f, (0,0), (pi,0), adaptive_tolerance=0.01)); n1
+        15
+        sage: n2 = len(adaptive_refinement(f, (0,0), (pi,0), adaptive_recursion=10, adaptive_tolerance=0.01)); n2
+        79
+        sage: n3 = len(adaptive_refinement(f, (0,0), (pi,0), adaptive_tolerance=0.001)); n3
+        26
+    """
+    if level >= adaptive_recursion:
+        return []
+    x = (p1[0] + p2[0])/2.0
+    try:
+        y = float(f(x))
+    except (ZeroDivisionError, TypeError, ValueError, OverflowError), msg:
+        sage.misc.misc.verbose("%s\nUnable to compute f(%s)"%(msg, x), 1)
+        # give up for this branch
+        return []
+
+    # this distance calculation is not perfect. 
+    if abs((p1[1] + p2[1])/2.0 - y) > adaptive_tolerance:
+        return adaptive_refinement(f, p1, (x, y),
+                    adaptive_tolerance=adaptive_tolerance,
+                    adaptive_recursion=adaptive_recursion,
+                    level=level+1) \
+                    + [(x, y)] + \
+            adaptive_refinement(f, (x, y), p2,
+                    adaptive_tolerance=adaptive_tolerance,
+                    adaptive_recursion=adaptive_recursion,
+                    level=level+1)
+    else:
+        return []