Euler Integration Step Size Change

4 ビュー (過去 30 日間)
Hamish Brown
Hamish Brown 2022 年 10 月 6 日
コメント済み: Hamish Brown 2022 年 10 月 6 日
I'm currently performing a Euler integration to solve the equations of motion for a rocket launch. I'm using a timestep of 1 at the moment, just for simplicity. When i try and change this to say, 0.1, the code breaks. Anyone know how to fix this?
for i = 1:0.1:500
%density calculator
if y<=11000
Temp(i+1) = 15.04-0.00649*y(i);
p(i+1) = 101.29*((Temp(i)+273.1)/288.08)^5.258;
elseif 11000 < y <=25000
Temp(i+1) = -56.46;
p(i+1) = 22.65*exp(1.73-0.000157*y(i));
else
Temp(i+1) = -131.21 + 0.00299*y(i);
p(i+1) = 2.488*((Temp(i)+273.1)/216.6)^-11.388;
end
rho(i+1) = p(i)/(0.2869*(Temp(i)+273.1));
% gravity model
g(i+1) = 9.81*(R_E/(R_E + y(i)))^2;
%drag force D
D(i+1) = 0.5*rho(i)*v(i)^2*Area*cd;
%speed magnitude v
v(i+1) = v(i) + (T/m0(i) - D(i)/m0(i) - g(i)*sin(gam(i)));
%path angle gam (from 90 to 0 degrees)
if i <= 25
gam(i+1) = gam(i); %vertical flight
elseif 25 < i <= 40
gam(i+1) = gam(i) - 0.008*gam(i); %kickover manouvre
else
gam(i+1) = gam(i) - ((g(i)/v(i) - v(i)/(R_E+y(i)))*cos(gam(i))); %gravity turn equation
end
%x distance downrange
x(i+1) = x(i) + (R_E/(R_E+y(i))*v(i)*cos(gam(i)));
%altitude y
y(i+1) = y(i) + (v(i)*sin(gam(i)));
VG(i) = g(i)*sin(gam(i));
VD(i) = D(i)/m0(i);
if i <= burntime
m0(i+1) = m0(i) - mdot;
else
m0(i+1) = m0(burntime);
T = 0;
end
q(i+1) = (rho(i) * v(i)^2)/2;
%Vertical Velocity
Vy(i+1) = v(i)*sin(gam(i));
%Horizontal Velocity
Vx(i+1) = v(i)*cos(gam(i));
end
It produces the error message:
Array indices must be positive integers or logical values.
Error in newtry (line 51)
Temp(i+1) = 15.04-0.00649*y(i);
Any help is appreciated.

採用された回答

Torsten
Torsten 2022 年 10 月 6 日
Only integer values can be used for array indexing. The loop must run from i=1 to i=5000 instead.
  3 件のコメント
Torsten
Torsten 2022 年 10 月 6 日
編集済み: Torsten 2022 年 10 月 6 日
500*0.1 = 50
So you would integrate up to 50 s.
The final value imax of the loop is not the end time. The end time is imax*dt.
Hamish Brown
Hamish Brown 2022 年 10 月 6 日
right, I get it now. Thank you for the help, Torsten.

サインインしてコメントする。

その他の回答 (0 件)

カテゴリ

Help Center および File ExchangeNumerical Integration and Differential Equations についてさらに検索

製品


リリース

R2021b

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by