MATLAB - Diffekvationsproblem

Hej,
Jag skulle behöva lite hjälp att lösa följande matlab-uppgift:

En långsträckt cylinder med radien R = 2 befinner sig i en inkompressibel vätska som strömmar i positiv x-riktning. Cylinderaxeln (x = 0) är vinkelrät mot flödesriktningen. Det hela kan betraktas som ett 2D-problem. Flödespartikelns läge (x(t), y(t)) vid tiden t bestäms av startpositionen och av differentialekvationerna

dx/dt=1-(R^2 (x^2-y^2 ))/(x^2+y^2 )^2
dy/dt=-(2xyR^2)/(x^2+y^2 )^2

Vid t=0 befinner sig fyra partiklar vid x=-4 med y-positionerna 0.2, 0.6, 1.0 och 1.4. Beräkna och rita deras strömningskurvor fram till tiden t=12. Notera läget för de fyra partiklarna vid denna tidpunkt. Den understa partikeln har hamnat på efterkälken. Beräkna med en effektiv algoritm hur lång tid som krävs för att den ska nå fram till samma x-position som den översta har vid t=12.

Jag lyckas lösa första delen genom att använda mig av function och ode45.

Jag har en funktionsfil, flode:

function dudt = flode(t,u)
x = u(1);
y = u(2);
R=2; %cylinderns radie
dxdt = 1 - (R^2*(x^2 - y^2))/(x^2 + y^2)^2;
dydt = -(2*x*y*R^2)/(x^2 + y^2)^2;
dudt=[dxdt; dydt];

Som jag anropar i filen:

for y0=[0.2 0.6 1.0 1.4];
[T,U]=ode45('flode',[0 12],[-4 y0]);
x=U(:,1);
y=U(:,2);

plot(x,y)
axis('equal')
hold on
xlabel('x')
ylabel('y')
end

Men sedan kommer jag inte på en bra lösning med hur jag snyggt plockar ut de olika x-värdena för partiklarnas slutpositioner. Jag plockade ut de här värdena ”för hand” (skrev ut x och y i for-slingan) i stället, men jag lyckas heller inte komma på en bra algoritm för att lösa sista uppgiften. Jag tycker att det borde gå att använda en while-slinga, började så här (men det här funkar ju inte...):

x=4.647214925231829; % sista x-värdet på partikel som kommit kortast.
y0=0.181156547069939; % sista y-värdet på partikeln som kommit kortast.
xend=7.851588396552807; % sista x-värdet för partikeln som kommit längst.
t=0;

while x<=xend
tend=t+0.2;
[T,U]=ode45('flode2',[t tend],[x y0]);
x=U(:,1);
y=U(:,2);
plot(x,y)
axis('equal')
hold on
xlabel('x')
ylabel('y')
end

Sedan undrar jag även hur jag kan göra för att lösa diffekvationsystemet en tidsetapp i taget (för en vidare uppgift)?

Jag är tacksam för hjälp!