Диссертация (1149675), страница 8

Просмтор этого файла доступен только зарегистрированным пользователям. Но у нас супер быстрая регистрация: достаточно только электронной почты!

Текст из файла (страница 8)

. , Rk (t), а начальным условием решения задачи Коши (под E понимается единичная матрица) соотношенияR0 (0) = 0, R1 (0) = E, R2 (0) = . . . = Rk (0) = 0.(3.4)48Данное соотношение очевидным образом следует из тождественности любого вектора состояния системы самому себе в начальный момент времени X(0) =X0 = EX0 . Для вывода системы дифференциальных уравнений, описывающихизменение матричного отображения во времени, необходимо продифферецировать решение (3.3) в силу системы (3.2). Производная по времени от (3.3) ивыражение (3.1) дают два равенства для производной вектора состоянияddddd[2][k]X = R0 (t) + R1 (t)X0 + R2 (t)X0 + .

. . + Rk (t)X0 ,dtdtdtdtdtdX =P0 (t) + P1 (t)X + P2 (t)X[2] + . . . + Pp (t)X[p] =dt=P0 (t)+()[2][k]+P1 (t) R0 (t) + R1 (t)X0 + R2 (t)X0 + . . . + Rk (t)X0 +()[2][2][k]+P2 (t) R0 (t) + R1 (t)X0 + R2 (t)X0 + . . . + Rk (t)X0+(3.5)+...+()[p][2][k].+Pp (t) R0 (t) + R1 (t)X0 + R2 (t)X0 + . . . + Rk (t)X0После раскрытия всех скобок и приведения подобных слагаемых получим систему обыкновенных дифференциальных уравнений, описывающих динамикуматриц R0 , R1 , R2 , .

. . , Rk∑d[i]R0 (t) =Pi (t)R0 ,dti=1p∑∂X[i]dRk (t) =Pi (t),[i] ∗dt∂(X)i=1(3.6)pk = 1, 2, . . .0Численно решая эту систему уравнений с начальным условием (3.4), можно получить искомое матричное отображение. Скобки в выражении (3.5) раскрываются по математическим правилам, применительно к кронекеровским операциям [3] аддитивности и коммутативности с учетом редуцирования размерности.Так соотношение AX[i] ⊗ BX[j] = (A ⊗ B)(X[i] ⊗ Y[j] ) следует записать в видеAX[i] ⊗ BX[j] = CX[i+j] , где матрица C получается из матрицы A ⊗ B сум-49мированием и редуцированием соответствующих элементов.

В данной реализации алгоритма подобные операции предлагается осуществлять алгоритмически[2][k]с учетом того, что решение X = R0 (t) + R1 (t)X0 + R2 (t)X0 + . . . + Rk (t)X0после перемножения всех матриц можно записать в виде вектораR1 0+(R11 , X0 )+ ... +[k](Rk1 , X0 )X=,...[k]nnnR0 + (R1 , X0 ) + . . . + (Rk , X0 )где под Rij понимается j-ая строка матрицы Ri , а под оператором (·, ·) — скалярное произведение векторов.Данный вектор уже легко может быть возведен в степень с учетом редуцирования размерности. После этой операции новый вектор, представляющий[k]искомую кронекеровскую степень, раскладывается по компонентам X0 , . .

. , X0 .3.2 Реализация алгоритмаВ данном параграфе приводится описание алгоритма на псевдокоде, а также рассмотрены его реализации на языках программирования Python и C#/C++.Язык Python выбран в виду своей простоты и используется для быстрого прототипирования алгоритма. На языках C#/C++ реализуются вычислительные библиотеки, которые впоследствии встраиваются в среду моделирования.3.2.1 Описание алгоритма на псевдокодеАлгоритм построения матричного отображения, приведенный в параграфе3.1 представляет собой метод пошагового интегрирования системы обыкновенных дифференциальных уравнений, применительно к системе (3.6).

Вычислениеправых частей заключается в перемножении соответствующих матриц и векторов с сохранением результата до заданного порядка нелинейности.50Алгоритм 1: численный метод вычисления матричного отображенияData: P0 , P1 , . . . , Pp , T, hResult: R0 (T ), R1 (T ), . . . , Rk (T )1t=02R0 = 0, R1 = E, .

. . , Rk = 03order = max(p, k)4X[1] , X[2] , . . . X[order] // генерируем список кронекервоских векторов5while t <= T do6P0 (t), P1 (t), . . . , Pp (t) −→ P0 , P1 , . . . , Pp // вычисляем матрицы7X1 = R0 + R1 X[1] + R2 X[2] + . . . + Rk X[k]8for i = 1 : order – 1 doXi+1 = X ⊗ Xi910end11for i = 0 : k do12DRi = 013for j = 0 : p do14R : Xp = RX[k] + . . . // вычисляем матрицу R15DRi = DRi + Pp Rend1617end18[R0 , R1 , . . . , Rk ] = integrate([DR0 , DR1 , . .

. , DRk ], h)19end20R0 (T ) = R0 , R1 (T ) = R1 , . . . , Rk (T ) = RkВ данном алгоритме под T понимается интервал, а под h — шаг интегриро-вания. Начальными данными являются матрицы Pi , выходными — Ri . Функцияintegrate осуществляет пересчет состояния на текущем шаге и может реализовываться любым пошаговым методом интегрирования, например, для простейшегометода Эйлера она будет выглядеть как Ri = Ri + hDRi .513.2.2 Реализация алгоритма на языках Python и C#/C++Описанный алгоритм реализован средствами таких языков программирования, как Python и C#/C++.

Python представляет собой гибкий интерпретируемыйязык, основными используемыми структурами данных в котором являются списки и словари, позволяющие быстро реализовать и протестировать разрабатываемый алгоритм. Построение вычислительных библиотек на языках C#/C++ позволяет достигнуть максимальной производительности вычислений. C# в данномслучае использован в небезопасном режим (unsafe mode), когда используютсяуказатели и явное управление памятью.Реализация на языке PythonОсновными структурами используемых данных являются вектор состояния(и его кронекеровские степени) и однородный полином.

Вектора состояния описываются списком, состоящим из строк. Каждая строка задает моном и представляет собой последовательностьзаписанных через пробел.переменных, x32x 2  x y    −→ [’x x x’ , ’x x y’, ’x y y’, ’y y y’].−→[’xx’,’xy’,’yy’],xy  2 xy  y2y3Для упорядочивания мономов, используется лексикографический порядок.Однородный полином описывается словарем, в котором ключем является строкамоном, а значением — вещественный коэффициент. Свободному члену соответствует пробел, например 3x2 + xy + 1 −→{’x x’:3.0, ’x y’:1.0, ’ ’:1.0}.

Описанный алгоритм реализован в виде библиотеки mode.py, которая опубликована воткрытом доступе. Программный скрипт содержит реализацию функций и класса, позволяющий задать произвольную систему дифференциальных уравненийв виде полиномиального разложения своих правых частей. Подробное описаниебиблиотеки с приведением модульных тестов для каждой функции представленов приложении B.52Реализация на языках C#/C++Основным элементом в логике работы программы является символьный полином. Эта структура данных представляет собой динамический массив, организованный по принципу работы хэш-таблицы, элементами которой являютсяслагаемые полинома. Хэш-коды определены на основе строкового представления слагаемых, которые упорядочиваются в лексикографическом порядке.

Каждое слагаемое реализовано как список переменных и числовой коэффициент.Для указанной структуры данных определены операции сложения, умножения,деления, подстановки полинома вместо определенной переменной и вычисление коэффициента при заданном списке параметров. Все операции выполняются до предопределенного порядка нелинейности. Вычисление кронекеровскойстепени вектора X[k] осуществляется на основе возведения полинома в степень(x1 + .

. . + xn )k , приведения подобных слагаемых и взятия частных производныхс целью преобразования результата к векторному представлению. Такой подходпозволяет унифицировать алгоритмы, связанные с кронекеровскими операциями, и производить необходимые вычисления автоматически для общего случая.Реализация структуры данных в виде хэш-таблицы позволяет увеличить скорость доступа к элементам полинома по ключу-моному. При реализации используется небезопасный код, адресация элементов массивов производится поуказателям.

В этом случае проверки выхода индекса за границы массива непроизводятся, что исключает потерю производительности, которая присутствует при использовании стандартных средств языка C#. Таким образом, код алгоритма написан в C-подобной нотации. Методы, описанные выше в данномпараграфе, реализованы в виде dll-библиотек, содержащих самодокументируемые XML-комментарии. Скомпилированные .NET-компоненты являются частьюразработанного программного комплекса и более подробно описаны в главе 4данной диссертации. Функционал этих библиотек полностью повторяет логикубиблиотеки mode.py, однако выигрывает у нее по производительность.533.3 Верификация алгоритма на модельных задачахВ качестве модельных выбраны некоторые наиболее хорошо изученные задачи нелинейной динамики.

Для численного анализа систем используются матричное отображение (на графиках представлены траекториями красного цвета)и метод Рунге – Кутта 4 порядка (траектории синего цвета).Скалярное однородное нелинейной уравнениеВ качестве первого примера рассмотрим применение алгоритма, описанногов этой главе, для уравнения видаdx= kx2 , x(0) = x0 ,dt(3.7)гдe k — вещественный коэффициент, x ∈ R1 . Будем искать решение уравнения ввиде матричного степенного ряда, который представляется в скалярном видеx = a + bx0 + cx20 + ex30 + . . .Подставляя это решение в уравнение и дифференцируя его получимdx da dbdcde=+ x0 + x20 + x30 + . .

. ,dtdt dtdtdtdx= k(a + bx0 + cx20 + ex30 + . . .)2 .dtПриводя теперь подобные слагаемые, можно построить систему уравнений (3.6)dbdcdeda= ka2 ,= 2kab,= k(2ac + b2 ),= k(2ad + 2bc), . . .dtdtdtdtУчитывая, что начальным условием для этой системы является тождественноепреобразование a = 0, b = 1, c = 0, d = 0, . . ., можно последовательно разрешить систему и получить решение в форме отображенияa(t) = 0, b(t) = 1, c(t) = kt, e(t) = k 2 t2 , .

Характеристики

Список файлов диссертации