Van der Pol Oscillator

Restricted second order information for the solution of optimal control problems using control vector parameterization. 2002, Eva Balsa Canto, Julio R. Banga, Antonio A. Alonso Vassilios S. Vassiliadis

Case Study 6.1: van der Pol oscillator

This case study has been studied by several authors, for example Morison, Gritsis, Vassiliadis and Tanartkit and Biegler.

Contents

Problem Description

The dynamic optimization problem is to minimize:

$$ J = x_3(t_f) $$

subject to:

$$ \frac{dx_1}{dt} = (1-x_2^2)*x_1-x_2+u $$

$$ \frac{dx_2}{dt} = x_1 $$

$$ \frac{dx_3}{dt} = x_1^2+x_2^2+u^2 $$

$$ -0.3 <= u <= 1.0 $$

$$ x(t_0) = [0 \ 1 \ 0]' $$

$$ t_f = 5 $$

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

Problem setup

toms t
p = tomPhase('p', t, 0, 5, 60);
setPhase(p);

tomStates x1 x2 x3
tomControls u

% Initial guess
x0 = {icollocate({x1 == 0; x2 == 1; x3 == 0})
    collocate(u == -0.01)};

% Box constraints
cbox = {-10  <= icollocate(x1) <= 10
    -10  <= icollocate(x2) <= 10
    -10  <= icollocate(x3) <= 10
    -0.3 <= collocate(u)   <= 1};

% Boundary constraints
cbnd = initial({x1 == 0; x2 == 1; x3 == 0});

% ODEs and path constraints
ceq = collocate({dot(x1) == (1-x2.^2).*x1-x2+u
    dot(x2) == x1; dot(x3) == x1.^2+x2.^2+u.^2});

% Objective
objective = final(x3);

Solve the problem

options = struct;
options.name = 'Van Der Pol';
solution = ezsolve(objective, {cbox, cbnd, ceq}, x0, options);
t  = subs(collocate(t),solution);
x1 = subs(collocate(x1),solution);
x2 = subs(collocate(x2),solution);
x3 = subs(collocate(x3),solution);
u  = subs(collocate(u),solution);

Problem type appears to be: lpcon
===== * * * =================================================================== * * *
TOMLAB - Tomlab Optimization Inc. Development license  999001. Valid to 2010-02-05
=====================================================================================
Problem: ---  1: Van Der Pol                    f_k       2.867259538084702400
                                       sum(|constr|)      0.000000020745141049
                              f(x_k) + sum(|constr|)      2.867259558829843500
                                              f(x_0)      0.000000000000000000

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

FuncEv    1 ConstrEv   26 ConJacEv   26 Iter   23 MinorIter  331
CPU time: 0.578125 sec. Elapsed time: 0.578000 sec.

Plot result

subplot(2,1,1)
plot(t,x1,'*-',t,x2,'*-',t,x3,'*-');
legend('x1','x2','x3');
title('vanDerPol state variables');

subplot(2,1,2)
plot(t,u,'+-');
legend('u');
title('vanDerPol control');