%
%Integrates two coupled  chua cisrcuits with synchronization
%


format compact

options= odeset('Reltol',1.e-5);
%options= odeset('Reltol',1.e-22,'OutputFcn','odephas3','Refine',222);
%options= odeset('Reltol',1.e-3,'OutputFcn','odephas3','Refine',2);
%options= odeset('Reltol',1.e-3);

 trans=0*100
 tf=100;
 fac=5;
 c3=25.58;
 DC = 0.1;
 CR = 0.48/2;

 N=3;
 S = 0.1;
 x0 = 0.1+S*(2*(rand(1,1)-0.5));
 y0 = 0.0+S*(2*(rand(1,1)-0.5));
 z0 = 1.0+S*(2*(rand(1,1)-0.5));
 v0 = [x0 y0 z0];
 x0 = 0.1+S*(2*(rand(1,1)-0.5));
 y0 = 0.0+S*(2*(rand(1,1)-0.5));
 z0 = 1.0+S*(2*(rand(1,1)-0.5));
 v0 = [v0 x0 y0 z0]
 tspan=[0:tf/2000:tf];

 v01=v0; v02=v0; v03=v0;
 if(trans)
  [T1,Y1]=ode45('chua_synchro',trans,v01,options,c3,CR-DC);
  [T2,Y2]=ode45('chua_synchro',trans,v02,options,c3,CR);
  [T3,Y3]=ode45('chua_synchro',trans,v03,options,c3,CR+DC);
  v01 = Y1(end,:); v02 = Y2(end,:); v03 = Y3(end,:);
 end

 [T1,Y1]=ode45('chua_synchro',tspan,v01,options,c3,CR-DC);
 [T2,Y2]=ode45('chua_synchro',tspan,v02,options,c3,CR);
 [T3,Y3]=ode45('chua_synchro',tspan,v03,options,c3,CR+DC);

 ifac=floor(length(tspan)/fac);
 Tfac = T1(1:ifac);
 
 figure(1);
 clf;
  subplot(3,3,1)
  plot(Y1(1:ifac,1),Y1(1:ifac,2),Y1(1:ifac,N+1),Y1(1:ifac,N+2))
  hold on
  plot(Y1(1,1),Y1(1,2),'*',Y1(1,N+1),Y1(1,N+2),'*')
  plot(Y1(ifac,1),Y1(ifac,2),'^',Y1(ifac,N+1),Y1(ifac,N+2),'^')
  axis tight; ylabel('y');
  subplot(3,3,4)
  plot(Y2(1:ifac,1),Y2(1:ifac,2),Y2(1:ifac,N+1),Y2(1:ifac,N+2))
  hold on
  plot(Y2(1,1),Y2(1,2),'*',Y2(1,N+1),Y2(1,N+2),'*')
  plot(Y2(ifac,1),Y2(ifac,2),'^',Y2(ifac,N+1),Y2(ifac,N+2),'^')
  axis tight; ylabel('y');
  subplot(3,3,7);
  plot(Y3(1:ifac,1),Y3(1:ifac,2),Y3(1:ifac,N+1),Y3(1:ifac,N+2))
  hold on
  plot(Y3(1,1),Y3(1,2),'*',Y3(1,N+1),Y3(1,N+2),'*')
  plot(Y3(ifac,1),Y3(ifac,2),'^',Y3(ifac,N+1),Y3(ifac,N+2),'^')
  axis tight; xlabel('x');ylabel('y');

  subplot(3,3,2)
  plot(Y1(1:ifac,1:N),Y1(1:ifac,N+1:2*N))
  hold on
  plot(Y1(1,1:N),Y1(1,N+1:2*N),'*')
  plot(Y1(ifac,1:N),Y1(ifac,N+1:2*N),'^')
  axis tight; ylabel('x2,y2,z2');
  subplot(3,3,5)
  plot(Y2(1:ifac,1:N),Y2(1:ifac,N+1:2*N))
  hold on
  plot(Y2(1,1:N),Y2(1,N+1:2*N),'*')
  plot(Y2(ifac,1:N),Y2(ifac,N+1:2*N),'^')
  axis tight; ylabel('x2,y2,z2');
  subplot(3,3,8)
  plot(Y3(1:ifac,1:N),Y3(1:ifac,N+1:2*N))
  hold on
  plot(Y3(1,1:N),Y3(1,N+1:2*N),'*')
  plot(Y3(ifac,1:N),Y3(ifac,N+1:2*N),'^')
  axis tight; ylabel('x2,y2,z2'); xlabel('x1,y1,z1');

  subplot(3,3,3)
  semilogy(T1,sqrt(sum((Y1(:,1:N)-Y1(:,N+1:2*N)).^2,2)))
  axis tight; ylabel('||x-y||');
  subplot(3,3,6)
  plot(T2,sqrt(sum((Y2(:,1:N)-Y2(:,N+1:2*N)).^2,2)))
  axis tight; ylabel('||x-y||');
  subplot(3,3,9)
  semilogy(T3,sqrt(sum((Y3(:,1:N)-Y3(:,N+1:2*N)).^2,2)))
  axis tight; ylabel('||x-y||');xlabel('t')

     
