Сер. 10. 2009. Вып. 3
ВЕСТНИК САНКТ-ПЕТЕРБУРГСКОГО УНИВЕРСИТЕТА
УДК 519.683.8 Е. И. Петренко
ОПТИМИЗАЦИЯ АЛГОРИТМОВ ПОСТРОЕНИЯ ИНВАРИАНТНЫХ МНОЖЕСТВ ДИНАМИЧЕСКИХ СИСТЕМ С ПОМОЩЬЮ ГЕНЕРАЦИИ КОДА
В программном комплексе исследования динамических систем с помощью символического образа - ориентированного графа, построенного по системе, - приходится использовать большое количество циклов по размерности, перебирать вершины гиперкуба этой размерности, хранить массивы соответствующего размера.
Цель нашей работы - применить генерацию кода на лету к задаче локализации цепно-рекуррентного множества с помощью символического образа. Мы применяем динамическую кодогенериацию, т. е. программа получает функцию системы, по размерности которой генерируется код некоторых компонент, компилирующийся и потом загружающийся в исходную программу.
Рассматриваемый подход представляет собой простейшую реализацию смешанных вычислений. В отличие от [1] мы не пытаемся построить общий смешанный вычислитель, а делаем акцент на то, что при фиксированной размерности можно сгенерировать реализацию многих алгоритмов и структур данных программного комплекса в виде явного перебора полей и/или набора простейших команд.
Методы построения и анализа символического образа даны в [2-5]. В частности, в работах [4, 5] описан способ, который позволяет повысить скорость построения символического образа. Для этого предлагается строить символический образ для большой степени исходного отображения. В [6] приводится подробное описание алгоритмов построения и анализа символического образа для общего случая. Следует отметить, что в реализации описанных подходов генерация кода не использовалась.
1. Применение генерации кода. Генерация кода используется во многих задачах. В компьютерной графике для построения проекции трехмерного пространства (сцены) на плоскость при помощи метода трассировки лучей оказывается выгодным специализировать алгоритм просчета пересечений для конкретной сцены [1]. В системах управления базами данных исполнитель запросов может генерировать план выполнения запроса на специальном внутреннем языке; в нейронных сетях для обучения сети оказывается удобным специализировать программу на конкретную топологию сети. В численных исследованиях, например в задаче n тел, для проведения расчетов эффективнее выполнять специализированную версию программы, в которой были подставлены значения всех параметров.
Петренко Евгений Игоревич — аспирант кафедры информатики математико-механического факультета Санкт-Петербургского государственного университета. Научный руководитель: Н. Б. Ампилова. Научные направления: символический образ динамической системы и приложения, алгоритмы исследования динамических систем, оптимизация программ, эффективная реализация, распределенные вычисления, проектирования сложных систем, разработка пользовательских программ в области исследования динамических систем. E-mail: [email protected].
© Е. И. Петренко, 2009
Созданный программный комплекс по исследованию динамических систем разработан на платформе Microsoft .NET Framework 3.5, использован язык C# 3.0.
Для реализации генерации кода мы используем встроенный в платформу .NET Framework компилятор языка C#.
1.1. Схема работы программного комплекса. В основной части программы определены интерфейсы компонент системы. По набору параметров программа генерирует код на языке C# для реализации специализированных версий некоторых интерфейсов. При помощи платформы .NET Framework этот код компилируется и загружается в память приложения. Для применения кода необходимо создать только классы, которые реализуют интерфейс компонент. Это можно сделать всего один раз, если использовать шаблон проектирования (design pattern) абстрактная фабрика (abstract factory) [7] и рефлексию (reflection) [8].
1.2. Достоинства подхода. В описанной системе можно подменять реализации только некоторых компонент. Это позволяет проводить разработку постепенно и по мере надобности. Вместе с тем в основной части программы можно разместить удобные базовые классы, которые будут расширены в сгенерированном коде, что позволит избежать генерации сложного кода.
Благодаря применению компонентной модели, в некоторых частях системы можно как включить, так и отключить генерацию кода, что дает возможность провести сравнительный анализ производительности системы: с использованием динамической генерации кода и без нее.
2. Определения. Пусть D - компактное подмножество m-мерного пространства Rm, f : D С Mm ^ D, f - диффеоморфизм.
Пусть C = {Di,..., Dk} - конечное покрытие компакта D замкнутыми множествами. Множества Di будем называть ячейками. По множеству D, диффеоморфизму f и множеству ячеек C можно построить ориентированный граф следующим образом [2] : каждой ячейке Di соответствует узел i графа. Между узлами i и j есть ориентированное ребро тогда и только тогда, когда f (Di) П Dj = 0. Построенный ориентированный граф будем называть символическим образом [2] динамической системы. Для непрерывных систем символический образ строится для отображения сдвига на фиксированный интервал по траекториям системы.
Приближение к инвариантному множеству динамической системы можно получить, если взять все ячейки, которые соответствуют вершинам, содержащимся в компонентах сильной связности символического образа [2]. Для построения инвариантного множества используется техника подразбиений [6].
3. Построение образа ячейки. Построение символического образа динамической системы опирается на методы построения образа ячейки. От их реализации зависят как точность, так и скорость построения символического образа.
Способы построения образа ячейки более подробно рассмотрены в [6, 9, 10]. Приведем описание некоторых из них, приведенных в [6]. Для каждого метода дано описание выполненной нами доработки.
3.1. Линейный метод [6]. В этом методе оцениваются возможные значения функции на заданной ячейке. Для этого производятся следующие действия:
1. Берем ячейку разбиения Di.
2. Вычисляем значения функции системы f в вершинах Di.
3. Строим минимальный параллелепипед E, ориентированный по осям координат, который содержит образы вершин Di.
4. Образом ячейки будем считать множество ячеек покрытия, которые пересекаются с построенным параллелепипедом D.
3.2. Ускорение работы линейного метода. Во всех методах построения образа ячейки требуется перебирать все вершины некоторого многомерного параллелепипеда.
Реализация такого алгоритма для общего случая, с одной стороны, будет содержать большое количество циклов по размерности и условных переходов. С другой - для фиксированной размерности можно явно выписать список координат всех вершин. Для того чтобы минимизировать количество присваиваний, используются коды Грея [11]. Это позволяет осуществить перебор вершин, делая только одно присваивание для вычисления координаты каждой следующей вершины.
Далее, если в программе нужно перебрать все вершины некоторого параллелепипеда, реализацию алгоритма перебора можно получить из абстрактной фабрики (abstract factory) [7], которая позволяет создавать сгенерированные для конкретной размерности реализации.
Отметим, что в программном комплексе всегда используется динамическая генерация специализированных реализаций этого алгоритма.
3.3. Адаптивный метод [6]. Пусть задано некоторое множество точек P С D и рефлексивное симметричное отношение ф на множестве P. Будем называть такое отношение ф отношением соседства. Пусть задана функция правой части системы f = (fi,f2,---,fm). Фиксируем £ = (е*)”=1, £i > 0. Тогда будем называть множество P е-множеством для функции f и отношения ф, если для любых двух соседних (ф(х,у)) точек x,y € P выполнено соотношение
\fi(x) - fi(y)\ < £i.
Будем обозначать такое множество через Pf ф.
Пусть задан набор m положительных чисел е = (ei)m=i. Тогда е-окрестностью точки х = (xi, Х2,..., xm) € Rm будем называть множество
[xi — £1, xi + £1] х [x2 — £2, x2 + £2] x [xm — ет, xm + em] С Rm.
Пусть заданы множество Pf ,ф и ячейка Dj, все вершины которой принадлежат множеству P. Тогда образ ячейки будет состоять из £-окрестностей образов точек, которые лежат в P П Dj. При этом, если элементы £ достаточно малы, т. е.
di
* 2 ’
к результату следует добавлять соседние ячейки для точек, которые попадают близко к границе ячейки.
Введенное отношение соседства ф можно представить как неориентированный граф, где вершины соответствуют точкам множества P, а ф^,у) эквивалентно существованию ребра между вершинами x и y или совпадению этих вершин x = y. В [6] приводится подробное описание алгоритма построения Pf ф.
Вершины графа разбиения соответствуют точкам фазового пространства исследуемой системы. Результат и скорость работы адаптивного метода зависит от выбора начального разбиения ячейки. Опишем графы разбиения ячейки для разных случаев:
1) для одномерного случая достаточно взять взять граф с одним ребром и двумя вершинами, которые соответствуют границам одномерной ячейки;
2) в двумерном случае будем рассматривать граф из 5 вершин - центр ячейки и ее четыре вершины и ребра, которые связывают соседние вершины;
3) в n-мерном случае рассматриваются все вершины ячейки и центр ячейки. Соединяем ребром соседние вершины.
3.4. Ускорение работы адаптивного метода. Для этого метода определенную сложность представляет алгоритм построения начального разбиения ячейки, так как его придется повторять для каждой ячейки (скажем, сотни тысяч раз) при построении символического образа динамической системы.
Как было указано выше, начальное разбиение ячейки включает в себя и все ребра ячейки как многомерного параллелепипеда, и ребра от каждой вершины к центру. Наша задача - построить граф соседних точек отношения PJ ^.
При помощи простейшей реализации алгоритма построения начального разбиения можно получить список всех ребер для графа отношения P^ ^. Далее генерируем реализацию алгоритма построения начального разбиения для фиксированной размерности, в коде которой явно перечислены все ребра.
Отметим, что в программном комплексе присутствует только одна, динамически генерируемая, реализация построения начального разбиения ячейки для адаптивного метода.
4. Представление ячейки. В реализации представления ячейки мы полагаем, что все ячейки - одинаковые параллелепипеды, ориентированные по осям координат, а область исследования - параллелепипед, ориентированный по осям координат. Целочисленная система координат [6] получается, если взять систему координат, в которой размер ячейки по любому из направлений будет равен 1 и оси параллельны осям координат. Преобразование координат ячейки из целочисленной системы и обратно - линейное.
4.1. Реализация системы координат. Генерируем код реализации объекта для координаты ячейки в виде объекта с m полями, где m - размерность системы. Это позволит заменить массивы на наборы полей и не тратить дополнительную память на хранение объекта-представления массива в .NET Framework [12]. Генерируемая реализация системы координат будет использовать поля сгенерированного объекта для координаты напрямую.
Методы построения образа ячейки так или иначе решают задачу разбиения того или иного пространственного объекта на ячейки. Чтобы избежать дублирования кода был выделен дополнительный уровень абстракции. Основные примитивы преобразования пространственных областей в наборы ячеек:
1. RadiusProcessor - позволяет представить в виде ячейки шар, или е-окрестность (см. п. 3.3);
2. RectProcessor - позволяет представить параллелепипед в виде набора ячеек;
3. PointProcessor - позволяет представить точку в виде ячейки или набора ячеек. Здесь есть два способа:
3.1) простой (point) - всегда возвращает не более одной ячейки;
3.2) с наложением (overlaped) - возвращает ячейку, в которую попадает точка. Если точка близка к границе, то возвращаются все ячейки с этой границей.
Реализация данных примитивов для общего случая описана в [6]. Для всех примитивов мы генерируем реализацию, которая напрямую работает со сгенерированным представлением координаты. Во всех примитивах, кроме PointProcessor, в генерируемом коде применяем вложенные циклы, количество которых равно размерности системы.
5. Анализ производительности программного комплекса в целом. В работе [5] приведено описание способа повышения скорости построения символического образа. Поскольку в ней не приведены какие-либо численные данные о скорости построения символического образа, мы не можем провести сравнение результатов. Будем сравнивать работу нашей системы в режиме с кодогенерацией и в обычном [6]. Как было упомянуто ранее, некоторые компоненты реализованы только с помощью динамической генерации кода. Потому сравниваем производительность системы при реализации системы координат с генерацией кода и без нее. В качестве примеров рассмотрим отображение Икеда [13], трехмерную систему Линдстрома - «хищник-жертва-суперхищник» [14], систему Лоренца [15] и биологическую модель с задержкой. Построение образа ячейки будем проводить двумя описанными методами: линейным и адаптивным. В качестве результата эксперимента будем выписывать ускорение метода с динамической кодогенерацией по сравнению с реализацией из [6].
5.1. Отображение «хищник-жертва-суперхищник» (система Линдстрома [14]). В работах [2, 5, 14] рассматривается следующая динамическая система:
х ______цохехр(-у)____
1 + ж max(exp( — y),k(z)k(y)) ’
y ^ xy exp(—z)k(4yz)k(y),
z ^ Ц2yz,
где k(t) = LjEELJ}^ і^0 = з; 4001; ці = 1 и Ц2 = 4.
На множестве D = [— 1,4] x [-1, 4] x [0,1.6] при заданных параметрах система имеет цепно-рекуррентное множество [4, 5, 15]. Проведем 8 и 10 шагов его построения.
5.2. Система Лоренца. В работах [2, 4, ,5, 15] рассматривается следующая динамическая система:
X = a(y — x),
y = x(r — z) — y,
z = by — bz,
в которой a = 10, b = 8/3 и r = 20.
На множестве D = [—35, 35] x [—35, 35] x [—35, 35] при заданных параметрах система имеет цепно-рекуррентное множество [4, ,5, 15]. Для этой системы используем метод Рунге-Кутта, 100 итераций с шагом 0.001, т. е. строим символический образ для отображения сдвига на t = 1. Проведем 8 шагов построения инвариантного множества.
5.3. Отображение Икеда. Оно имеет вид [13]
Сз
т(x, y) = C\ -
1 + x2 + y2 ’
x ^ d — C2(x cos t(x,y) — y sin т(x,y)), y ^ C2(xsinт(x,y)+y cosт(x,y)),
где d =2; Ci = 0,4; C2 = 0, 9; C3 = 6.
В области D = [—10,10] x [—10,10] и при заданных параметрах отображение имеет цепно-рекуррентное множество [13, 16]. Проведем по 8 и 10 шагов его построения для каждого из методов построения образа ячейки.
5.4. Биологическая система с задержкой. Ее можно представить следующим образом:
т(а, Ь) = 1.16 ^2 - 1+0.09(а+Ь)) ’ x ^ Т(x,y), y ^ т(У,т(x,y))-
Эта система была предложена Г. С. Осипенко. Она имеет инвариантное множество в области [0.01, 2] х [0, 2]. Проведем 8 и 10 шагов его построения.
5.5. Результаты эксперимента. В таблице приведено отношение скорости работы алгоритма построения инвариантного множества с использованием динамической кодогенериации к скорости и без нее.
Результаты численного эксперимента
Система Кол-во шагов Метод Ускорение
Линдстрома 8 Лин. 1.447
Адапт. 1.524
10 Адапт. 1.173
Лоренца 8 Лин. 1.080
Адап. 1.153
Икеда 8 Лин. 1.667
Адапт. 1.384
10 Лин. 1.151
Адапт. 1.242
С задержкой 8 Лин. 1.316
Адапт. 1.345
10 Лин. 1.234
Адапт. 1.251
Из данных таблицы видно, что реализация программного комплекса с использованием генерации кода на лету приводит к ускорению работы комплекса по исследованию динамических систем. В любом случае динамическая генерация кода дает возможность экономить на оперативной памяти, так как мы заменяем ссылочный тип (reference type) массива на набор полей - значимых типов (value type). В свою очередь, это позволяет проделать больше шагов построения символического образа.
Важно отметить, что генерируемый объект для представления ячейки занимает меньше оперативной памяти.
6. Заключение. В данной работе были рассмотрены особенности реализации системы построения и анализа динамических систем с помощью символического образа. Был разработан способ генерации специализированного кода в реальном времени на языке С#. Сравнение скорости работы сгенерированного кода для случая фиксированной размерности и кода, написанного для произвольной размерности, показало, что сгенерированный код требует меньше памяти и работает быстрее, чем реализация для общего случая.
Было проведено сравнение производительности на задаче построения инвариантного множества для систем Био, Икеда, Линдстрома и Лоренца. На всех примерах описанная реализация стала работать быстрее по сравнению с реализацией, в которой генерация кода не применялась. Описанный подход позволил уменьшить размер одной из структур данных, что привело к уменьшению требований к оперативной памяти.
Литература
1. Jones N. D., Gomard C. K., Sestoft P. Partial Evaluation and Automatic Program Generation (Prentice-Hall International Series in Computer Science). Prentice Hall, 1993. 400 p.
2. Осипенко Г. С., Ампилова Н. Б. Введение в символический анализ динамических систем. СПб.: Изд-во С.-Петерб. ун-та, 2005. 237 с.
3. Dellnitz M., Hohmann A. A subdivision algorithm for the computation of unstable manifolds and global attractors // Numerische Mathematik. 1997. Vol. 75, N 3. P. 293—317.
4. Fundinger D. Investigating dynamics by symbolic analysis: Tunings for an efficient computation of the symbolic image // Electronic J. of Differentional Equations and Control Processes. 2005. Vol. 3. P. 16—37. URL: http://www.neva.ru/journal.
5. Fundinger D. Investigating Dynamics by Multilevel Phase Space Discretization: Ph. D. thesis. Institut fur Parallele und Verteilte Systeme, Abteilung Bildverstehen. 2006. URL: http://elib.uni-stuttgart.de/opus/volltexte/2006/2614.
6. Петренко Е. И. Разработка и реализация алгоритмов построения символического образа // Электр. журн. Дифференц. уравнения и процессы управления. 2006. № 3. URL: http://www. neva.ru/journal.
7. Gamma E., Helm R., Johnson R., Vlissides J. M. Design Patterns: Elements of Reusable Object-Oriented Software (Addison-Wesley Professional Computing Series). Reading MA: Addison-Wesley Professional, 1994. 416 p.
8. .NET Framework Developer’s Guide. Reflection // MSDN. URL: http://msdn.microsoft.com/en-us/library/cxz4wk15.aspx.
9. Осипенко Г. С., Романовский И. В., Ампилова Н. Б., Петренко Е. И. О вычислении спектра
морса // Проблемы математического анализа. 2004. T. 27. C. 151—169.
10. Петренко Е. И. Разработка и реализация программного комплекса исследования динамических систем методами символической динамики // Технологии Microsoft в теории и практике программирования: Материалы межвуз. конкурса-конференции студентов, аспирантов и молодых ученых Северо-Запада. СПб.: Изд-во С.-Петерб. политехн. ун-та. 2005. С. 186—187.
11. Романовский И. В. Дискретный анализ: учеб. пособие для студентов, специализирующихся по прикладной математике и информатике. 3-е изд. СПб.: Невский Диалект; БХВ-Петербург. 2004. 320 с.
12. .NET Framework Class Library. Array Class // MSDN. URL: http://msdn.microsoft.com/en-us/lib-rary/system.array.aspx.
13. Ikeda K. Multiple-valued stationary state and its instability of the transmitted light by
a ring cavity system // Optics Communications. 1979. Vol. 30, N 2. P. 257—261. URL: http://cdsads.
u-strasbg.fr/abs/1979OptCo..30..257I.
14. Lindstrom T. Dependencies between competition and predation — and their consequences for initial value sensitivity // SIAM Journal on Applied Mathematics. 1999. Vol. 59, N 4. P. 1468-1486. URL: http://dx.doi.org/10.1137/S0036139997320858.
15. Lorenz E. N. Deterministic nonperiodic flow // J. of the Atmospheric Sciences. 1963. Vol. 20, N 2. P. 130-141.
16. Osipenko G. Numerical explorations of the Ikeda mapping dynamics // Electronic J. of Differentional Equations and Control Processes. 2004. Vol. 2. P 43-67. URL: http://www.math.spbu.ru/diffjournal/j/.
Статья рекомендована к печати член-кор. РАН, проф. Г. А. Леоновым.
Статья принята к печати 5 марта 2009 г.