PROPT Minimum Climb Time (English Units)

From TomWiki
Jump to navigationJump to search

Notice.png

This page is part of the PROPT Manual. See PROPT Manual.

Problem description

Example about climbing (increase altitude)

% Copyright (c) 2007-2008 by Tomlab Optimization Inc.

Problem setup

alt0     = 0;
altf     = 65600;
speed0   = 424.26;
speedf   = 968.148;
fpa0     = 0;
fpaf     = 0;
mass0    = 42000/32.208;

altmin   = 0;
altmax   = 69000;
speedmin = 10;
speedmax = 3000;
fpamin   = -40*pi/180;
fpamax   = -fpamin;
massmin  = 50/32.208;

toms t t_f
p = tomPhase('p',t,0,t_f,50);
setPhase(p);

% Altitude, speed, flight path angle, mass
tomStates h v fpa m

% Angle of Attack
tomControls aalpha

guess = {
    t_f == 300;
    icollocate({
    h == alt0 + t/t_f*(altf-alt0);
    v == speed0 + t/t_f*(speedf-speed0);
    fpa == 10*pi/180;
    m == mass0;
    })};

cbox = {
    100 <= t_f <= 800;
    collocate(-pi*20/180 <= aalpha <= pi/20*180)
    icollocate(altmin <= h <= altmax)
    icollocate(speedmin <= v <= speedmax)
    icollocate(fpamin <= fpa <= fpamax)
    icollocate(massmin <= m <= mass0)
    };

bnd = {
    initial(h) == alt0;
    initial(v) == speed0;
    initial(fpa) == fpa0;
    initial(m) == mass0;
    final(h) == altf;
    final(v) == speedf;
    final(fpa) == fpaf;
    };

% US1976 data
hTab = (-2000:2000:86000);
rhoTab = [1.478 1.225 1.007 0.8193 0.6601 0.5258 0.4135 0.3119 ...
    0.2279 0.1665 0.1216 0.08891 0.06451 0.04694 0.03426 0.02508 ...
    0.01841 0.01355 0.009887 0.007257 0.005366 0.003995 0.002995 ...
    0.002259 0.001714 0.001317 0.001027 0.0008055 0.0006389 0.0005044 ...
    0.0003962 0.0003096 0.0002407 0.000186 0.0001429 0.0001091 ...
    8.281e-005 6.236e-005 4.637e-005 3.43e-005 2.523e-005 1.845e-005 ...
    1.341e-005 9.69e-006 6.955e-006];
sosTab = [347.9 340.3 332.5 324.6 316.5 308.1 299.5 295.1 295.1 ...
    295.1 295.1 295.1 296.4 297.7 299.1 300.4 301.7 303 306.5 310.1 ...
    313.7 317.2 320.7 324.1 327.5 329.8 329.8 328.8 325.4 322 318.6 ...
    315.1 311.5 308 304.4 300.7 297.1 293.4 290.7 288 285.3 282.5 ...
    279.7 276.9 274.1];

Mtab   = [0; 0.2; 0.4; 0.6; 0.8; 1; 1.2; 1.4; 1.6; 1.8];
alttab = [0 5000 10000 15000 20000 25000 30000 40000 50000 70000];
Ttab = 1000*[24.2 24.0  20.3 17.3 14.5 12.2 10.2 5.7 3.4 0.1;
    28.0 24.6 21.1 18.1 15.2 12.8 10.7 6.5 3.9 0.2;
    28.3 25.2 21.9 18.7 15.9 13.4 11.2 7.3 4.4 0.4;
    30.8 27.2 23.8 20.5 17.3 14.7 12.3 8.1 4.9 0.8;
    34.5 30.3 26.6 23.2 19.8 16.8 14.1 9.4 5.6 1.1;
    37.9 34.3 30.4 26.8 23.3 19.8 16.8 11.2 6.8 1.4;
    36.1 38.0 34.9 31.3 27.3 23.6 20.1 13.4 8.3 1.7;
    36.1 36.6 38.5 36.1 31.6 28.1 24.2 16.2 10.0 2.2;
    36.1 35.2 42.1 38.7 35.7 32.0 28.1 19.3 11.9 2.9;
    36.1 33.8 45.7 41.3 39.8 34.6 31.1 21.7 13.3 3.1];
M2     = [0 0.4 0.8 0.9 1.0 1.2 1.4 1.6 1.8];
Clalphatab = [3.44 3.44 3.44 3.58 4.44 3.44 3.01 2.86 2.44];
CD0tab    = [0.013 0.013 0.013 0.014 0.031 0.041 0.039 0.036 0.035];
etatab    = [0.54 0.54 0.54 0.75 0.79 0.78 0.89 0.93 0.93];
M         = Mtab;
alt       = alttab;
Re        = 20902900;
mmu        = 0.14076539e17;
S         = 530;
g0        = 32.208;
ISP       = 1600;
H         = 23800;
rho0      = 0.002378;

rho = interp1(hTab,rhoTab,h*0.3048,'pchip')*0.001941;
sosI = interp1(hTab,sosTab,h*0.3048,'pchip')./0.3048;
Mach = v/sosI;

CD0       = interp1(M2,CD0tab,Mach,'pchip');
Clalpha   = interp1(M2,Clalphatab,Mach,'pchip');
eta       = interp1(M2,etatab,Mach,'pchip');

T = interp2(alttab, Mtab, Ttab, h, Mach, 'spline');
CD      = CD0 + eta.*Clalpha.*aalpha.^2;
CL      = Clalpha.*aalpha;
dynpres = 0.5.*rho.*v.^2;
D       = dynpres.*S.*CD;
L       = dynpres.*S.*CL;

equ = collocate({
    dot(h) == v.*sin(fpa);
    dot(v) == ((T.*cos(aalpha)-D)./m - mmu.*sin(fpa)./(Re+h).^2);
    dot(fpa) == (T.*sin(aalpha)+L)./(m.*v)+cos(fpa).*(v./(Re+h)-mmu./(v.*(Re+h).^2));
    dot(m) == -T./(g0.*ISP);
    });

options = struct;
options.name = 'Minimum Time to Climb (English)';
options.scale = 'auto';

Solve the problem

ezsolve(t_f,{cbox,bnd,equ},guess,options);
Problem type appears to be: lpcon
Auto-scaling
Time for symbolic processing: 0.79862 seconds
Starting numeric solver
===== * * * =================================================================== * * *
TOMLAB - TOMLAB Development license  999007. Valid to 2011-12-31
=====================================================================================
Problem: ---  1: Minimum Time to Climb (English)  f_k     318.348640855936940000
                                         sum(|constr|)      0.000000000187118413
                                f(x_k) + sum(|constr|)    318.348640856124060000
                                                f(x_0)    300.000000000000000000

Solver: snopt.  EXIT=0.  INFORM=1.
SNOPT 7.2-5 NLP code
Optimality conditions satisfied

FuncEv    1 ConstrEv   91 ConJacEv   90 Iter   49 MinorIter 2246
CPU time: 4.118426 sec. Elapsed time: 1.062000 sec. 

Plot result

subplot(2,1,1)
ezplot([1e-3*h, 1e-1*v, fpa*180/pi, 1e-1*m])
legend('1e-3*h', '1e-1*v', 'fpa*180/pi', '1e-1*m');
title('Minimum Time to Climb (English Units) state variables');

subplot(2,1,2)
ezplot(aalpha);
title('Minimum Time to Climb (English Units) control');

MinimumClimbEng 01.png