Матричный подход используется для моделирования и исследования устойчивости осесимметричных течений несжимаемой жидкости. Для решения системы нелинейных дискретизованных уравнений Навье-Стокса используется метод Ньютона в матричном виде. Для исследования устойчивости течения относительно трехмерных возмущений используется метод обратной итерации и метод Арнольди. Решение систем линейных уравнений реализуется прямым мультифронтальным методом, реализованным в пакете MUMPS. На основе применения матричного метода для численного бенчмарка для модели роста кристаллов по Чохральскому удалость получить решение на сверхподробной сетке с 125M неизвестных при использовании 65 процессоров с памятью 24 ГБ на СК «Ломоносов».
Целью данной работы является определение границ устойчивости течения расплавов в
представленных конфигурациях роста кристаллов [
На основе проведенных расчетов для бенчмарка по росту кристаллов по методу
Чохральского для больших чисел Прандтля [
Гидродинамическая модель роста кристаллов по Чохральскому представляет собой цилиндрический объем жидкости (расплава) высоты H и радиуса Rc (рис. 1). Боковая и нижняя стороны ограничены тиглем, с верхней стороны при r < Rx находится диск (кристалл), вне диска для Rx < r < Rc поверхность является свободной.
Поведение жидкости определяется уравнениями Навье-Стокса в приближении Буссинеска для ньютоновской жидкости, поверхностное натяжение которой предполагается линейной функцией температуры T :
0 (T T0 ) . Безразмерные уравнения неразрывности, движения и притока тепла в цилиндрических координатах (r, , z) с компонентами скорости и давлением P имеют вид: где Граничные на твердых поверхностях имеют вид: Граничные на свободной поверхности имеют вид: Безразмерные параметры включают два безразмерных геометрических параметра: , Ar= а также термокапиллярное число Рейнольдса Re , число Грасгофа Gr ,число Прандтля Pr и число Био Bi.
На первом этапе ищутся стационарные осесимметричные решения выписанных уравнений. На втором этапе найденное решение рассматривается как базовое и для него рассматривается устойчивость относительно малых возмущений, которые представляются в виде: получается обобщенная задача на соб3. Метод решения и тесты
Для численного решения уравнений используется метод контрольного объема с использованием разнесенных сеток. При пространственной дискретизации радиальная и осевая компоненты скорости определяются на гранях сеточных объемов, а скалярные величины, давление и температура, а также азимутальная компонента скорости определяются в центрах этих объемов. Конвективые производные аппроксимируются центральными разностями, обеспечивающими второй порядок аппроксимации. Введение, при необходимости, существенно неравномерной сетки позволяет воспроизвести тонкую структуру пограничных слоев вблизи точек соприкосновения твердых и свободных поверхностей без введения регуляризации.
Для решения стационарных осесимметричных разностных (дискретизованных) уравнений на первом шаге используется итерационный метод Ньютона в матричной форме для f (x) 0 : x[k1] x[k] f (x[k]) f (x[k])
, k 0,1,... для некоторого начального значения x[0] . Под вектором x подразумевается все множество неизвестных задачи, а под функцией f все множество уравнений и граничных условий. Так как якобиан f (x[k]) является матрицей, то для нахождения решения на следующей итерации требуется решать систему линейных алгебраических уравнений.
Использование матричных методов решения уравнений имеет ряд преимуществ по сравнению с традиционными явными и полунеявными методами: 1. возможность построения эффективных полностью неявных схем без расщепления и обременительных ограничений, связанных с шагом интегрирования по времени. 2. отсутствие необходимости разработки метода решения разностных уравнений; разностные уравнения решаются матричными методами. 3. возможность получения решений, неустойчивых для нестационарных уравнений; это необходимо при исследовании линейной устойчивости решений. 4. Высокая скорость сходимости при использовании метода Ньютона для решения нелинейных уравнений; решение уравнений находится за несколько итераций.
Недостатки матричных методов связаны с необходимостью решения систем линейных
уравнений высокого порядка. Матрица якобиана при использовании разностной
аппроксимации имеет блочно-диагональный характер и является сильно разреженной. Число элементов
матрицы примерно в 10 раз больше размерности вектора неизвестных. В данной работе для
решения систем линейных уравнений используются прямой метод решения уравнений,
основанный на на реализации мультифронтального метода [
Количество даже открытых пакетов для решения задач линейной алгебры весьма велико
[8]. Выбор конкретного пакета следует предпринимать из опубликованного опыта
предшественников, что не слишком просто из-за практического отсутствия опыта перекрестного
тестирования. Вопросы, связанные с выбором программных средств для моделирования задач
сеточными методами на суперкомпьютерах, являются очень важными, и они подробно
обсуждаются в [
Значительной экономии памяти теоретически можно было бы достичь используя итерационный метод решения градиентного типа для системы линейных уравнений, однако, построение эффективного предобуславливателя для полной системы уравнений конвекции является достаточно сложной задачей.
Для решения обобщенной задачи на собственные значения используется метод обратной итерации и метод Арнольди. Последний используется в режиме shift-and-invert пакета ARPACK. Для определения устойчивости течения достаточно найти собственное значение с наибольшей действительной частью. Если действительная части такого собственного значения положительно, то течение является неустойчивым, в противном случае оно устойчиво. Использование стандартных методов решения спектральной задачи затруднено вырожденным характером матрицы B. Для обоих подходов необходимо снова решение системы линейных уравнений, на этот раз над полем комплексных чисел.
Данный подход был разработан и тщательно тестирован для различных задач
вынужденной, термогравитационной и термокапиллярной конвекции [
Расчеты проводились на сетках до 5000×5000 точек (125 млн. неизвестных) для нахождения стационарного решения и до 3000×3000 точек (45 млн. неизвестных) для решения задачи устойчивости. Для задач использовалось до 65 процессоров с памятью 24 ГБ и до 17 процессоров с памятью 48 ГБ. Результатом решения задачи являются стационарные осесимметричные течения и нейтральные кривые, определяющие границы устойчивости, критические частоты и собственные вектора, определяющие наиболее опасные возмущения. (Будут приведены данные по временным затратам; в данный момент они отсутствуют так расчеты проводились для несколько другой постановки задачи).
В таблице 1 представлено сравнение результатов по бенчмарку для высоких чисел
Прандтля [
Таблица 1. Сравнение данных по бенчмарку [
(z) Vr max (r=0.5)
(z) Vz min (z=0.5H)
(r) Vz max (z=0.5H)
(z) Ekin Nu (дно) Nu (стенка) Nu (кристалл) Nu (поверхность) Pr=0.005 Pr=0.05 Pr=0.6 Pr=4 Pr=10 Pr=30 Рис. 2. Примеры собственных функций задачи (температура – слева, азимутальная компонента скорости - справа) для чисел Прандтля от 0.005 до 30. 5. Выводы
Для осесимметричных конвективных течений роста кристаллов из расплава по методу Чохральского разработана методика решения стационарных задач и исследования устойчивости течений по отношению к трехмерным возмущениям для кластерных систем. Для решения задач используется матричный подход, основанный на применении методов Ньютона, обратных итераций, Арнольди и линейного решателя (солвера) MUMPS для прямого решения разреженных систем.
Использование данного подхода позволило получить численное решение на сетках до
5000×5000 точек (125 млн. неизвестных) на 65 процессорах с памятью 24 ГБ или на 17
процессорах с памятью 48 ГБ. Результаты получены для международного бенчмарка по росту
кристаллов из расплава по Чохральскому для высоких чисел Прандтля [
Разработанные методы и коды могут быть также использованы для исследования устойчивости течений в задачах внешнего обтекания, крупномасштабной конвекции в мантии и атмосфере Земли и динамики плазмы в осесимметричных конфигурациях.
Работа выполнена с использованием ресурсов суперкомпьютерного комплекса МГУ имени
М.В. Ломоносова [
URL: http://www.netlib.org/utk/people/JackDongarra/la-sw.html. Matrix methods for fluid mechanics problems Michael Ermakov A matrix approach is used for modeling and stability analysis of axisymmetrical fluid flows for the Navier-Stokes equations in Boussinesq approximation. Advantages of the matrix approach consist in an unnecessity of equations splitting and an unnecessity of a partial method building for equations solution. For the solution of a system of discretized NavierStokes equations the matrix Newton method is employed. At each iterations a linear system for a vector of the full set of unknowns and for the Jacobian matrix is solved. For a flow stability analysis against 3-D disturbances, are being considered as a decomposition of normal modes in azimuthal direction, the inverse iterations and Arnoldi methods (ARPACK) are used. For a linear-system solution the direct method (MUMPS) is used. Due to application of matrix approach for numerical benchmark for Czochralski method of crystal growth, it was possible to get a solution for superfine mesh with 125M unknowns by use of 65 processors (24 GB) of the “Lomonosov” supercomputer.