Реализация точных методов анализа устойчивости нелинейных динамических систем

Исходные уравнения, особые точки, анализ устойчивости при малом начальном отклонении

В качестве объекта исследования необходимо рассмотреть некоторую абстрактную САР, математическая модель динамики которой описывается следующей системой нелинейных дифференциальных уравнений в форме Коши:

Особые точки находятся из системы при равных нулю левых частях уравнений динамики (условия стационара). Данная динамическая система имеет три особых точки:

1-я точка → (0, 0) → тривиальное решение;

2-я точка → (1, 0.5);

1-я точка → (1, – 0.5).

Для анализа типа особых точек (устойчивое или неустойчивое равновесие) обычно используют линеаризацию уравнений динамики в особой точке и рассматривают поведение линеаризованной системы в малой окрестности особой точки:

где Δx₁=[x₁(t)-x⁰₁]; Δx₂=[x₂(t)-x⁰₂] – малые отклонения от особой точки;

x⁰₁ и x⁰₂ – координаты особой точки, а коэффициенты aⁱ_j вычисляются по соотношениям

Оценку типа особых точек выполним на основании корней характеристического уравнения, которое для системы записывается в матричном виде:

det(A – λ·E) = 0, где

Преобразуя определитель, получаем характеристическое уравнение в виде:

λ²-(a₁₁+a₂₂)⋅λ+(a₁₁⋅a₂₂-a₁₂⋅a₂₁)=0.

Необходимо вычислить корни уравнения для каждой особой точки.

Для первой особой точки коэффициенты равны: а₁₁ = 0; а₁₂ = 0; а₂₁ = 0; а₂₂ = 1. Тогда характеристическое уравнение принимает вид λ² – λ = 0. По структуре это уравнение соответствует неустойчивому инерционно-интегрирующему звену. Корни уравнения равны «0» и «1» (первый корень – в начале координат; второй корень – в правой полуплоскости).

Положительный корень указывает, что первая точка может быть неустойчивой особой точкой при малом начальном отклонении. Данное утверждение необходимо проверить на этапе анализа результатов прямого расчёта фазовых траекторий.

Для второй особой точки коэффициенты равны: а₁₁ = -2; а₁₂ = 4; а₂₁ = -0.5; а₂₂ = 0. Характеристическое уравнения принимает вид λ² +2 · λ +2 = 0. Корни этого уравнения равны (-1 ± i), т.е. корни комплексно-сопряженные и лежат в левой полуплоскости.

Резюме: вторая точка является устойчивой при малом начальном отклонении (устойчивый фокус).

Для третьей особой точки коэффициенты равны: а₁₁ = -2; а₁₂ = -4; а₂₁ = 0.5; а₂₂ = 0. Характеристическое уравнения принимает такой же вид, что и для второй точки (λ² +2 · λ +2 = 0), поэтому корни уравнения равны (-1 ± i). Третья точка является устойчивой при малом начальном отклонении (устойчивый фокус).

Анализ движения автономной системы на фазовой плоскости

Используя освоенные в предыдущей лабораторной работе методы структурного моделирования для решения обыкновенных дифференциальных уравнений, выполнить в SimInTech решение системы уравнений.

Необходимо сформировать структурную схему для решения системы уравнений согласно рисунку (Рисунок 1).

Для построения фазовых портретов стоит воспользоваться блоком «Язык программирования», который позволяет реализовать численное решение системы дифференциальных уравнений динамики САР, записанной в форме Коши.

Открыть редактор блока «Язык программирования» и заполнить его согласно рисунку (Рисунок 2).

В свойствах блока с подписью «X1» задать свойство «Начальные условия» равным «-1», а в свойствах блока с подписью «X2» задать свойство «Начальные условия» равным «1». Таким образом была установлена начальная точка на фазовой плоскости с координатами (-1, 1) в начальный момент времени.

В параметрах расчёта проекта установить параметр «Конечное время расчёта» равным «10».

Запустить проект на расчёт и задать вид графиков согласно рисунку (Рисунок 3). Фазовая траектория асимптотически огибает вторую особую точку с координатами (1, 0.5), которая является устойчивым фокусом.

В свойствах блока с подписью «X2» задать свойство «Начальные условия» равным «1». Таким образом была установлена начальная точка на фазовой плоскости с координатами (-1, -1) в начальный момент времени.

Запустить проект на расчёт и задать вид графиков согласно рисунку (Рисунок 4). Фазовая траектория асимптотически огибает третью особую точку с координатами (1, -0.5), которая является устойчивым фокусом.

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

Для одновременного построения большого числа фазовых траекторий необходимо выбрать рассматриваемую часть фазовой плоскости. Пусть x1 ⊂ [-2,4], x2 ⊂ [-3,3].

В рассматриваемой области заданы двадцать начальных точек, координаты которых приведены в таблице (Таблица 1).

Заполнить редактор блока «Язык программирования» согласно рисунку (Рисунок 5). Такая запись означает, что на вход блока будет подан вектор из двадцати составляющих, сигнал на выходе также будет векторизованным, с той же размерностью.

Открыть окно «Скрипт страницы» и заполнить его согласно рисунку (Рисунок 6), сформировав вектор «х1» и «х2».

В блоке с подписью «Х1» задать свойства «Коэффициент усиления» равным «20#1», свойство «Начальные условия» равным «x1». В блоке с подписью «Х2» задать свойства «Коэффициент усиления» равным «20#1», свойство «Начальные условия» равным «x2».

Запустить проект на расчёт. Группа фазовых траекторий образовало фазовый портрет, вид которого должен быть близким к рисунку (Рисунок 7) и свидетельствовать, что фазовые траектории, начинающиеся строго в верхней полуплоскости стремятся ко второй особой точке с координатами (1, 0.5), а начинающиеся строго из нижней полуплоскости стремятся к третьей особой точке с координатами (1, – 0.5).

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

Для наблюдения за процессом построения фазового портрета необходимо в параметрах расчёта проекта на вкладке «Синхронизация» включить указатель «Синхронизация с реальным временем».

Необходимо увеличить рассматриваемую область фазовой плоскости. Для этого необходимо увеличить каждый элемент массивов «x1» и «x2» в десять раз. Запустить проект на расчёт и задать вид графиков согласно рисунку (Рисунок 8).

Практически вся фазовая плоскость (кроме отрицательной оси абсцисс) является областью устойчивых фазовых траекторий. Траектории стремятся обогнуть вторую или третью особые точки.

Резюме: только при отклонении системы в первой особой строго в отрицательную сторону по оси абсцисс нелинейная САР никогда не вернется в какое-то равновесное состояние.

Описание критерия абсолютной устойчивости В.М. Попова

Одним из точных критериев анализа устойчивости нелинейных САР, не утративших свою актуальность и в настоящее время, является критерий абсолютной устойчивости В.М. Попова.

В критерии абсолютной устойчивости В.М. Попова нелинейная САР условно разделена на линейную часть, обычно расположенную в прямой цепи, и нелинейную часть, обычно расположенную в цепи обратной связи (Рисунок 9).

В классическом варианте доказательства данного критерия принят ряд допущений:

1. Нелинейная часть – безинерционна.
2. Статическая характеристика нелинейной части является однозначной и вписывается в Гурвицев угол К (где «0 < К < ∞»).
3. Линейная часть должна быть устойчивой, или в особых случаях иметь не более двух полюсов, расположенных на мнимой оси, при всех остальных полюсах передаточной функции, расположенных в левой полуплоскости.
4. В особых случаях должна иметь место предельная устойчивость.
5. Видоизмененная АФЧХ обозначается W*(iω) и определяется соотношением: W*(iω)=u^*(ω)+i⋅ν^*(ω), где u^*(ω)≡u(ω); ν^*(ω)≡T⋅ω⋅ν(ω); T=1; u(w), n(w) – действительная и мнимая части АФЧХ линейной части, соответственно.

Существуют аналитическая и геометрическая формулировки абсолютной устойчивости по В.М. Попову. Более наглядной является геометрическая формулировка. Для того, чтобы имела место абсолютная устойчивость, в угле [0; К], в основном и в угле [eps; К] (где eps – бесконечно малое положительное число) в особых случаях, достаточно, чтобы в плоскости W^*(iω) можно было выбрать прямую, проходящую через точку действительной оси с абсциссой «–1/K» так, чтобы годограф W^*(iω) весь лежал строго справа от этой прямой и чтобы в особых случаях имела место предельная устойчивость.

На рисунке (Рисунок 10) представлена иллюстрация критерия Попова при анализе устойчивости нелинейной САР, где пунктирной линией представлен традиционный годограф Найквиста (годограф АФЧХ) для линейной части САР (W_лин), сплошной линией представлен видоизмененный годограф Попова, а точка на оси абсцисс с координатой «-1/K» (K – Гурвицев угол) расположена левее точки пересечения годографа Попова с осью абсцисс.

Очевидно, что через точку «-1/К» можно провести множество прямых.

На рисунке (Рисунок 10) одна из множества прямых проведена так, что видоизмененный годограф Попова лежит строго справа от этой прямой.

Преобразование линейной САР в нелинейную

Для анализа нелинейной САР с использованием критерия абсолютной устойчивости В.М. Попова необходимо воспользоваться структурной схемой из лабораторной работы №1, которая представлена на рисунке (Рисунок 11).

Объект управления представлен в виде блока «Колебательное звено» с подписью «W2(s)». В блоке «Колебательное звено» свойство «Коэффициент усиления» равно «1», свойство «Постоянная времени» равно «1», свойство «Коэффициент демпфирования» равно «0.5».

Локальная обратная связь представлена в виде блока «Апериодическое звено 1-го порядка» с подписью «W3(s)». В блоке «Апериодическое звено 1-го порядка» свойство «Коэффициент усиления» равно «0.6», свойство «Постоянная времени» равно «5».

Интегрирующий регулятор представленного блоком «Апериодическое звено 1-го порядка» с подписью «W₁(s)». В блоке «Апериодическое звено 1-го порядка» свойства «Коэффициент усиления» равно «0.35».

Добавить на схему блок «Константа» с вкладки «Источники», блок «Релейное с зоной нечувствительности» с вкладки «Нелинейные», блок «Построение частотных характеристик» с вкладки «Анализ и оптимизация». Придать вид структурной схемы и задать подписи блокам согласно рисунку (Рисунок 12).

В свойствах блока с подписью «Управляющее реле» задать свойство «Нижняя граница переключения» равным «-0.02», свойство «Нижняя граница зоны нечувствительности» равным «-0.01», свойство «Верхняя граница зоны нечувствительности» равным «0.01», свойство «Верхняя граница переключения» равным «0.02», свойство «Нижнее значение функции» равным «-1», свойство «Верхнее значение функции» равным «1».

В свойствах блока с подписью «Годографы Попова и Найквиста» задать свойство «Начальная круговая частота» равным «0.01», свойство «Конечная круговая частота» равным «1000», свойство «Количество точек вывода» равным «250», в свойство «Типы характеристик» установить «Годограф Попова; Годограф Найквиста».

Формулировка заданий к анализу устойчивости нелинейной САР

1. Используя критерий абсолютной устойчивости В.М. Попова, определить скоростную эффективность интегрирующего регулятора (блок с подписью «W₁(s)»), при которой нелинейная САР будет абсолютно устойчивой.
2. Определить тип устойчивости нелинейной САР, используя прямое моделирование переходного процесса в автономной системе при ненулевых начальных условиях.
3. Выполнить расчет переходного процесса в САР (нулевые начальные условия) при подаче управляющего воздействия, равного 0.8·1(t).

Необходимо проверить, удовлетворяет ли нелинейная часть нелинейной САР допущениям, рассмотренным выше.

Для нелинейной части системы:

• Типовое нелинейное звено, внесенное в структурную схему представляет из себя нелинейную части САР и является безынерционным.
• Статическая характеристика нелинейного звена с введенными свойствами является однозначной и ее статическая характеристика вписывается в Гурвицев угол [0; 50] (при 1/0.02 = 50 = К).

Для линейной части системы:

• Линейная часть САР соответствует варианту особого случая, так как она имеет один нулевой полюс, при всех остальных полюсах, расположенных в левой полуплоскости.
• Для линейной части с введенным значением скоростной эффективности привода (при «k»равном «0.35») существует предельная устойчивость (при замыкании линейной части отрицательной жесткой обратной связью с бесконечно малым коэффициентом усиления САР будет устойчивой).

При замыкании скорректированной линейной части даже единичной обратной связью САР является устойчивой, поэтому при меньшем коэффициенте усиления в блоке с подписью «W3» скорректированная линейная САР в замкнутом состоянии тем более будет устойчивой.

Анализ устойчивости с использованием критерия В.М. Попова

Необходимо привести созданную нелинейную САР к автономному виду.

В свойствах блока с подписью «Управляющее воздействие» задать свойство «Конечное состояние» равным «0».

Запустить проект на расчёт. Двойным нажатием по блоку «Построение частотных характеристик» открыть графики расчёта частотных характеристик. Привести графики к виду, согласно рисунку (Рисунок 13).

Одновременное отображение графиков годографа Найквиста и видоизмененного годографа Попова показывает, что при одном и том же значении частоты вещественные части у обоих годографов одинаковые, поэтому точки пересечения линий этих годографов с осью абсцисс совпадают. Это позволяет во многих случаях использовать для анализа устойчивости нелинейных САР не видоизмененный годограф, а годограф Найквиста.

При частоте, стремящейся, к нулю годограф Попова стремится в точку с координатами (b1 – а1,b1 – b0), где коэффициенты a1,b0и b1 – коэффициенты передаточной функции линейной части САР, определяемой выражением

W_лин(s)=W_экв(s)=W₁(s)⋅W₂(s)/(1+W₂(s)⋅W₃(s))=(b₁⋅s+b₀)/s⋅(a₃⋅s₃+a₂⋅s₂+a₁⋅s+1)

Необходимо уточнить, можно ли через точку с координатами (-1/K, 0), где «К» равно «50» – верхняя граница Гурвицева угла, провести прямую так, чтобы годограф Попова лежал строго справа от этой прямой.

Исходя из графиков годографов можно установить, что точка с абсциссой «-0.02» находится внутри годографа Попова и любая прямая, проведенная через эту точку, пересечет линию годографа Попова.

Этот результат свидетельствует о том, что рассматриваемая замкнутая автономная нелинейная САР (структурная схема которой получена вставкой дополнительного нелинейного звена в структурную схему устойчивой линейной САР) не будет абсолютно устойчивой.

Удалить линию связи между блоком «Константа» и блоком с подписью «Главное сравнивающее устройство». Провести линию связи между свободным портом блока с подписью «Главное сравнивающее устройство» и выходом блока с подписью «W2». В свойствах блока с подписью «W2» задать свойство «Начальные условия» равным «0.1».

Запустить проект на расчёт. График переходного процесса примет вид близкий к рисунку (Рисунок 14).

Вид переходного процесса показывает, что в автономной замкнутой нелинейной САР при ненулевых начальных условиях устанавливается режим приблизительно гармонических автоколебаний, амплитуда которых примерно в двадцать пять раз превышает зону нечувствительности в блоке с подписью «Управляющее реле».

Резюме: в данной нелинейной САР устанавливается режим относительно больших автоколебаний, поэтому САР не удерживает стационарное состояние с погрешностью и, следовательно, такая система должна считаться практически неустойчивой.

Коррекция САР и определение типа устойчивости

Из графика годографа Попова (Рисунок 13) следует: чтобы замкнутая нелинейная САР стала устойчивой, необходимо либо уменьшить приблизительно в двадцать раз коэффициент скоростной эффективности в интегрирующем регуляторе, либо в двадцать раз уменьшить значение управляющего воздействия в блоке с подписью «Управляющее реле».

В свойствах блока «Интегратор» в поле «Формула» задать свойство «Коэффициенты усиления» равным «0.35/20».

Запустить проект на расчёт и убедиться, что вид переходного процесса асимптотически устойчив.

Начальное отклонение относительно быстро стабилизируется и САР асимптотически возвращается в свое равновесное состояние (y_стационарное = 0).

Резюме: скорректированная автономная нелинейная САР асимптотически устойчива.

Анализ устойчивости скорректированной нелинейной САР

Придать вид структурной схемы согласно рисунку (Рисунок 12). В свойствах блока с подписью «W₂(s)» задать свойство «Начальные условия» равным «0».

Запустить проект на расчёт. Двойным нажатием по блоку «Построение частотных характеристик» открыть графики расчёта частотных характеристик. Привести графики к виду, согласно рисунку (Рисунок 16). Точка пересечения годографа Попова с осью абсцисс расположена правее точки с абсциссой -1/К = -0.02.

Нажать правой кнопкой мыши по графику и выбрать пункт «Анализ данных» –> «Табличный просмотр»: таким образом данные расчёта отобразаятся в виде таблицы. Прокруткой таблицы найти строку, в которой график с именем «годограф Попова» пересекает абцису с координатами «-0.02» и убедиться, что ордината меньше значения «0».

Из таблицы следует, что при смене знака мнимой части вещественная часть равна приблизительно «–0.0198», поэтому через точку пересечения можно провести прямую, относительно которой годограф Попова будет расположен строго справа.

Резюме: скорректированная нелинейная замкнутая САР абсолютно устойчива.

Расчет переходных процессов при подаче управляющего воздействия

В свойствах блока с подписью «Управляющее воздействие» задать свойство «Конечное состояние» равным «0.8».

Удалить линию связи между блоком «Константа» и блоком с подписью «Главное сравнивающее устройство». Провести линию связи между свободным портом блока с подписью «Главное сравнивающее устройство» и выходом блока с подписью «W2». В параметрах проекта задать параметр «Конечное время расчёта» равным «200».

Запустить проект на расчёт.

График переходного процесса (Рисунок 18) свидетельствует, что скорректированная нелинейная САР отработала управляющее воздействие, однако время переходного процесса значительно больше (около 80 с).

В свойствах блока «Интегратор» в поле «Формула» задать свойство «Коэффициент усиления» равным «0.35/10».

Запустить проект на расчёт.

График переходного процесса (Рисунок 19) свидетельствует, что в нелинейной САР установились высокочастотные автоколебания с амплитудой примерно «0.05», что в пять раз превышает ширину зону нечувствительности в блоке с подписью «Управляющем реле».

Резюме: нелинейная САР не отработала управляющее воздействие потому, что при значении «k» равным «0.35/10» не выполняются условия критерия абсолютной устойчивости В.М. Попова (точка «–1/K» расположена внутри годографа Попова).