Jakob Svenningsson
Jakob Svenningsson

Reputation: 4305

matlab - 2 second degree ODE plot with runge kutta - numerical methods

The following equation is given.

enter image description here

I know that i need to break the 2 second order ODE's into 4 first order ODE's. This i what i've got. Ive first introduced the new variable u and in the bottom of the picture i've written my matlab function that i use with ODE45.

enter image description here

Now the problem is that im supposed to get a parabola shaped figure(blue line) but that is not what I am getting.

enter image description here

I've gone through my code a thousand times with no results. Can detect any errors in my function?

main program

global g H R alfa 
alfa=pi/2;
g = 20.0;
R = 1;
H=2.3;
k = 0:0.01:2;
[T,Y] = ode45(@fspace,k,[H 0 0 0]);
plot(T,Y(:,1))
hold on 


fi= 0:2*pi/60:2*pi;
xx =R*cos(fi);
yy =R*sin(fi);
plot(xx,yy)

function f

function f = fspace(x,u)
global g R H alfa G
G=(g*R.^2)./((R+H).^2);
f = [u(2) G*cos(alfa)-g*((R.^2)/u(1).^2)+u(1)*u(4)^2 u(4) (G*sin(alfa)-2*u(2)*u(4))/u(1)];

Upvotes: 2

Views: 971

Answers (1)

am304
am304

Reputation: 13876

I think your problem lies with those lines:

fi= 0:2*pi/60:2*pi;
xx =R*cos(fi);
yy =R*sin(fi);
plot(xx,yy)

You are plotting a circle of radius R. phi is part of the solution from the ode solver, so you should have instead:

plot(R*cos(Y(:,3)),R*sin(Y(:,3)))

but that will always give you a circle of radius R, never a parabola. Or is the parabola meant to refer to plot(T,Y(:,1)).

The equations and the code appear to be correct as far as I can see. Replacing your definition of phi by Y(:,3) gives essentially the same plot, except that the resolution is less. As I said, you'll always get a circle by plotting yy vs xx. You need to clarify what the parabola should refer to.

enter image description here

Upvotes: 1

Related Questions