PROPT Nondifferentiable system
From TomWiki
Jump to navigationJump to search
This page is part of the PROPT Manual. See PROPT Manual. |
Global Optimization of Chemical Processes using Stochastic Algorithms 1996, Julio R. Banga, Warren D. Seider
Case Study IV: Optimal control of a nondifferentiable system
Problem description
This problem has been studied by Thomopoulos and Papadakis who report convergence difficulties using several optimization algorithms and by Luus using Iterative Dynamic Programming. The optimal control problem is:
Find u(t) to minimize:
Subject to:
with
where U(t-alpha) is the unit function such that U = 0 for t - alpha < 0 and U = 1 for t - alpha > 0. Hence d is a rectangular pulse of magnitude 100 from t=0.5 until t=0.6. These authors consider t_f = 2.0s and the initial conditions:
% Copyright (c) 2007-2008 by Tomlab Optimization Inc.
Problem setup
toms t1
p1 = tomPhase('p1', t1, 0, 0.5, 30);
toms t2
p2 = tomPhase('p2', t2, 0.5, 0.1, 20);
toms t3
p3 = tomPhase('p3', t3, 0.6, 2-0.6, 50);
x1p1 = tomState(p1,'x1p1');
x2p1 = tomState(p1,'x2p1');
x3p1 = tomState(p1,'x3p1');
up1 = tomControl(p1,'up1');
x1p2 = tomState(p2,'x1p2');
x2p2 = tomState(p2,'x2p2');
x3p2 = tomState(p2,'x3p2');
up2 = tomControl(p2,'up2');
x1p3 = tomState(p3,'x1p3');
x2p3 = tomState(p3,'x2p3');
x3p3 = tomState(p3,'x3p3');
up3 = tomControl(p3,'up3');
% Initial guess
x0 = {icollocate(p1,{x1p1 == 0;x2p1 == 0;x3p1 == 0})
icollocate(p2,{x1p2 == 0;x2p2 == 0;x3p2 == 0})
icollocate(p3,{x1p3 == 0;x2p3 == 0;x3p3 == 0})
collocate(p1,up1==-4-8*t1/0.5)
collocate(p2,up2==-12)
collocate(p3,up3==-12+14*t3/2)};
% Box constraints
cbox = {
icollocate(p1,{-100 <= x1p1 <= 100
-100 <= x2p1 <= 100
-100 <= x3p1 <= 100})
icollocate(p2,{-100 <= x1p2 <= 100
-100 <= x2p2 <= 100
-100 <= x3p2 <= 100})
icollocate(p3,{-100 <= x1p3 <= 100
-100 <= x2p3 <= 100
-100 <= x3p3 <= 100})
collocate(p1,-15 <= up1 <= 2)
collocate(p2,-15 <= up2 <= 2)
collocate(p3,-15 <= up3 <= 2)};
% Boundary constraints
cbnd = {initial(p1,{x1p1 == 0;x2p1 == 0;x3p1 == 0})
final(p3,x3p3) <= 60};
% ODEs and path constraints
ceq = {collocate(p1,{
dot(p1,x1p1) == x2p1
dot(p1,x2p1) == -x1p1-x2p1+up1
dot(p1,x3p1) == 5*x1p1.^2+2.5*x2p1.^2+0.5*up1.^2})
collocate(p2,{
dot(p2,x1p2) == x2p2
dot(p2,x2p2) == -x1p2-x2p2+up2+100
dot(p2,x3p2) == 5*x1p2.^2+2.5*x2p2.^2+0.5*up2.^2})
collocate(p3,{
dot(p3,x1p3) == x2p3
dot(p3,x2p3) == -x1p3-x2p3+up3
dot(p3,x3p3) == 5*x1p3.^2+2.5*x2p3.^2+0.5*up3.^2})};
% Objective
objective = final(p3,x3p3);
% Link phase
link = {final(p1,x1p1) == initial(p2,x1p2)
final(p1,x2p1) == initial(p2,x2p2)
final(p1,x3p1) == initial(p2,x3p2)
final(p2,x1p2) == initial(p3,x1p3)
final(p2,x2p2) == initial(p3,x2p3)
final(p2,x3p2) == initial(p3,x3p3)};
Solve the problem
options = struct;
options.name = 'Nondiff System';
constr = {cbox, cbnd, ceq, link};
solution = ezsolve(objective, constr, x0, options);
t = subs(collocate(p1,t1),solution);
t = [t;subs(collocate(p2,t2),solution)];
t = [t;subs(collocate(p3,t3),solution)];
x1 = subs(collocate(p1,x1p1),solution);
x1 = [x1;subs(collocate(p2,x1p2),solution)];
x1 = [x1;subs(collocate(p3,x1p3),solution)];
x2 = subs(collocate(p1,x2p1),solution);
x2 = [x2;subs(collocate(p2,x2p2),solution)];
x2 = [x2;subs(collocate(p3,x2p3),solution)];
x3 = subs(collocate(p1,x3p1),solution);
x3 = [x3;subs(collocate(p2,x3p2),solution)];
x3 = [x3;subs(collocate(p3,x3p3),solution)];
u = subs(collocate(p1,up1),solution);
u = [u;subs(collocate(p2,up2),solution)];
u = [u;subs(collocate(p3,up3),solution)];
Problem type appears to be: lpcon Time for symbolic processing: 0.34745 seconds Starting numeric solver ===== * * * =================================================================== * * * TOMLAB - TOMLAB Development license 999007. Valid to 2011-12-31 ===================================================================================== Problem: --- 1: Nondiff System f_k 58.065028469582899000 sum(|constr|) 0.000030760924864613 f(x_k) + sum(|constr|) 58.065059230507764000 f(x_0) 0.000000000000000000 Solver: snopt. EXIT=0. INFORM=1. SNOPT 7.2-5 NLP code Optimality conditions satisfied FuncEv 1 ConstrEv 45 ConJacEv 45 Iter 38 MinorIter 1056 CPU time: 0.530403 sec. Elapsed time: 0.534000 sec.
Plot result
subplot(2,1,1)
plot(t,x1,'*-',t,x2,'*-',t,x3,'*-');
legend('x1','x2','x3');
title('Nondiff System state variable');
subplot(2,1,2)
plot(t,u,'+-');
legend('u');
title('Nondiff System control');