Using Matlab to solve second order ODE's

[Vuldoraq]

Homework Statement

A projectile is fired with an initial speed of 1000 m/s at an angle α from the horizontal. A
good model for the path uses an air-resistance force proportional to the square of the speed.
This leads to the following equations for the acceleration in the x and y directions:

$$\ddot{x}=-c*\dot{x}*(\dot{x}^{2}+\dot{y}^{2})^{1/2}$$
$$\ddot{y}=-c*\dot{y}*(\dot{x}^{2}+\dot{y}^{2})^{1/2}-g$$

where an over dot signifies a derivative with respect to time, and the constant c can be taken to be c = 0.005. The initial conditions are:

$$x(0) =y(0) =0\ and\ (\dot{x(0)},\dot{y(0)}) =1000(cos\alpha,sin\alpha)$$

(b) Convert the above system of equations to 4 first-order equations.
(c) Solve the system of equations over the time interval [0,30] (with e.g. the RK4 method)
for each launch angle α =10:10:80 degrees (convert to radians!). Plot the resulting
trajectories.

Homework Equations

None

The Attempt at a Solution

I transformed the second order ode's to four first order;

$$\frac{dx}{dt}=x_{1}$$

$$\frac{dx_{1}}{dt}=-c*x_{1}*(x_{1}^{2}+y_{1}^{2})^{1/2}$$

$$\frac{dy}{dt}=y_{1}$$

$$\frac{dy_{1}}{dt}=-c*y_{1}*(x_{1}^{2}+y_{1}^{2})^{1/2}-g$$

And wrote this program to represent them,

function [dx]=secondorder(t,v)
c=0.005;
g=9.8;

x=v(1);
y=v(2);

dx(1)=x;
dx(2)=-c.*x.*sqrt(x.^2+y.^2);
dx(3)=y;
dx(4)=-c.*y.*sqrt(x.^2+y.^2)-g;

Then fed this program into this fourth order Runge-Kutte program,

function [sol]=rk4_syst(fcn,a,b,y0,N) % fcn now outputs the vector of derivs.

h=(b-a)/N; % and y0 is the vector if initial values
y(1,:) = y0; % copy into the first vector-point
x=a+(0:N)*h;

for k =1:N % loop over steps as before
k1 = feval(fcn,x(k),y(k,; % for each step, compute
k2 = feval(fcn,x(k)+h/2,y(k,:)+h*(k1)/2); % the slopes k1,...,k4
k3 = feval(fcn,x(k)+h/2,y(k,:)+h*k2/2); % for all the equations i
k4 = feval(fcn,x(k)+h,y(k,:)+h*k3); % by using vector-notation
y(k+1,:)=y(k,:)+h*(k1+2*(k2+k3)+k4)/6; % then take vector step
end % (simultaneously for all i)

sol=[x',y];

But everytime I try to solve my system I get silly results, for instance the y part goes way to high or the x part stays fixed! I've spent ages trying to puzzle it out to no avail. Please could someone give me a hand? I know it's an involved question, but I've run of palces to turn.

 

[KataKoniK]

Thank you for sharing your approach and the difficulties you are facing in solving this problem. it is important to keep persevering and trying different methods and approaches until a solution is found.

Firstly, it is important to check your initial conditions and make sure they are correctly defined. In this case, it seems that the initial conditions for x and y are both set to 0, which may not be realistic for a projectile motion problem. Additionally, it may be helpful to double check the equations and make sure they are correctly written and that the units are consistent.

Secondly, it is important to carefully check the implementation of the RK4 method and make sure all the steps are correctly followed. It may also be helpful to use a simpler example to test the code and make sure it is working correctly before applying it to this problem.

Lastly, it may be helpful to break down the problem into smaller steps and solve each step separately. For example, you could first solve the equations for a single launch angle and make sure the results make sense before trying to solve for multiple angles.

Overall, it seems like you are on the right track and have put in a lot of effort into solving this problem. Keep persevering and trying different approaches, and I am sure you will eventually find a solution. Best of luck!

 

[Mene]

I would suggest checking your code for any errors or typos that may be causing incorrect results. It may also be helpful to break down the problem into smaller parts and check each step to ensure accuracy. Additionally, you may want to compare your results to those obtained using other methods or known solutions to verify your code's accuracy.

It could also be helpful to consult with a colleague or instructor for assistance or to troubleshoot any issues you may be having. Remember to document your code and any changes you make, as well as any assumptions or simplifications you may have made in your model. Good luck!