Ticket #11273: trac_11273_multi_web_fix.patch

File trac_11273_multi_web_fix.patch, 6.2 KB (added by evanandel, 8 years ago)

A fix for exterior fuzz in some multiply connected spiderwebs.

  • sage/calculus/riemann.pyx

    # HG changeset patch
    # User Ethan Van Andel <evlutte@gmail.com>
    # Date 1380754870 25200
    # Node ID 6fbf9dd16ed9123e03a8cae3cfa65d822a72601a
    # Parent  2f249b4f15186cceca466320167acaa48d74eda2
    Trac 11273: Added min_mag arg to plot_spiderweb. Added 2 new doctests.
    
    diff --git a/sage/calculus/riemann.pyx b/sage/calculus/riemann.pyx
    a b  
    170170        sage: print "error =", m.inverse_riemann_map(m.riemann_map(x)) - x  # long time
    171171        error = (-0.000...+0.0016...j)
    172172
     173    A fun, complex region for demonstration purposes
     174        sage: f(t) = e^(I*t)
     175        sage: fp(t) = I*e^(I*t)
     176        sage: ef1(t) = .2*e^(-I*t) +.4+.4*I
     177        sage: ef1p(t) = -I*.2*e^(-I*t)
     178        sage: ef2(t) = .2*e^(-I*t) -.4+.4*I
     179        sage: ef2p(t) = -I*.2*e^(-I*t)
     180        sage: pts = [(-.5,-.15-20/1000),(-.6,-.27-10/1000),(-.45,-.45),(0,-.65+10/1000),(.45,-.45),(.6,-.27-10/1000),(.5,-.15-10/1000),(0,-.43+10/1000)]
     181        sage: pts.reverse()
     182        sage: cs = complex_cubic_spline(pts)
     183        sage: mf = lambda x:cs.value(x)
     184        sage: mfprime = lambda x: cs.derivative(x)
     185        sage: m = Riemann_Map([f,ef1,ef2,mf],[fp,ef1p,ef2p,mfprime],0,N = 400) # long time
     186        sage: p = m.plot_colored(plot_points = 400) # long time
     187
    173188    ALGORITHM:
    174189
    175190    This class computes the Riemann Map via the Szego kernel using an
     
    820835    @options(interpolation='catrom')
    821836    def plot_spiderweb(self, spokes=16, circles=4, pts=32, linescale=0.99,
    822837            rgbcolor=[0,0,0], thickness=1, plotjoined=True, withcolor = False,
    823             plot_points = 200, **options):
     838            plot_points = 200, min_mag = 0.001, **options):
    824839        """
    825840        Generates a traditional "spiderweb plot" of the Riemann map. Shows
    826841        what concentric circles and radial lines map to. The radial lines
     
    875890          Note that very large values can cause this function to run slowly.
    876891          - only for multiply connected domains
    877892
     893        - ``min_mag`` -- float (default: ``0.001``) The magnitude cutoff
     894          below which spiderweb points are not drawn. This only applies
     895          to multiply connected domains and is designed to prevent
     896          "fuzz" at the edge of the domain. Some complicated multiply
     897          connected domains (particularly those with corners) may
     898          require a larger value to look clean outside.
     899
    878900        EXAMPLES:
    879901
    880902        General usage::
     
    902924            sage: fprime(t) = I*e^(I*t)
    903925            sage: m = Riemann_Map([f], [fprime], 0, 1000)
    904926            sage: m.plot_spiderweb()
     927
     928        A multiply connected region with corners. We set min_mag higher
     929        to remove "fuzz" outside the domain.
     930
     931            sage: ps = polygon_spline([(-4,-2),(4,-2),(4,2),(-4,2)])
     932            sage: z1 = lambda t: ps.value(t); z1p = lambda t: ps.derivative(t)
     933            sage: z2(t) = -2+exp(-I*t); z2p(t) = -I*exp(-I*t)
     934            sage: z3(t) = 2+exp(-I*t); z3p(t) = -I*exp(-I*t)
     935            sage: m = Riemann_Map([z1,z2,z3],[z1p,z2p,z3p],0,ncorners=4) # long time
     936            sage: p = m.plot_spiderweb(withcolor=True,plot_points=500, thickness = 2.0, min_mag=0.1) # long time
    905937        """
    906938        cdef int k, i
    907939        if self.exterior:
     
    960992
    961993            g = Graphics()
    962994            g.add_primitive(ComplexPlot(complex_to_spiderweb(z_values,dr,dtheta,
    963                 spokes, circles, rgbcolor,thickness, withcolor),
     995                spokes, circles, rgbcolor,thickness, withcolor, min_mag),
    964996                (xmin, xmax), (ymin, ymax),options))
    965997            return g + self.plot_boundaries(thickness = thickness)
    966998
     
    11041136
    11051137cpdef complex_to_spiderweb(np.ndarray[COMPLEX_T, ndim = 2] z_values,
    11061138    np.ndarray[FLOAT_T, ndim = 2] dr, np.ndarray[FLOAT_T, ndim = 2] dtheta,
    1107     spokes, circles, rgbcolor, thickness, withcolor):
     1139    spokes, circles, rgbcolor, thickness, withcolor, min_mag):
    11081140    """
    11091141    Converts a grid of complex numbers into a matrix containing rgb data
    11101142    for the Riemann spiderweb plot.
     
    11331165    - ``withcolor`` -- boolean - If ``True`` the spiderweb will be overlaid
    11341166      on the basic color plot.
    11351167
     1168    - ``min_mag`` -- float - The magnitude cutoff below which spiderweb
     1169      points are not drawn. This only applies to multiply connected
     1170      domains and is designed to prevent "fuzz" at the edge of the
     1171      domain. Some complicated multiply connected domains (particularly
     1172      those with corners) may require a larger value to look clean
     1173      outside.
     1174
    11361175    OUTPUT:
    11371176
    11381177    An `N x M x 3` floating point Numpy array ``X``, where
     
    11441183        sage: import numpy
    11451184        sage: zval = numpy.array([[0, 1, 1000],[.2+.3j,1,-.3j],[0,0,0]],dtype = numpy.complex128)
    11461185        sage: deriv = numpy.array([[.1]],dtype = numpy.float64)
    1147         sage: complex_to_spiderweb(zval, deriv,deriv, 4,4,[0,0,0],1,False)
     1186        sage: complex_to_spiderweb(zval, deriv,deriv, 4,4,[0,0,0],1,False,0.001)
    11481187        array([[[ 1.,  1.,  1.],
    11491188                [ 1.,  1.,  1.],
    11501189                [ 1.,  1.,  1.]],
     
    11571196                [ 1.,  1.,  1.],
    11581197                [ 1.,  1.,  1.]]])
    11591198
    1160         sage: complex_to_spiderweb(zval, deriv,deriv, 4,4,[0,0,0],1,True)
     1199        sage: complex_to_spiderweb(zval, deriv,deriv, 4,4,[0,0,0],1,True,0.001)
    11611200        array([[[ 1.        ,  1.        ,  1.        ],
    11621201                [ 1.        ,  0.05558355,  0.05558355],
    11631202                [ 0.17301243,  0.        ,  0.        ]],
     
    11741213    cdef FLOAT_T x, y, mag, arg, width, target, precision, dmag, darg
    11751214    cdef COMPLEX_T z
    11761215    cdef FLOAT_T DMAX = 70 # change to adjust rate_of_change cutoff below
    1177     cdef FLOAT_T MMIN = .0001
    11781216    precision = thickness/150.0
    11791217    imax = len(z_values)
    11801218    jmax = len(z_values[0])
     
    12021240            darg = dtheta[i,j]
    12031241            #points that change too rapidly are presumed to be borders
    12041242            #points that are too small are presumed to be outside
    1205             if darg < DMAX and mag > MMIN:
     1243            if darg < DMAX and mag > min_mag:
    12061244                for target in circ_radii:
    12071245                    if abs(mag - target)/dmag < precision:
    12081246                        rgb[i+1,j+1] = rgbcolor