Модель сервопривода постоянного тока

Модель сервопривода постоянного тока представляет собой линейную динамическую систему, описанную в [1]. plant - непрерывная модель состояния и пространства двигателя. tau - максимально допустимый крутящий момент, используемый в качестве выходного ограничения.

[plant,tau] = mpcmotormodel;

Проектирование контроллера MPC

Установка имеет один вход, входное напряжение двигателя. Контроллер MPC использует этот вход в качестве управляемой переменной (MV). Установка имеет два выхода, угловое положение двигателя и крутящий момент вала. Угловое положение - измеряемый выходной сигнал (MO), и момент вращения вала не измеряется (UO).

plant = setmpcsignals(plant,'MV',1,'MO',1,'UO',2);

Наложение зависимости на управляемую переменную между +/- 220 вольт. Поскольку входы и выходы установки имеют разные порядки величин, для облегчения настройки используйте масштабные коэффициенты. Типичными вариантами масштабного коэффициента являются верхний/нижний предел или рабочий диапазон.

MV = struct('Min',-220,'Max',220,'ScaleFactor',440);

Угловое положение не ограничено. Укажите верхнюю и нижнюю границы крутящего момента вала в течение первых трех шагов горизонта прогнозирования. Чтобы определить эти границы, используйте tau.

OV = struct('Min',{-Inf, [-tau;-tau;-tau;-Inf]},...
    'Max',{Inf, [tau;tau;tau;Inf]},'ScaleFactor',{2*pi, 2*tau});

Задача управления заключается в достижении нулевой ошибки отслеживания углового положения. Так как имеется только одна управляемая переменная, допустите плавание крутящего момента вала в пределах его ограничения, установив его настроечный вес равным нулю.

Weights = struct('MV',0,'MVRate',0.1,'OV',[0.1 0]);

Укажите время и горизонты выборки и создайте контроллер MPC с помощью plant в качестве прогностической модели.

Ts = 0.1;           % Sample time
p = 10;             % Prediction horizon
m = 2;              % Control horizon
mpcobj = mpc(plant,Ts,p,m,Weights,MV,OV);

Моделирование в Simulink с помощью встроенного решателя QP

Для выполнения оставшегося примера требуется Simulink.

if ~mpcchecktoolboxinstalled('simulink')
    disp('Simulink is required to run this example.')
    return
end

Откройте модель Simulink, которая моделирует управление серводвигателем постоянного тока по замкнутому контуру с помощью контроллера MPC. По умолчанию MPC использует встроенный решатель QP, использующий алгоритм KWIK.

mdl = 'mpc_customQPcodegen';
open_system(mdl)

Запустить моделирование

sim(mdl)

-->Converting model to discrete time.
   Assuming no disturbance added to measured output channel #1.
-->The "Model.Noise" property of the "mpc" object is empty. Assuming white noise on each measured output channel.

Храните входные и выходные сигналы завода в рабочем пространстве MATLAB.

uKWIK = u;
yKWIK = y;

Моделирование в Simulink с помощью пользовательского решателя QP

Чтобы проверить поведение пользовательского решателя при тех же условиях, включите пользовательский решатель в контроллере MPC.

mpcobj.Optimizer.CustomSolver = true;

Необходимо также предоставить функцию MATLAB ®, удовлетворяющую следующим требованиям:

В этом примере используйте пользовательский решатель QP, определенный в файле шаблона mpcCustomSolverCodeGen_TemplateEML.txt, который реализует dantzig и подходит для генерации кода. Сохраните функцию в рабочей папке как mpcCustomSolver.m.

src = which('mpcCustomSolverCodeGen_TemplateEML.txt');
dest = fullfile(pwd,'mpcCustomSolver.m');
copyfile(src,dest,'f')

Смоделировать управление серводвигателем постоянного тока по замкнутому контуру и сохранить входные и выходные данные установки.

sim(mdl)
uDantzigSim = u;
yDantzigSim = y;

-->Converting model to discrete time.
   Assuming no disturbance added to measured output channel #1.
-->The "Model.Noise" property of the "mpc" object is empty. Assuming white noise on each measured output channel.

Создание кода с помощью пользовательского решателя QP

Для выполнения оставшегося примера требуется продукт Simulink Coder.

if ~mpcchecktoolboxinstalled('simulinkcoder')
    disp('Simulink(R) Coder(TM) is required to run this example.')
    return
end

Чтобы создать код из блока контроллера MPC, который использует пользовательский решатель QP, включите параметр пользовательского решателя для генерации кода в контроллере MPC.

mpcobj.Optimizer.CustomSolverCodeGen = true;

Необходимо также предоставить функцию MATLAB ®, которая удовлетворяет всем следующим требованиям:

В этом примере используйте тот же пользовательский решатель, что и в mpcCustomSolverCodeGen_TemplateEML.txt. Сохраните функцию в рабочей папке как mpcCustomSolverCodeGen.m.

src = which('mpcCustomSolverCodeGen_TemplateEML.txt');
dest = fullfile(pwd,'mpcCustomSolverCodeGen.m');
copyfile(src,dest,'f')

Просмотреть сохраненные mpcCustomSolverCodeGen.m файл.

function [x, status] = mpcCustomSolverCodeGen(H, f, A, b, x0)
%#codegen
% mpcCustomSolverCodeGen allows the user to specify a custom (QP) solver
% written in MATLAB to be used by MPC controller in code generation.
%
% Workflow:
%   (1) Copy this template file to your work folder and rename it to
%       "mpcCustomSolverCodeGen.m".  The work folder must be on the path. 
%   (2) Modify the "mpcCustomSolverCodeGen.m" to use your solver.
%       Note that your solver must use only fixed-size data.
%   (3) Set "mpcobj.Optimizer.CustomSolverCodeGen = true" to tell the MPC
%       controller to use the solver in code generation.
% To generate code:
%   In MATLAB, use "codegen" command with "mpcmoveCodeGeneration" (require MATLAB Coder) 
%   In Simulink, generate code with MPC and Adaptive MPC blocks
%
% To use this solver for simulation in MATLAB and Simulink, you need to:
%   (1) Copy "mpcCustomSolver.txt" template file to your work folder and 
%       rename it to "mpcCustomSolver.m".  The work folder must be on the path. 
%   (2) Modify the "mpcCustomSolver.m" to use your solver.
%   (3) Set "mpcobj.Optimizer.CustomSolver = true" to tell the MPC
%       controller to use the solver in simulation.
%
% The MPC QP problem is defined as follows: 
%
%       min J(x) = 0.5*x'*H*x + f'*x, s.t. A*x >= b.
%
% Inputs (provided by MPC controller at run-time):
%       H: a n-by-n Hessian matrix, which is symmetric and positive definite.
%       f: a n-by-1 column vector.
%       A: a m-by-n matrix of inequality constraint coefficients.
%       b: a m-by-1 vector of the right-hand side of inequality constraints.
%      x0: a n-by-1 vector of the initial guess of the optimal solution.
%
% Outputs (sent back to MPC controller at run-time):
%       x: must be a n-by-1 vector of optimal solution. 
%  status: must be an integer of:
%           positive value: number of iterations used in computation
%                        0: maximum number of iterations reached
%                       -1: QP is infeasible
%                       -2: Failed to find a solution due to other reasons
% Note that:
%   (1) When solver fails to find an optimal solution (status<=0), "x" 
%       still needs to be returned.
%   (2) To use sub-optimal QP solution in MPC, return the sub-optimal "x" 
%       with "status = 0".  In addition, you need to set
%       "mpcobj.Optimizer.UseSuboptimalSolution = true" in MPC controller. 
%
% DO NOT CHANGE LINES ABOVE

% This template implements a showcase QP solver using "Dantzig" algorithm
% by G. B. Dantzig, A. Orden, and P. Wolfe, "The generalized simplex method
% for minimizing a linear form under linear inequality constraints",
% Pacific J. of Mathematics, 5:183–195, 1955.
%
% User is expected to modify this template and plug in own custom QP solver
% that replaces the "Dantzig" algorithm.  

ZERO = zeros('like',H);
ONE = ones('like',H);
% xmin is a constant term that adds to the initial basis because "dantzig"
% requires positive optimization variables.  A fixed "xmin" does not work
% for all MPC problems.
xmin = -1e3*ones(size(f(:)))*ONE; 

maxiter = 200*ONE;
nvar = length(f);
ncon = length(b);
a = -H*xmin(:); 
H = H\eye(nvar);
rhsc = A*xmin(:) - b(:);
rhsa = a-f(:);
TAB = -[H H*A';A*H A*H*A'];
basisi = [H*rhsa; rhsc + A*H*rhsa];
ibi = -(1:nvar+ncon)'*ONE;
ili = -ibi*ONE;
%% Call EML function "qpdantzg"
[basis,ib,il,iter] = qpdantzg(TAB,basisi,ibi,ili,maxiter); %#ok<ASGLU>
%% status
if iter > maxiter
    status = ZERO;
elseif iter < ZERO
    status = -ONE;    
else
    status = iter;
end
%% optimal variable
x = zeros(nvar,1,'like',H);
for j = 1:nvar
    if il(j) <= ZERO
        x(j) = xmin(j);
    else
        x(j) = basis(il(j))+xmin(j);
    end
end

Создание исполняемого кода из модели Simulink с помощью slbuild от Simulink Coder.

slbuild(mdl)

### Starting build procedure for: mpc_customQPcodegen
-->Converting model to discrete time.
   Assuming no disturbance added to measured output channel #1.
-->The "Model.Noise" property of the "mpc" object is empty. Assuming white noise on each measured output channel.
### Successful completion of build procedure for: mpc_customQPcodegen

Build Summary

Top model targets built:

Model                Action                       Rebuild Reason                                    
====================================================================================================
mpc_customQPcodegen  Code generated and compiled  Code generation information file does not exist.  

1 of 1 models built (0 models already up to date)
Build duration: 0h 0m 22.665s

В системе Windows после завершения процесса построения программа добавляет исполняемый файл. mpc_customQPcodegen.exe в рабочую папку.

Запустите исполняемый файл. После успешного завершения исполняемого файла (status = 0), программное обеспечение добавляет файл данных mpc_customQPcodegen.mat в рабочую папку. Загрузите файл данных в рабочую область MATLAB и получите входные и выходные сигналы завода, сгенерированные исполняемым файлом.

if ispc
    status = system(mdl);
    load(mdl)
    uDantzigCodeGen = u;
    yDantzigCodeGen = y;
else
    disp('The example only runs the executable on Windows system.');
end

The example only runs the executable on Windows system.

Сравнение результатов моделирования

Сравните входные и выходные сигналы установки из всех моделей.

if ispc
    figure
    subplot(2,1,1)
    plot(u.time,uKWIK.signals.values,u.time,uDantzigSim.signals.values,...
        '+',u.time,uDantzigCodeGen.signals.values,'o')
    subplot(2,1,2)
    plot(y.time,yKWIK.signals.values,y.time,yDantzigSim.signals.values,...
        '+',y.time,yDantzigCodeGen.signals.values,'o')
    legend('KWIK','Dantzig Simu','Dantzig CodeGen','Location','northwest')
else
    figure
    subplot(2,1,1)
    plot(u.time,uKWIK.signals.values,u.time,uDantzigSim.signals.values,'+')
    subplot(2,1,2)
    plot(y.time,yKWIK.signals.values,y.time,yDantzigSim.signals.values,'+')
    legend('KWIK','Dantzig Simu','Location','northwest')
end

Сигналы от всех моделей идентичны.

Ссылки

[1] Бемпорад, А. и Моска, Е. «Выполнение жестких ограничений в неопределенных линейных системах путем управления ссылками». Automatica, том 34, номер 4, стр. 451-461, 1998.

bdclose(mdl)