Размещено на http://www.allbest.ru/

Постановка задачи

Описание модели

модель описывается системой дифференциальных уравнений 2-ого порядка;

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

Исследование стационарных точек

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

выполняется линеаризация системы уравнений в стационарных точках.

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

Обозначения:

- собственные значения якобиана ().

Устойчивый узел:

Собственные значения действительные (мнимая часть равна нулю), меньшие нуля:

Неустойчивый узел:

Собственные значения действительные, большие нуля:

Устойчивый фокус:

Собственные значения комплексно-сопряженные, действительная часть меньше нуля:

Неустойчивый фокус:

Собственные значения комплексно-сопряженные, действительная часть больше нуля:

"Центр”:

Собственные значения комплексно-сопряженные, действительная часть равна нулю:

"Седло”:

Собственные значения действительные, с разными знаками:

Построение фазовых траекторий

Фазовые траектории (портреты) строятся в окрестностях стационарных точек с помощью метода численного интегрирования.

Практическая реализация задачи

Задача реализована в среде системы Тополог.

Описание модели

Описание системы уравнений выполняется с использованием стандартных типов системы Тополог: функция и вектор. В виде функции описывается только правая часть уравнения.

Пример: пусть мы имеем систему из двух нелинейных уравнений вида:

тогда описание данной системы на языке системы Тополог будет выглядеть следующим образом:

данные

x: вещ

y: вещ

F1: функция (x, y) = x^2/100 + y^2/25 - 1

F2: функция (x, y) = y - 10 * Sin (X)

VF: вектор = (F1, F2) |

J: якобиан (VF)

Исследование стационарных точек

Решение системы уравнений в заданной прямоугольной области выполняется стандартным операционным блоком системы Найти_Все_Решения. Он имеет следующий синтаксис:

Найти_Все_Решения (<вектор функций>, <погрешность>,

<граница 1>, <граница 2>,<граница 3>, <граница 4>,

<шаг сетки>, <массив решений>)

Решения ищутся путем запуска метода Ньютона из узлов сетки, наложенной на прямоугольную область поиска, и накопления найденных различных решений. Прямоугольная область поиска задается параметрами: <граница 1>, <граница 2>,<граница 3>, <граница 4>. Шаг сетки равномерный по обеим координатам и задается параметром <шаг сетки>. После завершения работы блока все решения скапливаются в массиве точек, описанным как:

Решения: массив из Точка_вещ

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

Для получения собственных значений якобиана системы в найденных стационарных точках используется стандартный операционный блок системы Тополог Собственные_Значения. Он имеет следующий синтаксис:

Собственные_Значения (<матрица>,

<максимальное число итераций>, <погрешность>) - >

(<массив реальных частей>,

<массив мнимых частей>)

На выходе данный блок формирует два массива, содержащих реальные и мнимые части собственных значений квадратной матрицы (в нашем случае якобиана).

Численные методы нахождения собственных значений матриц делятся на два типа: прямые и итерационные. Прямые методы находят собственные значения путем решения характеристического уравнения матрицы с помощью одного из численных методов (метод Ньютона (касательных), метод секущих, метод хорд и др.). Итерационные методы при помощи преобразований, не меняющих собственных значений матрицы, приводят ее к диагональному (в случае действительных собственных значений) или блочно-диагональному (при наличии комплексно-сопряженных собственных значений) виду. Значения из диагонали полученной матрицы и являются собственными значениями исходной матрицы.

В основу блока Собственные_Значения заложен итерационный метод нахождения собственных значений по методу Якоби с понижением нормы для действительных матриц.

Основная суть метода Якоби:

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

Для реализации такого преобразования матрицы можно выбирать в виде произведения двух матриц , где матрица определяет так называемое вращение, а - сдвиг или комплексное вращение. Элементы этих двух матриц удовлетворяют соотношениям:

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

.

Параметр сдвига удовлетворяет соотношению

, где

; ;

;

.

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

Построение фазовых траекторий

Фазовые траектории образуются в процессе решения задачи интегрирования системы дифференциальных уравнений, описывающих исследуемую модель. Для этого использован метод численного интегрирования Рунге-Кутта четвертого порядка.

Алгоритм метода может быть представлен в следующем виде:

,,

,,

;

где -число уравнений, -номер шага,

,

,

,

.

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

Примеры.

Пример задачи из электротехники (линейная система)

Пусть необходимо исследовать переходные процессы и выявить стационарные состояния электротехнической цепи:

Якобиан имеет вид:

Варианты возможных переходных процессов представлены ниже.

Устойчивый узел

Начальные условия:

Собственные значения:

Параметры интегрирования:

нелинейная система анализ уравнение

Неустойчивый узел

Начальные условия:

Собственные значения:

Параметры интегрирования:

Устойчивый фокус

Начальные условия:

Собственные значения:

Параметры интегрирования:

Неустойчивый фокус

Начальные условия:

Собственные значения:

Параметры интегрирования:

"Центр”

Начальные условия:

Собственные значения:

Параметры интегрирования:

"Седло”

Реальной схема, описываемая исследуемой системой уравнений, не может иметь стационарную точку типа "седло”.

Приведенный пример был реализован в системе Тополог с помощью следующего описания модели:

ИМЯ "Решение системы дифференциальных уравнений"

ДАННЫЕ

{параметры модели}

E0: вещ R1: вещ R2: вещ L: вещ C: вещ k: вещ

{функциональные переменные и их рабочие значения}

Il: вещ Uc: вещ Iln: вещ Ucn: вещ

{параметры вычислительного эксперимента}

t: вещ tt: вещ tнач: вещ tкон: вещ h: вещ

{рабочие переменнные метода интегрирования}

di1: вещ di2: вещ di3: вещ di4: вещ

du1: вещ du2: вещ du3: вещ du4: вещ

{оптимизированный вариант}

dIldt: функция (Il,Uc) = (Il* (k - R1) - Uc + E0) / L

dUcdt: функция (Il,Uc) = (Il* (k + R2) - Uc + E0) / (C*R2)

{вектор производных}

V: вектор = (dIldt,dUcdt) |

{матрица Якоби}

J: якобиан (V)

{рабочие массивы}

Re: массив из вещ

Im: массив из вещ

{ Описание переменных и объектов для графического вывода }

i: цел n: цел

x: вещ y: вещ

mi: вещ mu: вещ

Точка: точка_вещ

Il_t: ломаная_вещ

Uc_t: ломаная_вещ

IlUc: ломаная_вещ

МОДЕЛЬ

{ параметры цепи }

E0 = 10.0

R1 = 1.0

R2 = 1.0

L = 1.0

C = 1.0

k = 1.5

{ начальные условия }

tнач = 0.0

Il = 0.0

Uc = 0.0

t = tнач

{ масштабные коэффициенты }

mi = 0.2

mu = 0.2

Собственные_Значения (J, 300, 1E-15) - > (Re, Im)

{ параметры вычислительного эксперимента }

h = 0.1

tкон = 10

{ установка начальной точки графиков Il (t) и Uc (t) }

Изменить_Абсциссу (t) - > (Точка)

Изменить_Ординату (Il*mi) - > (Точка)

Добавить_Элемент (Точка) - > (Il_t)

Изменить_Ординату (Uc*mu) - > (Точка)

Добавить_Элемент (Точка) - > (Uc_t)

Обновить (графика)

t = tнач

{ расчетный цикл интегрирования }

пока t < tкон

начало

{ расчет методом Рунге-Кутта 4-го порядка }

Iln = Il Ucn = Uc

{ 1-й подшаг }

tt = t

di1 = h * dIldt du1 = h * dUcdt

{ 2-й подшаг }

tt = t + h/2.

Il = Iln + di1/2. Uc = Ucn + du1/2.

di2 = h * dIldt du2 = h * dUcdt

{ 3-й подшаг }

tt = t + h/2.

Il = Iln + di2/2. Uc = Ucn + du2/2.

di3 = h * dIldt du3 = h * dUcdt

{ 4-й подшаг }

tt = t + h

Il = Iln + di3 Uc = Ucn + du3

di4 = h * dIldt du4 = h * dUcdt

Il = Iln + (di1 + 2. *di2 + 2. *di3 + di4) /6.0

Uc = Ucn + (du1 + 2. *du2 + 2. *du3 + du4) /6.0

t = t + h

Изменить_Абсциссу (t) - > (Точка)

Изменить_Ординату (Il*mi) - > (Точка)

Добавить_Элемент (Точка) - > (Il_t)

Изменить_Ординату (Uc*mu) - > (Точка)

Добавить_Элемент (Точка) - > (Uc_t)

конец

Число_Элементов (Il_t) - > (n)

цикл i = 1 до n

начало

Ордината (Il_t [i]) - > (y)

Ордината (Uc_t [i]) - > (x)

Добавить_точку (y,x) - > (IlUc)

конец

Обновить (ГРАФИКА)

КОНЕЦ

Пиримеры построения фазовых портретов для нелинейных систем

Дана система дифференциальных уравнеий:

Модель для системы Тополог:

ИМЯ "Решение системы дифференциальных уравнений"

ДАННЫЕ

x1: вещ

x2: вещ

x1s: вещ

x2s: вещ

x1_мин: вещ = - 1.0

x1_макс: вещ = 3.5

x2_мин: вещ = - 1.0

x2_макс: вещ = 1.5

Re_1: вещ

Re_2: вещ

R: вещ

Шаг_Сетки: вещ = 0.2

{ параметры вычислительного эксперимента }

t: вещ

tt: вещ

tнач: вещ

tкон: вещ

h: вещ

{ рабочие переменнные метода интегрирования }

X1n: вещ X2n: вещ

di1: вещ di2: вещ di3: вещ di4: вещ

du1: вещ du2: вещ du3: вещ du4: вещ

dX1: функция (x1, x2) = (3 - x1 - 5 * x2) * x1

dX2: функция (x1, x2) = (2 - 4 * x1 - 2 * x2) * x2

{ вектор производных }

V: вектор = (dX1, dX2) |

{ матрица Якоби }

J: якобиан (V)

{ рабочие массивы }

Re: массив из вещ

Im: массив из вещ

Решения: массив из точка_вещ

{ Описание переменных и объектов для графического вывода }

Число_Решений: цел

Число_Собственных_Значений: цел

Текущее_Решение: цел

Фазовый_Портрет: массив из ломаная_вещ

Текущая_Ломаная: ломаная_вещ

Устойчивые_Узлы: массив из точка_вещ

Неустойчивые_Узлы: массив из точка_вещ

Устойчивые_Фокусы: массив из точка_вещ

Неустойчивые_Фокусы: массив из точка_вещ

Центры: массив из точка_вещ

Седла: массив из точка_вещ

МОДЕЛЬ

{ начальные условия }

tнач = 0.0

{ параметры вычислительного эксперимента }

h = 0.05

tкон = 2

Найти_Все_Решения (V, 1E-8, x1_мин, x1_макс,

x2_мин, x2_макс, Шаг_Сетки, Решения)

Число_Элементов (Решения) - > (Число_Решений)

цикл Текущее_Решение = 1 до Число_Решений

начало

Абсцисса (Решения [Текущее_Решение]) - > (x1)

Ордината (Решения [Текущее_Решение]) - > (x2)

Собственные_Значения (J, 300, 1E-15) - > (Re, Im)

Число_Элементов (Re) - > (Число_Собственных_Значений)

если Число_Собственных_Значений = 2 то

начало

Re_1 = Re [1]

Re_2 = Re [2]

R = Re_1 * Re_2

если Im [1] = 0.0 то { вещественные СЗ }

если R < 0.0 то

Добавить_Элемент (Решения [Текущее_Решение]) - >

(Седла)

иначе

если Re_1 < 0.0 то

Добавить_Элемент (Решения [Текущее_Решение]) - >

(Устойчивые_Узлы)

иначе

Добавить_Элемент (Решения [Текущее_Решение]) - >

(Неустойчивые_Узлы)

иначе { комплексные СЗ }

если Re_1 < 0.0 то

Добавить_Элемент (Решения [Текущее_Решение]) - >

(Устойчивые_Фокусы)

иначе

если Re_1 = 0.0 то

Добавить_Элемент (Решения [Текущее_Решение]) - >

(Центры)

иначе

Добавить_Элемент (Решения [Текущее_Решение]) - >

(Неустойчивые_Фокусы)

конец

конец

x1s = 0

пока x1s < x1_макс

начало

x2s = 0

пока x2s < x2_макс

начало

гудок

x1 = x1s

x2 = x2s

Очистить_Массив (Текущая_Ломаная)

Добавить_Точку (x1, x2) - > (Текущая_Ломаная)

t = tнач

{ расчетный цикл интегрирования }

пока t < tкон

начало

{ расчет методом Рунге-Кутта 4-го порядка }

X1n = x1 X2n = x2

{ 1-й подшаг }

tt = t

di1 = h * dX1 du1 = h * dX2

{ 2-й подшаг }

tt = t + h / 2.0

x1 = X1n + di1/2.0 x2 = X2n + du1/2.0

di2 = h * dX1 du2 = h * dX2

{ 3-й подшаг }

tt = t + h / 2.0

x1 = X1n + di2/2.0 x2 = X2n + du2/2.0

di3 = h * dX1 du3 = h * dX2

{ 4-й подшаг }

tt = t + h

x1 = X1n + di3 x2 = X2n + du3

di4 = h * dX1 du4 = h * dX2

x1 = X1n + (di1 + 2.0*di2 + 2.0*di3 + di4) /6.0

x2 = X2n + (du1 + 2.0*du2 + 2.0*du3 + du4) /6.0

t = t + h

Добавить_Точку (x1, x2) - > (Текущая_Ломаная)

конец

Добавить_Элемент (Текущая_Ломаная) - >

(Фазовый_Портрет)

x2s = x2s + Шаг_Сетки

конец

x1s = x1s + Шаг_Сетки

конец

Обновить (ГРАФИКА)

КОНЕЦ

Фазовый портрет:

Дана система:

Модель для системы Тополог взята из предыдущего примера с параметрами:

x1_мин: вещ = - 1.0

x1_макс: вещ = 5.5

x2_мин: вещ = - 1.0

x2_макс: вещ = 2.5

dX1: функция (x1, x2) = (4 - 2.5 * x2) * x1

dX2: функция (x1, x2) = (-2 + x1) * x2

h = 0.01

tкон = 2

Получен фазовый портрет:

Размещено на Allbest.ru