Last active
February 1, 2019 14:26
-
-
Save isomorphisms/5a30e61fb305ee52bcff to your computer and use it in GitHub Desktop.
make complex-plane Wegert plots (plot only the argument, knowing that the length=modulus works the same as in a "regular" Cartesian plot)
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
#run convenience functions first.... | |
plat <- function(space,FUN,cex=2,chroma=45,a=66,b=4,c=3){ #a,b,c are tuning parameters | |
func=FUN(space) #eval for speedup | |
coleur=arg(func)%%360 | |
light=a+b*mod(func)+c*l(coleur) | |
plot(space, pch=46,cex=cex, col=hcl( h=coleur, c=chroma, l=light ) ) | |
} | |
#what this function does, is plot a boring series of points with nothing on them, | |
# moving the REAL ACTION to the colour parameter. | |
#So I'll draw for the space a square called Z, | |
# and then "pull back" the colours of what happens to the complex argument | |
# so all of the interesting computation happens in | |
# plot(...,col=hcl(HERE),...) | |
# instead of plot(HERE, ...). | |
#We could already see what happens to the modulus (size) with the usual Cartesian plots. | |
#The Wegert-plot shows us what's happening to the minus sign, | |
# with red being +, green being −, and other colours being "somewhere in between". | |
#This is why we like the complex plane: instead of +/− being an "on/off" switch (ℤ₂), | |
#it becomes either (in my %%360 way) a "circle" 𝕊₁ or | |
#(in Wegert's mod(x)/log(mod(x)) way), a spiraling parking-lot ramp. | |
#It's notable that, around a root (where the function equals/totals/evals to zero), | |
#you'll see a colour spiral. Try it with | |
plat(Z, function(x) (x-1)*(x-i)) | |
#for example. | |
#convenience | |
l <- function(x) x/100 - floor(x/100) | |
arg <- function(x) Arg(x)*360/2/pi | |
mod <- function(x) l(Mod(x)) | |
i <- complex(imaginary=1) #sqrt(-1) doesn't work … kinda sad … | |
#make a complex-plane template square | |
Z <- matrix(1:500,500,500) #you can do 501 if you want zero .. it's usually "bad" so I leave it out | |
Z + i*t(Z) -> Z | |
Z - complex(length.out=250, real=250,imaginary=250) -> Z | |
Z/100 -> Z | |
#fun things to plot in this way | |
plat( Z, zeta ) | |
plat(Z, function(x) zeta(x,n=15) ) | |
#airy | |
#bessel | |
plat( Z, log ) | |
plat( Z, exp ) | |
#(partial) geometric series | |
plat( Z, sinh ) | |
#polynomials | |
poly.1 <- function(x) x**9 - 3*x**5 + x**2 | |
poly.2 <- function(x) 5*x**4 + 2*i*x**3 - x**2/i | |
plat( Z, poly.1 ) | |
plat( Z, poly.2 ) | |
#ratios of polynomials | |
plat( Z, function(x) poly.1(x) / poly.2(x) ) | |
#random polynomials | |
2*pi -> tau | |
rcomplex <- function(n=1L, length=5) complex( argument=runif(n=n, min=0, max=tau), modulus=rpois(n=n, lambda=length) ) | |
rpoly <- function(x) rcomplex(1) * x**5 + rcomplex(1) * x**4 + rcomplex(1) * x**3 + rcomplex(1) * x**2 + rcomplex(1) * x | |
plat(Z, function(x) rpoly(x) - rpoly(x) ) #re-rolls the random's every time | |
plat(Z, function(x) rpoly(x) - rpoly(x) + rpoly(x) - rpoly(x) ) | |
plat(Z, function(x) rpoly(x)**2 - i*rpoly(x) + rpoly(x) - rpoly(x) ) | |
plat(Z, function(x) rpoly(x)**2 - i*rpoly(x) + rpoly(x)**3 - rpoly(x)**9 ) | |
#baseline reference | |
id <- function(x) x | |
plat( Z, id ) | |
#codes for special functions | |
zeta <- function(x,n=5) eval( parse(text= paste0( 1:n, '**-x', collapse=' + ') )) #tried it with sum(), but that discards imaginary parts *only* during plot( sum( ... )) | |
zeta <- function(x,M=1e9) (1-x)*( M**(1-x) - 1) | |
Lorentz <- function(x) 1/sqrt(1-x*x) | |
Airy <- function(x) exp(-2/3*x^1.5) / 2 / sqrt(pi) / x^.25 | |
Airy.plus <- function(x) Airy(x) + x | |
Airy.sinh <- function(x) Airy(x) + sinh(x) | |
Dilogarithm <- function(M) integrate( function(t) -log(1-t)/t, 0, M) | |
polynomial.1 <- function(x) x^3 + 5*x^2 - 2*x+i | |
polynomial.2 <- function(x) x^4/14 + x^3/13 | |
polynomial.1.2 <- function(x) polynomial.1(x) / polynomial.2(x) | |
polynomial.2.1 <- function(x) polynomial.2(x) / polynomial.1(x) | |
polynomial.222.11 <- function(x) polynomial.2(x)^3 / polynomial.1(x)^2 | |
Blaschke <- function(x,k=seed,const=rot) prod( (x-k) / (1-Conj(k)*x ) ) * rot | |
Bessel.1 <- function(x,a=1) integrate( function(t) cos(a*t - x*sin(t)), 0, pi) / pi | |
#complex digamma |
A polar plane:
U <- t( matrix( complex( modulus=rep(1:500,500), argument=1:500/500) ) )
Cool polynomials and a recursive way to generate them from graphs: en.wikipedia.org/wiki/Tutte_polynomial
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment
You can't copy-paste and run the code because some of the stuff is defined at the bottom. So either paste twice or make the
Z
and define the convenience functions first before running the whole code.