Требования по производительности

Рисунок ниже показов блок-схема системы с обратной связью. Схема включает номинальную модель самолета, контроллер K, а также элементы получая неопределенность модели и цели эффективности (см. следующие разделы для деталей).

Рисунок 1: устойчивая система управления для оси ответвления самолета

Цель проекта состоит в том, чтобы заставить самолет эффективно ответить на боковые входные параметры педали палки и руководящего принципа пилота. Требования по производительности включают:

HQ_p    = 5.0 * tf(2.0,[1 2.0]);
step(HQ_p), title('Desired response from lateral stick to roll rate (Handling Quality)')

Рисунок 2: Желаемый ответ от боковой палки, чтобы прокрутить уровень.

HQ_beta = -2.5 * tf(1.25^2,[1 2.5 1.25^2]);
step(HQ_beta), title('Desired response from rudder pedal to side-slip angle (Handling Quality)')

Рисунок 3: Желаемый ответ от педали руководящего принципа до угла заноса.

freq = 12.5 * (2*pi);  % 12.5 Hz
zeta = 0.5;
yaw_filt = tf(freq^2,[1 2*zeta*freq freq^2]);
lat_filt = tf(freq^2,[1 2*zeta*freq freq^2]);

freq = 4.1 * (2*pi);  % 4.1 Hz
zeta = 0.7;
roll_filt = tf(freq^2,[1 2*zeta*freq freq^2]);

AAFilters = append(roll_filt,yaw_filt,lat_filt);

От спецификаций до функций взвешивания

Алгоритмы проекта H-бесконечности стремятся минимизировать самое большое усиление с обратной связью через частоту (H-норма-по-бесконечности). Чтобы применить эти инструменты, мы должны сначала переделать технические требования проекта как ограничения на усиления с обратной связью. Мы используем функции взвешивания, чтобы "нормировать" технические требования через частоту и одинаково взвесить каждое требование.

Мы можем описать спецификации проекта в терминах функций взвешивания можно следующим образом:

W_act = ss(diag([1/50,1/20,1/60,1/30]));

W_n = append(0.025,tf(0.0125*[1 1],[1 100]),0.025);
clf, bodemag(W_n(2,2)), title('Sensor noise power as a function of frequency')

Рисунок 4: степень шума Датчика в зависимости от частоты

W_p = tf([0.05 2.9 105.93 6.17 0.16],[1 9.19 30.80 18.83 3.95]);
clf, bodemag(W_p), title('Weight on Handling Quality spec')

Рисунок 5: Вес при обработке качественной спецификации.

W_beta = 2*W_p;

Здесь мы масштабировали веса W_act, W_n\wp, и W_beta так усиление с обратной связью между всеми внешними входными параметрами и все взвесили выходные параметры, меньше 1 на всех частотах.

Номинальная модель самолета

Пилот может управлять боковым направленным ответом самолета с боковыми педалями палки и руководящего принципа. Самолет имеет следующие характеристики:

Номинальная боковая направленная модель LateralAxis имеет четыре состояния:

Эти переменные связаны уравнениями пространства состояний:

где x = [v; r; p; phi], u = [delta_stab; delta_rud], и y = [beta; p; r; yac].

load LateralAxisModel
LateralAxis

LateralAxis =
 
  A = 
               v         r         p       phi
   v      -0.116    -227.3     43.02     31.63
   r     0.00265    -0.259   -0.1445         0
   p    -0.02114    0.6703    -1.365         0
   phi         0    0.1853         1         0
 
  B = 
        delta_stab   delta_rud
   v        0.0622      0.1013
   r     -0.005252    -0.01121
   p      -0.04666    0.003644
   phi           0           0
 
  C = 
                 v          r          p        phi
   beta     0.2469          0          0          0
   p             0          0       57.3          0
   r             0       57.3          0          0
   yac   -0.002827  -0.007877    0.05106          0
 
  D = 
         delta_stab   delta_rud
   beta           0           0
   p              0           0
   r              0           0
   yac     0.002886    0.002273
 
Continuous-time state-space model.

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

A_S = [tf([25 0],[1 25]); tf(25,[1 25])];
A_S.OutputName = {'stab_rate','stab_angle'};

A_R = A_S;
A_R.OutputName = {'rud_rate','rud_angle'};

Составление моделирования ошибок

Номинальная модель только аппроксимирует истинное поведение самолета. С учетом несмоделированной динамики можно ввести относительный термин или мультипликативную неопределенность W_in*Delta_G во входе объекта, где ошибочная динамика Delta_G имейте усиление меньше чем 1 через частоты и функцию взвешивания W_in отражает частотные диапазоны, в которых модель более или менее точна. Там обычно больше моделируют ошибки на высоких частотах так W_in высокая передача.

% Normalized error dynamics
Delta_G = ultidyn('Delta_G',[2 2],'Bound',1.0);

% Frequency shaping of error dynamics
w_1 = tf(2.0*[1 4],[1 160]);
w_2 = tf(1.5*[1 20],[1 200]);
W_in = append(w_1,w_2);

bodemag(w_1,'-',w_2,'--')
title('Relative error on nominal model as a function of frequency')
legend('stabilizer','rudder','Location','NorthWest');

Рисунок 6: Относительная погрешность на номинальной модели самолета в зависимости от частоты.

Создавание неопределенной модели динамики самолета

Теперь, когда мы определили количество ошибок моделирования, мы можем создать неопределенную модель динамики самолета, соответствующей пунктирному, окружают рисунок 7 (то же самое как рисунок 1):

Рисунок 7: динамика Самолета.

Используйте connect функционируйте, чтобы объединить номинальную модель LateralAxis корпуса, модели A_S привода и A_R, и описание ошибки моделирования W_in*Delta_G в одну неопределенную модель Plant_unc отображение [delta_stab; delta_rud] на привод и объект выходные параметры:

% Actuator model with modeling uncertainty
Act_unc = append(A_S,A_R) * (eye(2) + W_in*Delta_G);
Act_unc.InputName = {'delta_stab','delta_rud'};

% Nominal aircraft dynamics
Plant_nom = LateralAxis;
Plant_nom.InputName = {'stab_angle','rud_angle'};

% Connect the two subsystems
Inputs = {'delta_stab','delta_rud'};
Outputs = [A_S.y ; A_R.y ; Plant_nom.y];
Plant_unc = connect(Plant_nom,Act_unc,Inputs,Outputs);

Это производит модель Plant_unc неопределенного пространства состояний (USS) из самолета:

Plant_unc

Plant_unc =

  Uncertain continuous-time state-space model with 8 outputs, 2 inputs, 8 states.
  The model uncertainty consists of the following blocks:
    Delta_G: Uncertain 2x2 LTI, peak gain = 1, 1 occurrences

Type "Plant_unc.NominalValue" to see the nominal value, "get(Plant_unc)" to see all properties, and "Plant_unc.Uncertainty" to interact with the uncertain elements.

Анализ, как ошибки моделирования влияют на ответы разомкнутого контура

Мы можем анализировать эффект моделирования неопределенности путем выбора случайных выборок несмоделированной динамики Delta_G и графический вывод номинальных и встревоженных ответов времени (анализ Монте-Карло). Например, для дифференциального канала стабилизатора, вес неопределенности w_1 подразумевает 5%, моделируя ошибку в низкой частоте, увеличиваясь до 100% после 93 рад/секунда, как подтверждено диаграммой Боде ниже.

% Pick 10 random samples
Plant_unc_sampl = usample(Plant_unc,10);

% Look at response from differential stabilizer to beta
figure('Position',[100,100,560,500]) 
subplot(211), step(Plant_unc.Nominal(5,1),'r+',Plant_unc_sampl(5,1),'b-',10)
legend('Nominal','Perturbed')

subplot(212), bodemag(Plant_unc.Nominal(5,1),'r+',Plant_unc_sampl(5,1),'b-',{0.001,1e3})
legend('Nominal','Perturbed')

Рисунок 8: Переходной процесс и диаграмма Боде.

Разработка контроллера Боковой Оси

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

Сначала мы создаем модель OLIC разомкнутого контура отображение внешних входных сигналов к связанным с эффективностью выходным параметрам как показано ниже.

Рисунок 9: модель разомкнутого контура отображение внешних входных сигналов к связанным с эффективностью выходным параметрам.

Чтобы создать эту модель, начните с блок-схемы системы с обратной связью, удалите блок K контроллера и используйте connect вычислить желаемую модель. Как прежде, возможность соединения задана путем маркировки вводов и выводов каждого блока.

Рисунок 10: Блок-схема для того, чтобы создать модель разомкнутого контура.

% Label block I/Os
AAFilters.u = {'p','r','yac'};    AAFilters.y = 'AAFilt';
W_n.u = 'noise';                  W_n.y = 'Wn';
HQ_p.u = 'p_cmd';                 HQ_p.y = 'HQ_p';
HQ_beta.u = 'beta_cmd';           HQ_beta.y = 'HQ_beta';
W_p.u = 'e_p';                    W_p.y = 'z_p';
W_beta.u = 'e_beta';              W_beta.y = 'z_beta';
W_act.u = [A_S.y ; A_R.y];        W_act.y = 'z_act';

% Specify summing junctions
Sum1 = sumblk('%meas = AAFilt + Wn',{'p_meas','r_meas','yac_meas'});
Sum2 = sumblk('e_p = HQ_p - p');
Sum3 = sumblk('e_beta = HQ_beta - beta');

% Connect everything
OLIC = connect(Plant_unc,AAFilters,W_n,HQ_p,HQ_beta,...
   W_p,W_beta,W_act,Sum1,Sum2,Sum3,...
   {'noise','p_cmd','beta_cmd','delta_stab','delta_rud'},...
   {'z_p','z_beta','z_act','p_cmd','beta_cmd','p_meas','r_meas','yac_meas'});

Это производит неопределенную модель в пространстве состояний

OLIC

OLIC =

  Uncertain continuous-time state-space model with 11 outputs, 7 inputs, 26 states.
  The model uncertainty consists of the following blocks:
    Delta_G: Uncertain 2x2 LTI, peak gain = 1, 1 occurrences

Type "OLIC.NominalValue" to see the nominal value, "get(OLIC)" to see all properties, and "OLIC.Uncertainty" to interact with the uncertain elements.

Вспомните, что конструкцией функций взвешивания, контроллер соответствует спецификациям каждый раз, когда усиление с обратной связью меньше 1 на всех частотах и для всех направлений ввода-вывода. Сначала спроектируйте контроллер H-бесконечности, который минимизирует усиление с обратной связью для номинальной модели самолета:

nmeas = 5;		% number of measurements
nctrls = 2;		% number of controls
[kinf,~,gamma_inf] = hinfsyn(OLIC.NominalValue,nmeas,nctrls);
gamma_inf

gamma_inf = 0.9700

Здесь hinfsyn вычисленный контроллер kinf это сохраняет усиление с обратной связью ниже 1, таким образом, спецификациям можно соответствовать для номинальной модели самолета.

Затем выполните Mu-Synthesis, чтобы видеть, можно ли спецификациям соответствовать надежно при принятии во внимание ошибок моделирования (неопределенность Delta_G). Используйте команду musyn выполнять синтез и использовать musynOptions установить сетку частоты, используемую для Mu-Analysis.

fmu = logspace(-2,2,60);
opt = musynOptions('FrequencyGrid',fmu);
[kmu,CLperf] = musyn(OLIC,nmeas,nctrls,opt);

D-K ITERATION SUMMARY:
-----------------------------------------------------------------
                       Robust performance               Fit order
-----------------------------------------------------------------
  Iter         K Step       Peak MU       D Fit             D
    1           5.097        3.487        3.488            12
    2            1.31        1.292        1.312            20
    3           1.242        1.242        1.693            12
    4           1.693        1.544        1.545            16
    5           1.223        1.223        1.551            12
    6           1.533        1.464        1.465            20
    7           1.278        1.277        1.302            12

Best achieved robust performance: 1.22

CLperf

CLperf = 1.2226

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

Сравнение частотного диапазона контроллеров

Сравните эффективность и робастность контроллера H-бесконечности kinf и контроллер mu kmu. Вспомните, что спецификации эффективности достигаются, когда усиление замкнутого цикла меньше 1 для каждой частоты. Используйте lft функция, чтобы замкнуть круг вокруг каждого контроллера:

clinf = lft(OLIC,kinf);
clmu = lft(OLIC,kmu);

Какова эффективность худшего случая (в терминах усиления с обратной связью) каждого контроллера для моделирования ошибок, ограниченных W_in? wcgain команда помогает вам ответить на этот трудный вопрос непосредственно без потребности в обширном gridding и симуляции.

% Compute worst-case gain as a function of frequency
opt = wcOptions('VaryFrequency','on');

% Compute worst-case gain (as a function of frequency) for kinf
[mginf,wcuinf,infoinf] = wcgain(clinf,opt);

% Compute worst-case gain for kmu
[mgmu,wcumu,infomu] = wcgain(clmu,opt);

Можно теперь сравнить номинал и эффективность худшего случая для каждого контроллера:

clf
subplot(211)
f = infoinf.Frequency;
gnom = sigma(clinf.NominalValue,f);
semilogx(f,gnom(1,:),'r',f,infoinf.Bounds(:,2),'b');
title('Performance analysis for kinf')
xlabel('Frequency (rad/sec)')
ylabel('Closed-loop gain');
xlim([1e-2 1e2])
legend('Nominal Plant','Worst-Case','Location','NorthWest');

subplot(212)
f = infomu.Frequency;
gnom = sigma(clmu.NominalValue,f);
semilogx(f,gnom(1,:),'r',f,infomu.Bounds(:,2),'b');
title('Performance analysis for kmu')
xlabel('Frequency (rad/sec)')
ylabel('Closed-loop gain');
xlim([1e-2 1e2])
legend('Nominal Plant','Worst-Case','Location','SouthWest');

Первый график показывает что в то время как контроллер H-бесконечности kinf соответствует спецификациям эффективности для номинальной модели объекта управления, ее эффективность может резко ухудшиться (пиковое усиление около 15) для некоторой встревоженной модели в наших ошибочных границах моделирования.

В отличие от этого mu контроллер kmu имеет немного худшую эффективность для номинального объекта когда по сравнению с kinf, но это последовательно обеспечивает эту эффективность для всех встревоженных моделей (усиление худшего случая около 1.25). mu контроллер поэтому более устойчив к моделированию ошибок.

Валидация временного интервала устойчивого контроллера

Далее протестировать робастность mu контроллера kmu во временном интервале можно сравнить ответы времени номинала и худшего случая модели с обратной связью с идеальной "Обработкой Качества" ответ. Для этого сначала создайте "истинную" модель CLSIM с обратной связью куда все функции взвешивания и образцы модели HQ были удалены:

kmu.u = {'p_cmd','beta_cmd','p_meas','r_meas','yac_meas'};
kmu.y = {'delta_stab','delta_rud'};

AAFilters.y = {'p_meas','r_meas','yac_meas'};

CLSIM = connect(Plant_unc(5:end,:),AAFilters,kmu,{'p_cmd','beta_cmd'},{'p','beta'});

Затем создайте тестовые сигналы u_stick и u_pedal показанный ниже

time = 0:0.02:15;
u_stick = (time>=9 & time<12);
u_pedal = (time>=1 & time<4) - (time>=4 & time<7);

clf
subplot(211), plot(time,u_stick), axis([0 14 -2 2]), title('Lateral stick command')
subplot(212), plot(time,u_pedal), axis([0 14 -2 2]), title('Rudder pedal command')

Можно теперь вычислить и построить идеал, номинал и ответы худшего случая на тестовые команды u_stick и u_pedal.

% Ideal behavior
IdealResp = append(HQ_p,HQ_beta);
IdealResp.y = {'p','beta'};

% Worst-case response
WCResp = usubs(CLSIM,wcumu);

% Compare responses
clf
lsim(IdealResp,'g',CLSIM.NominalValue,'r',WCResp,'b:',[u_stick ; u_pedal],time)
legend('ideal','nominal','perturbed','Location','SouthEast');
title('Closed-loop responses with mu controller KMU')

Ответ с обратной связью почти идентичен для номинала и худшего случая системы с обратной связью. Обратите внимание на то, что ответ уровня крена самолета отслеживает команду уровня крена хорошо первоначально и затем вылетает от этой команды. Это происходит из-за нуля правой полуплоскости в модели самолета в 0,024 рад/секунда.