Plotting an ODE with an event detection function

I am trying to get a plot of an ODE travelling over this surface (the contour plot and circle aren't important but necessary) and I am trying to find the periodic points.

When I am trying to plot though I just get the contour plot and the circle, I don't get the line the ODE creates, any help would be great.

**You should see the line start on the x-axis go up slightly and then come back down on itself, and stop again once it hits the x-axis

clc

tspan = [0:0.05:80]; %Range for independent var. i.e. time
r = 2;
centreX0 = 0;
centreY0 = 0;

y = 0;
x = -0.31417;
Px = 0;
%Total Energy = Kinetic Energy + Potential Energy
%KE = (1 ./ 2) .* (Px.^2 + Py.^2);
%V = (x.^3 - 3 .* x .* y^2);
%Total Energy = 1
Py = sqrt((1 - x.^3) .* 0.5);

y0 = [x y Px Py];

option = @isYnegative;
[x, y] = ode45(@ODEdbb,tspan,y0,option);

[m,n] = meshgrid(-3:.05:3,-3:.05:3);
contour(m,n,m.^3-3*n.^2.*m);

hold on

viscircles([centreX0,centreY0,],r);

axis square


hold off

function f = ODEdbb(A, B)
    x=B(1);
    y=B(2);
    Px=B(3);
    Py=B(4);
    
    %Differential Equations
    dPxdt = 3 * y^2 -3 * x^2;
    dPydt = 6 * x * y;
    dxdt = Px;
    dydt = Py;
    f = [dxdt; dydt; dPxdt; dPydt];
end


function value = isYnegative(Y)
    value = Y(2) == 0;
    isterminal = 1;
    direction = -1;
end

[Torsten]

Since you don't refer to the ODE solution in your plot, it cannot appear therein.

Answer 1

[Walter Roberson]

y0 = [x y Px Py];

[x, y] = ode45(@ODEdbb,tspan,y0,option);

The first output from ode45() is the time values. The second output from ode45() is an array that is the value of the boundary conditions corresponding to each time.

Your code appears to assume that because you are thinking of x and y as your first two input coordinates, that the first two outputs are x and y. That is not correct at all. You should be using something like

[t, out] = ode45(@ODEdbb,tspan,y0,option);
x = out(:,1);
y = out(:,2);

after which you can

plot(x, y)