function x = logmap(R, timesteps, xinit, color)

% logmap (R, timesteps, xinit) plots the x versus timesteps starting with 
% x(1) = xinit in the equation x(t+1) = R*x(t)*(1-x(t))
% Try logmap(2.5, 50, 0.2), logmap(3.5, 50, 0.2), logmap (3.9, 50, 0.2)
% logmap (3.9, 50, 0.20000001)
% returns x, a list of values of x at all timesteps

%Melanie Moses, 10/23/2012, based on Eqns from M. Mitchell Complexity,
%Oxford 2009

clc
x = zeros(timesteps,1);
x(1) = xinit;

for t = 1:length(x)-1
    x(t+1) = R*x(t)*(1-x(t));
end


plot(x, [color,'.-'])
xlabel('time (t)','fontsize', 18)
ylabel('x(t)', 'fontsize', 18)
title (['Logistic Map R = ',num2str(R)], 'fontsize',22)

%In class demos
%MEM 1-17-2013

%increase R
logmap(2,10,0.3, 'b')
logmap(2,10,0.7, 'b')

logmap(3,10,0.7, 'b')
logmap(3.2,20,0.4, 'b')

logmap(3.5,40,0.4, 'b') %period 4
logmap(3.565,100,0.4, 'b')%period 16

logmap(3.6,200,0.4, 'b')
logmap(4,200,0.4, 'b')

%sensitive dependence on initial conditions
logmap(3.2,60,0.7, 'b')
hold on
logmap(3.2,60,0.71, 'r')
hold off

%SDIC in the chaotic regime
logmap(3.8,60,0.7, 'b')
hold on
logmap(3.8,60,0.71, 'r')
hold off

%rounding errors

%TRY R = 4, x_0 = 0.3
%Round x_t+1 to 2 decimal points vs 1 decimal point
%how long before they diverge?



%Shaddowing Lemma
A = logmap(4,100,.0396,'r');
B = logmap(4,100,.039,'b');

hold on
plot ([1:100],A,'r.-')
plot ([1:100],B, 'b.-')
hold off

max(A)
max(B)
min(A)
min(B)
figure()
hist(A)
figure()
hist(B)

%Eventhough the trajectories diverge, statistical properties of 
%timeseries A & B are similar. 
%Despite SDIC, the logistic map characterizes the dynamics


%Fiegenbaum
R2 =3.44949  
R3 =3.54409  
R4 =3.564407

(R3-R2) ./ (R4-R3)

% the dynamics.