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

# Problem with the "Plot" function

Asked by Tae Yeong Kim 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?

## 1 Comment

Jan Simon on 28 Jan 2013

Please format the code properly. Thanks.

## Products

No products are associated with this question.

Answer by Image Analyst on 27 Jan 2013
Edited by Image Analyst 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 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 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 on 28 Jan 2013

Nope - I didn't even run it.

Answer by Youssef KHMOU 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')

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

## 1 Comment

Tae Yeong Kim 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!

Answer by Youssef KHMOU on 28 Jan 2013
Edited by Youssef KHMOU 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,

## 1 Comment

Tae Yeong Kim 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

Answer by Youssef KHMOU 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 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 on 29 Jan 2013

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

Tae Yeong Kim 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?