%
%   ResetOscillator.m
%
%  Reset oscillator
%  with chaotic behaviour
%
%  J. Lunze
%  16.6.2008
%  for ESAufl4: 1.8.2025
%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
close all
clear
%
%    Model:   xdot = A*x + e
%
delta=0.1;                 %    0 < delta < 1
A=[0 1; -1 2*delta];
e=[0; 1];
x0=[0.2; 0];
Tmax=10;
dt=0.001;                   % simulation time step
kmax=Tmax/dt;
%
xbar=-inv(A)*e;
dx0=x0-xbar;                % dx - deviation from xbar
%
%  Reset rule
%  if x_2=0 and x_1 < 0   then  x_1:=-x_1
eAdt=expm(A*dt);
dx(1,:)=dx0;
x(1,:)=x0;
for k=2:kmax
    dx(k,:)=(eAdt*dx(k-1,:)')';
    x(k,:)=dx(k,:)+xbar';
    if abs(x(k,2))<0.001
        if x(k,1) < 0
           x(k,1)=-x(k,1);
           dx(k,1)=x(k,1)-xbar(1);
        end
    end
end
plot(x(:,1), x(:,2), 'b')
axis([-1 2.5 -1.5 1]);
axis('square');
hold on
plot(x0(1,1), x0(2,1), 'bo');
text(0.5, -1.4, '\leftarrow');
plot(-0.504, 0, 'bd');
plot(0.504, 0, 'bd');
latexxlabel('$$x_1$$');
latexylabel('$$x_2$$');
rotateY;
yticks([-1 0 1]);
xticks([-1 0 1 2]);
rotateY;
%
%  longer trajectory
clear x dx;
figure(2)
Tmax=100;
dt=0.01;                   % simulation time step
eAdt=expm(A*dt);
kmax=Tmax/dt+1;
T=[0:dt:Tmax];
dx(1,:)=dx0;
x(1,:)=x0;
kjump=0;
for k=2:kmax
    dx(k,:)=(eAdt*dx(k-1,:)')';
    x(k,:)=dx(k,:)+xbar';
    if abs(x(k,2))<0.01
        if x(k,1) < 0
           x(k,1)=-x(k,1);
           dx(k,1)=x(k,1)-xbar(1);
           kjump=kjump+1;
           xjump(kjump)=x(k,1);
           Tjump(kjump)=(k-1)*dt;
        end
    end
end
plot(x(:,1), x(:,2), 'b')
axis([-1 2.5 -1.5 1]);
axis('square');
latexxlabel('$$x_1$$');
latexylabel('$$x_2$$');
rotateY;
yticks([-1 0 1]);
xticks([-1 0 1 2]);
rotateY;
%
figure(3)
subplot(6,1,[1,2])
plot(T, x(:,1), 'b');
axis([0 100 -1 2.3]);
xleer;
xticksempty([0 50 100 150]);
latexylabel('$$x_1$$');
rotateY;
yticks([-1 0 1 2]);
%
subplot(6,1,[3,4])
plot(Tjump, xjump,'bo');
axis([0 100 0 1]);
latexxlabel('$$t$$');
latexylabel('$$x_1$$');
rotateY;
yticks([-1 0 1]);
xticks([0 50 100 150]);
%
%
%
epsfigc('Resetoscillator');