'Program ORBIT3D.BAS plots orbit for 3D system of Smith, Thiffeault & Horton '(c) 2000 by J. C. Sprott 'Compile with PowerBASIC for DOS 3.5 '\$error all on SCREEN 12 line(20,0)-(639,459),,b for i%=1 to 9 line(20+619*i%/10,0)-(20+619*i%/10,10) line(20+619*i%/10,459)-(20+619*i%/10,449) line(20,459*i%/10)-(30,459*i%/10) line(639,459*i%/10)-(629,459*i%/10) next i% for xp=20 to 639 'Theory b=2*(xp-20)/619 yp=459-46*sqr(1-b/2) 'pset(xp,yp) next xp pset(639,459) defext a - z h = 0.01 'Iteration step size jmax&&=5e5 'Number of iterations for each value of b c=70 randomize timer xo=(log(1/c)-log(2-1/c))/2: yo=0: zo=0 x=xo: y=yo: z=zo for xp=638 to 21 step -1 b=2*(xp-20)/619 x=x+.1 rsum=0 for j&&=1 to jmax&& call rk4(x, y, z, h) if j&&>jmax&&/10 then rsum=rsum+sqr((x-xo)^2+y*y+z*z) end if next j&& rav=.9*rsum/jmax&& line-(xp,459-46*rav) if inkey\$=chr\$(27) then exit for next xp END def fnx(x,y,z) = y def fny(x,y,z) = z def fnz(x,y,z) = -b*z-y+1-c*(1+tanh(x)) SUB rk4 (x, y, z, h) 'Advances (x, y, z) with a fourth order Runge-Kutta method 'User must provide the functions fnx, fny, fnz, and a step size h k1x = h * fnx(x, y, z) k1y = h * fny(x, y, z) k1z = h * fnz(x, y, z) k2x = h * fnx(x + .5 * k1x, y + .5 * k1y, z + .5 * k1z) k2y = h * fny(x + .5 * k1x, y + .5 * k1y, z + .5 * k1z) k2z = h * fnz(x + .5 * k1x, y + .5 * k1y, z + .5 * k1z) k3x = h * fnx(x + .5 * k2x, y + .5 * k2y, z + .5 * k2z) k3y = h * fny(x + .5 * k2x, y + .5 * k2y, z + .5 * k2z) k3z = h * fnz(x + .5 * k2x, y + .5 * k2y, z + .5 * k2z) k4x = h * fnx(x + k3x, y + k3y, z + k3z) k4y = h * fny(x + k3x, y + k3y, z + k3z) k4z = h * fnz(x + k3x, y + k3y, z + k3z) x = x + (k1x + 2 * (k2x + k3x) + k4x) / 6 y = y + (k1y + 2 * (k2y + k3y) + k4y) / 6 z = z + (k1z + 2 * (k2z + k3z) + k4z) / 6 END SUB FUNCTION tanh (x) 'Returns hyperbolic tangent of x tanh = 1 - 2 / (EXP(2 * x) + 1) END FUNCTION