Discover MakerZone

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

Learn more

Discover what MATLAB® can do for your career.

Opportunities for recent engineering grads.

Apply Today

To resolve issues starting MATLAB on Mac OS X 10.10 (Yosemite) visit: http://www.mathworks.com/matlabcentral/answers/159016

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.

Tae Yeong Kim

Products

No products are associated with this question.

4 Answers

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.

8 Comments

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.

Image Analyst
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!

Youssef  KHMOU
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 .

Advice :

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

Youssef  KHMOU
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 ?

4 Comments

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?

Youssef  KHMOU

Contact us