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

15

Я пытаюсь решить двумерное уравнение Пуассона с помощью конечных разностей. В процессе, я получаю разреженную матрицу только с переменными в каждом уравнении. Например, если переменные были U , то дискретизация даст:5U

Ui1,j+Ui+1,j4Ui,j+Ui,j1+Ui,j+1=fi,J

Я знаю, что могу решить эту систему итеративным методом, но мне пришла в голову мысль, что если бы я упорядочил переменные надлежащим образом, я мог бы получить матрицу с полосами, которая могла бы быть решена прямым методом (т. Е. Исключение Гаусса w / o поворотный). Это возможно? Существуют ли стратегии для этого для других, возможно, менее структурированных разреженных систем?

Пол
источник
2
Что-то вроде Кутхилл-Макки?
JM
Интересно ... я никогда не слышал об алгоритме Cuthill-McKee раньше! :)
Пол
1
Также есть обратный Cuthill-McKee.
Джефф Оксберри
1
Я надеюсь, что это ясно из ответов, но вы не хотите использовать решатель с полосами для этой проблемы или выбирать порядок, который минимизирует пропускную способность. Возможно, вопрос или выбранный ответ можно отредактировать, чтобы прояснить это, иначе я боюсь, что этот миф будет увековечен. Я дал визуальное сравнение и сравнил заполнение scicomp.stackexchange.com/a/880/119 .
Джед Браун
@JedBrown: На самом деле, я не совсем работаю с проблемой Пуассона как таковой ... Моя проблема имеет структуру, аналогичную проблеме Пуассона ... Признаки переменных (i и j) точно такие же, и матрица доминирует по диагонали с недиагональными элементами (в пределах одной строки), добавляемыми к сумме диагональных элементов.
Павел

Ответы:

13

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

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

Джордж и Лю показали в своей классической книге , что 2d разреженных прямых решения требуют только работа и вывод ( п войти п ) память, в то время как 3d разреженных прямым требует вывод ( п 2 ) работа и О ( п 4О(N3/2)О(NжурналN)O(n2)память.O(n4/3)

Джек Полсон
источник
У вас есть ссылка на ссылку Джорджа и Лю?
Павел
Добавлен; Я собирался выйти из машины, когда впервые ее представил. Я знаю, что где-то существует свободно доступная версия книги онлайн (Джед знает, где она), но я не смог ее найти.
Джек Поулсон
Я обновил ссылку, чтобы она указала на PDF книги, а не на рецензию.
Джед Браун
@JedBrown Это была отличная ссылка! Спасибо! :)
Пол
1
@ Александр Все приписывают трехмерность, привязанную к Джорджу и Лю, хотя я не знаю, указывают ли они это явно в книге. Это очевидно из теории, однако. Минимальная вершина сепаратор для сетки п 2 / 3 = м × м . Плотная матрица , связанная с этой supernode имеет ( п 2 / 3 ) 2 = п 4 / 3 записи и требует ( п 2 / 3 ) 3Nзнак равном×м×мN2/3знак равном×м(N2/3)2знак равноN4/3(N2/3)3знак равноN2операции к фактору. Логарифмический термин в 2D-случае более тонкий и рассматривается в главе 8 «Вложенное рассечение», в которой достигается нижняя граница.
Джед Браун
5

Cuthill-McKee является стандартом де-факто для того, что вы хотите сделать. Если вы хотите поиграть с этим методом, есть простая в использовании реализация алгоритма (и его обратного) в библиотеке графов ускорения. (BGL), и документация содержит примеры, как его использовать.

Вольфганг Бангерт
источник
на самом деле обратный Cuhill-McKee; это обычно дает меньше заполнения. Но вложенный порядок разбиения намного превосходит порядок с низкой пропускной способностью.
Арнольд Ноймайер,
4

Говоря о мультифронтальных методах, Тим Дэвис , который работает над мультифронтальными методами для факторизации LU ( UMFPACK ), предлагает ряд процедур, которые будут переупорядочивать матрицы для минимизации заполнения. Вы можете найти их здесь как часть SuiteSparse . SuiteSparse использует MeTiS.

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

Джефф Оксберри
источник
Не за что, Пол. Если вам это нравится, проголосуйте.
Джефф Оксберри
3

В прикладных математических кругах есть алгоритм ADI (Alternating Direction Implicit) и оператор Split в кругах физики, который в основном выполняет то, что вы описываете. Это итеративный метод, и он следует этой основной процедуре:

  1. YИкс

  2. ИксY

  3. Повторите 1 и 2, пока ошибка не станет настолько маленькой, насколько вы хотите.

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

Дэн
источник
2
Если вы решите пойти по пути разделения операторов, вам следует быть осторожным в том, что разделение операторов, как известно, в некоторых случаях приводит к ошибкам в стационарных решениях. (Один из моих коллег разработал способ преодолеть эту трудность, но я не верю, что он его опубликовал.) Также известно, что разбиение операторов вызывает числовые ошибки. Существуют устоявшиеся способы апостериорной оценки этих ошибок ; Дон Эстеп проделал отличную работу в этой области.
Джефф Оксберри
@ GeoffOxberry Похоже, вы имеете в виду другое разделение. Вы можете использовать ADI в полностью неявной схеме, в которой нет ошибки расщепления, потому что это фактически решает систему. Существуют также методы IMEX, которые строго контролируют ошибки расщепления.
Джед Браун
ИксY
Я никогда не слышал о расщеплении Годунова и Странга. Я склонен делить своего оператора по формуле Бейкера-Кэмпбелла-Хаусдорфа. Это то же самое?
Дан