How to get the time vector from the function

Question

Mohamed Aburakhis 2016-3-13

0
链接

此问题的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/273329-how-to-get-the-time-vector-from-the-function

回答： Jan 2016-3-14

Dear all

I have a program the uses two function but I couldn't understand how to get the time to plot the system response. I saved t inside the function but it gives me a different length. So please any one have an idea The code is attached

Many thanks

close all;clear all; clc
global u_save t_save y_save  a_u eta_u segma_r eta_z eta_w a_w vhat;
tf=10;
x1_0=1;
x2_0=3;
uhat_0=2;
a_u=0.232;
eta_u=0.0866;
p11_0=3;
p21_0=0;
p12_0=0;
p22_0=3;
segma_r=0.000001;
eta_z=0.822;
zhat1_0=1;
zhat2_0=1;
eta_w=0.8;
a_w=0.2;
vhat=0;
time=[0:1:tf];       
x0=[x1_0;x2_0;uhat_0;p11_0;p21_0;p12_0;p22_0;zhat1_0;zhat2_0];
[t,x]=ode45('eqp1',time,x0);
figure(1)
plot(t_save,y_save);
xlabel('time (sec)');
ylabel('y');
%legend('y','')
grid on
figure(2)
plot(t_save,u_save);
xlabel('time (sec)');
ylabel('u_hat');
grid on
%%%%%%%%%%%%%%%this is the function
  function [x_dot] = eqp1(t, x)
t        
global u_save t_save y_save  a_u eta_u segma_r eta_z eta_w a_w vhat;
x1=x(1);
x2=x(2);
uhat=x(3);
p11=x(4);
p12=x(5);
p21=x(6);
p22=x(7);
zhat1=x(8);
zhat2=x(9);
u=uhat+a_w*sin(eta_w*t);
y_save=[y_save x2];
u_save=[u_save u];
t_save=[t_save t];
x_dot(1)=-x1+u;
x_dot(2)=-x2+0.5*x1^2;
x_dot(3)=-(a_u*eta_u*zhat2)/(eta_u+a_u*abs(zhat2));
x_dot(4)=eta_z*p11+1*(p12+p21)-eta_z*(p11^2+segma_r*p12*p21);
x_dot(5)=eta_z*p12+1*p22-eta_z*(p11*p12+segma_r*p12*p22);
x_dot(6)=eta_z*p21+1*p22-eta_z*(p11*p21+segma_r*p21*p22);
x_dot(7)=eta_z*p22-eta_z*(p12*p21^2+segma_r*p22^2);
x_dot(8)=(0-eta_z*segma_r*p12)*zhat2+eta_z*p11*(x2-zhat1-vhat);
x_dot(9)=(-eta_z*segma_r*p22)*zhat2+eta_z*p21*(x2-zhat1-vhat);
% x_dot(10)
x_dot=x_dot';

1 个评论
显示 -1更早的评论隐藏 -1更早的评论

Jan 2016-3-14

Please use the "{} Code" button top format your code. I've done this for you this time.

请先登录，再进行评论。

请先登录，再回答此问题。

Answer 1

Star Strider 2016-3-14

0
链接

此回答的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/273329-how-to-get-the-time-vector-from-the-function#answer_213481

在 MATLAB Online 中打开

I can’t understand what you want to do. It’s difficult to follow your code.

Also, I would not use global variables. They cause more problems than they solve. They used to be necessary to pass extra parameters to ODE functions, but that has not been true for at least 15 years (at least as I remember). You don’t say what version of MATLAB you’re using.

Looking through your code, if you want to plot ‘x(2)’ as a function of time (since that seems to be what ‘y_save’ is), I would just do:

... PREVIOUS CODE ...
[t,x]=ode45(@eqp1,time,x0);
figure(1)
plot(t, x(:,2));
... REMAINING CODE ...

2 个评论
显示无隐藏无

Mohamed Aburakhis 2016-3-14

编辑：Mohamed Aburakhis 2016-3-14

thanks very much, your command works but the length of "t" is not right, how can i make it as the tf of my codes? I use matlab 15a

Star Strider 2016-3-14

在 MATLAB Online 中打开

That value of ‘t’ is created by ode45. You gave ode45 a vector of times in your ‘time’ array, not a range, so if you plot all your variables as functions of the same ‘time’ vector, everything should work correctly. (The ‘t’ vector returned by ode45 will be identical to your ‘time’ vector.)

I cannot follow what you are doing in your code, so beyond that and what I wrote previously, I cannot help. It is extremely confusing with all the global statements — none of which are necessary. If you want to plot more of the ‘x’ array variables that ode45 solves, address them as I addressed ‘x(:,2)’ in the plot call, and plot them the same way.

Also, delete all these from your ‘eqp1’ ODE function:

y_save=[y_save x2]; 
u_save=[u_save u]; 
t_save=[t_save t];

If you want to create your ‘u’ array in your script workspace, you can do that simply by calculating it there after your ode45 call:

u = x(:,3) + a_w*sin(eta_w*t);

There should be no problems because both ‘x(:,3)’ and ‘t’ are column vectors, so ‘u’ will also be a column vector.

请先登录，再进行评论。

Answer 2

Jan 2016-3-14

1
链接

此回答的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/273329-how-to-get-the-time-vector-from-the-function#answer_213521

Storing the time inside the function to be integrated is not meaningful: The integrator can reject steps and calls the function several time to compute one step. So rely on the output "t" from ODE45 only.