Open in CoCalc with one click!
# Code for plotting complex-valued functions # # License: GPLv3 # Copyright: 2016 # Author: Gregory Bard <[email protected]> # Based in part on code written by Harald Schilly <[email protected]> def new_complex_plot3d( f, xmin, xmax, ymin, ymax, zmin, zmax, style, **kwargs ): """This plots a function whose sole input is a complex number, and whose sole output is a complex number. There are many ways to plot such a function, and they are chosen with the parameter style. The options for style are 'magnitude', 'real', 'imaginary', 'argument', or 'mixed'. Also 'phase' is a synonym for 'argument', 'rainbow' is a synonym for 'mixed', and 'norm' is a synonym for 'magnitude'. Consider the function to be plotted as f(s) = f(x+iy). To be clear, the input to the complex-function f is s=x+iy. The 3D plot has an x-axis, a y-axis, and a z-axis. The x-axis is always the real part of the input, and the y-axis is always the imaginary part of the input. The z-axis is determined by the choice of style. For 'magnitude' or 'norm' it is the magnitude/norm of the output of f(s) = f(x+iy). For 'real' it is the real part of the output of f(s) = f(x+iy). For 'imaginary' it is the imaginary part of the output of f(s) = f(x+iy). For 'argument' or 'phase' it is the argument/phase of the output of f(s) = f(x+iy). (In other words, if you think of the complex number in polar coordinates, it is the theta of the output, in the sense of x+iy = rho*(cos(theta) + i*sin(theta)). The rho is the magnitude. There is a known issue. If the z-coordinate for any point overflows the zmax, or underflows the zmin, then it will be plotted as if exactly equal to zmax or zmin respectively. This issue results in a bizarre "shelf" forming at z=zmax if the plot overflows for any non-trivial fraction of the plotting area, or a similar shelf forming at z=zmin if the plot underflows for any non-trivial fraction of the plotting area. """ g(x,y) = f( x + i*y ) assert xmin < xmax, "Error: xmax <= xmin." assert ymin < ymax, "Error: ymax <= ymin." assert zmin < zmax, "Error: zmax <= zmin." # the reason for the above asserts is to ensure that each of the following three # is always strictly positive, never zero nor never negative. xwide = xmax - xmin ywide = ymax - ymin zwide = zmax - zmin # This line of code will see if mesh has been specified by the calling # code. If it has *not* been specified, it will be set to True. kwargs['mesh'] = kwargs.get('mesh', True ) # This line of code will see if plot_points has been specified by the # calling code. If it has *not* been specified, it will be set to 25. kwargs['plot_points'] = kwargs.get('plot_points', 25) # This line of code will see if the aspect_ratio has been specified by # the calling code. If it has *not* been specified, it will be set # to [1/xwide, 1/ywide, 1/zwide ] kwargs['aspect_ratio'] = kwargs.get('apsect_ratio', [1/xwide, 1/ywide, 1/zwide] ) # This line of code will see if color has been specified by the # calling code. If it has *not* been specified, it will be set to seagreen. # Note: for 'mixed' or 'rainbow' plots, this will be over-ridden. kwargs['color'] = kwargs.get('color', 'seagreen') if ( (style=='argument') or (style=='Argument') or (style=='ARGUMENT') or (style=='phase') or (style=='Phase') or (style=='PHASE')): # since the theta always has -pi < theta <= pi, we can be pedantic and warn # the user if they choose sloppy bounds. if ((zmin < -pi)): print "Warning: zmin set to ", zmin, "but argument/phase is never less than -pi." print if ((zmax > pi)): print "Warning: zmax set to ", zmax, "but argument/phase is never greater than pi." print P=plot3d( lambda x, y : CC(g(x=x,y=y)).argument(), (x, xmin, xmax), (y, ymin, ymax), **kwargs ) elif ( (style=='magnitude') or (style=='Magnitude') or (style=='MAGNITUDE') or (style=='norm') or (style=='Norm') or (style=='NORM')): # since the magnitude is always positive, we can be pedantic and warn # the user if they choose sloppy bounds. if ((zmin < 0)): print "Warning: zmin set to ", zmin, "but magnitude/norm is never less than zero." print P=plot3d( lambda x, y : max(min(CC(g(x=x,y=y)).norm(), zmax), zmin), (x, xmin, xmax), (y, ymin, ymax), **kwargs ) elif ( (style=='rainbow') or (style=='Rainbow') or (style=='RAINBOW') or (style=='mixed') or (style=='Mixed') or (style=='MIXED')): # since the magnitude is always positive, we can be pedantic and warn # the user if they choose sloppy bounds. if ((zmin < 0)): print "Warning: zmin set to ", zmin, "but magnitude/norm is never less than zero." print T = lambda x,y : ( CC( g(x=x,y=y) ).argument() + N(pi) ) / (2*N(pi)) new_color_option = (T, colormaps.gist_rainbow) kwargs['color'] = new_color_option kwargs['viewer'] = 'tachyon' P=plot3d( lambda x, y : max(min(CC(g(x=x,y=y)).norm(), zmax), zmin), (x, xmin, xmax), (y, ymin, ymax), **kwargs ) elif ( (style=='real') or (style=='Real') or (style=='REAL')): P=plot3d( lambda x, y : max(min(CC(g(x=x,y=y)).real_part(), zmax), zmin), (x, xmin, xmax), (y, ymin, ymax), **kwargs ) elif ( (style=='imaginary') or (style=='Imaginary') or (style=='IMAGINARY')): P=plot3d( lambda x, y : max(min(CC(g(x=x,y=y)).imag_part(), zmax), zmin), (x, xmin, xmax), (y, ymin, ymax), **kwargs ) else: assert false, "Your choice of style is unknown. Use new_complex_plot3d? to get a list of known styles." # The follow cludge is hard to explain. # Essentially, plot3d will ignore zmin and zmax, and this ruins all sorts of situations. # Instead, plot3d uses the actual largest z and smallest z observed for all possible (x,y) in the plotting domain. # # For example, if plotting where some z's are large, like e^z, then the largest observed z value might be huge. # (This also ruins the aspect ratios, by the way.) # # so to force plot3d to behave, we add two tiny invisible dots # ... one is at the center of the coordinate plane, and z=zmax # ... the other is at the center of the coordinate plane, and z=zmin # # We also use the min() and max() commands to ensure that each plotted z-value is cut off, forcing zmin < z < zmax # # These measures forces the z bounds to be correct, and all behaves nicely (except the underflow/overflow shelf) P_top = point3d( ((xmin+xmax)/2, (ymin+ymax)/2, zmax), size=0 ) P_bottom = point3d( ((xmin+xmax)/2, (ymin+ymax)/2, zmin), size=0 ) return P + P_top + P_bottom
var("x y z") # new_complex_plot3d( exp(z), -3,3, -3,3, -pi,pi, 'argument') # working # new_complex_plot3d( exp(z), -3,3, -3,3, 0,12, 'magnitude') # working # new_complex_plot3d( exp(z), -3,3, -3,3, 0,12, 'mixed') # working # new_complex_plot3d( exp(z), -3,3, -3,3, -6,6, 'real') # working # new_complex_plot3d( exp(z), -3,3, -3,3, -6,6, 'imaginary') # working p(z) = z^3 - z^2 + z - 1 # new_complex_plot3d( p(z), -3,3, -3,3, -pi,pi, 'argument') # working # new_complex_plot3d( p(z), -3,3, -3,3, 0,12, 'magnitude') # working # new_complex_plot3d( p(z), -3,3, -3,3, 0,12, 'mixed') # working # new_complex_plot3d( p(z), -3,3, -3,3, -3,3, 'real') # working # new_complex_plot3d( p(z), -3,3, -3,3, -3,3, 'imaginary') # working # new_complex_plot3d( sqrt(z), -3,3, -3,3, -pi,pi, 'argument') # working # new_complex_plot3d( sqrt(z), 2.1,2.2, 2.1,2.2, -pi,pi, 'argument') # working # new_complex_plot3d( sqrt(z), -3,3, -3,3, 0,6, 'magnitude') # working # new_complex_plot3d( sqrt(z), -3,3, -3,3, 0,6, 'mixed') # working # new_complex_plot3d( sqrt(z), -3,3, -3,3, -3,3, 'real') # working # new_complex_plot3d( sqrt(z), -3,3, -3,3, -3,3, 'imaginary') # working # new_complex_plot3d( log(z), -3,3, -3,3, -pi,pi, 'argument') # working # new_complex_plot3d( log(z), -3,3, -3,3, 0,6, 'magnitude') # working # new_complex_plot3d( log(z), -3,3, -3,3, 0,6, 'mixed') # working # new_complex_plot3d( log(z), -3,3, -3,3, -3,3, 'real') # working # new_complex_plot3d( log(z), -3,3, -3,3, -3,3, 'imaginary') # working # new_complex_plot3d( sin(z), -2*pi,2*pi, -2*pi,2*pi, -pi,pi, 'argument') # working # new_complex_plot3d( sin(z), -pi,pi, -pi,pi, 0,2, 'magnitude') # working new_complex_plot3d( sin(z), -pi,pi, -pi,pi, 0,2, 'mixed') # working # new_complex_plot3d( sin(z), -pi,pi, -pi,pi, -2,2, 'real') # working # new_complex_plot3d( sin(z), -pi,pi, -pi,pi, -2,2, 'imaginary') # unknown
(x, y, z)
3D rendering not yet implemented