方程大概是这样的:
x_i (k+1)=f(x_i (k))+g(x_i (k-d(k)))+sum_{m=1}^{infinity}(u_m *h(x_i (k-m)))+sum_{j=1,j~=i}^{N}w_{ij} x_j (k)
是一个耦合的网络,不妨取N=3,即是3个网络耦合.下面是我编的一个画图程序,大家帮助看一下哪里错了,因为总是和手算的值不一样,这里的m的上界本应该是无穷,但是那样比较麻烦,我只取到1000,应该和实际比较接近才对.
function disc_nn()
tic; % 开始计时
clear;
k0=0; kfinal=50; h=1;
i=1;
for k=k0:h:kfinal % time range
d(i)=3+(1+(-1)^k)/2;
i=i+1;
end
u1(1,:)=fun1_ini(0);
u2(1,:)=fun2_ini(0);
u3(1,:)=fun3_ini(0);
i=2;
sumh11=0;
sumh12=0;
sumh21=0;
sumh22=0;
sumh31=0;
sumh32=0;
for k=k0+1:h:kfinal % time range
for m=1:h:1000
if (k-m)<=0 % 判断k-m的符号
u11(k,1)=u1(1,1);
u12(k,2)=u1(1,2);
u21(k,1)=u2(1,1);
u22(k,2)=u2(1,2);
u31(k,1)=u3(1,1);
u32(k,2)=u3(1,2);
else
u11(k,1)=u1((k-m),1);
u12(k,2)=u1((k-m),2);
u21(k,1)=u2((k-m),1);
u22(k,2)=u2((k-m),2);
u31(k,1)=u3((k-m),1);
u32(k,2)=u3((k-m),2);
end
sumh11=sumh11+(2^(-m-3))*h1(u11(k,1));
sumh12=sumh12+(2^(-m-3))*h2(u12(k,2));
sumh21=sumh21+(2^(-m-3))*h1(u21(k,1));
sumh22=sumh22+(2^(-m-3))*h2(u22(k,2));
sumh31=sumh31+(2^(-m-3))*h1(u31(k,1));
sumh32=sumh32+(2^(-m-3))*h2(u32(k,2));
end
if floor(k-d(k))<=0 % 判断k-d(k)近似整数的符号
u1(i,1)=0.05*u2(k,1)+0.05*u3(k,1)...
+f1(u1(k,1))+0.2*u1(k,2)+g1(u1(1,1))+sumh11;
u1(i,2)=0.05*u2(k,2)+0.05*u3(k,2)...
+f2(u1(k,2))+g2(u1(1,2))+sumh12;
u2(i,1)=0.05*u1(k,1)+0.05*u3(k,1)...
+f1(u2(k,1))+0.2*u2(k,2)+g1(u2(1,1))+sumh21;
u2(i,2)=0.05*u1(k,2)+0.05*u3(k,2)...
+f2(u2(k,2))+g2(u2(1,2))+sumh22;
u3(i,1)=0.05*u1(k,1)+0.05*u2(k,1)...
+f1(u3(k,1))+0.2*u3(k,2)+g1(u3(1,1))+sumh31;
u3(i,2)=0.05*u1(k,2)+0.05*u2(k,2)...
+f2(u3(k,2))+g2(u3(1,2))+sumh32;
else
u1(i,1)=0.05*u2(k,1)+0.05*u3(k,1)...
+f1(u1(k,1))+0.2*u1(k,2)+g1(u1(floor(k-d(k)),1))+sumh11;
u1(i,2)=0.05*u2(k,2)+0.05*u3(k,2)...
+f2(u1(k,2))+g2(u1(floor(k-d(k)),2))+sumh12;
u2(i,1)=0.05*u1(k,1)+0.05*u3(k,1)...
+f1(u2(k,1))+0.2*u2(k,2)+g1(u2(floor(k-d(k)),1))+sumh21;
u2(i,2)=0.05*u1(k,2)+0.05*u3(k,2)...
+f2(u2(k,2))+g2(u2(floor(k-d(k)),2))+sumh22;
u3(i,1)=0.05*u1(k,1)+0.05*u2(k,1)...
+f1(u3(k,1))+0.2*u3(k,2)+g1(u3(floor(k-d(k)),1))+sumh31;
u3(i,2)=0.05*u1(k,2)+0.05*u2(k,2)...
+f2(u3(k,2))+g2(u3(floor(k-d(k)),2))+sumh32;
end
i=i+1;
end
%%%%%%%%%%%%%%%%%%%%%%% plot %%%%%%%%%%%%%%%%%%%%%%%%%
n=[-1:h:kfinal];
x11=[55;u1(:,1)];
x12=[-80;u1(:,2)];
x21=[20;u2(:,1)];
x22=[145;u2(:,2)];
x31=[52;u3(:,1)];
x32=[40;u3(:,2)];
figure;
plot(n,x11,'r:.');ylabel('x1');%axis([-1 kfinal -30 60]);
hold on;
plot(n,x21,'b:.');
hold off;
figure;
plot(n,x12,'r:.');ylabel('x2');%axis([-1 kfinal -100 250]);
hold on;
plot(n,x22,'b:.');
hold off;
figure;
plot(n,x11,'r:.');ylabel('x1');%axis([-1 kfinal -30 60]);
hold on;
plot(n,x31,'b:.');
hold off;
figure;
plot(n,x12,'r:.');ylabel('x2');%axis([-1 kfinal -100 250]);
hold on;
plot(n,x32,'b:.');
hold off;
% u1(40:50,:)
% u2(40:50,:)
% u3(40:50,:)
elapse_time=toc % 返回程序的运行时间
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function xini_1=fun1_ini(s)
xini_1(1)=55;
xini_1(2)=-80;
return
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function xini_2=fun2_ini(s)
xini_2(1)=20;
xini_2(2)=145;
return
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function xini_3=fun3_ini(s)
xini_3(1)=52;
xini_3(2)=40;
return
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function f1=f1(x)
f1=tanh(0.2*x)-0.5*x;
return
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function f2=f2(x)
f2=0.95*x-tanh(0.75*x);
return
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function g1=g1(x)
g1=0.2*x-tanh(0.1*x);
return
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function g2=g2(x)
g2=0.1*x;
return
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function h1=h1(x)
h1=0.2*x-tanh(0.1*x);
return
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function h2=h2(x)
h2=0.1*x;
return