Goddard Rocket, Maximum Ascent

Benchmarking Optimization Software with COPS Elizabeth D. Dolan and Jorge J. More ARGONNE NATIONAL LABORATORY

Contents

Problem Formulation

Find u(t) over t in [0; T ] to minimize

$$ J = h(T) $$

subject to:

$$ \frac{dv}{dt} = \frac{1}{m}*(T-D)-g $$

$$ \frac{dh}{dt} = v $$

$$ \frac{dm}{dt} = -\frac{T}{c} $$

$$ D = D_0*v^2*exp^{-beta*\frac{h-h_0}{h_0}} $$

$$ g = g_0*\frac{h_0}{h}^{2} $$

$$ m(0) = 1 $$

$$ m(T) = 0.6 $$

$$ v(0) = 0 $$

$$ h(0) = 1 $$

$$ g_0  = 1 $$

$$ 0 <= u <= 3.5 $$

$$ D_0 = 0.5*620 $$

$$ c = 0.5 $$

$$ beta = 500 $$

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

Problem setup

toms t
toms tf

Solve the problem, using a successively larger number collocation points

for n=[20 50 100]
    p = tomPhase('p', t, 0, tf, n);
    setPhase(p);
    tomStates v h m
    tomControls T

    % Initial guess
    if n==20
        x0 = {tf == 1
            icollocate({v == 620; h == 1
            m == 1-0.4*t/tf})
            collocate(T == 0)};
    else
        x0 = {tf == tfopt
            icollocate({v == vopt; h == hopt
            m == mopt})
            collocate(T == Topt)};
    end

    % Box constraints
    cbox = {0.1 <= tf <= 1
        icollocate({
        0 <= v; 1 <= h
        0.6 <= m <= 1
        0   <= T <= 3.5})};

    % Boundary constraints
    cbnd = {initial({v == 0; h == 1; m == 1})
        final({m == 0.6})};

    beta = 500;
    D    = 0.5*620*v.^2.*exp(-beta*h);
    g    = 1./h.^2;
    c    = 0.5;

    % ODEs and path constraints
    ceq = collocate({dot(v) == (T-D)./m-g
        dot(h) == v; dot(m) == -T/c});

    % Objective
    objective = -final(h);

Solve the problem

    options = struct;
    options.name = 'Goddard Rocket';
    solution = ezsolve(objective, {cbox, cbnd, ceq}, x0, options);

    % Optimal v and more to use as starting guess
    vopt = subs(v, solution);
    hopt = subs(h, solution);
    mopt = subs(m, solution);
    Topt = subs(T, solution);
    tfopt = subs(tf, solution);
Problem type appears to be: lpcon
===== * * * =================================================================== * * *
TOMLAB - Tomlab Optimization Inc. Development license  999001. Valid to 2010-02-05
=====================================================================================
Problem: ---  1: Goddard Rocket                 f_k      -1.025133414041156300
                                       sum(|constr|)      0.000002519458375471
                              f(x_k) + sum(|constr|)     -1.025130894582780800
                                              f(x_0)     -0.999999999999998220

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

FuncEv    1 ConstrEv   41 ConJacEv   41 Iter   23 MinorIter 1141
CPU time: 0.171875 sec. Elapsed time: 0.172000 sec. 
Problem type appears to be: lpcon
===== * * * =================================================================== * * *
TOMLAB - Tomlab Optimization Inc. Development license  999001. Valid to 2010-02-05
=====================================================================================
Problem: ---  1: Goddard Rocket                 f_k      -1.025311927458321800
                                       sum(|constr|)      0.000016009288875488
                              f(x_k) + sum(|constr|)     -1.025295918169446300
                                              f(x_0)     -1.025133225224280400

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

FuncEv    1 ConstrEv   23 ConJacEv   23 Iter   14 MinorIter  435
CPU time: 0.250000 sec. Elapsed time: 0.250000 sec. 
Problem type appears to be: lpcon
===== * * * =================================================================== * * *
TOMLAB - Tomlab Optimization Inc. Development license  999001. Valid to 2010-02-05
=====================================================================================
Problem: ---  1: Goddard Rocket                 f_k      -1.025328777109888700
                                       sum(|constr|)      0.000000000010407547
                              f(x_k) + sum(|constr|)     -1.025328777099481200
                                              f(x_0)     -1.025311927458318500

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

FuncEv    1 ConstrEv   12 ConJacEv   12 Iter    7 MinorIter  533
CPU time: 0.671875 sec. Elapsed time: 0.703000 sec. 
end

t = subs(collocate(t),solution);
v = subs(collocate(vopt),solution);
h = subs(collocate(hopt),solution);
m = subs(collocate(mopt),solution);
T = subs(collocate(Topt),solution);

Plot result

subplot(2,1,1)
plot(t,v,'*-',t,h,'*-',t,m,'*-');
legend('v','h','m');
title('Goddard Rocket state variables');

subplot(2,1,2)
plot(t,T,'+-');
legend('T');
title('Goddard Rocket control');