Thursday, April 12, 2007

Lorenz attractor

#############
## LORENZ ##
#############

a=10; r=28; b=8/3; dt=0.01
X=0.01; Y=0.01; Z=0.01

n=60000
XYZ=array(0,dim=c(n,3))
t=0

for(i in 1:n)
{
X1=X; Y1=Y; Z1=Z

X=X1+(-a*X1+a*Y1)*dt
Y=Y1+(-X1*Z1+r*X1-Y1)*dt

Z=Z1+(X1*Y1-b*Z1)*dt
#t=t+dt
XYZ[i,]=c(X,Y,Z)
}


# install.packages("scatterplot3d") # Install if needed.
library(scatterplot3d)

x=XYZ[,1]
y=XYZ[,2]
z=XYZ[,3]
X11()

s3d<-scatterplot3d(x, y, z, highlight.3d=TRUE, xlab="x",ylab="y",zlab="z",type="p", col.axis="blue", angle=55, col.grid="lightblue", scale.y=0.7, pch=".", xlim=c(min(x),max(x)), ylim=c(min(y),max(y)), zlim=c(min(z),max(z)))

s3d.coords <- s3d$xyz.convert(x,y,z)
D3_coord=cbind(s3d.coords$x,s3d.coords$y)

title("Lorenz Attractor")




1 comment:

David said...

thanks dude, i was looking for exactly this! keep up the good work with R:)