Санкт-Петербургский государственный университет
Кафедра компьютерных технологий и систем
Мартынов Родион Сергеевич
Выпускная квалификационная работа бакалавра
Задача управления движением двухкорпусного судна
с малой площадью ватерлинии в условиях волнения
010400
Прикладная математика и информатика
Научный руководитель,
кандидат физ.-мат. наук,
доцент
Жабко Н.А.
Санкт-Петербург
2016
Содержание
Введение . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
Постановка задачи . . . . . . . . . . . . . . . . . . . . . . . . .
Обзор литературы . . . . . . . . . . . . . . . . . . . . . . . . .
Глава 1. Математическая модель . . . . . . . . . . . . . . . . .
1.1. Общие положения . . . . . . . . . . . . . . . . . . . .
1.2 Линейные уравнения движения . . . . . . . . . . . . .
1.3 Применение теории потенциалов . . . . . . . . . . . .
1.4 Моделирование неопределенностей в коэффициентах .
1.5 Моделирование волнения . . . . . . . . . . . . . . . .
Глава 2. Построение регулятора . . . . . . . . . . . . . . . . .
2.1 Анализ устойчивости . . . . . . . . . . . . . . . . . . .
2.2 H∞ синтез . . . . . . . . . . . . . . . . . . . . . . . . .
2.3 µ - синтез . . . . . . . . . . . . . . . . . . . . . . . . .
2.4 Понижение порядка системы. . . . . . . . . . . . . . .
2.5 Алгоритм решения задачи. . . . . . . . . . . . . . . .
Глава 3. Эксперимент и результаты . . . . . . . . . . . . . . .
Выводы . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
Заключение . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
Список литературы . . . . . . . . . . . . . . . . . . . . . . . .
Приложение . . . . . . . . . . . . . . . . . . . . . . . . . . . .
2
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
3
4
6
7
7
7
8
10
13
16
16
19
23
24
27
28
32
33
34
36
Введение
В данной работе рассматривается задача управления двухкорпусным
судном с малой площадью ватерлинии в условиях волнения. Подобные суда вызывают большой интерес в современном судостроительстве и морских исследованиях благодаря специфическим свойствам, обусловленным
их конструкцией. Двухкорпусное судно с малой площадью ватерлинии (далее СМПВ) представляет собой разновидность катамарана, однако в отличие от классической конструкции катамарана, корпуса СМПВ полностью
скрыты под водой и ватерлиния проходит по опорам, соединяющим корпуса
с надводной платформой. Это позволяет преодолевать волны без серьезной
бортовой качки, что выгодно выделяет данный тип судов среди остальных.
Однако, при определенных внешних возмущениях, СМПВ может терять устойчивость в вертикальной плоскости. В частности, при курсировании по морю при определенной скорости асимметричное распределение
давления на погруженные корпуса вызывает дестабилизирующий момент,
который может привести к тому, что судно может сильно сместиться по
дифференту. Этот момент имеет название «момент Мунка» и его величина
пропорциональна квадрату скорости судна.
В связи с этим было проведено исследование и доказано, что для
СМПВ момент Мунка может быть компенсирован с помощью горизонтальных рулей и была сформулирована задача построения стабилизирующего
регулятора для СМПВ. Однако, из-за высокой сложности построения математической модели и вычисления необходимых гидродинамических коэффициентов, оптимальное решение данной задачи до сих пор не было
найдено. Кроме того, позднее было показано, что модель вертикального
движения обладает высокой степенью неопределенности, сильно зависящей от скорости судна и характеристик моря.
Поэтому задача стабилизации СМПВ в условиях волнения не теряет
актуальности и по сей день, кроме того, она влечет за собой дальнейшие
работы по разработке управления для данного типа судов с использованием
более сложных моделей.
3
Постановка задачи
Математическая модель вертикального движения СМПВ может быть
представлена в виде линейной стационарной системы обыкновенных дифференциальных уравнений
(
ẋ = Ax + Bu + Hd
y = Cx + Du
(1)
x = [ḣ, h, θ̇, θ]
u = [Uf , Ua , Df , Dm ]
где h и θ обозначают смещение по глубине и по дифференту соответственно, Uf , Ua , Df , Dm - угол отклонения носового руля, угол отклонения кормового руля, силы и моменты вызываемые волнением соответственно, x вектор состояний, y - вектор измерений, y = [h, θ], u = [Uf , Ua ] - вектор
управлений и d = [Df , Dm ] - вектор внешних возмущений.
Матрицы [A, B, C, D] постоянные, однако из-за особенностей судна,
сделанных в процессе построения модели, реальная модель объекта отличается от номинальной:
(
ẋ = (A + δA)x + (B + δB)u
(2)
y = Cx + Du
где δA и δB - матрицы вариаций коэффициентов модели объекта, заданные
в определенных пределах.
Для этой модели необходимо построить регулятор, заданный передаточной матрицей K(s):
u = −K(s)y
(3)
такой, что замкнутая система (2), (3) будет удовлетворять требованиям по
робастной устойчивости и робастному качеству переходного процесса.
Обозначим эту систему M (P, K, ∆), где P - номинальная система
объекта, K - матрица регулятора, ∆ - неопределенность в коэффициентах
модели, заданная на некотором множестве: ∆ ∈ ∆. Тогда требования к
устойчивости можно выразить следующим образом:
Определение 1 Замкнутая система M (P, K, ∆) обладает свойством робастной устойчивости по отношению к неопределенности ∆, если для
всех ∆ ∈ ∆ корни характеристического полинома замкнутой системы
лежат в левой открытой комплексной полуплоскости C − . При этом будем говорить, что регулятор, выраженный передаточной матрицей K:
u = −Ky обеспечивает робастную устойчивость замкнутой системы.
Необходимо также дать определение робастного качества системы.
4
Согласно [12]:
Определение 2 Замкнутая система M (P, K, ∆) обладает определенным
робастным качеством для некоторого наперед заданного функционала J =
J(M (P, K, ∆)) если она робастно устойчива по отношению к ∆ и справедливо включение I ⊂ R ⊂ R1 , где I = J(M (P, K, ∆)) при фиксированном
K, и R - допустимое множество значений рассматриваемого функционала.
Таким образом, регулятор должен отвечать требованиям по робастному качеству переходного процесса в условиях волнения (ненулевого внешнего сигнала d) при нулевых начальных данных (x(0) = 0): амплитуда реакции на внешнее возмущение замкнутой системы должна быть значительно меньше амплитуды реакции на то же возмущение разомкнутой системы:
max yc → min
t∈[0,T ]
(4)
yo − yc → max
где yc и yo - модуль значений выходных переменных замкнутой и разомкнутой системы при одинаковом внешнем сигнале d.
5
Обзор литературы
Первые математические описания судов данного типа встречаются в
работах Lee C.M., в частности, в статье «Theoretical prediction of motion of
small waterplane area, twin-hull ships in waves», в которой автор предложил
полное описание гидродинамических коэффициентов математической модели движения СМПВ в условиях волнения. На его основе была построена
математическая модель движения СМПВ с применением двумерной теории потенциалов в [1], которая затем была развита в труде [2]. В каждой
из этих работ говорилось о необходимости разработки регулятора, обеспечивающего устойчивость при волнении, так как без этого невозможно
построение полноценного автопилота для СМПВ, а, следовательно, и эффективное использование судов данного типа.
Затем, в разное время были представлены решения задачи стабилизации вертикального движения СМВП-судна как без учета вариаций коэффициентов модели (нечеткий регулятор в [3], и PID-регулятор [4]), так
и с учетом вариаций коэффициентов (µ-синтез [5] и скользящий режим
управления [6]). Однако, несмотря на столь обширный список работ, оптимальное решение не было найдено из-за специфики судов такого типа и
сложности разработки управления, о чем отдельно говорится в [7].
Теория по математическому моделированию волнения и в целом морских объектов были подробно описаны в [8] и [9], включая задачи управления морскими объектами.
Теория робастного управления подробно описана в монографии [10],
а её прикладное применение в среде MATLAB описано в труде [11].
6
Глава 1. Математическая модель
Построение математической модели в виде системы обыкновенных
дифференциальных уравнений - необходимый шаг для начала любого анализа устойчивости системы и её различных свойств и качеств. В этой главе будет рассмотрено непосредственно построение модели СМВП-судна и
её линеаризация. Также будут рассмотрены причины возникновения параметрической неопределенности в системе и она будет должным образом
включена в модель для оценки робастных свойств.
1.1. Общие положения
Для исследования движения судна и, в частности, реакции судна
на морское волнение используются положения так называемой seakeeping
theory [8]. В данной теории формулируются методы получения математических моделей (линейных и нелинейных), с помощью которых воздействие
волн на судно может быть формализовано и изучено в терминах устойчивости и управляемости. Однако для её использования необходимо принять
некоторые допущения, а именно:
• движение судна, вызванное волнением, описывается линейными уравнениями;
• крутизна волн принимается невысокой, и сами волны равномерно распределены по морской поверхности;
• океан или море, в котором тестируется судно, полагается бесконечно
глубоким и не создающим никаких сил и моментов, способных вызвать
нелинейные реакции судна;
• скорость судна по курсу постоянна (что включает тривиальный случай
с нулевой скоростью);
• для погруженной части предполагается, что каждый корпус узок настолько, насколько это необходимо для того, чтобы гидродинамическое
давление на погруженные части судна было равномерно распределено
по каждой его части;
• система координат в данном случае связана с судном.
1.2 Линейные уравнения движения
Итак, используя второй закон Ньютона и разложение сил в ряд Тейлора для обеспечения линейности уравнений, получаем следующее пред-
7
ставление движений судна [1]:
6
X
(Aij ξ¨j + bij ξ˙j + cij ξj ) = Fi
j=1
где
• ξ1 - продольное смещение судна;
• ξ2 - поперечное смещение судна;
• ξ3 - вертикальное смещение судна;
• ξ4 - крен;
• ξ5 - дифферент;
• ξ6 - курс;
• Fi - возмущающие силы и моменты волн;
• Aij = Mij + aij - сумма матрицы массы и матрицы присоединенных
масс;
• bij - матрица демпфирования;
• cij - матрица восстанавливающих сил (restoring forces).
Присоединёнными массами называется толща воды, движущаяся благодаря движению погруженной части судна.
Для вычисления соответствующих коэффициентов модели могут использоваться различные подходы. Основным инструментом в таких случаях является теория потенциалов.
1.3 Применение теории потенциалов
Теория потенциалов разделяется на две основных категории: двумерная теория слоев и теория трехмерного потенциала. В данной работе разумно использовать двумерную теорию, так как для решения задачи необходимо рассматривать только вертикальное движение судна. Теория слоев
позволяет представить трехмерное погруженное тело как серию двумерных
слоев. Каждый слой имеет толщину δx, при этом полагающуюся малой, и
имеет собственные гидродинамические характеристики, которые соотносятся с коэффициентами всего судна в его уравнениях движения. Также
все испытываемые судном возмущения от морского волнения соотносятся
с каждым слоем по отдельности.
Теория слоев подразумевает, что локальные гидродинамические характеристики имеют тот же характер, что имели бы характеристики погруженного цилиндра бесконечной длины такого же поперечного сечения, что
8
и у каждого корпуса. Другими словами, все трехмерные эффекты, такие
как взаимодействие между слоями и рассеивание волн на концах судна, не
учитываются.
При применении теории слоев для СМПВ требуется учесть также
следующий ряд гидродинамических эффектов [1]:
1. Гидродинамические взаимодействия между двумя корпусами
2. Вязкое демпфирование
3. Эффекты стабилизирующих (горизонтальных) рулей
Таким образом, можно уточнить коэффициенты модели. В силу линейности системы разумно положить их выражение через суммы соответствующих гидродинамических эффектов:
A0i j = Ai j P + Ai f f ;
b0ij = bPij + bvij + bfij ;
c0i j = cvij + cij f f + cbij ;
Fi0 = FiP + Fiv + Fif ,
где верхние индексы b, P , v, и f обозначают коэффициенты, полученный от
эффектов подъемной силы воды (buoyancy), потенциала потока, вязкости
жидкости и эффекта рулей соответственно. Последнее стоит рассмотреть
подробнее, так как рули играют первостепенную роль в задаче стабилизации судна.
Силы и моменты, создаваемые движением рулей могут быть представлены в модели, аналогичной модели движения судна [4]:
X f
f
Fi =
[Aij + Mijf ]η̈j + Bijf η̇ + Cijf η − [Fife + Fik ],
j
где верхние индексы соответствуют отдельным рулям, e - возмущающая
сила волнения, k - управляющие силы, возникающие при активации рулей. В случае рассмотрения линейной модели, силы, стоящие в левой части
уравнения сводятся к подъемной силе (lift) Li :
1
L = ρVf2 Af CL (αe ),
2
где Af - площадь руля, CL (αe ) - наклон кривой подъемной силы, ρ и Vf плотность и скорость жидкости, действующей на руль.
Принимая во внимание тот факт, что влияние носовых рулей на кормовые рули незначительно, в рамках поставленной для СМПВ-судов задачи
9
(стабилизация вертикального движения) можно записать [3]:
1
L(t)h = ρU 2 A(f ) CLα α(t)
2
1
M (t)θ = ρU 2 A(f ) CLα α(t)Lθ
2
Здесь CLα есть эффективное соотношение сторон руля, которое может быть
вычислено следующим образом:
CLα = (KW (B) + KB(W ) )(CLα )W .
Индекс W обозначает руль, как самостоятельное тело, B(W ) - подъемную силу, возникающую от действия руля на тело, W (B) - подъемная
сила, действующая от тела на руль. Коэффициенты KW (B) и KB(W ) обозначают соответствующие взаимодействия и вычисляются как
1
a
a
KW (B) + KB(W ) = [1 + ]2 + 2
2
s
s
a - радиус корпуса, s - длина руля, включая a.
CLα =
1.8πAe
p
.
1.8 + A2e + 4
Ae - соотношение сторон руля.
Несмотря на то, что данное выражение кажется чересчур упрощенным, оно вписывается в систему условий, обозначенных при постановке
задачи и формулировании математической модели.
1.4 Моделирование неопределенностей в
коэффициентах
Приведенная выше модель не предназначена для использования современных методов анализа управляемости, устойчивости и непосредственно синтеза регулятора. В связи с этим преобразование в модели в пространство состояний является необходимым шагом.
Для начала следует выделить из модели уравнения, задающие движение в вертикальной плоскости. Это возможно сделать благодаря выбранному режиму движения, линейности модели и продольной симметрии судна.
Результатом будут следующие уравнения:
(M + A33 )ξ¨3 + B33 ξ˙3 + C33 ξ3 + A35 ξ¨5 + B35 ξ˙5 + C35 ξ5 = F3
A53 ξ¨3 + B53 ξ˙3 + C53 ξ3 + (I5 + A55 )ξ¨5 + B55 ξ˙5 + C55 ξ5 = F5
(5)
Теперь становится возможным переписать эти уравнения в матричной фор10
ме, добавив гидродинамические коэффициенты при стабилизирующих рулях:
ξ3
(M + A33 ) 0
A35
0
0
1
1
0
ξ˙3
; x=
T =
A53
0 (I5 + A55 ) 0
ξ5
0
0
0
1
ξ˙5
−B33 −C33 −B35 −C35
b11 b12 1 0
1
0
0
0
; B 0 = 0 0 0 0
A0 =
−B53 −C53 −B55 −C55
b31 b32 0 1
0
0
1
0
0 0 0 0
Uf
Ua
0
1
0
0
u=
Dξ3 ; c = 0 0 0 1
Dξ5
где Uf , UA , Dξ3 , Dξ5 - углы отклонения переднего и заднего стабилизирующего руля а также силы и моменты, задействованные в их активации и
влияющие на вертикальное движение судна.
В итоге данные уравнения запишутся в виде
T ẋ = A0 x + B 0 u
которые в свою очередь очевидными алгебраическими преобразованиями
приводятся к классической форме модели в пространстве состояний:
ẋ = Ax + Bu
y = Cx
Коэффициенты матрицы управления B определяются следующим образом:
f
f
b11 = −0.5ρAf CLc
cos(Cang
);
b21 = −b11 Lf
f
f
b12 = −0.5ρAa CLc
cos(Cang
)
b22 = −b12 La
Основная особенность данной модели состоит в том, что её коэффициенты заданы в частотной области и зависят от многих параметров,
большинство из которых невозможно учесть явно [7]. Поэтому при моделировании пользуются одним из основных частотных показателей судна в
условиях волнения: частотой встречи ωe (encounter frequency) [1]:
ω 2 U cos(β)
ωe = ω −
g
где U - собственная скорость судна по курсу, ω - частота волнения мо11
ря/океана, β - угол между курсом судна и направлением волнения (0◦ для
курса по волнам и 180◦ для курса строго против волн). Следовательно, для
использования данной модели должна быть выбрана определенная частота
ωe .
Наиболее частым подходом является выбор частоты, дающий наихудший случай движения (максимальную амплитуду движения). Однако
данный подход не гарантирует того, что построенный для такой модели регулятор будет отвечать поставленным требованиям, то есть обеспечивать
устойчивость замкнутой системы и необходимое качество переходного процесса. Поэтому разумным выходом в данном случае является формулирование модели с неопределенностями в коэффициентах, которая сможет
отразить возможные вариации движения судна.
Для построения подобной модели нужно предпринять следующие шаги: в первую очередь частота волнения дискретизируется в некотором диапазоне. Затем вычисляется частота встречи, основываясь на основных скоростных режимах рассматриваемого судна и интересующего угла β. Для
каждого значения ωe вычисляются коэффициенты модели в пространстве
состояний и составляется диапазон вариаций каждого коэффициента. Без
учета сил и моментов, создаваемых волнением, это выразится следующим
образом:
n
a11 + ar11 δ1 an12 + ar12 δ2 an13 + ar13 δ3 an14 + ar14 δ4
1
0
0
0
Au =
an31 + ar31 δ5 an32 + ar32 δ6 an33 + ar33 δ7 an34 + ar34 δ8
0
0
1
0
n
b11 + br11 δ9 bn12 + br12 δ10
0
0
Bu =
r
n
r
n
b31 + b11 δ11 b32 + b32 δ12
0
0
0 1 0 0
0 0
Cu =
Du =
0 0 0 1
0 0
Здесь коэффициенты с индексом n обозначают номинальные значения и
равны среднему значению среди всех вычисленных вариантов соответствующих коэффициентов, а коэффициенты с индексом r равны половине диапазона возможных значений. δi представляют собой изменения соответствующего коэффициента, с ограничением: |δi | ≤ 1(i = 1, 2, ..., 12).
Теперь, для применения современных методов теории робастности,
необходимо представить систему в виде номинальной модели, замкнутой
обратной связью с неопределенностями в коэффициентах (M − ∆ представление [11]). Это выполняется с помощью верхнего дробно-линейного
12
преобразования (рис. 1):
Gu (s) = Fu (M, ∆)
Рис. 1: M − ∆ конфигурация системы.
где
M11
M=
M21
δ1
δ2
M12
∆=
...
M22
δ12
Матрица M представляет собой некоторую комбинацию номинальной модели и вариаций её коэффициентов (зачастую, их включают в выходы системы). Тогда
Fu (M, ∆) = M22 + M21 ∆(I − M11 ∆)−1 M12
Все методы, описанные в следующих главах, применяются к замкнутой
системе Gu , полученной выше.
1.5 Моделирование волнения
Волнение создает давление на корпуса судна, которое в свою очередь создает силы и моменты, влияющие на движение. Возникающее в таком случае внешнее воздействие оказывает огромное влияние на движение
судна. Это внешнее воздействие выражается как передаточная функция
между амплитудой волнения и силами и моментами [9]. Соответствующая
амплитуда волнения Ai выражается следующим образом:
A2i = 2S(ωi )∆ω
13
где S(ωi ) – функция плотности спектра волны, ωi – частота i-го компонента, ∆ω - интервал между частотами. Добавив к этому величину ki = 2π
λi –
номер волны, где λi – длина i-ой волны и φi – случайный фазовый угол из
интервала [0; 2π], можно ввести функцию, моделирующую высоту волны
ξ(x, t) =
N
X
i=1
Ai cos(ωi t − ki x + φi ) +
N
X
1
i=1
2
ki A2i cos2(ωi t − ki x + φi ) + O(A3i )
Так как функция ξ(x, t) повторяет свои значения с периодом 2π/∆ω,
рекомендуется выбирать ωi случайно из интервала ∆ω для адекватного
моделирования случайных волн.
Для вычисления этой функции необходимо выбрать функцию S(ωi ).
Существует множество её вариантов, содержащее различное количество
параметров, однако все они должны быть ограничены при возрастании
частоты (ω 1, α > 0):
S(ω) −→ αg 2 ω −5
В данной работе рассматриваются два основных используемых спектра:
улучшенный спектр Пирсона-Московица и JONSWAP-спектр (рис. 2).
Для вычисления реакции морских объектов на волнение в открытом
море чаще всего используется улучшенный спектр Пирсона-Московица:
4π 3 Hs2
−16π 3
S(w) = 4 5 exp
Tz ω
Tz4 ω 4
Это двухпараметрический спектр, в котором Hs – значимая высота волны (обычно принимается 31 от высоты волны от поверхности моря),
Tz = 0.921T1 , где T1 – средний период волнения. Его следует использовать
только для морей, в которых влияние ветра на создание волнения близко к постоянному (состояние моря не зависит от расстояния, на котором
действует ветер и от длительности его действия). В случае, если это не
выполняется, рекомендуется использовать JONSWAP-спектр:
Hs2
−944
× (3.3)Y
S(ω) = 155 4 4 exp
4
4
T1 ω
T1 ω
"
2 #
0.191ωT1 − 1
√
Y = exp −
2σ
(
0.07, ω ≤ 5.24/T1
σ=
0.09, ω > 5.24/T1
В целом, алгоритм моделирования волнения в рамках данной задачи
описывается следующим образом:
14
1. Выбор функции плотности спектра S(ω)
2. Выделение N интервалов значений выбранной на шаге 1 функции с
длиной ∆ω
3. Выбрать случайную частоту ωi для каждого интервала
p
4. Вычислить амплитуду Ai = 2S(ωi )∆ω и номер волны ki = ωi2 /g
5. Вычислить функцию ξ(x, t).
Рис. 2: Сравнение спектров волнения для волнения в 3 балла.
15
Глава 2. Построение регулятора
В данной главе будет описан процесс анализа построенной математической модели в терминах теории управления, робастности и устойчивости,
а также приведены алгоритмы для синтеза регулятора, решающего поставленную задачу.
2.1 Анализ устойчивости
Устойчивость является основным предметом рассмотрения при синтезе управления.
Определение 3 Частное решение x = ξ(t) автономной системы ẋ =
f (t, x) называется устойчивым по Ляпунову, если для любого заданного
числа ε > 0 можно указать такое число δ > 0, зависящее от чисел и
t0 , что любое другое решение x = η(t) этой же системы, для которого
выполняется неравенство kη(t0 )−ξ(t0 )k < δ, удовлетворяет неравенству
kη(t) − ξ(t)k < ε для любого момента времени t > t0
Однако, для построения качественной модели необходимо, чтобы система обладала более сильным свойством - свойством асимптотической
устойчивости.
Определение 4 Частное решение x = ξ(t) этой же системы называется асимптотически устойчивым по Ляпунову, если оно является устойчивым в смысле определения 3 и если существует такое число H > 0,
что любое решение x = η(t) этой системы, для которого выполняется
неравенство kη(t0 ) − ξ(t0 )k < H, удовлетворяет соотношению
lim kη(t) − ξ(t)k = 0
t→∞
При этом для линейных стационарных систем вида ẋ = Ax, x ∈ E n
все решения являются одновременно либо асимптотически устойчивыми,
либо асимптотически неустойчивыми, поэтому в дальнейшем под устойчивостью системы будем пониматься асимптотическая устойчивость всех её
решений.
Неустойчивые системы не несут никакого практического смысла, что
делает анализ устойчивости как разомкнутой, так и замкнутой систем обязательным этапом при разработке управления.
Модель вертикального движения СМПВ, представленного системой
(5), можно свести к одному уравнению четвертого порядка, имеющему вид:
aλ4 + bλ3 + cλ2 + dλ + e = 0.
16
Где
a = (M + A33 )(I5 + A5 5) − A35 A53 ;
b = (M + A33 )B55 + (I5 + A5 5)B33 − A35 B53 − A53 B35 ;
c = (M + A33 )C55 + (I5 + A5 5)C33 + B33 B55 − A35 C53 − A53 C35 − B35 B53 ;
d = B33 C55 + B55 C33 − B53 C35 − B35 C53 ;
e = C33 C55 − C35 C53 .
Корни данного уравнения для СМПВ представляют собой две пары
комплексно-сопряженных значений, соответствующие вертикальной и килевой качке [15]. Если для определенности обозначить индексом h корни,
соответствующими вертикальной качке, а индексом θ корни, соответствующие килевой качке, то будут выполняться соотношения ( [15], [7]):
|Reλh | < |Reλθ |; |Imλh | > |Imλθ |
(6)
После приведения системы к представлению в пространстве состояний:
(
ẋ = Ax + Bu
(7)
y = Cx + Du
x ∈ E n, u ∈ E m, y ∈ E p
соотношение (6) между корнями (которые в данном случае будут являться
собственными числами матрицы A) сохранится.
Для подобных систем устойчивость гарантируется в случае, если (для
модели, представленной в пространстве состояний) все собственные числа
матрицы A располагаются в левой открытой комплексной полуплоскости
(имеют отрицательные вещественные части).
Однако, в случае, когда в модели присутствуют неопределенности
(как параметрические, т.е. колебания коэффициентов, так и неструктурированные, например, не моделируемая динамика), критерии устойчивости
необходимо расширить, чтобы учесть все возможные возмущения объекта.
Основой для такого подхода является теорема о малом возмущении (smallgain theorem) [11]:
Теорема 1 Рассмотрим систему, представленную на рисунке 3, где G1 и
G2 – передаточные матрицы соответствующих линейных стационарных
систем. Если G1 (s) и G2 (s) устойчивы, то замкнутая система устойчива тогда и только тогда, когда
||G1 G2 ||∞ < 1
и
||G2 G1 ||∞ < 1
17
Рис. 3: Замкнутая система.
||G||∞ = max σ̄(ω),
ω∈[0,∞)
где σ̄(ω) – максимальное сингулярное число матрицы G(iω) (корень квадратный из максимального собственного значения матрицы G0 (iω)G(iω)).
Данная теорема дает необходимое условие устойчивости замкнутой
системы. Замкнутая система называется робастно устойчивой, если как
система, представляющая объект, так и регулятор, выраженный передаточной матрицей K, замыкающий эту систему, являются устойчивыми для
любых заранее определенных возмущений и/или колебаний коэффициентов ∆. Это включает в себя также и то, что данный регулятор должен
стабилизировать и номинальную модель G. Тогда, объединяя эти условия,
получим следующую теорему [11]:
Теорема 2 Замкнутая система является робастно устойчивой для
устойчивого возмущения ∆, если K(s) стабилизирует Gn – номинальную
систему объекта G и выполняется:
||K(I + GK)−1 ||∞ <
1
||∆||∞
Теорема 2 дает практический критерий оценки качества построенного регулятора в терминах робастной устойчивости. Однако, для учета
неопределенности в конфигурации M − ∆ разумно ввести дополнительную
характеристику, учитывающую полную неопределенность ∆ (как параметрическую, так и неструктурированную).
Неопределенность, представленная в виде ∆-блока в M − ∆ конфигурации, имеющего специальную структуру, учитывается и анализируется
при помощи µ-теории [12]. Основным её понятием является структурированное сингулярное число. Пусть
∆ = { diag[ δ1 Ir1 , ..., δs Irs , ∆1 , ..., ∆f ] : δi ∈ C, ∆j ∈ C mj ×mj }
18
P
P
где si=1 ri + fj=1 mj = n, n – размерность блока ∆, ∆i , i = 1, f - блок,
содержащий диагональную матрицу с постоянными комплексными компонентами, представляющий часть параметрической неопределенности системы. Также отметим, что множество ∆ ограничено. Тогда структурированное сингулярное число µ∆ (M ) для системы M − ∆ определяется, как
µ∆ (M (s)) = min {σ̄(∆) : det(I − M ∆) = 0}
∆∈∆
Если не существует такого ∆ ∈ ∆, что det(I − M ∆) = 0, то µ∆ (M ) =
0. Обратная величина для структурированного сингулярного числа отражает границу устойчивости в частотной области. Тогда критерий робастной устойчивости с учетом структурированной неопределенности ∆ будет
определяться следующим образом:
Теорема 3 Пусть номинальная замкнутая система M (s) устойчива, а
возмущение ∆ ограничено: ||∆||∞ ≤ β ∀∆ ∈ ∆. Тогда возмущенная система M − ∆ робастно устойчива тогда и только тогда, когда
µ∆ (M (s)) <
1
β
(8)
Теорема 3 дает необходимое, и при этом практически эффективное
условие робастной устойчивости. Это позволяет использовать её при синтезе регуляторов, например, при H∞ синтезе. Однако, основной проблемой
при её применении является сложность вычисления структурированного
сингулярного числа, поэтому на практике точное значение часто заменяют
верхней границей его оценки.
Данная оценка получается из следующей теоремы [11]:
Теорема 4 Пусть
D = {D = diag[D1 , ..., Ds , d1 Im1 , ..., df Imf ] :
Di ∈ C ri ×ri , Di = Di∗ > 0, dj > 0}
т.е. совпадает по структуре с блоком неопределенностей ∆ системы.
Тогда имеет место следующее неравенство:
µ(M ) ≤ inf σ̄(DM D−1 )
D∈D
2.2 H∞ синтез
Одним из наиболее распространенных методов построения регулятора для робастных систем является H∞ синтез.
19
Рассмотрим модель в пространстве состояний следующего вида:
ẋ(t) = Ax(t) + B1 d(t) + B2 u(t)
e(t) = C1 x(t) + D11 d(t) + D12 u(t)
y(t) = C2 x(t) + D21 d(t)
где x(t) ∈ E n - вектор состояний, d(t) ∈ E l - вектор внешних возмущений, u(t) ∈ E m - вектор управлений, e(t) ∈ E p1 - вектор контролируемых
переменных и y(t) ∈ E p2 - вектор измерений. Данную систему можно представить в виде передаточной матрицы P (s):
e
d
P11 P12
=
P21 P22
y
u
B1 B2
A
A
B
C1 D11 D12 =
C D
C2 D21 0
Замкнутая система, с которой работает метод, представляет собой
результат нижнего дробно-линейного преобразования (рис. 4):
Ted (s) = Fl (P, K) = P11 + P12 K(I − P22 K)−1 P21
где K - передаточная матрица регулятора.
Рис. 4: Блок-схема замкнутой системы управления.
Общей идеей метода является минимизация нормы H∞ передаточной
функции от входа d к выходу e. При этом учитывается чувствительность
системы к величине входного сигнала и устойчивость относительно внешних возмущений. Однако, в этом подходе сложно учитывать параметрическую неопределенность системы, а также для сложных систем задание
требований к качеству переходного процесса является нетривиальной задачей.
Таким образом, постановка задачи для оптимального случая: найти
20
передаточную матрицу регулятора K, который по измеряемым выходам y
системы P будет создавать управляющее воздействие u, минимизируя норму замкнутой системы. Сам регулятор ищется в виде системы пространства
состояний с матрицами:
¯˙x = Ak x̄ + Bk y
u = Ck x̄ + Dk y
В дальнейшем, говоря о матрице регулятора K будет иметься в виду передаточная матрица K(s) = Ck (Is − Ak )−1 Bk + Dk .
Подобную задачу обычно рассматривают в паре с так называемой
«субоптимальной» задачей H∞ синтеза: построить матрицу регулятора K,
обеспечивающий выполнение неравенства
||Ted (s)||∞ = max σ̄(Ted (jω)) ≤ γ
ω
Тогда задача оптимального синтеза состоит в том, чтобы найти γmin =
minK γ.
Однако, в подобной постановке и классической конфигурации замкнутой системы сложно определить какие-либо требования к качеству
переходного процесса, такие как уменьшение перерегулирования, подавление внешних возмущений и т.д.. Для того, чтобы учесть их при синтезе, требуется ввести некоторые весовые функции в систему вместе с соответствующими им «виртуальными» выходами. На практике чаще всего
пользуются двумя специальными весовыми функциями We и Wu , которые
накладываются на контролируемый выход системы и выход регулятора соответственно (рис. 5).
Рис. 5: Конфигурация замкнутой системы с весовыми функциями.
В рамках задачи H∞ синтеза рассматривается функция чувствитель21
ности системы
S(s) =
1
1 + K(s)P (s)
(9)
и добавочной чувствительности
K(s)S(s) =
K(s)
1 + K(s)P (s)
и соответствующие веса задаются относительно этих функций [13] следующим образом:
1
|S(jω)| ≤
|We (jω)|
1
|K(jω)S(jω)| ≤
|Wu (jω)|
We (s) и Wu (s) обычно являются передаточными функциями первого
порядка. Обязательным условием является то, что эти функции должны
быть устойчивы, т.е. их полюса должны лежать в левой открытой полуплоскости. Для их поиска можно сформировать формулы с варьируемыми
параметрами, которые будет необходимо подбирать в ходе синтеза управления.
Функция We задает некоторый желаемый вид функции чувствительности системы (9) и позволяет задать требования по качеству переходного
процесса, подавлению шумов и внешних возмущений. Её можно искать в
следующем виде:
s + ωb
1
=
(10)
We (s) s/Ms + wb
где ' 0, Ms < 2, а основной варьируемый параметр wb . Его увеличение
ведет к лучшему переходному процессу, однако вызывает рост нагрузки на
управляющие элементы (рули).
В свою очередь, функция Wu позволяет задать требования для управляющих элементов системы - скорость работы приводов объекта, углы отклонения рулей и т.д., которые всегда ограничены при рассмотрении модели реального объекта. Одним из шаблонов для поиска Wu является формула
1
1 s + ωbc
=
(11)
Wu (s) s + wbc /Mu
где 1 ' 0, Mu выбирается исходя из ограничений, накладываемых на динамику управляющих элементов, ωbc - основной варьируемый элемент, уменьшение которого ведет к лучшему соблюдению наложенных ограничений,
однако ухудшает качество переходного процесса.
После построения весовых функций и включения их в систему задача
22
H∞ синтеза выглядит следующим образом:
!
WS
e
||Ted ||∞ =
Wu KS
≤γ
∞
"
Ted (s) =
#
"
#
We
−We P
+
K(I + P K)−1 I
0
Wu
Поиск значений коэффициентов указанных весовых функций - основная задача при использовании H∞ синтеза.
Для построенной модели СМПВ-судна соответствующие веса будут
являться передаточными матрицами [We ]4×4 и [Wu ]2×2 , при условии, что
контролируемые переменные e = (ξ˙3 , ξ3 , ξ˙5 , ξ5 ), а управляющий сигнал u =
(Uf , Ua ). Сама задача решается с помощью уравнения Риккати либо с помощью линейных матричных неравенств, что подробно описано в [10].
Регулятор, полученный при помощи H∞ синтеза, затем должен быть
проверен на соответствие выбранным критериям робастной устойчивости
замкнутой системы и качества переходного процесса. Однако, не всегда
удается достигнуть необходимой границы устойчивости, так как при H∞
синтезе сложно учесть параметрическую неопределенность, но его можно
использовать для получения регулятора, который будет выполнять требуемые критерии по качеству переходного процесса и затем использовать при
синтезе регулятора, который будет учитывать данный вид неопределенности в том числе. Одним из методов, позволяющих учесть это, является
µ-синтез.
2.3 µ - синтез
µ-синтез позволяет учесть параметрическую неопределенность в системе непосредственно, получив в итоге регулятор, который соответствует
требованиям как робастной устойчивости, так и робастного качества переходного процесса. Для этого необходимо воспользоваться понятием структурированного сингулярного числа системы. На его основе и построен основной алгоритм µ-синтеза - алгоритм D − K итераций.
Алгоритм опирается на теорему 4, согласно которой структурированное сингулярное число системы µ(M ) может быть приближено его верхней
границей, зависящей от некоторой матрицы D, повторяющей структуру
неопределенностей системы ∆. Рассмотрим замкнутую систему, состоящую
из трех блоков: M − ∆ конфигурация системы (номинальная модель P и
блок, содержащий неопределенности ∆) и матрица регулятора K. Обозначим её как M (P, K). Тогда задачу µ-синтеза можно сформулировать следующим образом: требуется найти стабилизирующий контроллер K, при
23
котором выполняется:
sup µ[M (P, K)(jω)] <
ω∈R
1
β
при условии k∆k∞ ≤ β. Для её решения был предложен метод D − K итераций [17], основанный на решении следующей оптимизационной задачи:
sup inf σ̄[DM (P, K)D−1 (jω)] <
ω∈R D∈D
1
β
.
Каждая итерация состоит из следующих шагов:
1. Выбор фиксированной матрицы D и решение задачи H∞ -синтеза:
K = arg inf kM (P, K)k∞
K
2. Для полученной на шаге 1 матрицы регулятора K решается оптимизационная задача вида:
D(jω) = arg inf σ̄[DM (P, K)D−1 (jω)]
D∈D
для заранее выбранного диапазона частот ω
3. Аппроксимация полученных на предыдущем шаге матриц D устойчивой передаточной матрицей D(s). Затем перейти на шаг 1.
Для начала работы алгоритма нужно либо выбрать начальную матрицу D - в таком случае её зачастую принимают единичной, - либо начальную матрицу регулятора K - тогда первая итерация начнется со второго шага. Таким образом, если сначала использовать H∞ - синтез для
получения «начального» регулятора, удовлетворяющего критериям качества переходного процесса, то с помощью D − K-итераций можно добиться
робастной устойчивости замкнутой системы.
2.4 Понижение порядка системы.
Регуляторы, полученные методом D − K итераций, зачастую имеют
очень высокий порядок, что может привести к высокой стоимости расчетов и усложнению их использования. Поэтому, после µ-синтеза разумно
рассмотреть задачу понижения порядка регулятора.
24
Рассмотрим систему, заданную в пространстве состояний:
ẋ(t) = Ax(t) + B1 d(t) + B2 u(t)
e(t) = C1 x(t) + D11 d(t) + D12 u(t)
y(t) = C2 x(t) + D21 d(t) + D22 u(t)
где x(t) ∈ E n - вектор состояний, d(t) ∈ E l - вектор внешних возмущений,
u(t) ∈ E m - вектор управлений, e(t) ∈ E p1 - вектор контролируемых переменных и y(t) ∈ E p2 - вектор измерений. Передаточная матрица от входа
u к выходу e будет выглядеть следующим образом:
M (s) = C(sI − A)−1 B + D
где A ∈ E n × n, B = [B1 B2 ] : n × (m + l), C = [C1 C2 ]T : (p1 + p2 ) × n,
D11 D12
D=
: (p1 + p2 ) × (m + l).
D21 D22
Тогда система пониженного порядка представится как
Mr (s) = Cr (sI − Ar )−1 Br + Dr
где размерности матриц равны: Ar : r×r, B : r×(m+l), C : (p1 +p2 )×r, D :
(p1 +p2 )×(m+l) и r < n. Суть задачи состоит в том, чтобы получить такую
систему Mr (s), для которой некоторая норма ошибки была бы малой, т.е.
минимизировать некоторую норму (например, норму H∞ ):
kM (s) − Mr (s)k
Для решения этой задачи необходимо, чтобы система M (s) была
устойчива (либо было возможно разложение на устойчивую и неустойчивую составляющие) и представлена в сбалансированной реализации:
Определение 5 Устойчивая система M (s) называется сбалансированной, если решения P и Q матричных уравнений Ляпунова:
AP + P AT + BB T = 0
AT Q + QA + C T C = 0
такие, что P = Q = diag(σ1 , σ2 , ..., σn ) = Σ, где σ1 ≥ σ2 ≥ ...σn > 0
P и Q называются грамианом управляемости и грамианом наблюдаемости соответственно, σi , i = 1, ..., n - i-ое сингулярное число Ганкеля.
Для системы в несбалансированной реализации существует следующий алгоритм приведения её в нужную форму [11]:
1. Вычислить грамианы P и Q.
2. Для Q вычислить разложение Холецкого: Q = RT R.
25
3. Сформировать положительно определенную матрицу RP RT и привести её к диагональному виду:
RP RT = U Σ2 U T
где U - ортонормированная матрица, U T U = I, Σ = diag(σ1 , σ2 , ..., σn ),
σ1 ≥ σ2 ≥ ...σn > 0.
1
4. Найти матрицу балансирующего преобразования T = Σ− 2 U T R. Тогда
[T AT −1 , T B, CT −1 ] - сбалансированная реализация системы M (s).
Одним из методов для решения задачи уменьшения порядка системы является метод минимизации нормы Ганкеля, которая определяется как
максимальное сингулярное число Ганкеля системы. Суть метода заключается в следующем: пусть M (s) - устойчивая система в сбалансированной
реализации с матрицами [A, B, C, D]. Грамианы управляемости и наблюдаемости этой системы равны: P = Q = diag(Σ1 , σIl ), где σ - минимальное
сингулярное число Ганкеля кратности l и каждый диагональный элемент
матрицы Σ1 больше σ.
Пусть теперь матрицы [A, B, C] могут быть представлены в виде:
A11 A12
B1
A=
B=
C = [C1 C2 ]
A21 A22
B2
Тогда система Mr (S) порядка (n − l) строится как:
Ar = Γ−1 (σ 2 AT11 + ΣT1 A11 Σ1 − σC1T U B1T )
Br = Γ−1 (Σ1 B1 + σC1T U )
Cr = C1 Σ1 + σU B1T
Dr = D − σU
где U - ортонормированная матрица, удовлетворяющая уравнению
B2 = −C2T U
и
Γ = Σ21 − σ 2 I
Полученная система устойчива и удовлетворяет следующим неравенствам
[14]:
kM (s) − Mr (s)kH = σ
(12)
kM (s) − Mr (s)k∞ = σ
К полученной системе можно применять описанный выше метод итеративно, на каждой итерации получая систему меньшего порядка, но с
большим значением ошибки.
26
Итак, в данной главе были рассмотрены основные понятия и методы, используемые при решении поставленной задачи. На их основе можно
сформировать следующий алгоритм
2.5 Алгоритм решения задачи.
Основываясь на выводах из предыдущих пунктов, можно составить
следующий алгоритм поиска стабилизирующего регулятора для СМПВсудна:
1. Для данной модели СМПВ-судна (7) сформировать M − ∆ - конфигурацию.
2. Провести анализ робастной устойчивости и робастного качества переходного процесса полученной модели.
3. Полученные на шаге 2 результаты отразить в конфигурации для H∞ синтеза, т.е. решить задачу поиска весовых функций We и Wu , в соответствии с (10) и (11), таким образом, чтобы удовлетворять требованиям по качеству переходного процесса для номинальной модели,
используя варьирование параметров так, как это указано в параграфе
2.2.
4. Для полученных весов провести процедуру H∞ -синтеза и оценить робастные качества полученного регулятора.
5. При несоответствии полученных робастных качеств требуемым, применить метод D − K итераций для улучшения обеспечиваемой робастной
устойчивости и робастного качества регулятора.
6. При необходимости, воспользоваться методами редукции полученного
регулятора, например, методом аппроксимации нормы Ганкеля.
Следует отметить, что в целом данный алгоритм применим для поиска стабилизирующих регуляторов, обладающих робастными свойствами,
для самых различных объектов, не ограничиваясь только лишь СМПВ.
В следующей главе будет описан пример применения данного алгоритма для модели СМПВ-судна (CHINA 2001).
27
Глава 3. Эксперимент и результаты
Для подтверждения сделанных в предыдущих главах теоретических
выводов в рамках данной работы была рассмотрена модель СМПВ судна
водоизмещением 400 тонн (CHINA 2001) [3] и произведено моделирование
в среде MATLAB с использованием пакетов Control System Toolbox, Robust
Control Toolbox и Marine Systems Simulator. Они содержат необходимые
функции для построения систем в пространстве состояний, синтеза обратных связей а также анализа робастных качеств системы. Marine Systems
Simulator использовался для моделирования волнения, указанного в пункте 1.4.
Номинальные коэффициенты модели следующие:
−0.3129 −8.9758 −0.5743 42.196
1
0
0
0
An =
0.0021
0.8279 −0.1077 −4.8088
0
0
1
0
1.3918 1.3918 0.5945 0.5945
0
0
0
0
Bn =
−0.0568 −0.0568 0.5945 0.5945
0
0
0
0
0 1 0 0
0 0 0 0
Cu =
Du =
0 0 0 1
0 0 0 0
Вариации соответствующих коэффициентов указаны в таблице 1.
ar11 0.0626
ar12 1.7952
ar13 0.1149
ar14 8.4392
ar31 0.0004
ar32 0.1656
ar33 -0.0215
ar34 -0.9618
br11 0.2784
br12 0.2784
br31 0.0114
br32 0.0114
Таблица 1: Вариации коэффициентов
Следуя алгоритму, описанному в предыдущей главе, необходимо построить M − ∆ конфигурацию системы. Блок M состоит из номинальной
28
модели и вариаций её коэффициентов:
An Bn
Cn Dn L
S
0
где S = [S1 . . . S12 ], Si - вектор-строка размерности 1 × 12. Если обозначить
ненулевые элементы S в виде Si,j , то все элементы матрицы S равны нулю,
кроме S1,1 = ar11 , S2,2 = ar12 , S3,3 = ar13 , S4,4 = ar14 , S5,1 = ar31 , S6,2 = ar32 ,
S7,3 = ar33 , S8,4 = ar34 , S9,5 = br11 , S10,6 = br12 , S11,5 = br31 , S11,6 = br32 . Блок L =
[L1 . . . L12 ], где Li = [1 0 0 0 0 0]T для i = 1, 2, 3, 4, 9, 10 и Li = [0 0 1 0 0 0]T
для i = 5, 6, 7, 8, 11, 12. ∆ = diag(δ1 . . . δ1 2), |δi | ≤ 1, k∆k ≤ 1.
Теперь необходимо провести анализ устойчивости модели. Номинальная модель устойчива (все её собственные значения лежат в левой открытой полуплоскости), однако существуют возмущения, способные сделать
модель неустойчивой (рис. 6). Отсюда следует, что необходим синтез регу-
Рис. 6: Собственные числа номинальной и возмущенной системы.
лятора, обеспечивающего робастную устойчивость и качество переходного
процесса замкнутой системы. Разумно применить H∞ - синтез.
Для этого сначала необходимо задать весовые функции We и Wu , которые будут определять качество переходного процесса и ограничения на
управляющие элементы. В данном случае, угол отклонения рулей не должен превышать 20◦ . We представляет собой диагональную передаточную
29
матрицу размерности 2 × 2:
wh 0
We =
0 wθ
где соответствующие элементы wh и wθ отвечают контролируемым выходам eh и eθ (смещение по глубине и смещение по дифференту соответственно). Параметры подбирались таким образом, чтобы максимизировать
подавление внешних возмущений, следуя сформулированному условию (4).
Соответствующие функции были получены с применением формулы (10).
График амплитудно-частотной характеристики данных весовых функций
представлен в приложении А (рис. 7 и рис. 8).
Следует отметить, что зачастую от замкнутой системы также требуется соблюдение качественных характеристик переходного процесса в режиме при заданных начальных данных x(0) = x0 и в отсутствие внешних
возмущений (d = 0):
Tp = inf{tm : ρ(t) ∈ D(ρ0 , δ), ∀t ≥ tm } ≤ t∗
Jp = (ρm − ρ0 )/ρ0 ∗ 100% ≤ ρ∗
где Tp - длительность переходного процесса, ρ(t) = |ỹ(t)|, ρ0 = |ye | - действительное и желаемое значения выхода y системы (2) соответственно,
D(ρ0 , δ) - заданная δ-окрестность числа ρ0 . Jp - величина перерегулирования в процентах, где ρm = supt∈[0,t) - точная верхняя граница функции,
причем ρm ≥ ρ0 > 0. Однако, учитывая высокую степень неопределенности
системы (все значимые коэффициенты могут варьироваться), решение задачи о качестве переходного процесса в двух режимах при помощи одного
регулятора не находится тривиальным путем. Одним из возможных методов является построение регуляторов специальной структуры (скоростные
регуляторы [16]), либо построение нескольких регуляторов и переключение
между ними в зависимости от режима работы. В случае же одного регулятора и при учете неопределенностей в коэффициентах модели, обычно
данную задачу опускают.
В свою очередь, Wu - диагональная передаточная матрица 2 × 2,
c постоянными элементами, которые отражают масштабирование сигнала u = −K(s)y в соответствии с наложенными ограничениями. В данном
случае применение формул (11) оказалось необязательным, так как ограничения, накладываемые на управляющие элементы, возможно выразить
соответствующими константами.
Полученный регулятор стабилизирует номинальную модель и отвечает требованиям качества переходного процесса: амплитуда в режиме с
ненулевым волнением выходного сигнала замкнутой системы в несколько
раз ниже, чем у разомкнутой системы при соблюдении ограничений на ди30
намику рулей (Приложение B, рис. 11), однако система не отвечает требованиям робастной устойчивости. Так как k∆k ≤ 1, то должно выполняться
(8):
µ∆ (M (s)) < 1
однако верхняя граница структурированного сингулярного числа больше
1 на некоторых частотах (Приложение C, рис. 12), что означает, что существуют такие возмущения ∆ ∈ ∆, которые могут сделать систему неустойчивой. Следовательно, нужно расширить границы робастной устойчивости
с помощью µ-синтеза методом D − K итераций, используя регулятор, полученный при H∞ -синтезе в качестве начального для алгоритма.
В итоге был выбран регулятор, полученный на 10 итерации с минимальной достигнутой границей структурированного сингулярного числа
равной 0.8951 (Приложение D, рис. 13). Полученный регулятор имеет двадцатый порядок, поэтому следующим шагом стоит рассмотреть вопрос о
возможности снижения порядка замкнутой системы.
Для этого воспользуемся методом аппроксимации нормы Ганкеля.
Сингулярные числа Ганкеля замкнутой системы представлены на рисунке 14. Применяя метод, описанный в параграфе 2.4, с уровнем ошибки (12)
σ = 0.05, получаем регулятор, выраженный матрицей Kr 9 порядка при
ошибке kK − Kr k∞ = 0.0385, который при этом сохранил робастные свойства регулятора полного порядка (рис. 15).
31
Выводы
В рамках данной работы были достигнуты следующие результаты:
1. Проведен анализ математической модели вертикального движения СМПВ
в пространстве состояний с учетом параметрической неопределенности,
зависящей от гидродинамических коэффициентов;
2. Предложен алгоритм синтеза регулятора на основе методов H∞ -синтеза
и D − K итераций, которые обеспечивают выполнение необходимых
требований по робастной устойчивости и качеству переходного процесса в тестовом режиме;
3. Произведено моделирование движения судна в различных условиях в
среде MATLAB, сделаны соответствующие выводы о допустимости использования предложенного алгоритма для решения поставленной задачи.
32
Заключение
В данной работе была рассмотрена задача управления СМПВ в условиях волнения, в частности, проблема потери устойчивости при определенном режиме движения и состоянии моря. Для её решения была построена
математическая модель вертикального движения СМПВ в системы линейных стационарных дифференциальных уравнений и определены границы
вариаций её коэффициентов, возникающие в связи с зависимостью гидродинамических коэффициентов от частоты моря и скорости судна. Затем
данная модель была приведена к M − ∆ конфигурации и был проведен
робастный анализ устойчивости и качества переходного процесса.
Для построения регулятора, решающего поставленную задачу, использовались метод H∞ синтеза, с помощью которого обеспечивались требования к качеству переходного процесса и ограничениям на управляющие
элементы, а затем метод µ - синтеза, с помощью которого обеспечивались
робастные свойства замкнутой системы. В результате поставленная задача
была решена - был получен регулятор, обеспечивающий робастную устойчивость и приемлемое качество переходного процесса в режиме волнения.
В дальнейшем приведенный алгоритм решения задачи и непосредственно регулятор, построенный для модели судна «DQ0 SWATH», могут
быть применены при работе с более сложными математическими моделями, задающими движение других СМПВ в 6 степенях свободы. Также регулятор может быть использован как часть стабилизирующего регулятора
при построении автопилота и системы управления курсом. В целом, данная
работа представляет собой первый шаг в разработке полноценной системы
управления СМПВ, решающей различные задачи мореходства судов данного типа.
33
Список литературы
[1] Chinn N. L., Roberts G. N., Scrace R. G., Owens D.H.. Mathematical
Modelling of a small waterplane area twin hulled (SWATH) vessel. //
Control’94 International Conference on (Volume: 2) 1994. С. 1560 – 1565.
[2] Dubrovsky V. A., Lyakhoviyskiy A.. Multi-Hull ships. Backbone Publishing
Company, 2001. 495 c.
[3] Lihua Liang, Baohua Wang, Ming Ji. Adaptive Fuzzy Control for
SWATH Ship Seakeeping Characteristics. // International Conference on
Mechatronics and Automation 2012. С. 440 – 445.
[4] Omar Bin Yaakob. Development of a semi-SWATH craft for Malaysian
waters. University of technology, Malaysia, 2006. 198 с.
[5] Qiang Liu. Research on Robust Control of Longitudinal Motion of
SWATH.// Journal of computers 2011. Т. 6. В. 12. С. 2706 – 2710.
[6] Jing-yun Zhao, Yue-feng Liao, Xing Pu, Xiao Han. On the
Active Disturbance Rejection of SWATH Ship based on Terminal
Sliding Mode Control. 2016. International Journal of Control and
Automation. Т. 9. В. 1. С. 399 – 408.
[7] Dubrovsky V. A., Matveev K., Sutulo S.. Small waterplane area
ships. Backbone Publishing Company , 2007. 255 c.
[8] Thor I. Fossen. Handbook of marine craft hydrodynamics and
motion control. Norwegian University of Science and Technology,
Norway, 2011. 575 с.
[9] Thor I. Fossen. Guidance and control of Ocean Vehicles. University of
Trondheim, Norway, 1994. 480 с.
[10] Поляк Б. Т., Щербаков П.С.. Робастная устойчивость и управление , 2002. 303 c.
[11] Gu D .W., Petkov P. Hr., Konstantinov M. M.. Robust Control Design
with MATLAB. Springer, 2005. 389 с.
[12] Веремей Е. И.. Линейные системы с обратной связью. 2013. 445 с.
[13] Olivier S.. Robust control of MIMO systems. 2007 138 c. 2015.
[14] Glover K.. All optimal Hankel-norm approximations of linear
multivariable systems and their l∞ -error bounds.// International Journal
of Control, с. 39:1115– 1193, 1984.
34
[15] Бородай И. К., Мореншильдт В.А., Виленский Г.В., Дубицкий В.М.,
Смирнов Б.Н.. Прикладные задачи динамики судов на волнении., Судостроение, с. 137-152, 1989.
[16] Веремей Е. И., Корчанов В. М., Коровкин М. В., Погожев С. В. Компьютерное моделирование систем управления движением морских подвижных объектов. СПб.: НИИ Химии СПбГУ, 2002. 370 с.
[17] Doyle J.C.. Structured uncertainty in control system design.//In
Proceedings of the 24 IEEE Conference on Decision and
Control. 1985. 260–265 c.
35
Приложение
Приложение A
На графиках представлены соответствующие амплитудно-частотные
характеристики весовых функций wh и wθ , использованные в практической
части данной работы.
Рис. 7: Амплитудно-частотная характеристика wh .
Рис. 8: Амплитудно-частотная характеристика wθ .
36
Приложение B
Реакция СМПВ в двух режимах работы: при ненулевом внешнем возмущении и нулевых начальных условиях и при нулевом внешнем возмущении и ненулевых начальных условиях для системы, замкнутой регулятором, получившимся в результате H∞ синтеза.
Рис. 9: Выход системы и отклонение рулей при волнении в 3 балла.
37
Рис. 10: Выход системы и отклонение рулей при волнении в 4 балла.
38
Рис. 11: Выход системы и отклонение рулей при волнении в 5 баллов.
39
Приложение C
График верхней границы структурированного сингулярного числа
системы, замкнутой регулятором, выраженным матрицей K, получившимся в результате H∞ синтеза.
Рис. 12: Граница структурированного сингулярного числа замкнутой системы.
40
Приложение D
Сравнение характеристик робастной устойчивости регуляторов, полученных H∞ и µ-синтезом, как сравнение верхней границы структурированного сингулярного числа замкнутой системы.
Рис. 13: Границы структурированного сингулярного числа.
41
Приложение E
Рис. 14: Сингулярные числа Ганкеля замкнутой системы.
Рис. 15: Граница сингулярного числа для полного и редуцированного регуляторов.
42
Отзывы:
Авторизуйтесь, чтобы оставить отзыв