trylorenz.m

The programme ahout scientific computation

function trylorenz
%  交互实验trylorenz，  模拟 Lorenz 方程的吸引子
%   Lorenz方程由3个非线性的微分方程组成，形式很简单。该方程是MIT的
%   Edward Lorenz 提出的。在所有的常微分方程中，被研究的最多的毫无
%   疑问就是Lorenz方程。
%   在图形界面上你看到的两个绿色圈，是该方程的平衡解。方程解的曲线
%   在3维空间上是有界（环绕它的两个平衡解）的，它不是周期的但也不
%   收敛。对充分大的时间，解的曲线定义出的集合就是著名的Lorenz吸引子
%
%   用鼠标可以改变3维相空间的透视角度，界面上的按钮提供如下功能：
%   "暂停pause", "继续resume", "停止stop", "重启restart", "清除clear", 
%   "关闭close".
%
%   界面右下角提供了参数rho的选择值。其中第一个值 28是在各种文献中
%   最常用的，它生成Lorenz混沌吸引子。
%   其他的四个值吸引子是周期的。要改变参数，首先要执行"stop"，参数
%   调整后再执行"restart"。
%   
%   注意，吸引子是指时间充分大时解的性态。为了看到吸引子，应当在
%   程序运行一段时间后，执行"暂停pause"和"清除clear"，这样将时间
%   小的时候的解去掉，然后再"继续resume"
%

% If this is first entry, just initialize the figure window.

   if ~isequal(get(gcf,'name'),'Lorenz吸引子')
      lorenzplot([],[],'init',[]);
      return
   end

% The differential equation is ydot = A(y)*y

   rhopick = findobj('tag','rhopick');
   rhos = get(rhopick,'userdata');
   rho = rhos(get(rhopick,'value'));
   sigma = 10;
   beta = 8/3;
   eta = sqrt(beta*(rho-1));
   A = [ -beta    0     eta
            0  -sigma   sigma 
         -eta   rho    -1  ];

% With this value of eta, A is singular.
% The critical points are the null vectors of A.

   yc = [rho-1; -eta; -eta];
   yc = [rho-1; eta; eta];

% The initial value of y(t) is near one of the critical points.
   
   y0 = yc + [0; 0; 3];

% Integrate forever, or until the stop button is toggled.
% The eta's in A will be replaced by y(2) during the integration.

   tspan = [0 Inf];
   opts = odeset('reltol',1.e-6,'outputfcn',@lorenzplot,'refine',4);
   ode45(@lorenzeqn, tspan, y0, opts, A);