function SIRpp
global d p 

set(0,                           ...
   'defaultaxesfontsize', 20,   ...
   'defaultaxeslinewidth', 1.0, ...
   'defaultlinelinewidth', 1.2, ...
   'defaultpatchlinewidth', 0.7);

cmykflg=0;
%if cmykflg=1, this code will make cmyk figure files
%cmykflg = 1;

st = 0.2;
p = (st-1)/log(st)  % p is the single parameter of the problem
 
dlist = [0.8;1.8];
for jj = 1:2
    d=dlist(jj); % d = 1/c;
 
s = [0.1:.01:1.1];
eta = 1+p*log(s) - s;

s0= [st;0.01];
tstep = 0.01;
t_end = 200;
tspan = [0:tstep:t_end];

[T,S] = ode23s(@deRHS,tspan,s0);   
se = 1.4;
figure(1)

plot( S(:,1), S(:,2),'linewidth',2)
hold on
axis([0 se -0.05 0.2])
xlabel('S','fontsize',20)
annotation('textarrow',[0.47,0.57],[0.83 0.75],'linewidth',2)
ylabel('U','fontsize',20)
end
plot(s,eta,'--',[st,1],[0,0],'*',[0,se],[0,0],'--','linewidth',2)

legend('1/c = 0.8','1/c = 1.8','fontsize',20)
text(0.22,-0.01,'S^*','fontsize',20)
hold off

if(cmykflg==1)
    print('../../figs_c/chapt_6/sir_wave_pp','-deps','-cmyk')
end
 
S = [0.001:0.001:0.9999];
eta = (-S)./log(1-S);
figure(2)
plot(eta, S, 'linewidth',2)
xlabel('\eta=\alpha/\beta s_0','fontsize',20)
ylabel('I^* = fraction of Affected Susceptibles','fontsize',20)
function s_prime=deRHS(t,s)
global d p 
v = s(1);
w = s(2);
 
Fv = v*w;
Fw = (1-v+p*log(v)-w)/d;

s_prime = [Fv Fw]';

