matlab编写的Lyapunov指数计算程序汇总
申明：以下各程序为个人在网络上收集的Lyapunov指数计算程序，未经过验证，不保证程序的正确性和计算结果的正确性，请大家见谅，也欢迎大家探讨！计算连续方程Lyapunov指数的程序，比较好用的
其中给出了两个例子：
计算Rossler吸引子的Lyapunov指数
计算超混沌Rossler吸引子的Lyapunov指数
http://www.pudn.com/downloads39/sourcecode/math/detail132231.html 基于Volterra滤波器混沌时间序列多步预测
----------------------------------------------
参考文献：
1、张家树.混沌时间序列的Volterra自适应预测.物理学报.2000.03
2、Scott C.Douglas, Teresa H.-Y. Meng, Normalized Data Nonlinearities for LMS Adaptation. IEEE Trans.Sign.Proc. Vol.42 1994
----------------------------------------------
文件说明：
1、original_MultiStepPred_main.m 程序主文件,直接运行此文件即可
2、original_train.m 训练函数
3、original_test.m 测试函数
4、LorenzData.dll 产生Lorenz离散序列
5、normalize_1.m 归一化
6、PhaSpaRecon.m 相空间重构
7、PhaSpa2VoltCoef.m 构造 Volterra 自适应 FIR 滤波器的输入信号矢量 Un
8、TrainTestSample_2.m 将特征矩阵前 train_num 个为训练样本，其余为测试样本
9、FIR_NLMS.dll NLMS自适应算法
http://www.pudn.com/downloads45/sourcecode/math/detail150408.html recnstitution重构相空间，在非线性系统分析中有重要的作用
function =lyapunov(n,tstart,stept,tend,ystart,ioutp);
global DS;
global P;
global calculation_progress first_call;
global driver_window;
global TRJ_bufer Time_bufer bufer_i;
%
% Lyapunov exponent calcullation for ODE-system.
%
% The alogrithm employed in this m-file for determining Lyapunov
% exponents was proposed in
%
% A. Wolf, J. B. Swift, H. L. Swinney, and J. A. Vastano,
% "Determining Lyapunov Exponents from a Time Series," Physica D,
% Vol. 16, pp. 285-317, 1985.
%
% For integrating ODE system can be used any MATLAB ODE-suite methods.
% This function is a part of MATDS program - toolbox for dynamical system investigation
% See: http://www.math.rsu.ru/mexmat/kvm/matds/
%
% Input parameters:
% n - number of equation
% rhs_ext_fcn - handle of function with right hand side of extended ODE-system.
% This function must include RHS of ODE-system coupled with
% variational equation (n items of linearized systems, see Example).
% fcn_integrator - handle of ODE integrator function, for example: @ode45
% tstart - start values of independent value (time t)
% stept - step on t-variable for Gram-Schmidt renormalization procedure.
% tend - finish value of time
% ystart - start point of trajectory of ODE system.
% ioutp - step of print to MATLAB main window. ioutp==0 - no print,
% if ioutp>0 then each ioutp-th point will be print.
%
% Output parameters:
% Texp - time values
% Lexp - Lyapunov exponents to each time value.
%
% Users have to write their own ODE functions for their specified
% systems and use handle of this function as rhs_ext_fcn - parameter.
%
% Example. Lorenz system:
% dx/dt = sigma*(y - x) = f1
% dy/dt = r*x - y - x*z = f2
% dz/dt = x*y - b*z = f3
%
% The Jacobian of system:
% | -sigmasigma0 |
% J = | r-z -1 -x |
% | y x -b |
%
% Then, the variational equation has a form:
%
% F = J*Y
% where Y is a square matrix with the same dimension as J.
% Corresponding m-file:
% function f=lorenz_ext(t,X)
% SIGMA = 10; R = 28; BETA = 8/3;
% x=X(1); y=X(2); z=X(3);
%
% Y= [X(4), X(7), X(10);
% X(5), X(8), X(11);
% X(6), X(9), X(12)];
% f=zeros(9,1);
% f(1)=SIGMA*(y-x); f(2)=-x*z+R*x-y; f(3)=x*y-BETA*z;
%
% Jac=[-SIGMA,SIGMA,0; R-z,-1,-x; y, x,-BETA];
%
% f(4:12)=Jac*Y;
%
% Run Lyapunov exponent calculation:
%
% =lyapunov(3,@lorenz_ext,@ode45,0,0.5,200,,10);
%
% See files: lorenz_ext, run_lyap.
%
% --------------------------------------------------------------------
% Copyright (C) 2004, Govorukhin V.N.
% This file is intended for use with MATLAB and was produced for MATDS-program
% http://www.math.rsu.ru/mexmat/kvm/matds/
% lyapunov.m is free software. lyapunov.m is distributed in the hope that it
% will be useful, but WITHOUT ANY WARRANTY.
%
%
% n=number of nonlinear odes
% n2=n*(n+1)=total number of odes
%
options = odeset('RelTol',DS(1).rel_error,'AbsTol',DS(1).abs_error,'MaxStep',DS(1).max_step,...
'OutputFcn',@odeoutp,'Refine',0,'InitialStep',0.001);
n_exp = DS(1).n_lyapunov;
n1=n; n2=n1*(n_exp+1);
neq=n2;
%Number of steps
nit = round((tend-tstart)/stept)+1;
% Memory allocation
y=zeros(n2,1);
cum=zeros(n2,1);
y0=y;
gsc=cum;
znorm=cum;
% Initial values
y(1:n)=ystart(:);
for i=1:n_exp y((n1+1)*i)=1.0; end;
t=tstart;
Fig_Lyap = figure;
set(Fig_Lyap,'Name','Lyapunov exponents','NumberTitle','off');
set(Fig_Lyap,'CloseRequestFcn','');
hold on;
box on;
timeplot = tstart+(tend-tstart)/10;
axis();
title('Dynamics of Lyapunov exponents');
xlabel('t');
ylabel('Lyapunov exponents');
Fig_Lyap_Axes = findobj(Fig_Lyap,'Type','axes');
for i=1:n_exp
PlotLyap{i}=plot(tstart,0);
end;
uu=findobj(Fig_Lyap,'Type','line');
for i=1:size(uu,1)
set(uu(i),'EraseMode','none') ;
set(uu(i),'XData',[],'YData',[]);
set(uu(i),'Color',);
end
ITERLYAP = 0;
% Main loop
calculation_progress = 1;
while t<tend
tt = t + stept;
ITERLYAP = ITERLYAP+1;
if tt>tend, tt = tend; end;
% Solutuion of extended ODE system
% = feval(fcn_integrator,rhs_ext_fcn,,y);
while calculation_progress == 1
= integrator(DS(1).method_int,@ode_lin,,y,options,P,n,neq,n_exp);
first_call = 0;
if calculation_progress == 99, break; end;
if ( T(size(T,1))<tt ) & (calculation_progress~=0)
y=Y(size(Y,1),:);
y(1,1:n)=TRJ_bufer(bufer_i,1:n);
t = Time_bufer(bufer_i);
calculation_progress = 1;
else
calculation_progress = 0;
end;
end;
if (calculation_progress == 99)
break;
else
calculation_progress = 1;
end;
t=tt;
y=Y(size(Y,1),:);
first_call = 0;
%
% construct new orthonormal basis by gram-schmidt
%
znorm(1)=0.0;
for j=1:n1 znorm(1)=znorm(1)+y(n1+j)^2; end;
znorm(1)=sqrt(znorm(1));
for j=1:n1 y(n1+j)=y(n1+j)/znorm(1); end;
for j=2:n_exp
for k=1:(j-1)
gsc(k)=0.0;
for l=1:n1 gsc(k)=gsc(k)+y(n1*j+l)*y(n1*k+l); end;
end;
for k=1:n1
for l=1:(j-1)
y(n1*j+k)=y(n1*j+k)-gsc(l)*y(n1*l+k);
end;
end;
znorm(j)=0.0;
for k=1:n1 znorm(j)=znorm(j)+y(n1*j+k)^2; end;
znorm(j)=sqrt(znorm(j));
for k=1:n1 y(n1*j+k)=y(n1*j+k)/znorm(j); end;
end;
%
% update running vector magnitudes
%
for k=1:n_exp cum(k)=cum(k)+log(znorm(k)); end;
%
% normalize exponent
%
rescale = 0;
u1 =1.e10;
u2 =-1.e10;
for k=1:n_exp
lp(k)=cum(k)/(t-tstart);
%Plot
Xd=get(PlotLyap{k},'Xdata');
Yd=get(PlotLyap{k},'Ydata');
if timeplot<t
u1=min(u1,min(Yd));
u2=max(u2,max(Yd));
end;
Xd=; Yd=;
set(PlotLyap{k},'Xdata',Xd,'Ydata',Yd);
end;
if timeplot<t
timeplot = timeplot+(tend-tstart)/20;
figure(Fig_Lyap);
axis(); end;
drawnow;
% Output modification
if ITERLYAP==1
Lexp=lp;
Texp=t;
else
Lexp=;
Texp=;
end;
if (mod(ITERLYAP,ioutp)==0)
for k=1:n_exp
txtstring{k}=['\lambda_' int2str(k) '=' num2str(lp(k))];
end
legend(Fig_Lyap_Axes,txtstring,3);
end;
end;
ss=warndlg('Attention! Plot of lyapunov exponents will be closed!','Press OK to continue!');
uiwait(ss);
delete(Fig_Lyap);
fprintf('\n \n Results of Lyapunov exponents calculation: \n t=%6.4f',t);
for k=1:n_exp fprintf(' L%d=%f; ',k,lp(k)); end;
fprintf('\n');
if ~isempty(driver_window)
if ishandle(driver_window)
delete(driver_window);
driver_window = [];
end;
end;
calculation_progress = 0;
update_ds;
wolf 方法计算李雅普诺夫指数
http://www.pudn.com/downloads52/sourcecode/math/detail178303.html 给出了分形计算的源代码的matlab程序，可以迅速帮助大家进行分形的分析与计算，参数容易设置
function =new1lyapunov(n,rhs_ext_fcn,fcn_integrator,tstart,stept,tend,ystart,ioutp,d);
%
% Lyapunov exponent calcullation for ODE-system.
%
% The alogrithm employed in this m-file for determining Lyapunov
% exponents was proposed in
%
% A. Wolf, J. B. Swift, H. L. Swinney, and J. A. Vastano,
% "Determining Lyapunov Exponents from a Time Series," Physica D,
% Vol. 16, pp. 285-317, 1985.
%
% For integrating ODE system can be used any MATLAB ODE-suite methods.
% This function is a part of MATDS program - toolbox for dynamical system investigation
% See: http://www.math.rsu.ru/mexmat/kvm/matds/
%
% Input parameters:
% n - number of equation
% rhs_ext_fcn - handle of function with right hand side of extended ODE-system.
% This function must include RHS of ODE-system coupled with
% variational equation (n items of linearized systems, see Example).
% fcn_integrator - handle of ODE integrator function, for example: @ode45
% tstart - start values of independent value (time t)
% stept - step on t-variable for Gram-Schmidt renormalization procedure.
% tend - finish value of time
% ystart - start point of trajectory of ODE system.
% ioutp - step of print to MATLAB main window. ioutp==0 - no print,
% if ioutp>0 then each ioutp-th point will be print.
%
% Output parameters:
% Texp - time values
% Lexp - Lyapunov exponents to each time value.
%
% Users have to write their own ODE functions for their specified
% systems and use handle of this function as rhs_ext_fcn - parameter.
%
% Example. Lorenz system:
% dx/dt = sigma*(y - x) = f1
% dy/dt = r*x - y - x*z = f2
% dz/dt = x*y - b*z = f3
%
% The Jacobian of system:
% | -sigmasigma0 |
% J = | r-z -1 -x |
% | y x -b |
%
% Then, the variational equation has a form:
%
% F = J*Y
% where Y is a square matrix with the same dimension as J.
% Corresponding m-file:
% function f=lorenz_ext(t,X)
% SIGMA = 10; R = 28; BETA = 8/3;
% x=X(1); y=X(2); z=X(3);
%
% Y= [X(4), X(7), X(10);
% X(5), X(8), X(11);
% X(6), X(9), X(12)];
% f=zeros(9,1);
% f(1)=SIGMA*(y-x); f(2)=-x*z+R*x-y; f(3)=x*y-BETA*z;
%
% Jac=[-SIGMA,SIGMA,0; R-z,-1,-x; y, x,-BETA];
%
% f(4:12)=Jac*Y;
%
% Run Lyapunov exponent calculation:
%
% =lyapunov(3,@lorenz_ext,@ode45,0,0.5,200,,10);
%
% See files: lorenz_ext, run_lyap.
%
% --------------------------------------------------------------------
% Copyright (C) 2004, Govorukhin V.N.
% This file is intended for use with MATLAB and was produced for MATDS-program
% http://www.math.rsu.ru/mexmat/kvm/matds/
% lyapunov.m is free software. lyapunov.m is distributed in the hope that it
% will be useful, but WITHOUT ANY WARRANTY.
%
%
% n=number of nonlinear odes
% n2=n*(n+1)=total number of odes
%
n1=n; n2=n1*(n1+1);
%Number of steps
nit = round((tend-tstart)/stept);
% Memory allocation
y=zeros(n2,1); cum=zeros(n1,1); y0=y;
gsc=cum; znorm=cum;
% Initial values
y(1:n)=ystart(:);
for i=1:n1 y((n1+1)*i)=1.0; end;
t=tstart;
% Main loop
for ITERLYAP=1:nit
% Solutuion of extended ODE system
= unit(t,stept,y,d);
t=t+stept;
y=Y(size(Y,1),:);
for i=1:n1
for j=1:n1 y0(n1*i+j)=y(n1*j+i); end;
end;
%
% construct new orthonormal basis by gram-schmidt
%
znorm(1)=0.0;
for j=1:n1 znorm(1)=znorm(1)+y0(n1*j+1)^2; end;
znorm(1)=sqrt(znorm(1));
for j=1:n1 y0(n1*j+1)=y0(n1*j+1)/znorm(1); end;
for j=2:n1
for k=1:(j-1)
gsc(k)=0.0;
for l=1:n1 gsc(k)=gsc(k)+y0(n1*l+j)*y0(n1*l+k); end;
end;
for k=1:n1
for l=1:(j-1)
y0(n1*k+j)=y0(n1*k+j)-gsc(l)*y0(n1*k+l);
end;
end;
znorm(j)=0.0;
for k=1:n1 znorm(j)=znorm(j)+y0(n1*k+j)^2; end;
znorm(j)=sqrt(znorm(j));
for k=1:n1 y0(n1*k+j)=y0(n1*k+j)/znorm(j); end;
end;
%
% update running vector magnitudes
%
for k=1:n1 cum(k)=cum(k)+log(znorm(k)); end;
%
% normalize exponent
%
for k=1:n1
lp(k)=cum(k)/(t-tstart);
end;
% Output modification
if ITERLYAP==1
Lexp=lp;
Texp=t;
else
Lexp=lp;
Texp=t;
end;
for i=1:n1
for j=1:n1
y(n1*j+i)=y0(n1*i+j);
end;
end;
end; 小数据量法计算 Lyapunov 指数的 Matlab 程序 - (mex 函数，超快)
----------------------------------------------
参考文献：
Michael T.Rosenstein,
"A practical method for calculating largest lyapunov exponents from small sets",
Physica D,1993,65: 117-134
文件说明：
----------------------------------------------
1、LargestLyapunov_example1.m 程序主文件1，直接运行此文件即可
2、LargestLyapunov_example2.m 程序主文件2，直接运行此文件即可
3、LorenzData.dll 产生 Lorenz 离散数据
4、PhaSpaRecon.m 相空间重构
5、Lyapunov_luzhenbo.dll Lyapunov 计算主函数
6、lyapunov_buffer.dll Lyapunov 计算缓存
http://www.pudn.com/downloads63/sourcecode/math/detail221870.html 计算lyapunov指数的程序
program lylorenz
parameter(n=3,m=12,st=100)
integer::i,j,k
real y(m),z(n),s(n),yc(m),h,y1(m),a,b,r,f(m),k1,k2,k3
y(1)=10.
y(2)=1.
y(3)=0.
a=10.
b=8./3.
r=28.
t=0.
h=0.01
!!!!!initial conditions
do i=n+1,m
y(i)=0.
end do
doi=1,n
y((n+1)*i)=1.
s(i)=0
end do
open(1,file='lorenz1.dat')
open(2,file='ly lorenz.dat')
do 100 k=1,st !!!!!!!!st iterations
call rgkt(m,h,t,y,f,yc,y1)
!!!!normarize vector 123
z=0.
do i=1,n
do j=1,n
z(i)=z(i)+y(n*j+i)**2
enddo
if(z(i)>0.)z(i)=sqrt(z(i))
do j=1,n
y(n*j+i)=y(n*j+i)/z(i)
enddo
end do
!!!!generate gsr coefficient
k1=0.
k2=0.
k3=0.
do i=1,n
k1=k1+y(3*i+1)*y(3*i+2)
k2=k2+y(3*i+3)*y(3*i+2)
k3=k3+y(3*i+1)*y(3*i+3)
end do
!!!!conduct new vector2 and 3
do i=1,n
y(3*i+2)=y(3*i+2)-k1*y(3*i+1)
y(3*i+3)=y(3*i+3)-k2*y(3*i+2)-k3*y(3*i+1)
end do
!!!generate new vector2 and 3,normarize them
do i=2,n
z(i)=0.
do j=2,n
z(i)=z(i)+y(n*j+i)**2
enddo
if(z(i)>0.)z(i)=sqrt(z(i))
do j=2,n
y(n*j+i)=y(n*j+i)/z(i)
end do
end do
!!!!!!!update lyapunov exponent
do i=1,n
if(z(i)>0)s(i)=s(i)+log(z(i))
enddo
100 continue
do i=1,n
s(i)=s(i)/(1.*st*h*1000.)
write(2,*)s(i)
enddo
end
!!!!!!!!!!!!!!!!!!!!!!!!!
subroutine rgkt(m,h,t,y,f,yc,y1)
real y(m),f(m),y1(m),yc(m),a,b,r
integer::i,j
do j=1,1000
call df(m,t,y,f)
t=t+h/2.0
do i=1,m
yc(i)=y(i)+h*f(i)/2.0
y1(i)=y(i)+h*f(i)/6.0
end do
call df(m,t,yc,f)
do i=1,m
yc(i)=y(i)+h*f(i)/2.0
y1(i)=y1(i)+h*f(i)/3.0
end do
call df(m,t,yc,f)
t=t+h/2.0
do i=1,m
yc(i)=y(i)+h*f(i)
y1(i)=y1(i)+h*f(i)/3.0
end do
call df(m,t,yc,f)
do i=1,m
y(i)=y1(i)+h*f(i)/6.0
end do
if(j>500)write(1,*)t,y(1),y(2),y(3)
end do
return
end
!!!!!!!!!!!!!!!!!!!!!!!!
subroutine df(m,t,y,f)
real y(m),a,b,r,f(m)
common a,b,r
a=10.
b=8./3.
r=28.
f(1)=a*(y(2)-y(1))
f(2)=y(1)*(r-y(3))-y(2)
f(3)=y(1)*y(2)-b*y(3)
doi=0,2
f(4+i)=a*y(7+i)-y(4+i)
f(7+i)=y(4+i)*(r-y(3))-y(7+i)-y(1)*y(10+i)
f(10+i)=y(2)*y(4+i)-b*y(10+i)+y(1)*y(7+i)
enddo
return
end C-C方法计算时间延迟和嵌入维数计算Lyapunov指数计算关联维数混沌时间序列预测
C-C方法计算时间延迟和嵌入维数
主程序：C_CMethod.m, C_CMethod_independent.m
子函数：correlation_integral.m(计算关联积分)
disjoint.m(将时间序列拆分成t个不相关的子序列)
heaviside.m（计算时间序列的海维赛函数值）
参考文献Nonlinear dynamics, delay times, and embedding windows。
计算Lyapunov指数：
largest_lyapunov_exponent.m（用吕金虎的方法计算最大Lyapunov指数）
参考文献：基于Lyapunov指数改进算法的边坡位移预测。
lyapunov_wolf.m(用wolf方法计算最大Lyapunov指数)
计算关联维数：G_P.m(G-P算法)
混沌时间序列预测
主函数
MainPre_by_jiaquanyijie_1.m(该程序用加权一阶局域法对数据进行进行一步预测)
MainPre_by_jiaquanyijie_n.m(该程序用加权一阶局域法对数据进行进行n步预测)
MainPre_by_Lya_1.m(基于最大Lyapunov指数的一步预测)
MainPre_by_Lya_n.m(基于最大Lyapunov指数的n步预测)
nearest_point.m(计算最后一个相点的最近相点的位置及最短距离)
子函数
jiaquanyijie.m(该函数用加权一阶局域法(xx)、零级近似(yy)和基于零级近似的加权一阶局域法(zz)对时间数据进行一步预测)
pre_by_lya.m(基于最大Lyapunov指数的预测方法)
pre_by_lya_new.m(改进的基于最大Lyapunov指数的预测方法)
如果需要该工具包，大家到http://www.pudn.com/downloads25/sourcecode/others/detail82624.html下载 求取lyapunov指数的小数据量方法，采用混合编程
lylorenz.f90的程序如下：
program lylorenz
parameter(n=3,m=12,st=100)
integer::i,j,k
real y(m),z(n),s(n),yc(m),h,y1(m),a,b,r,f(m),k1,k2,k3
y(1)=10.
y(2)=1.
y(3)=0.
a=10.
b=8./3.
r=28.
t=0.
h=0.01
!!!!!initial conditions
do i=n+1,m
y(i)=0.
end do
doi=1,n
y((n+1)*i)=1.
s(i)=0
end do
open(1,file='lorenz1.dat')
open(2,file='ly lorenz.dat')
do 100 k=1,st !!!!!!!!st iterations
call rgkt(m,h,t,y,f,yc,y1)
!!!!normarize vector 123
z=0.
do i=1,n
do j=1,n
z(i)=z(i)+y(n*j+i)**2
enddo
if(z(i)>0.)z(i)=sqrt(z(i))
do j=1,n
y(n*j+i)=y(n*j+i)/z(i)
enddo
end do
!!!!generate gsr coefficient
k1=0.
k2=0.
k3=0.
do i=1,n
k1=k1+y(3*i+1)*y(3*i+2)
k2=k2+y(3*i+3)*y(3*i+2)
k3=k3+y(3*i+1)*y(3*i+3)
end do
!!!!conduct new vector2 and 3
do i=1,n
y(3*i+2)=y(3*i+2)-k1*y(3*i+1)
y(3*i+3)=y(3*i+3)-k2*y(3*i+2)-k3*y(3*i+1)
end do
!!!generate new vector2 and 3,normarize them
do i=2,n
z(i)=0.
do j=2,n
z(i)=z(i)+y(n*j+i)**2
enddo
if(z(i)>0.)z(i)=sqrt(z(i))
do j=2,n
y(n*j+i)=y(n*j+i)/z(i)
end do
end do
!!!!!!!update lyapunov exponent
do i=1,n
if(z(i)>0)s(i)=s(i)+log(z(i))
enddo
100 continue
do i=1,n
s(i)=s(i)/(1.*st*h*1000.)
write(2,*)s(i)
enddo
end
!!!!!!!!!!!!!!!!!!!!!!!!!
subroutine rgkt(m,h,t,y,f,yc,y1)
real y(m),f(m),y1(m),yc(m),a,b,r
integer::i,j
do j=1,1000
call df(m,t,y,f)
t=t+h/2.0
do i=1,m
yc(i)=y(i)+h*f(i)/2.0
y1(i)=y(i)+h*f(i)/6.0
end do
call df(m,t,yc,f)
do i=1,m
yc(i)=y(i)+h*f(i)/2.0
y1(i)=y1(i)+h*f(i)/3.0
end do
call df(m,t,yc,f)
t=t+h/2.0
do i=1,m
yc(i)=y(i)+h*f(i)
y1(i)=y1(i)+h*f(i)/3.0
end do
call df(m,t,yc,f)
do i=1,m
y(i)=y1(i)+h*f(i)/6.0
end do
if(j>500)write(1,*)t,y(1),y(2),y(3)
end do
return
end
!!!!!!!!!!!!!!!!!!!!!!!!
subroutine df(m,t,y,f)
real y(m),a,b,r,f(m)
common a,b,r
a=10.
b=8./3.
r=28.
f(1)=a*(y(2)-y(1))
f(2)=y(1)*(r-y(3))-y(2)
f(3)=y(1)*y(2)-b*y(3)
doi=0,2
f(4+i)=a*y(7+i)-y(4+i)
f(7+i)=y(4+i)*(r-y(3))-y(7+i)-y(1)*y(10+i)
f(10+i)=y(2)*y(4+i)-b*y(10+i)+y(1)*y(7+i)
enddo
return
end 计算各种混沌系统李雅普洛夫指数的MATLAB 源程序
http://www.pudn.com/downloads50/sourcecode/math/detail172856.html 正好帮上忙了
让我看看程序的正确性怎么样
先表示感谢
这个工作量是很大的，辛苦了 想问一下对于不连续的非线性方程,有没有计算程序,
我知道有文章讨论这个问题
P.Muller,Chaos Solitions Fractals 5(1995),1671
但还找不到程序 离散问题的论坛也有好几个程序，自己搜索一下吧
[ 本帖最后由 ChaChing 于 2010-5-11 10:36 编辑 ]
应用中出现问题，向大侠求救ing
各位好，我在用楼主所提供的C-C方法求时间延迟时遇到了下面的错误提示，请各位帮小弟诊断一下，谢谢！t =
1
t =
2
??? Index exceeds matrix dimensions.
有朋友试过这个C-C这个程序吗？怎么我试了下出现问题的 ，
[ 本帖最后由 ChaChing 于 2010-5-11 10:35 编辑 ] 不知道你给的data是什么样的？我这里没有这样的序列
从错误提示上看估计是因为提供序列有问题造成的
[ 本帖最后由 ChaChing 于 2010-5-11 10:36 编辑 ]