## Problem with the "Plot" function

### Tae Yeong Kim (view profile)

on 27 Jan 2013

Hi, I need some help from those of you who are good at Matlab sinc I am not a profession.. I appreciate the time and effort.

I need to plot my time in x axis and my roll angle(x1) as y axis

----This is my parameters------- function dx = wheelchaircorrectversion (t,x) global m g d F h Ix w m=10; g=-9.81; F=5; d=10; h=20; Ix=30; w=20;

if t>0 && t>=0.01 F=5 elseif t>0.01 F=0 end

z=(sqrt((.5*w).^2+(.5*h).^2)).*(sin(atan((h/2)/w/2))); d=(sqrt((.5*w).^2+(.5*h).^2)).*(cos(atan((h/2)/w/2)));

dx=zeros(2,1); dx(1)=x(2); dx(2)=(-(m*g*d)+(F+z))/Ix; end

---- my main function---

clear all clc global m g d F h Ix w m=10; g=-9.81; F=5; d=10; h=20; Ix=30; w=20;

x1_ini = [0 0]'; t_dim = 1; i=1; dt = 0.01; duration = 0.05;

for t = 0.0 : dt : duration t0 = t; i=i+1; tf = (t+dt)*(1+eps); tspan = [t0 tf];

```    [t1,x1] = ode45('wheelchaircorrectversion',tspan,x1_ini);
n1 = length(x1); %length of a vector
tmp1 = x1(n1,:); %the last value
tmp2 = t1(t_dim)```

x1_ini = tmp1; % update initial conditions

```    traj_p0(i,:) = x1_ini;
traj_t(i)= tmp2
end```

subplot(211) %Creates axes in tiled positions plot(traj_t,x1(:,1)) xlabel('Time') ylabel('Roll Angle')

subplot(212) plot(traj_t,x1(:,2)) xlabel('Time') ylabel('Roll Velocity')

----

Whenever I try to plot these functions, it says that vectors needs to be in same length for plot(traj_t, x1(;,1)) and for plot(traj_t,x1(:,2))

I have tried everything to make the number of vectors equal to each other, but it didnt work..

could someone help?

Jan Simon

### Jan Simon (view profile)

on 28 Jan 2013

Please format the code properly. Thanks.

## Products

No products are associated with this question.

### Image Analyst (view profile)

on 27 Jan 2013
Edited by Image Analyst

### Image Analyst (view profile)

on 27 Jan 2013

More accurately, x1(:,1) needs to be the same number of elements at traj_t. x1(:,2) needs to be that length also. Do this before you call plot

```size(traj_t)
size(x1)
```

Tell us what it says.

Image Analyst

### Image Analyst (view profile)

on 27 Jan 2013

Set a breakpoint and type this into the command window:

```t = 0.0 : dt : duration
```

See how many elements it is. That is the loop where you're setting traj_t so that is how many elements traj_t will have. If you don't like it, you have to adjust dt or duration.

Tae Yeong Kim

### Tae Yeong Kim (view profile)

on 27 Jan 2013

Wow, thanks for your help. I finalized it, and it is working. However, I'm not sure if it's physically making sense. Do you have any thoughts about my code?

Image Analyst

### Image Analyst (view profile)

on 28 Jan 2013

Nope - I didn't even run it.

### Youssef Khmou (view profile)

on 27 Jan 2013

Try this :

1.Delete the traj_t variable that exists inside the loop.

2. create a vector time, outside the loop, of length 41 such as : time=linspace(0,duration,41).

3. replace,in the plot section, traj_t wit time'.

The code now works well, but are the results reasonable In Physic viewpoint ?

As you posted non organized code , i post the new version of the code here for further discussions by other users :

--------------------------------------------------------------------------

function dx = wheelchaircorrectversion (t,x)

global m g d F h Ix w

m=10; g=-9.81; F=5; d=10; h=20; Ix=30; w=20;

if t>0 && t>=0.01

`    F=5;`

elseif t>0.01

`    F=0; `

end

z=(sqrt((.5*w).^2+(.5*h).^2)).*(sin(atan((h/2)/w/2)));

d=(sqrt((.5*w).^2+(.5*h).^2)).*(cos(atan((h/2)/w/2)));

dx=zeros(2,1); dx(1)=x(2); dx(2)=(-(m*g*d)+(F+z))/Ix;

end -------------------------------------------------------------------

And the main program :

-------------------------------------------------------------------

clear all ,clc

global m g d F h Ix w

m=10;g=-9.81; F=5; d=10; h=20; Ix=30; w=20;

x1_ini = [0 0]'; t_dim = 1; i=1; dt = 0.01; duration = 0.05;

for t = 0 : dt : duration

`    t0 = t;i=i+1; tf = (t+dt)*(1+eps); tspan = [t0 tf];`
`    [t1,x1] = ode45('wheelchaircorrectversion',tspan,x1_ini);`

n1 = length(x1); %length of a vector

tmp1 = x1(n1,:); %the last value

`    tmp2 = t1(t_dim);`
`    x1_ini = tmp1; % update initial conditions`
`    traj_p0(i,:) = x1_ini;`
`    %traj_t(i)= tmp2;`

end

time=linspace(0,duration,41);time=time';

subplot(211), plot(time,x1(:,1)), xlabel('Time'), ylabel('Roll Angle')

subplot(212), plot(time,x1(:,2)), xlabel('Time'), ylabel('Roll Velocity')

--------------------------------------------------------------------------

Tae Yeong Kim

### Tae Yeong Kim (view profile)

on 27 Jan 2013

Thank you for your response. I have one more question about my code: for the last part,subplot(211), plot(time,x1(:,1)), xlabel('Time'), ylabel('Roll Angle'), I am trying to plot the function before it went through ODE45 so that I get un-integrated (original) function to get roll velocity on subplot 212.

I am not sure how I have to write out the function..

Thanks a lot!

### Youssef Khmou (view profile)

on 28 Jan 2013
Edited by Youssef Khmou

### Youssef Khmou (view profile)

on 28 Jan 2013

Tae Yeong Kim,

so its about ordinary differential equation ! THE TIME IS ALREADY THERE !!

[t1,x1]=ode('wheelchaircorrectversion',tspan,x1_ini);

Position and velocity are in the 41x2 matrix x1 and the time axis is t1(41,1)

you can replace the vector "time" that we created with t1, ITS THE SAME .

About your last question i do not understand what you mean : you have a position X and Velocity V , you said you want to integrate the x1(:,1) to get the Velocity, i think its wrong :

YOU HAVE POSITION X , VELOCITY V and ACCELERATION A

DX/Dt= V(t) , DV/Dt= A, INTEG(V(t))=X(t), INTEG(A)= V(t)

so lets see :

---------------------------------------------------------------------------

plot(t1,x1(:,2)) % This is the velocity

hold on

Velocity2=diff(x(:,1))./diff(t1); % Diff function truncates to n-1

Velocity2(end+1)=Velocity2(end); % add the last element

plot(t1,Veclocity2,'r'),

-----------------------------------------------------------------------

Numerically its the same .

1.do not declare the variables again the testing script .

2. try to comment every variable to explain to the reader , for example :

------------------------------------------------------------------------

m=10; % mass in Kg

g=-9.81; % Gravity acc in m/s²

F=5; % Force in Newton

d=10; % Distance between ..and .. in meteres....

I hope that helps

KHMOU Youssef,

Tae Yeong Kim

### Tae Yeong Kim (view profile)

on 28 Jan 2013

Thank you so much. It answered everything I needed.

If you could check out my equations and see if they make sense, it would be awesome.

I am not sure if variable z and d are correctly derived.. Thanks

### Youssef Khmou (view profile)

on 28 Jan 2013

Hi

i think there is a problem, based on what you wrote :

Z/d=h/w. and based on numerical values u gave h=w=20. in that case z=d which is wrong ,

what is the point in the center of the rectangle? what is its dynamic ?

Tae Yeong Kim

### Tae Yeong Kim (view profile)

on 29 Jan 2013

Oh, and I forgot to tell you that z and d changes as it the object is moving from the effect of the force.

So they are variables, but I do not know the right equations for them if w and h are not equal to each other.

Jan Simon

### Jan Simon (view profile)

on 29 Jan 2013

Please, Tae Yeong Kim, take the time to format your code. Currently it looks very messy.

Tae Yeong Kim

### Tae Yeong Kim (view profile)

on 2 Feb 2013
```clear all
clc
global m g F h Ix w
m=10;
g=9.81;
F=5;
h=20;
Ix=30;
w=20;
```
```x1_ini = [0 0]';
i=1;
dt = 0.01;
duration = 0.05;
```
```for t = 0.0 : dt : duration
t0 = t;
i=i+1;
tf = (t+dt)*(1+eps);
tspan = [t0 tf];
```
```      [t1,x1] = ode45('wheelchaircorrectversion',tspan,x1_ini);
n1 = length(x1); %length of a vector
tmp1 = x1(n1,:); %the last value```
```      n2 = length(t1);
tmp2 = t1(n2);```
`      x1_ini = tmp1; % update initial conditions`
```      traj_p0(i,:) = x1_ini;
traj_t(i)= tmp2;
end```
```timetraj = traj_t';
timetraj = linspace(0,duration,41)% Generates "duration" points between X1 and X2
```
```subplot(211) %Creates axes in tiled positions
plot(timetraj,x1(:,1))
xlabel('Time')
ylabel('Roll Angle')
```
```subplot(212)
plot(timetraj,x1(:,2))
xlabel('Time')
ylabel('Roll Velocity')
----------------------------------------------------------------------
```

Does this look any better?

MATLAB and Simulink resources for Arduino, LEGO, and Raspberry Pi test