# Two second order differential equation

1. Oct 12, 2009

### people

hi, everybody
i got for homework a pendulum on a cart.

i solved the system and got two equations.

(m+M)x'' + M L theta'' = F

x'' + L theta'' + g* theta = 0

F = 1000 N, m = 500 kg, M = 1250 kg, l = 10m

i know how to put them in state space and solve them with SS block in
simulink.but i don't know how to solve the system with ode23. there are
no initial conditions stated in a task, i guess there are zero.

i understand this sintax
function xdot = vdpol(t,x)
xdot = [x(1).*(1-x(2).^2)-x(2); x(1)]

t0 = 0; tf = 20;
x0 = [0 0.25]'; % Initial conditions
[t,x] = ode23('vdpol',t0,tf,x0);
plot(t,x)

but dont know how to use it on those eq.

2. Oct 13, 2009

### matematikawan

Do you know how to solve higher order differential equation using Runge-Kutta method?
Change your system of equations to the following form first before you use ode23

$$\frac{d\vec{y}}{dt} = \vec{f}(t,\vec{y})$$

where
$$\vec{y} = ( x , \dot{x} , \theta , \dot{\theta})^t$$

3. Oct 14, 2009

### people

i dont know if that is the method i am using but i know how to transform
this second order to first order eq.
and i got.
$$\ddot{\theta} +\frac{F}{m*l}+ \frac{(m+M)*g}{l*m}\theta=0$$
$$\ddot{x}-\frac{F}{m}-\frac{g*M}{m}\theta=0$$

$$x_{1} = \dot{x} -> \dot{x}_{1}= \ddot{x}$$
$$x_{2} = x -> \dot{x}_{2} = \dot{x} = x_{1}$$
$$\theta_{1} = \dot{\theta} -> \dot{\theta}_{1}= \ddot{\theta}$$
$$x_{2} = \theta -> \dot{\theta}_{2} = \dot{\theta} = \theta_{1}$$

and now i transform in
$$\dot{\theta}_{1} = -\frac{F}{m*l}-\frac{(m+M)*g}{l*m}\theta_{2}$$
$$\dot{\theta}_{2} = \theta_{1}$$
$$\dot{x}_{1} = \frac{F}{m}+\frac{g*M}{m}\theta_{2}$$
$$\dot{x}_{2} = x_{1}$$

and now i know how to put that in matlab
thetadot = [-0.2 - theta(2).*3.4335; theta(1)] % here are numbers because matlab
% cant "see" the variables F,m,M...
t0 = 0;
tf = 10;
theta0 = [0 0]';
[t,theta] = ode23('thetam',t0,tf,theta0);
plot(t,theta)
plot(t, theta(84:1:166))

and i get two graphs. one of them is the on i need.
can you explain me what is the second one, please.

but i dont know how to write in matlab solution for x
becase it has a theta in it eq??

4. Oct 14, 2009

### matematikawan

try checking the size for theta. I think it's size is n by 2. That's why you have two graphs I think. It is just like using the command plot(t,theta(:,1),t,theta(:,2)).

The equation that you solve with matlab actually can be solve analytically because it is an inhomogeneous linear DE with constant coefficients.

5. Oct 14, 2009

### people

yes, i know it has two columns. what is in a first? also solution to the eq?
problem is i dont know how to solve it for x... can you say exactly how
because i have no clue

6. Oct 14, 2009

### matematikawan

My guess is that your second graph is t against theta. You must have known better. You wrote the program.

Try using the substitution that I wrote in my first post.

$$y_1= x , y_2=\dot{x},y_3=\theta , y_4=\dot{\theta}$$

7. Oct 14, 2009

### people

i am sorry, but i cant see how is your's supstitution different form
mine when i wrote what i did to that point

8. Oct 14, 2009

### matematikawan

Not much different actually . But you use subscripts for both x and theta that confuse me. It looks like you want to solve the two equations separately. What I have in mind is solving the two equations simultaneously using the four variables y1, y2, y3, and y4.

Then the plots for your solution are plot(t,y(:,1)) for x and plot(t,y(:,3)) for theta.

9. Oct 21, 2009

### people

i forgot to write the solution... if anyone was interested or have similar problem
$$y_{1} = \dot{x}$$
$$y_{2} = x$$
$$y_{3} = \dot{\theta}$$
$$y_{4} = \theta$$
$$\dot{y}_{1} = \ddot{x}$$
$$\dot{y}_{2} = \dot{x} = y_{1}$$
$$\dot{y}_{3} = \ddot{\theta}$$
$$\dot{y}_{4} = \dot{\theta} = y_{3}$$

$$\dot{y}_{1} = \frac{F}{m} + \frac{gM}{m}y_{4}$$
$$\dot{y}_{3} = -\frac{F}{ml} - \frac{(m+M)g}{ml}y_{4}$$
$$\dot{y}_{4} = y_{3}$$
$$\dot{y}_{2} = y_{1}$$

and when we put it in matlab
_________________________
function ydot = ydot(t, y)
ydot = [2 + y(4).*24.525; y(1); -0.2 - y(4).*3.4335; y(3)]
%this is ydot.m file
_________________________
M = 1250;
m=500;
l = 10;
g = 9.81;
F = 1000;
t0 = 0;
tf = 10;
poc0 = [0 0 0 0]';
[t,rj]= ode23('xm', t0, tf, poc0)
plot(t,rj)

problem is similar to http://www.myphysicslab.com/pendulum_cart.html#indirect", except force is constant in my problem.

Last edited by a moderator: Apr 24, 2017