CC BY
74
ИНФОРМАЦИОННЫЕ ТЕХНОЛОГИИ
УДК 519.852.2
ПАРАЛЛЕЛЬНЫЙ АЛГОРИТМ НАХОЖДЕНИЯ ОБЩЕГО РЕШЕНИЯ СИСТЕМЫ ЛИНЕЙНЫХ НЕРАВЕНСТВ*
© 2009 г. Н.Ю. Золотых, С.С. Лялин
Нижегородский госуниверситет им. Н.И. Лобачевского {Nikolai.Zolotykh,Sergey.Lyalin}@gmail.com
Поступила в редакцию 25.05.2009
Предлагается параллельная модификация метода двойного описания многогранного конуса. Приводятся результаты вычислительного эксперимента на многопроцессорной машине с общей памятью и сравнение с другими реализациями алгоритмов, решающих данную задачу. Эксперименты показали близкие к линейным результаты по масштабируемости.
Ключевые слова: полиэдр, многогранный конус, параллельный алгоритм, метод двойного описания.
Введение
Во многих прикладных задачах возникает такой объект, как полиэдр, или многогранное множество [1], - множество точек многомерного пространства, удовлетворяющих системе линейных неравенств. Например, в биологии статический анализ метаболических сетей приводит к решению системы линейных уравнений и неравенств [2]. В программной оптимизации полиэдры возникают при рассмотрении линейного класса программ и его расширений [3].
Работа с полиэдрами, как правило, заключается в выполнении следующих операций: пересечение полиэдров, выпуклое объединение полиэдров, проекция на линейные многообразия, обход всех граней некоторой заданной размерности и т.д. В качестве одной из базовых задач, которую приходится решать при осуществлении перечисленных операций, выступает задача перехода между двумя представлениями полиэдра.
В данной статье даётся схематичное описание алгоритма, реализованного в программе Skeleton 3 - параллельной модификации алгоритма двойного описания многогранного конуса, который решает задачу перехода между двумя представлениями полиэдра. Сообщаются последние результаты по вычислительной эф-
* Статья рекомендована к печати программным комитетом Международной научной конференции «Параллельные вычислительные технологии 2009» (http:// agora.guru.ru/pavt).
фективности в сравнении с программами других авторов, решающими ту же задачу, а также масштабируемость Skeleton 3 при вычислении с вещественной арифметикой и арифметикой многократной точности.
Необходимые определения
Для краткости изложения будем рассматривать постановку задачи в однородной форме. Пусть дана матрица A є Rmxd . Многогранным выпуклым конусом C в пространстве Rd называют множество решений системы линейных неравенств Ax > 0 :
C = {x є Rd | Ax > 0}. (1)
Кроме представления (1) конус C может
быть записан в параметрическом виде, а именно для некоторой матрицы E є Rdxn:
C = {x є Rd | x = ВЛ, Лє R + }, (2)
где R+ - неотрицательные действительные
числа. Алгоритм двойного описания [4], известный также как алгоритм Моцкина-Бургера [1], по одному из описаний, (1) или (2), строит другое. При этом метод находит неприводимую систему векторов в матрице В (неприводимую систему неравенств в матрице A соответственно). Везде далее мы будем считать, что рассматривается задача перехода от (1) к (2). Для простоты изложения будем считать, что ранг матрицы A равен d. В этом случае неприводи-
Рис. 1. Зависимости между стадиями алгоритма
мая система векторов в описании (2) называется остовом конуса. Векторы остова (лучи) определяются единственным образом с точностью до умножения на положительные константы. Так как порядок строк в А и столбцов в В , очевидно, не имеет значения, то далее для упрощения изложения символы А и В используются также для обозначения множества строк и столбцов соответствующих матриц.
Существующие алгоритмы и реализации
Авторам известны несколько других алгоритмов (и соответствующих реализаций), решающих поставленную задачу. Родственным алгоритмом, который построен на тех же идеях, что и описываемый алгоритм, является еще одна модификация [5] метода двойного описания, реализация которой, cdd и cdd+, доступна в Интернете [6]. Хорошо известны ещё два алгоритма: 1ге [7] и qhu11 [8]. Все перечисленные реализации доступ-
ны в последовательной форме. В [9] описана параллельная реализация алгоритма [4].
Алгоритм
Рассматриваемый алгоритм является параллельной модификацией метода двойного описания [4], построенной на основе модификации [10] и реализации [11]. Предварительные результаты по улучшению версии [10] анонсировались ранее в [12-14].
Опишем схему алгоритма. Оригинальный метод двойного описания имеет итеративный характер: на каждой итерации рассматривается очередное неравенство из матрицы А и проводятся все необходимые вычисления, после чего алгоритм переходит к следующему неравенству. При переходе с итерации на итерацию перестраиваются две структуры данных: текущий остов конуса В (множество лучей) и текущий подграф графа смежности лучей Е (множество пар лучей). Перечислим
Поток k1 Поток k2
Итерация i Шаг 1 1 Шаг 2 1
Шаг 7 1
Итерация i + 1 Шаг 1 1 Шаг 2 1 г
Шаг 7 1
Шаг 1
Шаг 2
Итерация i
Шаг 7
Шаг 1 Итерация
Шаг 2
i + 1
Шаг 7
Условные обозначения:
Блокировка
Рис. 2. Пример реализации параллельного потока управления при синхронизации в конце каждой стадии. Указаны номера шагов алгоритма 1
основные стадии, из которых состоит одна итерация алгоритма:
Алгоритм 1. Последовательный алгоритм построения остова многогранного конуса.
Шаг 0. [Начальный симплициальный конус.] Выбрать произвольно невырожденную квадратную подматрицу А0 размера й матрицы А. Положить В = А- и Е равным всевозможным упорядоченным парам из В.
Шаг 1. [Выбор неравенства.] Выбирается очередное неравенство из А . Хотя неравенства могут выбираться не последовательно, в данной работе рассматривается только последовательный выбор (в том порядке, в котором идут строки в А ).
Шаг 2. [Классификация.] Существующие лучи конуса В разделяются на три группы: (а) удовлетворяющие выбранному неравенству строго, (б) удовлетворяющие неравенству как равенству, (в) не удовлетворяющие неравенству.
Шаг 3. [Получение новых лучей.] Перебором всех элементов из Е комбинируются сущест-
вующие лучи, входящие в пары, и определённым образом порождаются новые лучи. Они включаются в В и добавляются в группу (б), созданную на стадии 2.
Шаг 4. [06pa3oeaoue nap.] Подготавливаются пары из лучей группы (б), которые будут проверены на потенциальную смежность в следующей стадии.
Шаг 5. [Пomepцuaлъpaя cмeжpocmъ.] Лучи, попавшие в группу (б), попарно проверяются на потенциальную смежность в конусе.
Шаг 6. [Смeжpocmъ.] Для каждого луча группы (б) из числа потенциально смежных лучей определённым образом выбираются смежные лучи; все такие пары лучей добавляются в E .
Шаг 7. \ydaMeoue pepужpнlх лучєй.] Все лучи из группы (в) удаляются.
При завершении алгоритма, когда все неравенства A рассмотрены, В содержит остов конуса, т.е. описание (2).
На рис. 1 дана схема зависимостей между стадиями одной итерации алгоритма. Нисходящие зависимости (т.е. идущие от стадий с меньшими номерами к стадиям с большими номерами) указывают на зависимости стадий
Поток k1 Поток к2
Итерация i Шаг 1 Шаг 2 Шаг 7
Итерация i + 1 Шаг 1 Шаг 2 Шаг 7
Шаг 1 Итерация i
Шаг 2
Шаг 7
Шаг 1 Итерация i + 1
Шаг 2
Шаг 7
Шаг 6
Условные обозначения:
Ожидание данных
Рис. 3. Пример реализации параллельного потока управления алгоритма 2. Указаны номера шагов алгоритма 1
Число потоков
Число потоков
Рис. 4. Масштабируемость Skeleton 3 на вещественной арифметике
внутри одной итерации. Восходящие зависимости (т.е. идущие от стадий с большими номерами к стадиям с меньшими номерами) указывают на зависимость следующей итерации от текущей итерации. Зависимость по данным, которая обозначена сплошной стрелкой на рис. 1 , можно рассматривать как очередь из данных от одной стадии к другой.
Рис. 5. Масштабируемость Skeleton 3 на целочисленной арифметике
Анализ показал, что при наличии такого рода зависимостей можно так реализовать алгоритм, чтобы избежать обязательного упорядочивания работы по основным стадиям. Другим словами, часть работы стадии с большим номером может быть выполнена до окончания работы стадии с меньшим номером, если нет дополнительных зависимостей.
Таблица
Время работы Skeleton 2 и Skeleton 3 на нескольких наборах данных
Набор данных Параметры задачи Время работы Skeleton 2 (вещественная арифметика), сек Время работы Skeleton 3 на одном потоке (вещественная арифметика), сек Время работы Skeleton 3 на одном потоке (целочисленная арифметика GMP), сек
d m n
cube16.ine 17 32 65536 5.959 4.548 6.929
mit729-9.ine 9 729 4862 54.147 15.368 731.327
prodmT5.ine 20 711 76 845.146 269.169 863.840
prodst62.ext 25 3461 168 1195.584 320.192 2352.040
ccp7.ext 22 64 116764 10257.774 2269.150 2539.730
Параллельный алгоритм и его реализация
Работа каждой из стадий 1-7 в отдельности может быть организована параллельно. При этом если не пользоваться свойством, указанным выше, в конце каждой стадии требуется точка синхронизации всех исполнительных устройств (барьер). Это ухудшает балансировку загрузки между вычислительными ресурсами: время работы одной стадии трудно предсказать и трудно (если вообще возможно) разделить данные между вычислительными устройствами так, чтобы добиться идеальной балансировки. Неизбежна ситуация когда в конце работы стадии многие вычислительные устройства будут простаивать (рис. 2).
Использование свойства отсутствия обязательного упорядочивания между стадиями позволяет построить более эффективный параллельный алгоритм, так как устраняет обязательные барьерные синхронизации в конце каждой стадии (рис. 3).
Skeleton 3 [15] - это реализация описанного алгоритма с использованием указанного свойства. Программа работает на многопроцессорной машине с общей памятью и использует дешёвый способ синхронизации через атомарно модифицируемые целочисленные переменные (interlocked-операции). Каждый из буферов (овалов на рис. 1) реализован в виде односвязного списка с неблокируемым доступом. Все потоки равнозначны и синхронизируются через небольшой набор глобальных переменных. Синхронизации типа «барьер» отсутствуют. Такое понятие, как «текущая итерация» (итерация алгоритма, обрабатываемая в данный момент времени), для всех потоков исполнения теряет свой смысл, так как в текущий момент времени разные потоки могут обрабатывать данные с разных итераций. Понятие «текущая итерация для алгоритма» заменяется «текущей итерацией для потока». Далее для краткости мы
опускаем слово «текущий» в словосочетаниях «текущая итерация» и «текущая стадия», имея в виду ту итерацию и ту стадию, которые исполняет данный поток.
Итак, в описываемой реализации все потоки равнозначны и каждый исполняет представленный алгоритм следующим образом.
Алгоритм 2. Исполнение алгоритма 1 каждым потоком параллельной реализации.
Шаг 1. [Выбор стадии.] Опираясь на метки-номера итераций для буферов данных, выбрать стадию алгоритма 1, которую можно исполнить и для которой сейчас есть данные. Если выбор сделан, то перейти на шаг 2, если готовых к исполнению данных нет, то перейти к шагу 4.
Шаг 2. [Обработка.] Обработать одну порцию данных на выбранной стадии, после чего проверить условия завершения стадии для итерации.
Шаг 3. [Завершение стадии.] Если условия завершения стадии на текущей итерации, проверенные на шаге 2, удовлетворяются, то поток завершает стадию локально (в контексте этого потока), и если он последний завершает данную стадию (т.е. все остальные потоки уже завершили стадию), то стадия помечается завершённой для итерации глобально (для всех потоков). Если это не завершающая стадия на последней итерации, то переход на шаг 1, иначе - завершение работы.
Шаг 4. [Подглядывание.] Проверить, есть ли необрабатываемые порции данных у других потоков в каких-либо стадиях, которые сейчас можно обработать. Если есть, то перейти на шаг
5. Если нет - на шаг 6.
Шаг 5. [Заимствование.] Переместить порцию данных из другого потока в соответствующий буфер текущего потока и перейти к шагу 2.
Шаг 6. [Ожидание.] Сообщить операционной системе об отказе от текущего выделенного кванта времени (выполнить yield) или средст-
вами операционной системы пассивно подождать период времени, равный среднему времени обработки порции данных (требует сбор статистики по времени исполнения). Затем перейти на шаг 1.
Как уже было описано в [15], эффективность реализации по сравнению с первыми результатами Skeleton 3 [12] улучшилась. Более эффективное использование кэша процессора, некоторые низкоуровневые оптимизации и использование инструкций Intel® SSE и Intel® SSE 2 позволили ускорить работу Skeleton 3 по сравнению со своим предшественником Skeleton 2 [10, 11] в несколько раз при исполнении на одном потоке (см. табл.). При этом нельзя сказать об однозначно хорошей масштабируемости программы, так как есть такие наборы данных, на которых Skeleton 3 показывает плохие результаты по этому параметру (например, cube16.ine; см. рис. 4 и рис. 5).
Для сравнения производительности и оценки масштабируемости использовались несколько наборов данных, которые позаимствованы из дистрибутива cdd [б] и фигурируют в различных сравнительных исследованиях программ построения остова многогранного конуса (см. табл.). В экспериментах использовались две различные реализации арифметики: 1) вещественная арифметика двойной точности, выраженная типом double в С++, и 2) целочисленная «длинная» арифметика GMP [16], подключённая к библиотеке Arageli [17].
Вычисления производились на машине с двумя четырёхъядерными процессорами Intel® Xeon® X5355 2.66 ГГц, 8 Гбайт оперативной памяти. Реализация использует потоки Win32 API операционной системы Microsoft Windows и библиотеку Arageli [7]. На рис. 4 и рис. 5 показаны графики масштабируемости Skeleton 3 при запуске на различном числе ядер для вещественной (рис. 4) и целочисленной (рис. 5) арифметики.
При сравнении с программой cdd [9] было получено, что на одних наборах данных быстрее работает Skeleton 3 (например, cube16.ine и mit729-9.ine), а на других - cdd (например, prodst62.ext). Это связано с одним нюансом в работе рассматриваемого алгоритма и cdd, который связан с порядком проверки критериев отбора подходящих пар лучей на стадиях 4 и 5 основного алгоритма 1. В настоящий момент авторам известна модификация, которая улучшит работу Skeleton 3 на наборах, подобных prodst62.ext, что позволит приблизиться или обогнать cdd на таких примерах.
Skeleton 3 доступен on-line через сервис, разработанный С.В. Лобановым [18].
Реализация алгоритма для распределенной памяти
Продолжается работа над параллельной версией алгоритма для машин с распределённой памятью. Используется та же схема, что приведена на рис. 1, но появляются дополнительные стадии по асинхронной пересылке данных с узла на узел. Проводятся эксперименты по спекулятивной (упреждающей) посылке данных, чтобы минимизировать потерю производительности из-за большого времени ожидания коммуникационной сети.
С этой же целью проверяется возможность спекулятивных вычислений, когда один узел, обрабатывая какую-то порцию данных, ещё не знает, что какой-то другой узел начал обрабатывать ту же порцию данных. Эта ситуация возникает как следствие спекулятивной посылки данных, а также в случае применения разных эвристик на шагах 1, 4, 5 и 6: при избытке вычислительных ресурсов на недостаточно больших задачах можно использовать разные группы вычислительных узлов для различных эвристик. Для продолжения вычислений берется результат, вычисленный быстрее всех одной из групп, а на всех остальных группах вычисления отменяются. Данный подход позволяет значительно увеличить скорость работы для задач, для которых неизвестно поведение Skeleton, следовательно, неизвестно, какая из эвристик даст результат за минимальное время.
Список литературы
1. Черников С.Н. Линейные неравенства. М.: Наука, 1968.
2. Stelling J., Klamt S., Bettenbrock K. et al. Metabolic network structure determines key aspects of functionality and regulation // Nature. 2002. V. 420. P. 190193.
3. Cousot P., Halbwachs N. Automatic discovery of linear restraints among variables of a program // Conference Record of the Fifth Annual ACM SIGPLAN-SIGACT Symposium on Principles of Programming Languages, Tucson, Arizona. New York: ACM Press, 1978. P. 84-97.
4. Motzkin T.S., Raiffa H., Thompson G.L., Thrall R.M. The double description method // Contributions to the Theory of Games. Princeton: Princeton University Press, 1953.
5. Fukuda K., Prodon A. Double description method revisited // Lecture Notes in Computer Science. V. 1120. Springer-Verlag, 1996. P. 91-111.
6. Cdd Home Page: http://www.ifor.math.ethz.ch/ ~fukuda/cdd home/cdd.html.
7. Avis D. lrs: A Revised Implementation of the Reverse Search Vertex Enumeration Algorithm // Polytopes - Combinatorics and Computation. DMV Seminar Band 29, Birkhauser-Verlag, 2000. P. 177-198.
8. Barber C.B., Dobkin D.P., Huhdanpaa H.T. The Quickhull algorithm for convex hulls // ACM Trans. on Mathematical Software. 1996. V. 22. N. 4. P. 469-483.
9. Агафонов Е.А., Земскова Е.Л., Золотых Н.Ю. Параллельный алгоритм построения остова многогранного конуса и его программная реализация // Информационный бюллетень Ассоциации математического программирования № 10. Конференция «Математическое программирование и приложения» (тезисы докладов). Екатеринбург: УрО РАН, 2003. С. 15-16.
10. Золотых Н.Ю. Новая модификация метода двойного описания для построения остова многогранного конуса // III Всероссийская конференция «Проблемы оптимизации и экономические приложения». Тезисы докладов. Омск, 2006. С. 108.
11. Skeleton 2 Home Page: http://www.uic.
nnov.ru/~zny/skeleton.
12. Золотых Н.Ю., Лялин С.С. Оптимизация одной модификации алгоритма двойного описания для
многогранного конуса // Технологии Microsoft в теории и практике программирования. Материалы конференции. Нижний Новгород: Изд-во Нижегородского госуниверситета, 2007. С. 333-336.
13. Лялин С.С. Параллельная модификация алгоритма двойного описания для многогранного конуса // Высокопроизводительные параллельные вычисления на кластерных системах. Материалы Седьмой Международной конференции-семинара. Нижний Новгород: Изд-во Нижегородского госуниверситета, 2007. С. 240-243.
14. Лялин С.С. Skeleton 3 - реализация параллельного алгоритма построения остова многогранного конуса // Высокопроизводительные параллельные вычисления на кластерных системах. Материалы Восьмой Международной конференции-семинара. Казань: Изд-во КГТУ, 2008. С. 167-171.
15. Skeleton 3 Home Page: http://www.
arageli.org/skeleton.
16. GNU Multiple Precision Arithmetic Library Home Page: http://gmplib.org.
17. Arageli Home Page: http://www.arageli.org.
18. Skeleton 3 on-line demonstration: http://www. arageli.org/skeletondemo.
A PARALLEL ALGORITHM FOR FINDING THE GENERAL SOLUTION OF THE SYSTEM OF LINEAR INEQUALITIES
N. Yu. Zolotykh, S. S. Lyalin
A parallel modification of the double description method for a polyhedral cone is proposed. The results of the computer experiment on a shared-memory multiprocessor machine are given as well as their comparison with other algorithm realizations of this problem. The experiments have shown nearly linear results in terms of scalability.
Keywords: polyhedron, polyhedral cone, parallel algorithm, double description method.
