Ускорение граничных элементов: Быстрый мультипольный метод (FMBEM)

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

Метод граничных элементов (МГЭ, или BEM) обладает колоссальным преимуществом перед объемными методами: он снижает размерность физической задачи на единицу. Для расчета распределения электромагнитного поля вокруг огромного самолета инженеру нужно покрыть сеткой только поверхность обшивки фюзеляжа, а не кубические километры окружающего пространства. Однако за это математическое изящество приходится платить суровую цену. Переход к поверхностным интегральным уравнениям порождает систему линейных алгебраических уравнений (СЛАУ), матрица которой является плотной (dense). Она заполнена ненулевыми элементами на 100%, потому что каждая точка (заряд) на поверхности самолета непрерывно взаимодействует со всеми остальными точками посредством закона Кулона.

Для сетки всего из 100 000 элементов матрица потребует около 80 гигабайт оперативной памяти для хранения и колоссального количества процессорного времени (порядка O(N^2) операций для умножения на вектор и O(N^3) для прямого решения методом Гаусса). Это означает, что классический МГЭ упирается в «стеклянный потолок» размерности: попытка решить задачу для детализированной сетки современного подводного корабля или антенной решетки из миллионов элементов становится математически невозможной даже для самых мощных суперкомпьютеров министерства обороны.

Интеграция с Быстрым мультипольным методом (FMM)

Спасением для метода граничных элементов стала его интеграция с революционным Быстрым преобразованием мультиполей (Fast Multipole Method, FMM), разработанным Владимиром Рохлиным. Этот синтез получил название Быстрого мультипольного метода граничных элементов (Fast Multipole Boundary Element Method, FMBEM). Идея алгоритма основана на строгом математическом факте: влияние группы далеких источников (точек на поверхности) на приемник можно сгруппировать и заменить одним эквивалентным источником, используя разложение фундаментального решения (функции Грина) в ряды по сферическим гармоникам.

С вычислительной точки зрения алгоритм FMBEM полностью отказывается от сборки и хранения в памяти плотной глобальной матрицы. Алгоритм строит иерархическое дерево кубических ячеек (Octree), заключающее в себе всю поверхность модели. Если два элемента сетки находятся в соседних, близких ячейках, их взаимное влияние вычисляется точно, прямым интегрированием функций формы, и эти матричные коэффициенты образуют сильно разреженную локальную матрицу (ближнее поле). Если же элементы находятся в далеких ячейках (дальнее поле), матричные элементы не вычисляются вообще!

Мультипольные разложения и ускорение до O(N log N)

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

Использование FMM в качестве способа быстрого умножения плотной матрицы на вектор позволяет применять итерационные методы подпространств Крылова (например, GMRES). В результате вычислительная сложность и затраты оперативной памяти радикально падают с O(N^2) до O(N log N) или даже O(N). Это технологический прорыв колоссального масштаба: сегодня коммерческие FMBEM солверы позволяют моделировать рассеяние акустических волн на турбинах подводных лодок, используя сетки из 10-50 миллионов узлов, решая задачи, которые еще 20 лет назад считались абсолютно нерешаемыми в вычислительной акустике.