Численное моделирование фазовых переходов: задача Стефана и метод энтальпии

Математика таяния льда и кристаллизации металлов

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

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

Методы явного отслеживания фронта (Front Tracking)

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

Этот подход (например, с использованием ALE сеток) дает высочайшую точность определения формы кристалла и распределения температур. Однако он терпит сокрушительное фиаско, если топология фронта начинает меняться. Если в процессе литья в форме образуется несколько независимых кусков льда, которые затем сливаются (или наоборот, одна жидкая ванна распадается на части), алгоритм отслеживания сетки запутывается в узлах, порождает математические ошибки и часто останавливает программу. Для сложной 3D-металлургии этот подход оказался слишком хрупким.

Метод энтальпии: сквозной счет без выделения границы

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

Математически скрытая теплота плавления не выделяется как сингулярный источник на границе, а искусственно «размазывается» (сглаживается) по небольшому температурному интервалу (например, между 0 и 1 градусом Цельсия). В этом искусственном интервале (называемом двухфазной зоной или mushy zone) эффективная теплоемкость материала становится просто гигантской. Теперь мы решаем одно единое, нелинейное уравнение энергии для всей расчетной области (и для твердой, и для жидкой фазы одновременно), вообще не задумываясь о том, где находится фронт! Позиция фронта легко восстанавливается после расчетов просто путем поиска изотермы, соответствующей температуре плавления. Метод энтальпии легко справляется с любыми топологическими изменениями (слиянием и распадом зон кристаллизации), легко программируется на фиксированных декартовых сетках и сегодня является стандартом во всех литейных CAE-модулях (например, в ProCAST или MAGMASOFT).