关于duffing系统求Lyapunov指数的程序，请各位高手指点
function dX = DF(t,X,flag,alf)global alf;a=1;b=0.25;c=1;omega=1;x=X(1);y=X(2);z=X(3);Y = ; dX = zeros(12,1);dX(1)=y;dX(2)=alf*cos(z)+a*x-c*x^3-b*y;dX(3)=omega;J=;dX(4:12) = J*Y;endclear all;
clc;
% Z1=[ ];
% Z2=[ ];
% Z3=[ ];
global alf;
options=odeset('RelTol',1e-8,'AbsTol',1e-9);
for mm=1:100;
alf=(mm-1)*0.05
a0=0.1;b0=0.1
yinit = ;
orthmatrix = eye(3);
IC(1:3) = yinit;
IC(4:12) = orthmatrix;
InitialTime = 0; % 时间初始值
TimeStep = pi/180; % 时间步长
wholetimes = 300*pi; % 总的循环时间
UpdateStepNum =360; % 每次演化的步数
iteratetimes = wholetimes/(UpdateStepNum*TimeStep); % 演化的次数
DiscardItr=10; % 抛弃的不稳定迭代次数
Iteration=UpdateStepNum*TimeStep; % 迭代次数
DiscardTime=DiscardItr*Iteration+InitialTime; % 抛弃的时间段
mod = zeros(1,3);
d=3; % 维数
Sum=zeros(1,d);
n=0; % 总的迭代计数
k=0; % 对结果有作用的迭代计数
xData=[];
yData=[];
T1=InitialTime;
T2=T1+Iteration;
TSpan=;
for i=1:iteratetimes
n=n+1
=ode45('DF', TSpan, IC,options);
=size(X); %提取行数和列数
L=cX-d*d; % 取积分得到的最后一个时刻的值
for j=1:d
m1=L+1+(j-1)*d;
m2=m1+d-1;
A(:,j)=(X(rX,m1:m2))';
end
y0 = ThreeGS(A); %正交归一化
mod(1) = sqrt(y0(:,1)'*y0(:,1)); % 取三个向量的模
mod(2) = sqrt(y0(:,2)'*y0(:,2));
mod(3) = sqrt(y0(:,3)'*y0(:,3));
y0(:,1) = y0(:,1)/mod(1);
y0(:,2) = y0(:,2)/mod(2);
y0(:,3) = y0(:,3)/mod(3);
y(4:12) = y0';
if (T2>DiscardTime)
k=k+1;
T=k*Iteration;
TT=n*Iteration+InitialTime;
%计算lyapunov指数
Sum=Sum+log(abs(mod));
lambda=Sum/T;
%按降序排列指数
Lambda=fliplr(sort(lambda));
xData=;
yData=;
end
% 重新定义起始时刻
ic=X(rX,1:L);
T1=T2;
T2=T2+Iteration;
TSpan=;
IC=;
end
YY(mm,:)=Lambda'
end
本程序是从论坛下载，修改的，DF函数我还是不太明白另外结果好像不对，请赐教
DF函数应该就是你所要计算的duffing方程 请问你解决了吗我也在做这些可以想你请教下吗 1713573225 发表于 2015-11-12 21:13
请问你解决了吗我也在做这些可以想你请教下吗
有问题建议直接在论坛提出来，大家一起讨论
因为楼主未必能够看到你发的内容
页:
[1]