Riemann Mapping
AUTHORS:
Development supported by NSF award No. 0702939.
Bases: object
The Riemann_Map class computes a Riemann or Ahlfors map from supplied data. It also has various methods to provide information about the map. A Riemann map conformally maps a simply connected region in the complex plane to the unit disc. The Ahlfors map does the same thing for multiply connected regions.
Note that all the methods are numeric rather than analytic, for unusual regions or insufficient collocation points may give very inaccurate results.
INPUT:
The following inputs must all be passed in as named parameters:
EXAMPLES:
The unit circle identity map:
sage: f(t) = e^(I*t)
sage: fprime(t) = I*e^(I*t)
sage: m = Riemann_Map([f], [fprime], 0) # long time (4 sec)
The unit circle with a small hole:
sage: f(t) = e^(I*t)
sage: fprime(t) = I*e^(I*t)
sage: hf(t) = 0.5*e^(-I*t)
sage: hfprime(t) = 0.5*-I*e^(-I*t)
sage: m = Riemann_Map([f, hf], [fprime, hfprime], 0.5 + 0.5*I)
A square:
sage: ps = polygon_spline([(-1, -1), (1, -1), (1, 1), (-1, 1)]) # long time
sage: f = lambda t: ps.value(real(t)) # long time
sage: fprime = lambda t: ps.derivative(real(t)) # long time
sage: m = Riemann_Map([f], [fprime], 0.25, ncorners=4) # long time
sage: m.plot_colored() + m.plot_spiderweb() # long time
Compute rough error for this map:
sage: x = 0.75 # long time
sage: print "error =", m.inverse_riemann_map(m.riemann_map(x)) - x # long time
error = (-0.000...+0.0016...j)
ALGORITHM:
This class computes the Riemann Map via the Szego kernel using an adaptation of the method described by [KT].
REFERENCES:
[KT] | N. Kerzman and M. R. Trummer. “Numerical Conformal Mapping via the Szego kernel”. Journal of Computational and Applied Mathematics, 14(1-2): 111–123, 1986. |
Returns a discretized version of the Szego kernel for each boundary function.
INPUT:
The following inputs must all be passed in as named parameters:
OUTPUT:
A list of points of the form [t value, value of the Szego kernel at that t].
EXAMPLES:
Generic use:
sage: f(t) = e^(I*t) - 0.5*e^(-I*t)
sage: fprime(t) = I*e^(I*t) + 0.5*I*e^(-I*t)
sage: m = Riemann_Map([f], [fprime], 0)
sage: sz = m.get_szego(boundary=0)
sage: points = m.get_szego(absolute_value=True)
sage: list_plot(points)
Extending the points by a spline:
sage: s = spline(points)
sage: s(3*pi / 4)
0.00121587378429...
The unit circle with a small hole:
sage: f(t) = e^(I*t)
sage: fprime(t) = I*e^(I*t)
sage: hf(t) = 0.5*e^(-I*t)
sage: hfprime(t) = 0.5*-I*e^(-I*t)
sage: m = Riemann_Map([f, hf], [fprime, hfprime], 0.5 + 0.5*I)
Getting the szego for a specifc boundary:
sage: sz0 = m.get_szego(boundary=0)
sage: sz1 = m.get_szego(boundary=1)
Returns an array of points of the form [t value, theta in e^(I*theta)], that is, a discretized version of the theta/boundary correspondence function. For multiply connected domains, get_theta_points will list the points for each boundary in the order that they were supplied.
INPUT:
The following input must all be passed in as named parameters:
OUTPUT:
A list of points of the form [t value, theta in e^(I*theta)].
EXAMPLES:
Getting the list of points, extending it via a spline, getting the points for only the outside of a multiply connected domain:
sage: f(t) = e^(I*t) - 0.5*e^(-I*t)
sage: fprime(t) = I*e^(I*t) + 0.5*I*e^(-I*t)
sage: m = Riemann_Map([f], [fprime], 0)
sage: points = m.get_theta_points()
sage: list_plot(points)
Extending the points by a spline:
sage: s = spline(points)
sage: s(3*pi / 4)
1.62766037996...
The unit circle with a small hole:
sage: f(t) = e^(I*t)
sage: fprime(t) = I*e^(I*t)
sage: hf(t) = 0.5*e^(-I*t)
sage: hfprime(t) = 0.5*-I*e^(-I*t)
sage: m = Riemann_Map([f, hf], [hf, hfprime], 0.5 + 0.5*I)
Getting the szego for a specifc boundary:
sage: tp0 = m.get_theta_points(boundary=0)
sage: tp1 = m.get_theta_points(boundary=1)
Returns the inverse Riemann mapping of a point. That is, given pt on the interior of the unit disc, inverse_reimann_map() will return the point on the original region that would be Riemann mapped to pt.
Note
This method does not work for multiply connected domains.
INPUT:
OUTPUT:
The point on the region that Riemann maps to the input point.
EXAMPLES:
Can work for different types of complex numbers:
sage: f(t) = e^(I*t) - 0.5*e^(-I*t)
sage: fprime(t) = I*e^(I*t) + 0.5*I*e^(-I*t)
sage: m = Riemann_Map([f], [fprime], 0)
sage: m.inverse_riemann_map(0.5 + sqrt(-0.5))
(0.406880548363...+0.361470279816...j)
sage: m.inverse_riemann_map(0.95)
(0.486319431795...-4.90019052...j)
sage: m.inverse_riemann_map(0.25 - 0.3*I)
(0.165324498558...-0.180936785500...j)
sage: import numpy as np
sage: m.inverse_riemann_map(np.complex(-0.2, 0.5))
(-0.156280570579...+0.321819151891...j)
Plots the boundaries of the region for the Riemann map. Note that this method DOES work for multiply connected domains.
INPUT:
The following inputs must all be passed in as named parameters:
EXAMPLES:
General usage:
sage: f(t) = e^(I*t) - 0.5*e^(-I*t)
sage: fprime(t) = I*e^(I*t) + 0.5*I*e^(-I*t)
sage: m = Riemann_Map([f], [fprime], 0)
Default plot:
sage: m.plot_boundaries()
Big blue collocation points:
sage: m.plot_boundaries(plotjoined=False, rgbcolor=[0,0,1], thickness=6)
Draws a colored plot of the Riemann map. A red point on the colored plot corresponds to a red point on the unit disc. Note that this method DOES work for multiply connected domains.
INPUT:
The following inputs must all be passed in as named parameters:
EXAMPLES:
Given a Riemann map m, general usage:
sage: f(t) = e^(I*t) - 0.5*e^(-I*t)
sage: fprime(t) = I*e^(I*t) + 0.5*I*e^(-I*t)
sage: m = Riemann_Map([f], [fprime], 0)
sage: m.plot_colored()
Plot zoomed in on a specific spot:
sage: m.plot_colored(plot_range=[0,1,.25,.75])
High resolution plot:
sage: m.plot_colored(plot_points=1000) # long time (29s on sage.math, 2012)
To generate the unit circle map, it’s helpful to see what the colors correspond to:
sage: f(t) = e^(I*t)
sage: fprime(t) = I*e^(I*t)
sage: m = Riemann_Map([f], [fprime], 0, 1000)
sage: m.plot_colored()
Generates a traditional “spiderweb plot” of the Riemann map. Shows what concentric circles and radial lines map to. Note that this method DOES NOT work for multiply connected domains.
INPUT:
The following inputs must all be passed in as named parameters:
EXAMPLES:
General usage:
sage: f(t) = e^(I*t) - 0.5*e^(-I*t)
sage: fprime(t) = I*e^(I*t) + 0.5*I*e^(-I*t)
sage: m = Riemann_Map([f], [fprime], 0)
Default plot:
sage: m.plot_spiderweb()
Simplified plot with many discrete points:
sage: m.plot_spiderweb(spokes=4, circles=1, pts=400, linescale=0.95, plotjoined=False)
Plot with thick, red lines:
sage: m.plot_spiderweb(rgbcolor=[1,0,0], thickness=3)
To generate the unit circle map, it’s helpful to see what the original spiderweb looks like:
sage: f(t) = e^(I*t)
sage: fprime(t) = I*e^(I*t)
sage: m = Riemann_Map([f], [fprime], 0, 1000)
sage: m.plot_spiderweb()
Returns the Riemann mapping of a point. That is, given pt on the interior of the mapped region, riemann_map() will return the point on the unit disk that pt maps to. Note that this method only works for interior points; it breaks down very close to the boundary. To get boundary corrospondance, use get_theta_points().
INPUT:
OUTPUT:
A complex number representing the point on the unit circle that the input point maps to.
EXAMPLES:
Can work for different types of complex numbers:
sage: f(t) = e^(I*t) - 0.5*e^(-I*t)
sage: fprime(t) = I*e^(I*t) + 0.5*I*e^(-I*t)
sage: m = Riemann_Map([f], [fprime], 0)
sage: m.riemann_map(0.25 + sqrt(-0.5))
(0.137514...+0.87669602...j)
sage: m.riemann_map(1.3*I)
(-1.56...e-05+0.989694...j)
sage: I = CDF.gen()
sage: m.riemann_map(0.4)
(0.733242677...+3.2...e-06j)
sage: import numpy as np
sage: m.riemann_map(np.complex(-3, 0.0001))
(1.405757...e-05+8.06...e-10j)
Convert from a (Numpy) array of complex numbers to its corresponding matrix of RGB values. For internal use of plot_colored() only.
INPUT:
OUTPUT:
An floating point Numpy array X, where
X[i,j] is an (r,g,b) tuple.
EXAMPLES:
sage: from sage.calculus.riemann import complex_to_rgb
sage: import numpy
sage: complex_to_rgb(numpy.array([[0, 1, 1000]],dtype = numpy.complex128))
array([[[ 1., 1., 1.],
[ 1., 0., 0.],
[-998., 0., 0.]]])
sage: complex_to_rgb(numpy.array([[0, 1j, 1000j]],dtype = numpy.complex128))
array([[[ 1.00000000e+00, 1.00000000e+00, 1.00000000e+00],
[ 5.00000000e-01, 1.00000000e+00, 0.00000000e+00],
[ -4.99000000e+02, -9.98000000e+02, 0.00000000e+00]]])
TESTS:
sage: complex_to_rgb([[0, 1, 10]])
Traceback (most recent call last):
...
TypeError: Argument 'z_values' has incorrect type (expected numpy.ndarray, got list)