x0=2/3;
xf=1;
n=1000;
h=(xf-x0)/n;
x=x0;
len=length(3.9:0.001:4);
q=zeros(1,len);
j=1;
for u=3.9:0.001:4 
i=1;
y=zeros(1,len);
while i<=n
z=u*x*(1-x);
x=x+h;
z1=u*x*(1-x);
f=log(abs(z1-z)/h);
y(i)=f;
i=i+1;
x=z;
end
i=1:n
s=sum(y(i));
p=s/n;
q(j)=p;
j=j+1;
plot(u,q) 
axis([3.9 4 -1 1])
end