В этом примере показано, как преобразовать самолет в линейную модель пространства состояний инварианта времени (LTI) для линейного анализа.
Этот пример описывает:
Импорт и заполнение данных из файла DATCOM.
Построение самолета из данных DATCOM.
Вычисление статической устойчивости самолета.
Линеаризация самолета вокруг начального состояния.
Валидация статического анализа устойчивости с динамической характеристикой.
Изоляция передаточной функции «лифт на тангаж» и разработка контроллера обратной связи для лифта.
Этот пример использует файл DATCOM, созданный для самолета Sky Hogg.
Во-первых, импортируйте выходной файл DATCOM с помощью datcomimport.
allData = datcomimport('astSkyHoggDatcom.out', false, 0);
skyHoggData = allData{1}skyHoggData = struct with fields:
case: 'SKYHOGG BODY-WING-HORIZONTAL TAIL-VERTICAL TAIL CONFIG'
mach: [0.1000 0.2000 0.3000 0.3500]
alt: [1000 3000 5000 7000 9000 11000 13000 15000]
alpha: [-16 -12 -8 -4 -2 0 2 4 8 12]
nmach: 4
nalt: 8
nalpha: 10
rnnub: []
hypers: 0
loop: 2
sref: 225.8000
cbar: 5.7500
blref: 41.1500
dim: 'ft'
deriv: 'deg'
stmach: 0.6000
tsmach: 1.4000
save: 0
stype: []
trim: 0
damp: 1
build: 1
part: 0
highsym: 1
highasy: 0
highcon: 0
tjet: 0
hypeff: 0
lb: 0
pwr: 0
grnd: 0
wsspn: 18.7000
hsspn: 5.7000
ndelta: 5
delta: [-20 -10 0 10 20]
deltal: []
deltar: []
ngh: 0
grndht: []
config: [1x1 struct]
version: 1976
cd: [10x4x8 double]
cl: [10x4x8 double]
cm: [10x4x8 double]
cn: [10x4x8 double]
ca: [10x4x8 double]
xcp: [10x4x8 double]
cma: [10x4x8 double]
cyb: [10x4x8 double]
cnb: [10x4x8 double]
clb: [10x4x8 double]
cla: [10x4x8 double]
qqinf: [10x4x8 double]
eps: [10x4x8 double]
depsdalp: [10x4x8 double]
clq: [10x4x8 double]
cmq: [10x4x8 double]
clad: [10x4x8 double]
cmad: [10x4x8 double]
clp: [10x4x8 double]
cyp: [10x4x8 double]
cnp: [10x4x8 double]
cnr: [10x4x8 double]
clr: [10x4x8 double]
dcl_sym: [5x4x8 double]
dcm_sym: [5x4x8 double]
dclmax_sym: [5x4x8 double]
dcdmin_sym: [5x4x8 double]
clad_sym: [5x4x8 double]
cha_sym: [5x4x8 double]
chd_sym: [5x4x8 double]
dcdi_sym: [10x5x4x8 double]
Затем подготовьте интерполяционные таблицы DATCOM.
Интерполяционные таблицы DATCOM могут иметь отсутствующие значения из-за того, что таблицы заполняют только одно значение для всего столбца.
Эти отсутствующие данные представлены в интерполяционных таблицах как значение 99999 и могут быть заполнены с помощью «предыдущего» метода заполнения.
В этом примере, , , , и иметь отсутствующие данные.
skyHoggData.cyb = fillmissing(skyHoggData.cyb, "previous", "MissingLocations", skyHoggData.cyb == 99999); skyHoggData.cnb = fillmissing(skyHoggData.cnb, "previous", "MissingLocations", skyHoggData.cnb == 99999); skyHoggData.clq = fillmissing(skyHoggData.clq, "previous", "MissingLocations", skyHoggData.clq == 99999); skyHoggData.cmq = fillmissing(skyHoggData.cmq, "previous", "MissingLocations", skyHoggData.cmq == 99999);
При заполненных недостающих данных могут быть построены самолеты.
Во-первых, самолет готовят с желаемым именем самолета.
Опционально имя самолета может быть извлечено из поля «case» на struct DATCOM путем передачи пустого объекта самолета.
skyHogg = Aero.FixedWing(); skyHogg.Properties.Name = "Sky_Hogg"; skyHogg.DegreesOfFreedom = "3DOF"; [skyHogg, cruiseState] = datcomToFixedWing(skyHogg, skyHoggData);
datcomToFixedWing преобразует все совместимые данные из struct datcom в объект fixwing и его состояние. Однако возвращенное состояние все еще нуждается в обработке, чтобы получить желаемые начальные условия самолета.
В этом примере окружение, масса, инерция, воздушная скорость и центр давления нуждаются в регулировании.
h = 2000; [T, a, p, rho] = atmoscoesa(convlength(h, "ft", "m")); cruiseState.AltitudeMSL = h; cruiseState.Environment.Temperature = convtemp(T,"K", "R"); cruiseState.Environment.SpeedOfSound = convvel(a, "m/s", "ft/s"); cruiseState.Environment.Pressure = convpres(p, "pa", "psf"); cruiseState.Environment.Density = convdensity(rho, "kg/m^3", "slug/ft^3"); cruiseState.U = 169.42; cruiseState.Mass = 1299.214; cruiseState.Inertia.Variables = [5787.969 0 117.64;0 6928.93 0;-117.64 0 11578.329]; cruiseState.CenterOfPressure = [0.183, 0, 0];
Выполнение статического анализа устойчивости помогает определить динамическую устойчивость системы без вычисления динамической системы отклика.
stability = staticStability(skyHogg, cruiseState)
stability=6×8 table
U V W Alpha Beta P Q R
________ _________ ________ ________ _________ _________ ________ _________
FX "Stable" "" "" "" "" "" "" ""
FY "" "Neutral" "" "" "" "" "" ""
FZ "" "" "Stable" "" "" "" "" ""
L "" "" "" "" "Neutral" "Neutral" "" ""
M "Stable" "" "" "Stable" "" "" "Stable" ""
N "" "" "" "" "Neutral" "" "" "Neutral"
При статически устойчивых силах и моментах, с возмущениями к прямым и вертикальным скоростям и возмущениями к углу атаки, динамическая устойчивость системы имеет тенденцию к колебанию установившегося состояния при возмущении передней и вертикальной скоростей.
Чтобы проверить это поведение, используйте Control System Toolbox™.
Чтобы использовать инструменты в Control System Toolbox™, линеаризируйте самолет вокруг состояния.
Это делается с помощью метода линеаризации с таким же круизным состоянием, как и раньше.
linSys = linearize(skyHogg, cruiseState)
linSys =
A =
XN XD U W Q
XN 0 0 1 0 0
XD 0 0 0 1 0
U 0 1.319e-06 -0.001714 -0.0007608 0
W 0 0 -0.002705 -0.3319 2.557
Q 0 0 0.02379 -2.495 -26.41
Theta 0 0 0 0 1
Theta
XN -2.586e-07
XD -2.957
U -0.5617
W -4.867e-08
Q 0
Theta 0
B =
Delta
XN 0
XD 0
U -0.0004314
W -0.02084
Q -2.321
Theta 0
C =
XN XD U W Q Theta
XN 1 0 0 0 0 0
XD 0 1 0 0 0 0
U 0 0 1 0 0 0
W 0 0 0 1 0 0
Q 0 0 0 0 1 0
Theta 0 0 0 0 0 1
D =
Delta
XN 0
XD 0
U 0
W 0
Q 0
Theta 0
Continuous-time state-space model.
С помощью построенной линейной модели пространства состояний можно построить график динамического поведения системы.
Чтобы проверить результаты статической устойчивости с динамическим поведением системы, постройте график модели пространства состояний с учетом начальных условий.
Чтобы вызвать возмущение в системе, к сигналу лифта добавляется 5 степень шаг в течение 1 секунды.
x0 = getState(cruiseState, linSys.OutputName); t = linspace(0, 50, 500); u = zeros(size(t)); u(t > 1 & t < 2) = 5; lsim(linSys,u,t,x0)
Как ожидалось из статического анализа устойчивости, скорость воздуха и скорость тангажа стабильны при реагировании на небольшие возмущения в лифте.
В дополнение к статической верификации устойчивости, изоляция управляющих поверхностей к их предполагаемой динамической реакции может помочь спроектировать контроллеры, характерные для отдельных поверхностей.
В этом случае имеется только одна управляющая поверхность - лифт.
Лифт управляет тангажом ответом самолета. Чтобы показать ответ тангажа, изолируйте вход лифта от угла тангажа до передаточной функции лифта.
linSysElevatorTF = tf(linSys(6,1))
linSysElevatorTF =
From input "Delta" to output "Theta":
-2.321 s^3 - 0.7222 s^2 - 0.001232 s - 8.935e-09
------------------------------------------------------------------
s^5 + 26.75 s^4 + 15.19 s^3 + 0.03932 s^2 + 0.008227 s + 5.712e-08
Continuous-time transfer function.
step(linSysElevatorTF)
Как видно на графике шага, ответ тангажа на вход лифта имеет нежелательную колебательную природу и большую установившуюся ошибку.
Путем добавления ПИД обратной связи контроллера к входу лифта может быть достигнута гораздо более желательная характеристика тангажа.
C = pidtune(linSysElevatorTF, "PID")C =
1
Ki * ---
s
with Ki = -0.000524
Continuous-time I-only controller.
elevatorFeedback = feedback(linSysElevatorTF * C, 1)
elevatorFeedback =
From input to output "Theta":
0.001217 s^3 + 0.0003786 s^2 + 6.459e-07 s + 4.684e-12
------------------------------------------------------------------------
s^6 + 26.75 s^5 + 15.19 s^4 + 0.04053 s^3 + 0.008605 s^2 + 7.03e-07 s
+ 4.684e-12
Continuous-time transfer function.
step(elevatorFeedback)
