[t,y]=ode45(@duffing3,[0 1000],[1 1]);
[t,z]=ode45(@duffing4,[0 1000],[1 1]);
figure(1);
subplot(1,2,1),plot(y(:,1),y(:,2)),xlabel('x'),ylabel('y');title('duffing r=0.8 混沌 tspan=1000');
subplot(1,2,2),plot(z(:,1),z(:,2));xlabel('x'),ylabel('y'),title('duffing r=0.82 臨界混沌');