Как эффективно реализовать граничные условия Дирихле в глобальных разреженных матрицах жесткости конечных элементов

9

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

Кзнак равно[520-102410001632-1037000203]и правый векторбзнак равно[б1б2б3б4б5]

Затем применить условие Дирихле к первому узлу (Икс1знак равнос) мы обнуляем первый ряд, ставим 1 на К11и вычтите первый столбец с правой стороны. Например, наша система станет:

Кзнак равно[1000004100016320037000203]и правый векторбзнак равно[сб2-2×сб3-0×сб4+1×сб5-0×с]

Это все хорошо в теории, но если наша K-матрица хранится в формате сжатых строк (CRS), то перемещение столбцов в правую часть становится дорогим для больших систем (многие узлы являются дирихлетами). Альтернативой было бы не перемещать столбцы, соответствующие условию Дирихле, в правую часть, то есть наша система стала бы:

Кзнак равно[100002410001632-1037000203]и правый векторбзнак равно[сб2б3б4б5]

Это, однако, имеет существенный недостаток в том, что система больше не является симметричной, и поэтому мы больше не можем использовать предобусловленный сопряженный градиент (или другие симметричные решатели). Одним интересным решением, с которым я столкнулся, является «Метод больших чисел», который я нашел в книге Геннадий Никишкова «Программирование конечных элементов на Java». Этот метод использует тот факт, что двойная точность содержит только около 16 цифр точности. Вместо того, чтобы поставить 1 вК11Положение мы ставим большое количество. Например, наша система становится:

Кзнак равно[1,0е6420-102410001632-1037000203]и правый векторбзнак равно[с×1,0е64б2б3б4б5]

Преимущества этого метода в том, что он поддерживает симметрию матрицы, а также очень эффективен для разреженных форматов хранения. Мои вопросы тогда следующие:

Как граничные условия Дирихле обычно реализуются в кодах конечных элементов для тепла / жидкости? Люди обычно используют метод больших чисел или делают что-то еще? Есть ли недостаток в методе больших чисел, который кто-то может увидеть? Я предполагаю, что, вероятно, существует какой-то стандартный эффективный метод, используемый в большинстве коммерческих и некоммерческих кодов, который решает эту проблему (очевидно, я не ожидаю, что люди будут знать все внутренние механизмы каждого коммерческого решателя конечных элементов, но эта проблема кажется базовой / фундаментальной Достаточно того, что кто-то, вероятно, работал над такими проектами и мог обеспечить руководство).

Джеймс
источник
2
У вас есть доказательства того, что это существенно замедляет вас?
Билл Барт
@BillBarth Да, хотя всегда есть вероятность, что я делаю что-то неэффективно. Сам Геннадий пишет, что, хотя явный метод прост для полных 2d-массивов, «… не всегда легко получить доступ к строкам и столбцам матрицы, когда матрица находится в компактном формате». предполагая, что явный метод может быть более сложным для эффективной реализации. Поскольку мой код в настоящее время написан, явный метод может занять больше времени, чем фактическое решение.
Джеймс
1
сделайте это, как говорит Вольфганг, и примените граничные условия к элементным матрицам перед сборкой.
Билл Барт
@BillBarth Да, я думаю, что сделаю это. Его видео потрясающие! Я просто оставил ему комментарий / вопрос о том, нужно ли собирать глобальные матрицы на каждом временном шаге, после чего, я думаю, я приму его ответ.
Джеймс

Ответы:

11

В сделке II ( http://www.dealii.org - : я один из основных авторов этой библиотеки) мы исключаем целые строки и столбцы, и это не слишком дорого в целом. Хитрость заключается в том, чтобы использовать тот факт, что шаблон разреженности, как правило, симметричен, поэтому вы знаете, какие строки нужно просматривать при удалении целого столбца.

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

Я никогда не слышал о методе «больших чисел» и не буду его использовать, потому что, несомненно, он приведет к ужасно плохо обусловленным проблемам.

Для справки, алгоритмы, которые мы используем в deal.II, концептуально описаны в лекциях 21.6 и 21.65 на http://www.math.colostate.edu/~bangerth/videos.html . Они точно соответствуют вашему описанию.

Вольфганг Бангерт
источник
2
В случае нестационарной задачи (скажем, уравнения теплопроводности) вы собираете глобальную матрицу на каждом временном шаге? Причина, по которой я спрашиваю, состоит в том, что в случае ненулевых условий Дирихле при изменении правой части вам нужна информация из исходной глобальной матрицы, но если вы обнулили эти столбцы во время предыдущего временного шага, то эта информация будет потеряна (если вы не сохраните ее в дополнительных массивах). Это не было бы проблемой, если бы глобальная матрица перекомпоновывалась каждый раз, хотя это то, что я собираюсь сделать, и что должно быть сделано в любом случае, если использовать адаптивную сетку.
Джеймс
1
Это зависит от приложения. Все «большие» коды решают нелинейные задачи, зависящие от времени, и для них очевидно, что вам нужно пересобрать так или иначе. Для линейных кодов вы можете просто сохранить исходную матрицу и на каждом временном шаге копировать ее в другое место, применять граничные условия, а затем использовать ее в решателе. Это просто требует больше памяти, но в остальном дешево.
Вольфганг Бангерт
1
Ах, я вижу, это то, что я подозревал. Я буду реализовывать, как вы предложили. Хорошо, это за вашу помощь. Эти видео уроки deallii действительно хороши!
Джеймс
2

Ноль БК легки. Для ненулевых БК вы также можете использовать множители Лагранжа. Например, смотрите здесь . Одним из преимуществ LM является то, что вы можете использовать любое уравнение ограничения, хотя система становится неопределенной, поэтому вам нужен соответствующий решатель.

Stali
источник