本文主要基于以下参考:

[1] John T. Betts. Survey of Numerical Methods for Trajectory Optimization.

[2] Anil V. Rao. A Survey of Numerical Methods For Optimal Control.

[3] John T. Betts. Practical Methods for Optimal Control and Estimation Using Nonlinear Programming 2nd.

[4] A E. Bryson. Applied Optimal Control.

[5] KIRK. Optimal Control Theory: An Introduction.

[6] Matthew Kelly. An Introduction to Trajectory Optimization: How to Do Your Own Direct Collocation.

上一篇文章中,最优控制问题的终端状态是完全指定的,而没有用到关于协态终端的关系式。

在[4] 2.4节,讨论了只对一部分终端状态带约束的最优控制问题。

考虑引入乘子函数 λ(t) ,通过该函数将系统微分方程引入到优化指标中,得到如下:

定义标量函数 H 为:

此时,公式(1)变为:

对上式最后一项作分部积分,可以得到:

考虑由于控制变量的变分在固定初末时间条件下导致的优化指标函数的变分如下:

上式中,由于状态变量与控制变量满足系统微分方程,因此控制变量的变化会导致状态的变化。由于上述原因,控制变量的变分与状态变分之间满足某些关系,确定这种关系比较复杂。

下面给出一个例子:

最优控制问题:

边界条件满足:

首先类似之前,利用变分法写出必要最优条件:

边界条件为:

%% bvp solve optimal control

% given the initial guess
solinit = bvpinit(linspace(0, 5,10),[1 0 0 1]);
% solve the problem
sol = bvp4c(@sysdyn,@sysbc,solinit);

% eval points
t = linspace(0, 5, 50);
% eval
y = deval(sol,t);
% plot states
figure(1);
plot(t, y(1,:), 'o-', t, y(2,:),'o-');
grid;
xlabel('time/s');
legend('state 1', 'state 2');
% plot costates
figure(2);
plot(t, y(3,:), 'o-',t, y(4,:),'o-');
grid;
xlabel('time/s');
legend('costate 1', 'costate 2');
% plot control
figure(3);
plot(t, -y(3,:)./y(4,:),'o-');
grid;
xlabel('time/s');
legend('control');


% plot H funcs
u = -y(3,:) ./ y(4,:);
x1 = y(1,:);    x2 = y(2,:);
cx1 = y(3,:);   cx2 = y(4,:);

H = sum([cx1;cx2] .* [0.5*x1+u; x1.^2+0.5*u.^2]);
figure(4);
plot(t, H,'o-');
grid;
xlabel('time/s');
legend('H funcs');

% 状态依次: x1,x2,lambda1,lambda2
function dydx = sysdyn(t,y)
    dydx = [ 0.5*y(1)+(-y(3)/y(4)); y(1)^2+0.5*(-y(3)/y(4))^2; -(0.5*y(3)+2*y(1)*y(4)); 0];
end

function res = sysbc(y0,y1)
    res = [ y0(1)-1;y0(2);y1(1)-0.5; y1(4)-1 ];
end