For loop broken ?
1 visualización (últimos 30 días)
Mostrar comentarios más antiguos
BioZ
el 12 de Oct. de 2021
Comentada: BioZ
el 12 de Oct. de 2021
Hi all, I have a for loop in my matlab code which uses forward euler method to calculate ceratain parameters. The first parameter is rho (line 47) which is suppose to be dependent on h(n). However, it seems as the for loop skips the rho as it is not being updated when the loop iterates, althought, below parameter rhotest (line64) with manual h(55) in it gives a correct answer. Any help would be appreciated :) thank you.
clc; clear all; close all
%Aircraft parameters
g=9.80665;
m=41000;
W = m*g;
Cd0 = 0.02
S=120%m^2
b=34 %m
AR=b^2/S
ef=0.82; % Efficiency Factor
K=0.03 % combined induced drag factor, = 𝑘 ∕ (𝜋 𝐴𝑅)
Ta = 110000; %Thrust Available
Tp=1; %Thrust Percentage
Tu=Ta*Tp; % (Thrust Used)
Cl=sqrt((pi()/3)*ef*AR*Cd0);
Cd = K+Cd0*Cl^2
%Cd=Cd0+((Cl^2)/(pi()*ef*AR));
tf = 60; %Final time
dt = 1; %Time step
t = 0:dt:tf;
V0 = 85; %Initial Velocity
X0 = 0; %Initial Displacement
h0 = 0; %Initial Height
gamma0 = 0.0872665 % 5deg initial climb.
%Pre-Fill with zeros in order to avoid buildup in the loop.
V = zeros(1,length(t))
V(1) = V0; % Initial velocity in dt frame
X = zeros(1,length(t))
X(1) = X0; % Initial X displacement in dt frame
h = zeros(1,length(t));
h(1) = h0; %Initial h displacement in dt frame
gamma = zeros(1,length(t));
gamma(1) = gamma0; %Initial gamma in dt frame
for n=2:length(t)
rho = (20-h(n)/1000)/(20+h(n)/1000)*1.225;
D = Cd*S*0.5*rho*V(n-1)^2;
V(n) = V(n-1)+((Ta-D-W*sin(gamma(n-1)))/m)*dt;
ROC = ((Ta*V(n)-D*V(n))/W);
gamma(n)= asin(ROC/V(n));
X(n) = X(n-1)+V(n)*dt*cos(gamma(n));
h(n) = h(n-1)+V(n)*dt*sin(gamma(n));
end
%% equation checks
rhotest= (20-h(55)/1000)/(20+h(55)/1000)*1.225;
tst2=((Ta-D-W*sin(gamma0))/m)*dt;
%V(n) = V(n-1)+(1/m)*(Ta-D-W*sin(gamma(n-1)));
tst=(1/m)*(Ta-D-W*sin(gamma(1)));
%% Plots
plot(X,h)
0 comentarios
Respuesta aceptada
David Hill
el 12 de Oct. de 2021
rho = (20-h(n)/1000)/(20+h(n)/1000)*1.225;%h(n) is always going to be zero
rho = (20-h(n-1)/1000)/(20+h(n-1)/1000)*1.225;%did you want h(n-1)?
Más respuestas (1)
Cris LaPierre
el 12 de Oct. de 2021
rho is calculated using the current value of h(n), which your code has defined as zero. Therefore, rho is always going to have the same value.
h = zeros(1,length(t));
Perhaps you meant to use h(n-1) instead?
0 comentarios
Ver también
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!