i have problems of my ode45 codes

Question

kyu hong lee 2016-6-10

0
链接

此问题的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/289135-i-have-problems-of-my-ode45-codes

回答： Torsten 2016-6-10

采纳的回答： Torsten

在 MATLAB Online 中打开

(Sorry about my english skill)

i want to calculate height of rocket. so i wrote some codes.

(maybe some values are not correct.. because i want to know my codes are work without any error)

.m file

function dH=rocket(t,H)
global km mildo press temp;
m = 6100;   % kg
a = 0.8;     % m^2
ve = 2060.1; % m/s
thrust = 13000*9.81; % thrust N
pe = thrust/a;  % pressure at engine
dmdt = -(thrust)/ve;
dpae=polyfit(km,mildo,15);  % density per altitude equation
ppae=polyfit(km,press,15);  % pressure per altitude equation
tpae=polyfit(km,temp,15);   % temperature per altitude equation
rho = polyval(dpae,H);  % density per altitude 
pa =  polyval(ppae,H);  % pressure per altitude 
tem = polyval(tpae,H);  % temperature per altitude 
%M = (331.5+(0.6*tem));      % mach
veeff = ve-a*(pe-pa)/dmdt;  % ve.eff 
%v = ; % velocitiy of rocket
cd = 1;
dH(1) = H(2);
dH(2) = cd*rho(1)*a*H(1)^2/(2*m)+dmdt*veeff(1);

and

rocket_test.m file

[t H] = ode45(@rocket,[0 15],[0 1]);

i saw some examples and i replaced example variations to my own variations.

the error contents are

i
Warning: Polynomial is not unique; degree >= number of data points. 
> In polyfit (line 70)
In rocket (line 12)
In odearguments (line 87)
In ode45 (line 115)
In rocket_test (line 1) 
Warning: Polynomial is not unique; degree >= number of data points. 
> In polyfit (line 70)
In rocket (line 13)
In odearguments (line 87)
In ode45 (line 115)
In rocket_test (line 1) 
Warning: Polynomial is not unique; degree >= number of data points. 
> In polyfit (line 70)
In rocket (line 14)
In odearguments (line 87)
In ode45 (line 115)
In rocket_test (line 1) 
Error using odearguments (line 90)
ROCKET must return a column vector.
Error in ode45 (line 115)
  odearguments(FcnHandlesUsed, solver_name, ode, tspan, y0, options, varargin);
Error in rocket_test (line 1)
[t H] = ode45(@rocket,[0 15],[0 1]);

tell me what is the problem (T_T)

0 个评论
显示 -2更早的评论隐藏 -2更早的评论

请先登录，再进行评论。

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

Answer 1

Torsten 2016-6-10

0
链接

此回答的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/289135-i-have-problems-of-my-ode45-codes#answer_225041

在 MATLAB Online 中打开

1. Call polyfit before calling ODE45 and pass dpae, ppae and tpae to rocket:

dpae=polyfit(km,mildo,15);  % density per altitude equation
ppae=polyfit(km,press,15);  % pressure per altitude equation
tpae=polyfit(km,temp,15);   % temperature per altitude equation
[t H] = ode45(@(t,H)rocket(t,H,dpae,ppae,tpae),[0 15],[0 1]);

2. Use a polynomial of degree much smaller than the number of data points (i.e. reduce 15).

3. Use a column vector for dH:

function dH=rocket(t,H,dpae,ppae,tpae)
m = 6100;   % kg
a = 0.8;     % m^2
ve = 2060.1; % m/s
thrust = 13000*9.81; % thrust N
pe = thrust/a;  % pressure at engine
dmdt = -(thrust)/ve;
rho = polyval(dpae,H);  % density per altitude 
pa =  polyval(ppae,H);  % pressure per altitude 
tem = polyval(tpae,H);  % temperature per altitude 
%M = (331.5+(0.6*tem));      % mach
veeff = ve-a*(pe-pa)/dmdt;  % ve.eff 
%v = ; % velocitiy of rocket
cd = 1;
dH(1,1) = H(2);
dH(2,1) = cd*rho(1)*a*H(1)^2/(2*m)+dmdt*veeff(1);