Lorenz 吸引子

編輯
 
Lorenz 吸引子

> restart; with(plots); > f1 := diff(x(t), t) = sigma*(y(t)-x(t));

> f2 := diff(y(t), t) = r*x(t)-y(t)-x(t)*z(t);

> f3 := diff(z(t), t) = x(t)*y(t)-b*z(t);

> DE := f1, f2, f3; > > > params := sigma = 10, r = 28, b = 8/3;

> initv := x(0) = 4, y(0) = 5, z(0) = 6;

> dvars := [x(t), y(t), z(t)];

> chaodisplay := proc (chartname) EQ := [op(subs(params, [DE])), initv];

EQ1 := dsolve(EQ, numeric);

odeplot(EQ1, dvars, t = 0 .. 30, axes = frame, numpoints = 5000, color = blue, orientation = [-30, 100], title = chartname)

end proc;

chaodisplay("Lorentz Attractor")

若斯叻吸引子

編輯
 
Rossler 吸引子

> restart; with(plots); > DE := diff(x(t), t) = -y(t)-z(t), diff(y(t), t) = x(t)+a*y(t), diff(z(t), t) = b-c*z(t)+x(t)*z(t); > params := a = .3, b = .2, c = 5;

> initv := x(0) = 1, y(0) = 2, z(0) = .3;

> dvars := [x(t), y(t), z(t)]; > chaodisplay := proc (chartname)

EQ := [op(subs(params, [DE])), initv];
EQ1 := dsolve(EQ, numeric); 

odeplot(EQ1, dvars, t = 0 .. 200, axes = frame, numpoints = 5000, color = blue, orientation = [-30, 100], title = chartname) end proc; ` > > chaodisplay("Rossler Attractor");

陳氏吸引子

編輯
 
陳氏吸引子

> restart; with(plots); > DE := diff(x(t), t) = a*(y(t)-x(t)), diff(y(t), t) = (c-a)*x(t)-x(t)*z(t)+c*y(t), diff(z(t), t) = x(t)*y(t)-b*z(t);

      > params := a = 40, c = 28, b = 3;
                    a = 40, c = 28, b = 3

> initv := x(0) = -.1, y(0) = .5, z(0) = -.6;

             x(0) = -0.1, y(0) = 0.5, z(0) = -0.6

> > > > dvars := [x(t), y(t), z(t)];

> chaodisplay := proc (chartname) EQ := [op(subs(params, [DE])), initv]; EQ1 := dsolve(EQ, numeric); odeplot(EQ1, dvars, t = 0 .. 30, axes = frame, numpoints = 5000, orientation = [-30, 100], title = chartname) end proc;

> chaodisplay("Chen Attractor");

蔡氏電路

編輯
 
蔡氏電路 吸引子

> with(plots); with(DEtools); > f := c*x(t)+(.5*(d-c))*(abs(x(t)+1)-abs(x(t)-1));

> f1 := diff(x(t), t) = a*(-x(t)+y(t)-f);


> f2 := diff(y(t), t) = x(t)-y(t)+z(t); > f3 := diff(z(t), t) = -b*y(t);

> Chua3 := (f1, f2, f3);

> params := a = 15.6, b = 25.58, c = -5/7, d = -8/7;

> vars := {x(t), y(t), z(t)};

> inits := [[x(0) = 1.8, y(0) = -.70, z(0) = -2.85], [x(0) = -2, y(0) = .5, z(0) = 3.6]];

        > 

> > Chua6 := subs(params, {Chua3});

> p1t := DEplot3d(Chua6, vars, t = 0 .. 100, inits, scene = [x(t), y(t), z(t)], stepsize = 0.5e-1, thickness = 1, linecolor = blue);

> display(p1t, font = [Helvetica, 15], title = 'Chua*Attractor', titlefont = [Courier, 15]);

壓縮子

編輯
 
壓縮子

> f := (1/4)*(4*lambda*(1/3))*cos(x(t)-lambda*t)^2;

               1                             2
               - lambda cos(-x(t) + lambda t) 
               3                              

> with(plots); > animate3d(f, x = -10 .. 10, t = 0 .. 100, lambda = 0 .. 2);

Logistic 差分方程

編輯

  V[i+1]=beta*V[i]*(1-V[i])

logistic := proc (beta, size)

local X, V, i, h, s;

X := Vector(size, 1);

for i to size do X[i] := i end do;

for i to beta do V[i] := Vector(size, i); V[i][1] := .3 end do;

for h to size-1 do for s to beta do

V[s][h+1] := (1/50)*s*V[s][h]*(1-V[s][h]) end do

end do;

plots:-display([seq(plot(X, V[i]), i = 1 .. beta)], insequence = true)

end proc

logistic(200,200);

變奏

編輯
 

logistic := proc (beta, size) local X, V, i, h, s; X := Vector(size, 1); for i to size do X[i] := i end do; for i to beta do V[i] := Vector(size, i); V[i][1] := .3 end do; for h to size-1 do for s to beta do V[s][h+1] := (1/51)*s*V[s][h]*(1-V[s][h]^2) end do end do; plots:-display([seq(plot(X, V[i]), i = 1 .. beta)], insequence = true) end proc

logistic(153, 200)

變化

編輯