Hi,
I was asking some time a go about (pseudo) vector fields for discrete
dynamical systems build for example on complex function f(z,c):= z^2 +c;
Below is a code using vectors based on a program by Mario Rodr?guez
Riotorto
http://riotorto.users.sourceforge.net/gnuplot/vectors/index.html
It works, but image is not as I expected.
For simplest case
c :0
I was expected image similar to this :
plotdf([v,-k*z/m], [z,v], [parameters,"m=2,k=2"],
[sliders,"m=1:5"], [trajectory_at,6,0])$
Here are good images
http://www.clarku.edu/~djoyce/complex/powers.htmlhttp://facstaff.unca.edu/mcmcclur/Mathematica/ComplexDynamics/HTMLLinks/ComplexQuadAndCubicIteration_1.html
The problem is that vector draws a straight line, but I need a curve (
complex square is a rotation ).
I'm working to get better result.
Any help is welcome.
Adam
/* ---------- functions -----------------------*/
f(z):=z*z+c$
/* compute vectors at the given points vector ([x,y], [dx,dy])
vector([x,y], [dx,dy]) plots vector [dx,dy] with origin in [x,y].
*/
GiveVector(z0x,z0y):=block(
[z0,z1,dz],
z0:z0x+z0y*%i,
z1:f(z0),
dz:z1-z0,
vector([z0x,z0y],[realpart(dz)/10,imagpart(dz)/10])
)$
compile(all);
/*---------------------------- constants ----------------------*/
c:0;
HalfSide:1$
ZxMin:-HalfSide;
ZxMax:HalfSide;
dSide:ZxMax-ZxMin;
iXmin:0$
iXmax:5;
iLength : iXmax-iXmin-1$
PixelSize:dSide/iLength$
/* ----------- compute ------------------------------------ */
coord: setify(makelist(float(ZxMin+n*PixelSize),n,iXmin,iXmax))$
points2d: listify(cartesian_product(coord,coord))$ /* 2D grid */
/* make list of vectors for draw package
makelist (expr, x, list) returns a list, the j'th element of which is
equal to ev (expr, x=list[j]) for j equal to 1 through length (list).
*/
vectors: makelist(GiveVector(k[1],k[2]),k,points2d)$
/* ----------draw the vector field ----------------------------*/
draw2d(
title="Vectors ",
user_preamble = "set size square",
[head_length=0.04, color=blue],
vectors , /* Since Maxima 5.22, draw2d also accepts lists as arguments */
color=red,
polar(1,theta,0,2*%pi) /* unit circle = Julia set */
)$