如何调试此 Matlab 代码以模拟滑翔机下降?

How to debug this Matlab code to model glider descent?

我有这个 Matlab 项目,但出于某种原因我无法停止思考它,因为我无法让它工作。

Objective:

这是一个 MATLAB 脚本,可以计算从 10000 英尺下降的滑翔机的压力、温度和密度的变化。随着它的下降,我们想要使用计算出的这些新值,然后插入一个函数,该函数具有 4 个方程,需要使用 ode45 在每个点进行微分,以及 P T 和 Rho 的新值。

下面是主要代码:

% HouseKeeping:

clc
clear all
close all

% Constants:
S = 232;                                         % ft^2
Cd0 = 0.02;
K = 0.07;
W = 11000;                                       % lbf
Cl_max = sqrt(Cd0/K);
Cd_max = 2*K*Cl_max^2;
Rho_10000 = .001756;                             % slugs/ ft^3

%Initial conditions:
t = 0;                                           % Sec
x = 0;                                           % ft
h = 10000;                                       % ft
v = sqrt((2*W)/(Rho_10000*S*Cl_max));            %ft/s
gamma = -Cd_max/Cl_max; 

% Find Endurance:
V_RD= sqrt((2*W)/(S* Rho_10000* sqrt(3*Cd0/K)));
RD= V_RD/((sqrt(3*Cd0/K))/(2*Cd0))  ;            % ft/s
Endurcance= h/RD;                                % 958.3515 sec

% Sea Level values:
TSL = 518.69;                                    % Rankine
PSL = 2116.199414;                               % lbs/ft^2
RhoSL = 0.0023769;                          % slugs/ft^3

while h > 0

    tspan = [t t+1];
    i=1;
    X = [x;h;v;gamma;Rho_10000];
    Time(i)= t;

    % Calling ODE45:
    [F] = ode45(@ D,tspan, X)

    % Hight Varying Parameters:
    T = TSL - 0.00356616*h;
    P = (1.137193514E-11)*(T)^5.2560613;
    Rho = (RhoSL * TSL / PASL)*(P / T);
    a = 49.0214 * (T)^.5;

    H_Del(i) = (-Cd_max/Cl_max)*(plotted_x(i))+10000;
    V_Del(i) = sqrt((2*W)/(Rho*S*Cl_max));
    Gamma_Del(i) = -Cd_max/Cl_max;
    i= i+1;

    X = [ x ; H_Del(i); V_Del(i); Gamma_Del(i); Rho];

end



% Plots:

%1
figure (1)
plot(F(:,1),'-r',F(:,3),'-b')
title('Velocity vs Distance');
xlabel('x (ft)');
ylabel('v (ft/s)');


%2 
Figure (2)
plot(F(:,1),'-r',F(:,2),'-b')
title('Altitude vs Distance ');
xlabel('x (ft)');
ylabel('h (ft)');


%3
figure (3)
plot(F(:,1),'-r',F(:,4),'-b')
title('Gamma vs Distance');
xlabel('x (ft)');
ylabel('Gamma (rad)');

%4
figure (4)
plot(t,'-r',F(:,3),'-b')
title('Velocity vs Time');
xlabel(' t (s)');
ylabel('v  (ft/s)');

%5
figure (5)
plot(t,'-r',F(:,1),'-b')
title(' Distance vs Time ');
xlabel('t (s)');
ylabel('x (ft)');


%6
figure (6)
plot(t,'-r',F(:,4),'-b')
title('Gamma vs Time ');
xlabel('t (s)');
ylabel('Gamma (rad)');

%7 
figure (7)
plot(t,'-r',F(:,3),'-b')
title('Velocity vs Time');
xlabel('t (s)');
ylabel('V (ft/s)');

这是函数

function [F] = D(X)

% Constants:
S = 232;                                         % ft^2
Cd0 = 0.02;
K = 0.07;
W = 11000;                                       % lbf
Cl_max = sqrt(Cd0/K);
Cd_max = 2*K*Cl_max^2;
Rho_10000 = .001756;                             % slugs/ ft^3

% Initial conditions:
t = 0;                                           % Sec
x = 0;                                           % ft
h = 10000;                                       % ft
v = sqrt((2*W)/(Rho_10000*S*Cl_max));            % ft/s
gamma = -Cd_max/Cl_max; 
g= 32.2;                                         % ft/s^2

% Equations:
X_Pr = X(3)*cos(X(4));
H_Pr = X(3)*sin(X(4));
V_Pr = (-32.2./W)*(0.5*X(5)*Cd_max*S*X(3)^2 + W*sin(X(4)));
G_Pr = (32.2./(W*X(3)))*(0.5*X(5)*Cl_max*S*X(3)^2 - W*cos(X(4)));

F = [X_Pr;H_Pr;V_Pr;G_Pr];

我的 Matlab 不是很好,但我已经尽力了!我去找我的教授寻求帮助,但他们说他们太忙了。我什至跟踪了我认识的每一位学长,他们都说他们不知道该怎么做。我的下一个项目马上就要分配了,我想如果我不能做这个,那么我就不能做下一个。

您的代码产生以下错误:

Error using main>D

Too many input arguments.

这意味着 ode45 会尝试使用过多的输入参数调用您提供的函数 D。您应该在 ode45 文档中检查所需的 odefun 格式:dydt = odefun(t,y)

因此,您应该将 D 的函数声明更改为

function [F] = D(t, X)

这应该可以解决您的第一个问题,但弹出以下错误:

D returns a vector of length 4, but the length of initial conditions vector is 5. The vector returned by D and the initial conditions vector must have the same number of elements.

同样,您应该查看 ode45 文档。您的函数应该 return 所有输入变量的导数 XF = dX/dt。您还应该 return 第五个元素 Rho_10000 的导数。

接下来,我得到一些关于未定义变量的错误,例如 PASL。可能是因为您没有 post 完整代码。

除了错误之外,您真的应该再次检查您的代码。您已经编写了一个无限 while 循环 while h > 0。您永远不会在循环中更改 h,也不会在循环中使用 ode45 的输出。此外,您总是在循环开始时覆盖 iX 值,这可能不是您想要的。

这不是您问题的完整答案,但我希望您能够继续 post 更小、定义明确的问题,而不是一个非常困难的大问题完整回答。