Help with formulating conditional logic.

조회 수: 6(최근 30일)
Jeffrey Lewis
Jeffrey Lewis 2022년 11월 4일
댓글: Sam Chak 2022년 11월 4일
I am writing a stability and control code (for a class assignment) that takes into account 3 different airplane conditions; Climb data, Cruise data, and Approach data. I don't want to change variables for each condition, rather I have captured all of the data for each condition in the code with differing data for each condition.
Ultimately, I would like to output the transfer functions for each condition. I am not very savy with MATLAB and would love some help in finding the best way to formulate this logic.
The code:
clc
clear all
g = 32.174
theta1 = 0
%Cessna design parameters
S = 174 %Planform Area (ft^2)
c = 4.9 %Chord length (ft)
b = 36.0 %Wing span (ft)
%%%%%%%%%%%%%%%%Flight Conditions%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%Climb Flight Condition
altitude = 0
M = 0.120
U_1 = 133.5
q_bar = 21.2
c_frac = 26.2
a_1 = 5.4
%Cruise Flight Condition
h = 5000 %Altitude (ft)
M = 0.201 %Mach Number
U_1 = 220.1 %TAS (ft/s)
q_bar = 49.6 %Dynamic Pressure (lbs/ft^2)
c_frac = 26.4 %CG Location
a_1 = 0 %Angle of Attack (Deg)
%Approach Flight Condition
h = 0
M = 0.096
U_1 = 107.1
q_bar = 13.6
c_frac = 26.4
a_1 = 4
%%%%%%%%%%%%%%%%Mass Data For All Flight Conditions%%%%%%%%%%%%%%%%%%%%%%
%Mass Data in Climb, Cruise, and Approach Conditions
W = 2650 %Weight (lbs)
m = W/g %Mass
I_xx_B = 948 %Moment of Inertia along the X-axis (slug-ft^2)
I_yy_B = 1346 %Moment of Inertia along the Y-axis (slug-ft^2)
I_zz_B = 1967 %Moment of Inertia along the Z-axis (slug-ft^2)
I_xz_B = 0 %Moment of Inertia along the XZ plane (slug-ft^2)
%%%%%%%%%%%%%%%%%%%%%%%Longitudinal Control and Hinge Momoment $$$$$$$$$
%Climb Longitudinal Control and Hinge Moment Derivatives
Ch_a = -0.0545
Ch_del_e = -0.594
%Cruise Longitudinal Control and Hinge Moment Derivatives
Ch_a = -0.0584
Ch_del_e = -0.585
%Approach Longitudinal Control and Hinge Moment Derivatives
Ch_a = -0.0549
Ch_del_e = -0.594
%%%%%%%%%%%%%%%%%%%Lateral-Directional Stability%%%%%%%%%%%%%%%%%%%%
%Climb Lateral-Directional Stabillity Derivatives
Cl_B = -0.0895
Cl_p = -0.487
Cl_r = 0.1869
Cy_B = -0.404
Cy_p = -0.145
Cy_r = 0.267
Cn_B = 0.0907
Cn_TB = 0
Cn_p = -0.0649
Cn_r = -0.1199
%Cruise Lateral-Directional Stability Derivatives
Cl_B = -0.0923
Cl_p = -0.484
Cl_r = 0.0798
Cy_B = -0.393
Cy_p = -0.075
Cy_r = 0.214
Cn_B = 0.0587
Cn_TB = 0
Cn_p = -0.0278
Cn_r = -0.0937
%Approach Lateral-Directional Stability Derivatives
Cl_B = -0.0969
Cl_p = -0.494
Cl_r = 0.2039
Cy_B = -0.303
Cy_p = -0.213
Cy_r = 0.201
Cn_B = 0.0701
Cn_TB = 0
Cn_p = -0.0960
Cn_r = -0.1151
%%%%%%%%%%%%%%%%%%%%Lateral-Directional Control and Hing Moment%%%%%%%%%
%Climb Lateral-Directional Control and Hinge Moment Derivatives
Cl_dela = 0.229
Cl_delr = 0.0147
Cy_dela = 0
Cy_delr = 0.187
Cn_dela = -0.0504
Cn_delr = -0.0805
Ch_aa = 0
Ch_dela = -0.369
Ch_Br = 0.0819
Ch_delr = -0.579
%Cruise Lateral-Directional Control and Hinge Moment Derivatives
Cl_dela = 0.229
Cl_delr = 0.0147
Cy_dela = 0
Cy_delr = 0.187
Cn_dela = -0.0216
Cn_delr = -0.0645
Ch_aa = 0
Ch_dela = -0.363
Ch_Br = 0.0819
Ch_delr = -0.567
%Approach Lateral-Directional Control and Hinge moment Derivatives
Cl_dela = 0.229
Cl_delr = 0.0147
Cy_dela = 0
Cy_delr = 0.187
Cn_dela = -0.0504
Cn_delr = -0.0805
Ch_aa = 0
Ch_dela = -0.369
Ch_Br = 0.0819
Ch_delr = -0.579
%Lateral-Directional Dimensional Stability Derivatives
YB = (q_bar*S*Cy_B)/m
Yp = (S*b*Cy_p)/(2*m*U_1)
Yr = (q_bar*S*b*Cy_r)/(2*m*U_1)
Y_dela = (q_bar*S*Cy_dela)/m
Y_delr = (q_bar*S*Cy_delr)/m
LB = (q_bar*S*b*Cl_B)/I_xx_B
Lp = (q_bar*S*(b^2)*Cl_p)/(2*I_xx_B*U_1)
Lr = (q_bar*S*(b^2)*Cl_r)/(2*I_xx_B*U_1)
L_dela = (q_bar*S*b*Cl_dela)/I_xx_B
L_delr = (q_bar*S*b*Cl_delr)/I_xx_B
NB = (q_bar*S*b*Cn_B)/I_zz_B
NT_B = (q_bar*S*b*Cn_TB)/I_zz_B
Np = (q_bar*S*(b^2)*Cn_p)/(2*I_zz_B*U_1)
Nr = (q_bar*S*(b^2)*Cn_r)/(2*I_zz_B*U_1)
N_dela = (q_bar*S*b*Cn_dela)/I_zz_B
N_delr = (q_bar*S*b*Cn_delr)/I_zz_B
%Perturbed Lateral-Directional Equations of Motion
A_1 = I_xz/I_xx_B
B_1 = I_xz/I_zz_B
%Denominator Polynomial
A2 = U_1(1-A_1*B_1)
B2 = -YB*(1-A_1*B_1)-U_1*(Lp+Nr+A_1*Np+B_1*Lr)
C2 = U_1*(Lp*Nr-Lr*Np)+YB*(Nr+Lp+A_1*Np+B_1*Lr)-Yp*(LB+NB*A_1+NT_B*A_1)...
+U_1*(LB*B_1+NB+NT_B)-Yr*(LB*B_1+NB+NT_B)
D2 = -YB*(Lp*Nr-Lr*Np)+Yp*(LB*Nr-NB*Lr-NT_B*Lr)-g*cos(theta1)...
*(LB+NB*A_1+NT_B*A_1)+U_1*(LB*Np-NB*Lp-NT_B*Lp)-Yr*(LB*Np-NB*Lp-NT_B*Lp)
E2 = g*cos(theta1)*(LB*Nr-NB*Lr-NT_B*Lr)
D_2 = [A2 B2 C2 D2 E2]
%Numerator Polynomial
AB = Y_dela*(1-A_1*B_1)+Y_delr*(1-A_1*B_1)
BB = -Y_dela*(Nr+Lp+A_1*Np+B_1*Lr)-Y_delr*(Nr+Lp+A_1*Np+B_1*Lr)+Yp*(L_dela+N_dela*A_1)+Yp*(L_delr+N_delr*A_1)...
+Yr*(L_dela*B_1+N_dela)+Yr*(L_delr*B_1+N_delr)-U_1*(L_dela*B_1+N_dela)...
+U_1*(L_delr*B_1+N_delr)
CB = Y_dela*(Lp*Nr-Np*Lr)+Y_delr*(Lp*Nr-Np*Lr)+Yp*(N_dela*Lr-L_dela*Nr)+Yp*(N_delr*Lr-L_delr*Nr)...
+g*cos(theta1)*(L_dela+N_dela*A_1)+g*cos(theta1)*(L_delr+N_delr*A_1)...
+Yr*(L_dela*Np-N_dela*Lp)+Yr*(L_delr*Np-N_delr*Lp)+U_1*(L_dela*Np-N_dela*Lp)+U_1*(L_delr*Np-N_delr*Lp)
DB = g*cos(theta1)*(N_dela*Lr-L_dela*Nr)+g*cos(theta1)*(N_delr*Lr-L_delr*Nr)
N_B = [AB BB CB DB]
%Transfer Function
sys = tf(N_B,D_2)
damp(sys)

채택된 답변

Sam Chak
Sam Chak 2022년 11월 4일
There are a few ways to do it. But since you have written the data into the code, I think probably the simplest way is to write a function file for eash flight condition, and then call the specific function. The example is shown below:
1. Calling climbsys
climbsys
Climb_sys = 8.375 s^3 - 200.8 s^2 - 1.251e04 s + 195.1 ------------------------------------------------ 133.5 s^4 + 1392 s^3 + 2612 s^2 + 9028 s - 255.2 Continuous-time transfer function. Pole Damping Frequency Time Constant (rad/seconds) (seconds) 2.80e-02 -1.00e+00 2.80e-02 -3.57e+01 -6.80e-01 + 2.65e+00i 2.48e-01 2.74e+00 1.47e+00 -6.80e-01 - 2.65e+00i 2.48e-01 2.74e+00 1.47e+00 -9.09e+00 1.00e+00 9.09e+00 1.10e-01
2. Calling cruisesys
cruisesys
Cruise_sys = 19.59 s^3 - 1240 s^2 - 4.263e04 s + 2174 -------------------------------------------------- 220.1 s^4 + 3163 s^3 + 6232 s^2 + 3.028e04 s + 540 Continuous-time transfer function. Pole Damping Frequency Time Constant (rad/seconds) (seconds) -1.79e-02 1.00e+00 1.79e-02 5.59e+01 -6.74e-01 + 3.18e+00i 2.08e-01 3.25e+00 1.48e+00 -6.74e-01 - 3.18e+00i 2.08e-01 3.25e+00 1.48e+00 -1.30e+01 1.00e+00 1.30e+01 7.69e-02
3. Calling approachsys
approachsys
Approach_sys = 5.373 s^3 - 102.2 s^2 - 5475 s + 28.61 ------------------------------------------------- 107.1 s^4 + 897.5 s^3 + 1294 s^2 + 3403 s - 66.09 Continuous-time transfer function. Pole Damping Frequency Time Constant (rad/seconds) (seconds) 1.93e-02 -1.00e+00 1.93e-02 -5.19e+01 -5.37e-01 + 2.02e+00i 2.57e-01 2.09e+00 1.86e+00 -5.37e-01 - 2.02e+00i 2.57e-01 2.09e+00 1.86e+00 -7.32e+00 1.00e+00 7.32e+00 1.37e-01
function climbsys
g = 32.174;
theta1 = 0;
% Cessna design parameters
S = 174; % Planform Area (ft^2)
c = 4.9; % Chord length (ft)
b = 36.0; % Wing span (ft)
%%%%%%%%%%%%%%%%Flight Conditions%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Climb Flight Condition
altitude = 0;
M = 0.120;
U_1 = 133.5;
q_bar = 21.2;
c_frac = 26.2;
a_1 = 5.4;
% % Cruise Flight Condition
% h = 5000; % Altitude (ft)
% M = 0.201; % Mach Number
% U_1 = 220.1; % TAS (ft/s)
% q_bar = 49.6; % Dynamic Pressure (lbs/ft^2)
% c_frac = 26.4; % CG Location
% a_1 = 0; % Angle of Attack (Deg)
%
% % Approach Flight Condition
% h = 0;
% M = 0.096;
% U_1 = 107.1;
% q_bar = 13.6;
% c_frac = 26.4;
% a_1 = 4;
%%%%%%%%%%%%%%%%Mass Data For All Flight Conditions%%%%%%%%%%%%%%%%%%%%%%
% Mass Data in Climb, Cruise, and Approach Conditions
W = 2650; % Weight (lbs)
m = W/g; % Mass
I_xx_B = 948; % Moment of Inertia along the X-axis (slug-ft^2)
I_yy_B = 1346; % Moment of Inertia along the Y-axis (slug-ft^2)
I_zz_B = 1967; % Moment of Inertia along the Z-axis (slug-ft^2)
I_xz_B = 0; % Moment of Inertia along the XZ plane (slug-ft^2)
%%%%%%%%%%%%%%%%%%%%%%%Longitudinal Control and Hinge Momoment $$$$$$$$$
% Climb Longitudinal Control and Hinge Moment Derivatives
Ch_a = -0.0545;
Ch_del_e = -0.594;
% % Cruise Longitudinal Control and Hinge Moment Derivatives
% Ch_a = -0.0584;
% Ch_del_e = -0.585;
%
% % Approach Longitudinal Control and Hinge Moment Derivatives
% Ch_a = -0.0549;
% Ch_del_e = -0.594;
%%%%%%%%%%%%%%%%%%%Lateral-Directional Stability%%%%%%%%%%%%%%%%%%%%
% Climb Lateral-Directional Stabillity Derivatives
Cl_B = -0.0895;
Cl_p = -0.487;
Cl_r = 0.1869;
Cy_B = -0.404;
Cy_p = -0.145;
Cy_r = 0.267;
Cn_B = 0.0907;
Cn_TB = 0;
Cn_p = -0.0649;
Cn_r = -0.1199;
% % Cruise Lateral-Directional Stability Derivatives
% Cl_B = -0.0923;
% Cl_p = -0.484;
% Cl_r = 0.0798;
% Cy_B = -0.393;
% Cy_p = -0.075;
% Cy_r = 0.214;
% Cn_B = 0.0587;
% Cn_TB = 0;
% Cn_p = -0.0278;
% Cn_r = -0.0937;
%
% % Approach Lateral-Directional Stability Derivatives
% Cl_B = -0.0969;
% Cl_p = -0.494;
% Cl_r = 0.2039;
% Cy_B = -0.303;
% Cy_p = -0.213;
% Cy_r = 0.201;
% Cn_B = 0.0701;
% Cn_TB = 0;
% Cn_p = -0.0960;
% Cn_r = -0.1151;
%%%%%%%%%%%%%%%%%%%%Lateral-Directional Control and Hing Moment%%%%%%%%%
% Climb Lateral-Directional Control and Hinge Moment Derivatives
Cl_dela = 0.229;
Cl_delr = 0.0147;
Cy_dela = 0;
Cy_delr = 0.187;
Cn_dela = -0.0504;
Cn_delr = -0.0805;
Ch_aa = 0;
Ch_dela = -0.369;
Ch_Br = 0.0819;
Ch_delr = -0.579;
% % Cruise Lateral-Directional Control and Hinge Moment Derivatives
% Cl_dela = 0.229;
% Cl_delr = 0.0147;
% Cy_dela = 0;
% Cy_delr = 0.187;
% Cn_dela = -0.0216;
% Cn_delr = -0.0645;
% Ch_aa = 0;
% Ch_dela = -0.363;
% Ch_Br = 0.0819;
% Ch_delr = -0.567;
%
% % Approach Lateral-Directional Control and Hinge moment Derivatives
% Cl_dela = 0.229;
% Cl_delr = 0.0147;
% Cy_dela = 0;
% Cy_delr = 0.187;
% Cn_dela = -0.0504;
% Cn_delr = -0.0805;
% Ch_aa = 0;
% Ch_dela = -0.369;
% Ch_Br = 0.0819;
% Ch_delr = -0.579;
% Lateral-Directional Dimensional Stability Derivatives
YB = (q_bar*S*Cy_B)/m;
Yp = (S*b*Cy_p)/(2*m*U_1);
Yr = (q_bar*S*b*Cy_r)/(2*m*U_1);
Y_dela = (q_bar*S*Cy_dela)/m;
Y_delr = (q_bar*S*Cy_delr)/m;
LB = (q_bar*S*b*Cl_B)/I_xx_B;
Lp = (q_bar*S*(b^2)*Cl_p)/(2*I_xx_B*U_1);
Lr = (q_bar*S*(b^2)*Cl_r)/(2*I_xx_B*U_1);
L_dela = (q_bar*S*b*Cl_dela)/I_xx_B;
L_delr = (q_bar*S*b*Cl_delr)/I_xx_B;
NB = (q_bar*S*b*Cn_B)/I_zz_B;
NT_B = (q_bar*S*b*Cn_TB)/I_zz_B;
Np = (q_bar*S*(b^2)*Cn_p)/(2*I_zz_B*U_1);
Nr = (q_bar*S*(b^2)*Cn_r)/(2*I_zz_B*U_1);
N_dela = (q_bar*S*b*Cn_dela)/I_zz_B;
N_delr = (q_bar*S*b*Cn_delr)/I_zz_B;
% Perturbed Lateral-Directional Equations of Motion
A_1 = I_xz_B/I_xx_B;
B_1 = I_xz_B/I_zz_B;
% Numerator Polynomial
AB = Y_dela*(1 - A_1*B_1) + Y_delr*(1 - A_1*B_1);
BB = - Y_dela*(Nr + Lp + A_1*Np + B_1*Lr) - Y_delr*(Nr + Lp + A_1*Np + B_1*Lr) + Yp*(L_dela + N_dela*A_1) + Yp*(L_delr + N_delr*A_1) + Yr*(L_dela*B_1 + N_dela) + Yr*(L_delr*B_1 + N_delr) - U_1*(L_dela*B_1 + N_dela) + U_1*(L_delr*B_1 + N_delr);
CB = Y_dela*(Lp*Nr - Np*Lr) + Y_delr*(Lp*Nr - Np*Lr) + Yp*(N_dela*Lr - L_dela*Nr) + Yp*(N_delr*Lr - L_delr*Nr) + g*cos(theta1)*(L_dela + N_dela*A_1) + g*cos(theta1)*(L_delr + N_delr*A_1) + Yr*(L_dela*Np - N_dela*Lp) + Yr*(L_delr*Np - N_delr*Lp) + U_1*(L_dela*Np - N_dela*Lp) + U_1*(L_delr*Np - N_delr*Lp);
DB = g*cos(theta1)*(N_dela*Lr-L_dela*Nr)+g*cos(theta1)*(N_delr*Lr-L_delr*Nr);
N_B = [AB BB CB DB];
% Denominator Polynomial
A2 = U_1*(1 - A_1*B_1);
B2 = - YB*(1 - A_1*B_1) - U_1*(Lp + Nr + A_1*Np + B_1*Lr);
C2 = U_1*(Lp*Nr - Lr*Np) + YB*(Nr + Lp + A_1*Np + B_1*Lr) - Yp*(LB + NB*A_1 + NT_B*A_1) + U_1*(LB*B_1 + NB + NT_B) - Yr*(LB*B_1 + NB + NT_B);
D2 = - YB*(Lp*Nr - Lr*Np) + Yp*(LB*Nr - NB*Lr - NT_B*Lr) - g*cos(theta1)*(LB + NB*A_1 + NT_B*A_1) + U_1*(LB*Np - NB*Lp - NT_B*Lp) - Yr*(LB*Np - NB*Lp - NT_B*Lp);
E2 = g*cos(theta1)*(LB*Nr - NB*Lr - NT_B*Lr);
D_2 = [A2 B2 C2 D2 E2];
% Transfer Function
Climb_sys = tf(N_B, D_2)
damp(Climb_sys)
step(Climb_sys)
end
function cruisesys
g = 32.174;
theta1 = 0;
% Cessna design parameters
S = 174; % Planform Area (ft^2)
c = 4.9; % Chord length (ft)
b = 36.0; % Wing span (ft)
%%%%%%%%%%%%%%%%Flight Conditions%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% % Climb Flight Condition
% altitude = 0;
% M = 0.120;
% U_1 = 133.5;
% q_bar = 21.2;
% c_frac = 26.2;
% a_1 = 5.4;
% Cruise Flight Condition
h = 5000; % Altitude (ft)
M = 0.201; % Mach Number
U_1 = 220.1; % TAS (ft/s)
q_bar = 49.6; % Dynamic Pressure (lbs/ft^2)
c_frac = 26.4; % CG Location
a_1 = 0; % Angle of Attack (Deg)
% % Approach Flight Condition
% h = 0;
% M = 0.096;
% U_1 = 107.1;
% q_bar = 13.6;
% c_frac = 26.4;
% a_1 = 4;
%%%%%%%%%%%%%%%%Mass Data For All Flight Conditions%%%%%%%%%%%%%%%%%%%%%%
% Mass Data in Climb, Cruise, and Approach Conditions
W = 2650; % Weight (lbs)
m = W/g; % Mass
I_xx_B = 948; % Moment of Inertia along the X-axis (slug-ft^2)
I_yy_B = 1346; % Moment of Inertia along the Y-axis (slug-ft^2)
I_zz_B = 1967; % Moment of Inertia along the Z-axis (slug-ft^2)
I_xz_B = 0; % Moment of Inertia along the XZ plane (slug-ft^2)
%%%%%%%%%%%%%%%%%%%%%%%Longitudinal Control and Hinge Momoment $$$$$$$$$
% % Climb Longitudinal Control and Hinge Moment Derivatives
% Ch_a = -0.0545;
% Ch_del_e = -0.594;
% Cruise Longitudinal Control and Hinge Moment Derivatives
Ch_a = -0.0584;
Ch_del_e = -0.585;
% % Approach Longitudinal Control and Hinge Moment Derivatives
% Ch_a = -0.0549;
% Ch_del_e = -0.594;
%%%%%%%%%%%%%%%%%%%Lateral-Directional Stability%%%%%%%%%%%%%%%%%%%%
% % Climb Lateral-Directional Stabillity Derivatives
% Cl_B = -0.0895;
% Cl_p = -0.487;
% Cl_r = 0.1869;
% Cy_B = -0.404;
% Cy_p = -0.145;
% Cy_r = 0.267;
% Cn_B = 0.0907;
% Cn_TB = 0;
% Cn_p = -0.0649;
% Cn_r = -0.1199;
% Cruise Lateral-Directional Stability Derivatives
Cl_B = -0.0923;
Cl_p = -0.484;
Cl_r = 0.0798;
Cy_B = -0.393;
Cy_p = -0.075;
Cy_r = 0.214;
Cn_B = 0.0587;
Cn_TB = 0;
Cn_p = -0.0278;
Cn_r = -0.0937;
% % Approach Lateral-Directional Stability Derivatives
% Cl_B = -0.0969;
% Cl_p = -0.494;
% Cl_r = 0.2039;
% Cy_B = -0.303;
% Cy_p = -0.213;
% Cy_r = 0.201;
% Cn_B = 0.0701;
% Cn_TB = 0;
% Cn_p = -0.0960;
% Cn_r = -0.1151;
%%%%%%%%%%%%%%%%%%%%Lateral-Directional Control and Hing Moment%%%%%%%%%
% % Climb Lateral-Directional Control and Hinge Moment Derivatives
% Cl_dela = 0.229;
% Cl_delr = 0.0147;
% Cy_dela = 0;
% Cy_delr = 0.187;
% Cn_dela = -0.0504;
% Cn_delr = -0.0805;
% Ch_aa = 0;
% Ch_dela = -0.369;
% Ch_Br = 0.0819;
% Ch_delr = -0.579;
% Cruise Lateral-Directional Control and Hinge Moment Derivatives
Cl_dela = 0.229;
Cl_delr = 0.0147;
Cy_dela = 0;
Cy_delr = 0.187;
Cn_dela = -0.0216;
Cn_delr = -0.0645;
Ch_aa = 0;
Ch_dela = -0.363;
Ch_Br = 0.0819;
Ch_delr = -0.567;
% % Approach Lateral-Directional Control and Hinge moment Derivatives
% Cl_dela = 0.229;
% Cl_delr = 0.0147;
% Cy_dela = 0;
% Cy_delr = 0.187;
% Cn_dela = -0.0504;
% Cn_delr = -0.0805;
% Ch_aa = 0;
% Ch_dela = -0.369;
% Ch_Br = 0.0819;
% Ch_delr = -0.579;
% Lateral-Directional Dimensional Stability Derivatives
YB = (q_bar*S*Cy_B)/m;
Yp = (S*b*Cy_p)/(2*m*U_1);
Yr = (q_bar*S*b*Cy_r)/(2*m*U_1);
Y_dela = (q_bar*S*Cy_dela)/m;
Y_delr = (q_bar*S*Cy_delr)/m;
LB = (q_bar*S*b*Cl_B)/I_xx_B;
Lp = (q_bar*S*(b^2)*Cl_p)/(2*I_xx_B*U_1);
Lr = (q_bar*S*(b^2)*Cl_r)/(2*I_xx_B*U_1);
L_dela = (q_bar*S*b*Cl_dela)/I_xx_B;
L_delr = (q_bar*S*b*Cl_delr)/I_xx_B;
NB = (q_bar*S*b*Cn_B)/I_zz_B;
NT_B = (q_bar*S*b*Cn_TB)/I_zz_B;
Np = (q_bar*S*(b^2)*Cn_p)/(2*I_zz_B*U_1);
Nr = (q_bar*S*(b^2)*Cn_r)/(2*I_zz_B*U_1);
N_dela = (q_bar*S*b*Cn_dela)/I_zz_B;
N_delr = (q_bar*S*b*Cn_delr)/I_zz_B;
% Perturbed Lateral-Directional Equations of Motion
A_1 = I_xz_B/I_xx_B;
B_1 = I_xz_B/I_zz_B;
% Numerator Polynomial
AB = Y_dela*(1 - A_1*B_1) + Y_delr*(1 - A_1*B_1);
BB = - Y_dela*(Nr + Lp + A_1*Np + B_1*Lr) - Y_delr*(Nr + Lp + A_1*Np + B_1*Lr) + Yp*(L_dela + N_dela*A_1) + Yp*(L_delr + N_delr*A_1) + Yr*(L_dela*B_1 + N_dela) + Yr*(L_delr*B_1 + N_delr) - U_1*(L_dela*B_1 + N_dela) + U_1*(L_delr*B_1 + N_delr);
CB = Y_dela*(Lp*Nr - Np*Lr) + Y_delr*(Lp*Nr - Np*Lr) + Yp*(N_dela*Lr - L_dela*Nr) + Yp*(N_delr*Lr - L_delr*Nr) + g*cos(theta1)*(L_dela + N_dela*A_1) + g*cos(theta1)*(L_delr + N_delr*A_1) + Yr*(L_dela*Np - N_dela*Lp) + Yr*(L_delr*Np - N_delr*Lp) + U_1*(L_dela*Np - N_dela*Lp) + U_1*(L_delr*Np - N_delr*Lp);
DB = g*cos(theta1)*(N_dela*Lr-L_dela*Nr)+g*cos(theta1)*(N_delr*Lr-L_delr*Nr);
N_B = [AB BB CB DB];
% Denominator Polynomial
A2 = U_1*(1 - A_1*B_1);
B2 = - YB*(1 - A_1*B_1) - U_1*(Lp + Nr + A_1*Np + B_1*Lr);
C2 = U_1*(Lp*Nr - Lr*Np) + YB*(Nr + Lp + A_1*Np + B_1*Lr) - Yp*(LB + NB*A_1 + NT_B*A_1) + U_1*(LB*B_1 + NB + NT_B) - Yr*(LB*B_1 + NB + NT_B);
D2 = - YB*(Lp*Nr - Lr*Np) + Yp*(LB*Nr - NB*Lr - NT_B*Lr) - g*cos(theta1)*(LB + NB*A_1 + NT_B*A_1) + U_1*(LB*Np - NB*Lp - NT_B*Lp) - Yr*(LB*Np - NB*Lp - NT_B*Lp);
E2 = g*cos(theta1)*(LB*Nr - NB*Lr - NT_B*Lr);
D_2 = [A2 B2 C2 D2 E2];
% Transfer Function
Cruise_sys = tf(N_B, D_2)
damp(Cruise_sys)
step(Cruise_sys)
end
function approachsys
g = 32.174;
theta1 = 0;
% Cessna design parameters
S = 174; % Planform Area (ft^2)
c = 4.9; % Chord length (ft)
b = 36.0; % Wing span (ft)
%%%%%%%%%%%%%%%%Flight Conditions%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% % Climb Flight Condition
% altitude = 0;
% M = 0.120;
% U_1 = 133.5;
% q_bar = 21.2;
% c_frac = 26.2;
% a_1 = 5.4;
%
% % Cruise Flight Condition
% h = 5000; % Altitude (ft)
% M = 0.201; % Mach Number
% U_1 = 220.1; % TAS (ft/s)
% q_bar = 49.6; % Dynamic Pressure (lbs/ft^2)
% c_frac = 26.4; % CG Location
% a_1 = 0; % Angle of Attack (Deg)
% Approach Flight Condition
h = 0;
M = 0.096;
U_1 = 107.1;
q_bar = 13.6;
c_frac = 26.4;
a_1 = 4;
%%%%%%%%%%%%%%%%Mass Data For All Flight Conditions%%%%%%%%%%%%%%%%%%%%%%
% Mass Data in Climb, Cruise, and Approach Conditions
W = 2650; % Weight (lbs)
m = W/g; % Mass
I_xx_B = 948; % Moment of Inertia along the X-axis (slug-ft^2)
I_yy_B = 1346; % Moment of Inertia along the Y-axis (slug-ft^2)
I_zz_B = 1967; % Moment of Inertia along the Z-axis (slug-ft^2)
I_xz_B = 0; % Moment of Inertia along the XZ plane (slug-ft^2)
%%%%%%%%%%%%%%%%%%%%%%%Longitudinal Control and Hinge Momoment $$$$$$$$$
% % Climb Longitudinal Control and Hinge Moment Derivatives
% Ch_a = -0.0545;
% Ch_del_e = -0.594;
%
% % Cruise Longitudinal Control and Hinge Moment Derivatives
% Ch_a = -0.0584;
% Ch_del_e = -0.585;
% Approach Longitudinal Control and Hinge Moment Derivatives
Ch_a = -0.0549;
Ch_del_e = -0.594;
%%%%%%%%%%%%%%%%%%%Lateral-Directional Stability%%%%%%%%%%%%%%%%%%%%
% % Climb Lateral-Directional Stabillity Derivatives
% Cl_B = -0.0895;
% Cl_p = -0.487;
% Cl_r = 0.1869;
% Cy_B = -0.404;
% Cy_p = -0.145;
% Cy_r = 0.267;
% Cn_B = 0.0907;
% Cn_TB = 0;
% Cn_p = -0.0649;
% Cn_r = -0.1199;
%
% % Cruise Lateral-Directional Stability Derivatives
% Cl_B = -0.0923;
% Cl_p = -0.484;
% Cl_r = 0.0798;
% Cy_B = -0.393;
% Cy_p = -0.075;
% Cy_r = 0.214;
% Cn_B = 0.0587;
% Cn_TB = 0;
% Cn_p = -0.0278;
% Cn_r = -0.0937;
% Approach Lateral-Directional Stability Derivatives
Cl_B = -0.0969;
Cl_p = -0.494;
Cl_r = 0.2039;
Cy_B = -0.303;
Cy_p = -0.213;
Cy_r = 0.201;
Cn_B = 0.0701;
Cn_TB = 0;
Cn_p = -0.0960;
Cn_r = -0.1151;
%%%%%%%%%%%%%%%%%%%%Lateral-Directional Control and Hing Moment%%%%%%%%%
% % Climb Lateral-Directional Control and Hinge Moment Derivatives
% Cl_dela = 0.229;
% Cl_delr = 0.0147;
% Cy_dela = 0;
% Cy_delr = 0.187;
% Cn_dela = -0.0504;
% Cn_delr = -0.0805;
% Ch_aa = 0;
% Ch_dela = -0.369;
% Ch_Br = 0.0819;
% Ch_delr = -0.579;
%
% % Cruise Lateral-Directional Control and Hinge Moment Derivatives
% Cl_dela = 0.229;
% Cl_delr = 0.0147;
% Cy_dela = 0;
% Cy_delr = 0.187;
% Cn_dela = -0.0216;
% Cn_delr = -0.0645;
% Ch_aa = 0;
% Ch_dela = -0.363;
% Ch_Br = 0.0819;
% Ch_delr = -0.567;
% Approach Lateral-Directional Control and Hinge moment Derivatives
Cl_dela = 0.229;
Cl_delr = 0.0147;
Cy_dela = 0;
Cy_delr = 0.187;
Cn_dela = -0.0504;
Cn_delr = -0.0805;
Ch_aa = 0;
Ch_dela = -0.369;
Ch_Br = 0.0819;
Ch_delr = -0.579;
% Lateral-Directional Dimensional Stability Derivatives
YB = (q_bar*S*Cy_B)/m;
Yp = (S*b*Cy_p)/(2*m*U_1);
Yr = (q_bar*S*b*Cy_r)/(2*m*U_1);
Y_dela = (q_bar*S*Cy_dela)/m;
Y_delr = (q_bar*S*Cy_delr)/m;
LB = (q_bar*S*b*Cl_B)/I_xx_B;
Lp = (q_bar*S*(b^2)*Cl_p)/(2*I_xx_B*U_1);
Lr = (q_bar*S*(b^2)*Cl_r)/(2*I_xx_B*U_1);
L_dela = (q_bar*S*b*Cl_dela)/I_xx_B;
L_delr = (q_bar*S*b*Cl_delr)/I_xx_B;
NB = (q_bar*S*b*Cn_B)/I_zz_B;
NT_B = (q_bar*S*b*Cn_TB)/I_zz_B;
Np = (q_bar*S*(b^2)*Cn_p)/(2*I_zz_B*U_1);
Nr = (q_bar*S*(b^2)*Cn_r)/(2*I_zz_B*U_1);
N_dela = (q_bar*S*b*Cn_dela)/I_zz_B;
N_delr = (q_bar*S*b*Cn_delr)/I_zz_B;
% Perturbed Lateral-Directional Equations of Motion
A_1 = I_xz_B/I_xx_B;
B_1 = I_xz_B/I_zz_B;
% Numerator Polynomial
AB = Y_dela*(1 - A_1*B_1) + Y_delr*(1 - A_1*B_1);
BB = - Y_dela*(Nr + Lp + A_1*Np + B_1*Lr) - Y_delr*(Nr + Lp + A_1*Np + B_1*Lr) + Yp*(L_dela + N_dela*A_1) + Yp*(L_delr + N_delr*A_1) + Yr*(L_dela*B_1 + N_dela) + Yr*(L_delr*B_1 + N_delr) - U_1*(L_dela*B_1 + N_dela) + U_1*(L_delr*B_1 + N_delr);
CB = Y_dela*(Lp*Nr - Np*Lr) + Y_delr*(Lp*Nr - Np*Lr) + Yp*(N_dela*Lr - L_dela*Nr) + Yp*(N_delr*Lr - L_delr*Nr) + g*cos(theta1)*(L_dela + N_dela*A_1) + g*cos(theta1)*(L_delr + N_delr*A_1) + Yr*(L_dela*Np - N_dela*Lp) + Yr*(L_delr*Np - N_delr*Lp) + U_1*(L_dela*Np - N_dela*Lp) + U_1*(L_delr*Np - N_delr*Lp);
DB = g*cos(theta1)*(N_dela*Lr-L_dela*Nr)+g*cos(theta1)*(N_delr*Lr-L_delr*Nr);
N_B = [AB BB CB DB];
% Denominator Polynomial
A2 = U_1*(1 - A_1*B_1);
B2 = - YB*(1 - A_1*B_1) - U_1*(Lp + Nr + A_1*Np + B_1*Lr);
C2 = U_1*(Lp*Nr - Lr*Np) + YB*(Nr + Lp + A_1*Np + B_1*Lr) - Yp*(LB + NB*A_1 + NT_B*A_1) + U_1*(LB*B_1 + NB + NT_B) - Yr*(LB*B_1 + NB + NT_B);
D2 = - YB*(Lp*Nr - Lr*Np) + Yp*(LB*Nr - NB*Lr - NT_B*Lr) - g*cos(theta1)*(LB + NB*A_1 + NT_B*A_1) + U_1*(LB*Np - NB*Lp - NT_B*Lp) - Yr*(LB*Np - NB*Lp - NT_B*Lp);
E2 = g*cos(theta1)*(LB*Nr - NB*Lr - NT_B*Lr);
D_2 = [A2 B2 C2 D2 E2];
% Transfer Function
Approach_sys = tf(N_B, D_2)
damp(Approach_sys)
step(Approach_sys)
end
  댓글 수: 2
Sam Chak
Sam Chak 2022년 11월 4일
I actually added the step function in the code to generate step response for each flight condition. Step response always start
You can also try the lsim(sys, u, t, x0) function for a more realistic simulated response data because you can design the the input signal u, and specify the initial condition x0. However, you must convert sys from transfer function into a state-space model (preferably in controllable canonical form)

댓글을 달려면 로그인하십시오.

추가 답변(0개)

제품


릴리스

R2021a

Community Treasure Hunt

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

Start Hunting!

Translated by