Этот пример показывает, как использовать встроенный решатель QP, чтобы реализовать пользовательский алгоритм MPC управления, который поддерживает генерацию кода C в MATLAB.
Модель объекта управления является системой пространства состояний дискретного времени, и это - нестабильный разомкнутый цикл. Мы принимаем, что все состояния объекта измеримы. Поэтому мы избегаем потребности в разработке средства оценки состояния, которое выходит за рамки этого примера.
A = [1.1 2; 0 0.95];
B = [0; 0.0787];
C = [-1 1];
D = 0;
Ts = 1;
sys = ss(A,B,C,D,Ts);
x0 = [0.5;-0.5]; % initial states at [0.5 -0.5]
Разработайте неограниченный LQR с выходным взвешиванием. Этот контроллер служит базовой линией, чтобы соответствовать пользовательскому алгоритму MPC. Закон о надзоре LQ является u (k) =-K_lqr*x (k).
Qy = 1; R = 0.01; K_lqr = lqry(sys,Qy,R);
Запустите симуляцию с начальными состояниями в [0.5 - 0.5]. Ответ с обратной связью стабилен.
t_unconstrained = 0:1:10; u_unconstrained = zeros(size(t_unconstrained)); Unconstrained_LQR = tf([-1 1])*feedback(ss(A,B,eye(2),0,Ts),K_lqr); lsim(Unconstrained_LQR,'-',u_unconstrained,t_unconstrained,x0); hold on;
Разработайте пользовательский контроллер MPC с терминальным весом, примененным на последнем шаге прогноза.
Предсказанные последовательности состояния, X (k), сгенерированный линейной моделью и входной последовательностью, U (k), могут быть сформулированы как: X (k) = M*x (k) + CONV*U (k). В этом примере используйте четыре шага прогноза (N = 4).
M = [A;A^2;A^3;A^4]; CONV = [B zeros(2,1) zeros(2,1) zeros(2,1);... A*B B zeros(2,1) zeros(2,1);... A^2*B A*B B zeros(2,1);... A^3*B A^2*B A*B B];
Целевая функция MPC является J (k) = сумма (x (k) '*Q*x (k) + u (k) '*R*u (k) + x (k+N) '*Q_bar*x (k+N)). Чтобы гарантировать, что целевая функция MPC имеет ту же квадратичную стоимость как бесконечный горизонт квадратичная стоимость, используемая LQR, терминальный вес, Q_bar получен путем решения следующего уравнения Ляпонова:
Q = C'*C; Q_bar = dlyap((A-B*K_lqr)', Q+K_lqr'*R*K_lqr);
Преобразуйте проблему MPC в стандартную проблему QP, которая имеет целевую функцию J (k) = U (k) '*H*U (k) + 2*x (К) '*F '*U (k).
Q_hat = blkdiag(Q,Q,Q,Q_bar); R_hat = blkdiag(R,R,R,R); H = CONV'*Q_hat*CONV + R_hat; F = CONV'*Q_hat*M;
Когда нет никаких ограничений, оптимальная предсказанная входная последовательность U (k) сгенерированный контроллером MPC является-K*x, где K = inv (H) *F.
K = H\F;
На практике только первое перемещение управления u (k) =-K_mpc*x (k) применяется к объекту (отступающий управление горизонтом).
K_mpc = K(1,:);
Запустите симуляцию с начальными состояниями в [0.5 - 0.5]. Ответ с обратной связью стабилен.
Unconstrained_MPC = tf([-1 1])*feedback(ss(A,B,eye(2),0,Ts),K_mpc); lsim(Unconstrained_MPC,'*',u_unconstrained,t_unconstrained,x0) legend show
LQR и контроллеры MPC приводят к тому же результату, потому что законы о надзоре являются тем же самым.
K_lqr K_mpc
K_lqr =
4.3608 18.7401
K_mpc =
4.3608 18.7401
Ограничьте контроллер вывод, u (k), чтобы быть между-1 и 1. Контроллер LQR генерирует медленный и колебательный ответ с обратной связью из-за насыщения.
x = x0; t_constrained = 0:40; for ct = t_constrained uLQR(ct+1) = -K_lqr*x; uLQR(ct+1) = max(-1,min(1,uLQR(ct+1))); x = A*x+B*uLQR(ct+1); yLQR(ct+1) = C*x; end figure subplot(2,1,1) plot(t_constrained,uLQR) xlabel('time') ylabel('u') subplot(2,1,2) plot(t_constrained,yLQR) xlabel('time') ylabel('y') legend('Constrained LQR')
Одно из главных преимуществ использования контроллера MPC - то, что оно обрабатывает ограничения ввода и вывода явным образом путем решения задачи оптимизации в каждом интервале управления.
Используйте встроенный решатель KWIK QP, mpcqpsolver, чтобы реализовать пользовательский контроллер MPC, разработанный выше. Ограничительные матрицы заданы как Ac*x> =b0.
Ac = -[1 0 0 0;... -1 0 0 0;... 0 1 0 0;... 0 -1 0 0;... 0 0 1 0;... 0 0 -1 0;... 0 0 0 1;... 0 0 0 -1]; b0 = -[1;1;1;1;1;1;1;1];
Функция mpcqpsolver требует, чтобы первый вход был инверсией нижнего треугольного разложения Холесского матрицы H Гессиана.
L = chol(H,'lower');
Linv = L\eye(size(H,1));
Запустите симуляцию путем вызова mpcqpsolver на каждом шаге симуляции. Первоначально все неравенства неактивны (холодный запуск).
x = x0; iA = false(size(b0)); opt = mpcqpsolverOptions; opt.IntegrityChecks = false; for ct = t_constrained [u, status, iA] = mpcqpsolver(Linv,F*x,Ac,b0,[],zeros(0,1),iA,opt); uMPC(ct+1) = u(1); x = A*x+B*uMPC(ct+1); yMPC(ct+1) = C*x; end figure subplot(2,1,1) plot(t_constrained,uMPC) xlabel('time') ylabel('u') subplot(2,1,2) plot(t_constrained,yMPC) xlabel('time') ylabel('y') legend('Constrained MPC')
Контроллер MPC производит ответ с обратной связью с более быстрым временем установления и меньше колебания.
mpcqpsolver может использоваться в блоке MATLAB function, чтобы обеспечить симуляцию и генерацию кода в окружении Simulink.
mdl = 'mpc_customqp';
open_system(mdl)
Блоком Custom MPC Controller является блок MATLAB function. Дважды кликните его, чтобы исследовать код MATLAB. Начиная с Linv F, Ac, матрицы b0 и структура opt являются постоянными, они передаются в блок MATLAB function как параметры.
Запустите симуляцию в Simulink. Закрытые ответы LQR и контроллеров MPC идентичны своим дубликатам в симуляции MATLAB.
open_system([mdl '/u_lqr']) open_system([mdl '/y_lqr']) open_system([mdl '/u_mpc']) open_system([mdl '/y_mpc']) sim(mdl)
mpcqpsolver поддерживает генерацию кода C с MATLAB Coder. Примите, что у вас есть функциональный mycode, который совместим со стандартом генерации кода.
function [x,iter,iA1,lam] = mycode() %#codegen n = 5; m = 10; q = 2; H = diag(10*rand(n,1)); f = randn(n,1); A = randn(m,n); b = randn(m,1); Aeq = randn(q,n); beq = randn(q,1); Linv = chol(H,'lower')\eye(n); iA = false(m,1); Opt = mpcqpsolverOptions(); [x,iter,iA1,lam] = mpcqpsolver(Linv,f,-A,-b,Aeq,beq,iA,Opt);
Можно использовать следующую команду, чтобы сгенерировать код С с MATLAB Coder:
fun = 'mycode'; Cfg = coder.config('mex'); % or 'lib', 'dll', etc. codegen('-config',Cfg,fun,'-o',fun);
Этот пример вдохновлен примечаниями лекции профессора Марка Кэннона для Образцового Прогнозирующего класса Управления в Оксфордском университете. Модель объекта управления - тот же используемый в Примере 2.1 в разделе "Prediction and optimization".
bdclose(mdl)
mpcqpsolver | mpcqpsolverOptions