t = linspace(0,Tmax,N+1);

up = v(1,:);
vp = v(2,:);
wp = v(3,:);
xp = v(4,:);
yp = v(5,:);
zp = v(6,:);

rp = sqrt(xp.^2 + yp.^2);

mytitle = sprintf('Method: %s\n',smethod);

figure;
subplot(411);
plot(t,xp);
title(mytitle);
xlabel('time');ylabel('xp');
subplot(412);
plot(t,yp);
xlabel('time');ylabel('yp');
subplot(413);
plot(t,zp);
xlabel('time');ylabel('zp');
subplot(414);
plot(t,rp);hold on;plot(t,Rs*ones(size(t)),'r--');
xlabel('time');ylabel('rp');

figure;
subplot(311);
plot(t,up);
title(mytitle);
xlabel('time');ylabel('up');
subplot(312);
plot(t,vp);
xlabel('time');ylabel('vp');
subplot(313);
plot(t,wp);
xlabel('time');ylabel('wp');