Main menu

Вейвлет-Галеркинские методы: объединение МКЭ и спектрального анализа

Ограничения базисных функций в вычислительной математике

При численном решении дифференциальных уравнений в частных производных (УЧП) методом Галеркина выбор базисных функций определяет успех всего алгоритма. В классическом методе конечных элементов (МКЭ) используются локальные кусочные полиномы. Они отлично подстраиваются под геометрию, но дают медленную алгебраическую сходимость. В спектральных методах используются глобальные синусоиды (ряды Фурье) или ортогональные полиномы. Они дают феноменальную экспоненциальную сходимость для гладких функций, но полностью разваливаются (порождают глобальные осцилляции Гиббса), если в решении появляется локальный разрыв, ударная волна или острый пик. Глобальный базис не умеет «фокусироваться» на локальных особенностях.

Инженерам был необходим базис, который обладал бы лучшими свойствами обоих миров: строгой ортогональностью и спектральной точностью глобальных функций, а также локальностью и компактным носителем конечных элементов. В начале 1990-х годов, на волне бурного развития теории обработки сигналов, вычислительная математика нашла идеальный инструмент — вейвлеты (маленькие волны). Использование вейвлетов в качестве пробных функций породило мощный класс алгоритмов: вейвлет-Галеркинские методы.

Многомасштабный анализ и адаптивные сетки

Фундаментальным преимуществом вейвлетов (например, вейвлетов Добеши) является их способность к многомасштабному анализу (Multiresolution Analysis, MRA). Пространство решений разбивается на иерархию подпространств с разным разрешением. Решение строится как сумма грубого (низкочастотного) приближения с помощью функций масштабирования и тонких (высокочастотных) корректировок с помощью самих вейвлетов.

Это дает алгоритму потрясающее свойство: встроенную, естественную адаптивность сетки. Алгоритм вычисляет вейвлет-коэффициенты по всей области. В тех местах, где функция гладкая (например, равномерный поток воздуха), высокочастотные вейвлет-коэффициенты оказываются практически равными нулю, и алгоритм просто удаляет их из расчета (сжимает оператор), экономя гигантские объемы памяти. Но в области ударной волны или пограничного слоя высокочастотные коэффициенты велики, и алгоритм автоматически сгущает расчет только в этой узкой зоне, не затрагивая глобальную структуру данных. При этом вейвлеты остаются строго ортогональными, что обеспечивает идеальную обусловленность матриц (в отличие от адаптивного МКЭ).

Проблема коэффициентов связи (Connection Coefficients)

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

Гладкие ортогональные вейвлеты с компактным носителем (вейвлеты Добеши) не имеют аналитических формул — они задаются рекурсивными фильтрами. Соответственно, вычислить аналитически интеграл от их производной невозможно. Для решения этой проблемы была разработана теория вычисления «коэффициентов связи» (Connection coefficients) — точных значений этих интегралов через решение специальных систем алгебраических уравнений, вытекающих из свойств масштабирования. Сегодня вейвлет-методы активно применяются для моделирования турбулентности (DNS и LES), где необходимость отслеживать каскад энергии вихрей разных масштабов идеально ложится на многомасштабную физику самих вейвлетов.

Подробнее

Предобуславливание (Preconditioning): ключ к решению гигантских СЛАУ

Узкое место итерационных методов Крылова

Ранее мы обсуждали, что гигантские разреженные системы линейных алгебраических уравнений (СЛАУ), возникающие из метода конечных элементов или разностей, невозможно решить прямым методом Гаусса из-за проблемы заполнения нулей. Индустрия использует итерационные методы подпространств Крылова, такие как Метод сопряженных градиентов (CG) для симметричных матриц и GMRES для несимметричных.

Однако эффективность этих методов напрямую зависит от одного фундаментального математического свойства матрицы — числа обусловленности (Condition Number). Число обусловленности характеризует «растянутость» спектра матрицы (отношение максимального собственного значения к минимальному). В реальных физических задачах (особенно при сильных перепадах свойств материалов, например, металл и воздух, или при сильном измельчении сетки) число обусловленности стремится к бесконечности. Матрица становится «плохо обусловленной». В такой ситуации метод сопряженных градиентов начинает топтаться на месте: невязка колеблется, и алгоритму требуются сотни тысяч итераций для сходимости, что делает расчет недопустимо долгим.

Математическая трансформация системы

Для преодоления этого вычислительного тупика применяется предобуславливание (Preconditioning). Идея заключается в преобразовании исходной СЛАУ (Ax = b) в другую, математически эквивалентную систему, которая имеет точно такое же решение, но гораздо лучшее число обусловленности. Для этого вводится матрица предобуславливателя (обозначим ее M).

Мы умножаем левую и правую части системы на обратную матрицу предобуславливателя: M^(-1)Ax = M^(-1)b. Теперь алгоритм Крылова решает новую систему с матрицей M^(-1)A. Идеальным предобуславливателем была бы сама матрица A (тогда M^(-1)A превратилась бы в единичную матрицу, и алгоритм сошелся бы за 1 итерацию). Но обращение A — это именно то, чего мы избегаем. Следовательно, матрица M должна удовлетворять двум противоречивым требованиям: она должна быть «максимально похожа» на матрицу A (чтобы сгруппировать собственные значения), и при этом решение системы My = z должно вычисляться алгоритмически очень быстро и дешево (так как это придется делать на каждой итерации).

Популярные алгоритмы: ILU и диагональное масштабирование

Самым простым является диагональный предобуславливатель (метод Якоби), где M — это просто матрица, состоящая только из диагональных элементов исходной матрицы A. Обратить такую матрицу тривиально (просто поделить на диагональные числа), но ее эффективность для сложных задач крайне мала.

Золотым стандартом в индустрии стало Неполное LU-разложение (Incomplete LU factorization, ILU). Алгоритм пытается сделать стандартное разложение матрицы A на нижнюю и верхнюю треугольные матрицы, но принудительно отбрасывает (заменяет нулями) те новые элементы (fill-ins), которые возникают в процессе исключения. В результате мы получаем приближенные матрицы L и U, которые сохраняют идеальную разреженность исходной системы. Решение системы с такими треугольными матрицами выполняется мгновенно прямой и обратной подстановкой. Баланс между точностью ILU (задаваемой порогом отбрасывания) и скоростью одной итерации является главным искусством настройки современных коммерческих решателей. Для самых тяжелых задач применяются многосеточные предобуславливатели (AMG), которые обеспечивают теоретическую линейную масштабируемость при решении миллиардных систем.

Подробнее

Метод сглаженных частиц (SPH): бессеточная вычислительная гидродинамика

Проблема подвижных границ и сильных деформаций в сеточных методах

Традиционная вычислительная гидродинамика (CFD) на основе метода конечных объемов или конечных элементов опирается на эйлеров подход: мы строим неподвижную сетку в пространстве и смотрим, как жидкость протекает сквозь ее ячейки. Этот подход идеален для расчета аэродинамики самолета или течения в трубе. Однако он сталкивается с колоссальными трудностями при моделировании процессов с сильной деформацией свободной поверхности: например, при ударе волны о волнолом, падении метеорита в океан, разбрызгивании металла при сварке или фрагментации планет в астрофизике.

Попытка отслеживать свободную поверхность (поверхность раздела газ-жидкость) на неподвижной сетке требует использования сложных алгоритмов типа Volume of Fluid (VOF) или Level Set, которые размазывают границу и подвержены численной диффузии. Использование же деформируемых (лагранжевых) сеток приводит к их быстрому скручиванию и вырождению ячеек, требуя постоянного, вычислительно дорогого перестроения (ремешинга). Инженеры и физики остро нуждались в методе, который вообще не зависит от геометрической сетки.

Лагранжев подход: жидкость как ансамбль частиц

В 1977 году Люси, а также Гинголд и Монаган, независимо друг от друга предложили метод сглаженных частиц (Smoothed Particle Hydrodynamics, SPH) для решения задач астрофизики. SPH — это чисто лагранжев бессеточный метод. Среда (жидкость, газ или деформируемое твердое тело) представляется в виде набора дискретных макроскопических частиц. Каждая частица обладает фиксированной массой и несет в себе информацию о плотности, давлении, внутренней энергии и скорости.

Эти частицы не привязаны к какой-либо сетке — они свободно движутся в пространстве в соответствии с уравнениями Навье-Стокса (законами сохранения импульса и массы). Поскольку масса жестко привязана к частице, закон сохранения массы выполняется автоматически и абсолютно точно. Свободная граница жидкости определяется естественным образом — это просто внешняя граница облака частиц, что делает моделирование брызг и капель тривиальной задачей.

Ядро сглаживания и аппроксимация производных

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

Весовые коэффициенты определяются специальной колоколообразной функцией — ядром сглаживания (Smoothing kernel). Ядро имеет радиус действия (smoothing length). Частицы, находящиеся за пределами этого радиуса, не влияют на расчет, что позволяет использовать алгоритмы быстрого поиска соседей (например, KD-деревья). Самое важное свойство SPH заключается в том, что пространственная производная от физической величины аналитически переносится на саму функцию ядра! То есть, чтобы найти градиент давления для уравнения импульса, компьютеру не нужно брать разности между частицами; он просто берет аналитическую производную от известной функции ядра. Благодаря своей гибкости метод SPH сегодня является стандартом не только в моделировании баллистики и цунами, но и в компьютерной графике для создания реалистичной воды в голливудских фильмах.

Подробнее

ГДЗ по математике

1 класс 2 класс 3 класс
4 класс 5 класс 6 класс

Книги

Фильмы

Старый сайт

Соц. сети