ЖУРНАЛ ВЫЧИСЛИТЕЛЬНОЙ МАТЕМАТИКИ И МАТЕМАТИЧЕСКОЙ ФИЗИКИ, 2015, том 55, № 3, с. 502-511
УДК 519.634
АНАЛИЗ РЕЖИМОВ МАГНИТОГАЗОДИНАМИЧЕСКОГО ВЗАИМОДЕЙСТВИЯ ТОКОВОГО СЛОЯ С ПОТОКОМ АРГОНА
В рамках нестационарной трехмерной модели магнитной газодинамики исследована динамика токового слоя, взаимодействующего с поперечным магнитным полем в сверхзвуковом потоке аргона в канале постоянного сечения. По результатам вычислительного моделирования проведен анализ режимов магнитогазодинамического взаимодействия и особенностей формирования структуры токового слоя при различных значениях сопротивления нагрузки и ширины канала. Библ. 11. Фиг. 8.
Ключевые слова: вычислительное моделирование, токовый слой, магнитогазодинамическое взаимодействие, неустойчивость Рэлея—Тейлора.
Взаимодействие токового слоя (ТС) с потоком газа является сложным газодинамическим процессом, зависящим от многих факторов. На характеристики такого процесса значительное влияние оказывают параметры течения на входе в канал, теплофизические и радиационные свойства газа, геометрические размеры канала, величина индукции внешнего магнитного поля, параметры электрической цепи. От этих факторов зависит пространственная форма и размеры ТС, которые, в свою очередь, определяют режим магнитогазодинамического (МГД) взаимодействия, характеристики и структуру течения в канале. Режимы с максимальной эффективностью взаимодействия реализуются при условии, когда поперечное сечение канала полностью перекрывается ТС, который в таком случае выполняет роль плазменного поршня в газовом потоке. В режимах с обтекаемым ТС, который может иметь один или несколько разрядных шнуров, эффективность взаимодействия снижается и зависит от степени перекрытия канала. В целом динамику формирования ТС и характер течения в канале определяет совместное действие энергетических механизмов и баланс электромагнитных и газодинамических сил.
Результаты ранее проведенных исследований пространственной формы ТС и эффективности МГД-взаимодействия в воздухе приведены в [1], [2]. Зависимость характеристик МГД-взаимо-действия ТС в аргоне от сопротивления внешней нагрузки изучалась в [3]. Данная работа является продолжением исследования режимов взаимодействия ТС со сверхзвуковым потоком аргона при различных значениях сопротивления внешней нагрузки и ширины канала.
Рассматривается динамика формирования ТС, взаимодействующего с газовым потоком и внешним поперечным магнитным полем в канале постоянного прямоугольного сечения (см. фиг. 1). Моделирование процесса проведено на основе следующей нестационарной трехмерной математической модели магнитной газодинамики:
© 2015 г. Е. Н. Васильев, Д. А. Нестеров
(660036Красноярск, Академгородок, Ин-т вычисл. моделирования СО РАН)
e-mail: ven@icm.krasn.ru Поступила в редакцию 01.07.2014 г. Переработанный вариант 12.09.2014 г.
DOI: 10.7868/S0044466915030205
1. ВВЕДЕНИЕ
2. МАТЕМАТИЧЕСКАЯ МОДЕЛЬ
dU + dE + 5F + dG = S dt dx dy dz
U = (p,pu,pv,pw, Et),
= S,
(2.1)
(2.2)
2
E = (pu,pu + p,puv,puw,(Et + p)u + qx), (2.3)
F = (pv,pvu,pv + p,pvw,(Et + p)v + qy), (2.4)
G = (pw, pwu, pwv, pw + p,(Et + p)w + qz), (2.5)
S = (0,fx,fy,fz,Qj - Qr + fxu + fyv + fw), (2.6)
p = RpT, Et = p^e +1 + ^ + W] (2.7)
Здесь t — время, p,p, T — плотность, давление и температура газа; u, v, w — компоненты вектора скорости газа v; Et — полная энергия единицы объема газа; e = RT/(y — 1) — внутренняя энергия единицы массы газа; fx, fy, fz — компоненты вектора силы f, действующей на газ вдоль координатных осей x, y и z, соответственно; QJ — объемная мощность джоулевой диссипации, QR — мощность радиационных потерь энергии, определяемая в приближении объемного излучателя (см. [3]); qx = XdT/dx, qy = XdT/dy и qz = XdT/dz — компоненты вектора плотности теплового потока, обусловленного теплопроводностью газа, R — индивидуальная газовая постоянная, у — показатель адиабаты, X — коэффициент теплопроводности газа.
Для расчета электродинамических параметров использовано неоднородное эллиптическое уравнение с переменными коэффициентами
div(a - Уф) = div(a • v х B). (2.8)
Формула (2.8) получена из уравнения непрерывности для электрического тока div j = 0 (токи смещения пренебрежимо малы) с учетом закона Ома j = a(E + v х B) и введения скалярной функции потенциала электрического поля ф, удовлетворяющей соотношению E = -Уф. Решение уравнения (2.8) позволяет найти распределение потенциала ф^, y, z), по которому с помощью приведенных выше соотношений определяются распределения напряженности электрического поля E(x, y, z), плотности тока j(x, y, z), а также джоулевой диссипации QJ = j2/a и силы Лоренца f=j х B. Вектор индукции магнитного поля B = (0, B, 0) учитывает только внешнее поле. Для характерных значений электропроводности газа a = 3 х 103 Ом-1 м-1, скорости токового слоя uT = = 103 м/с и радиуса разрядной области r = 0.075 м величина индуцированного магнитного поля является малой, Bind/B = ^,0au-f/4 ~ 0.03 [1], здесь ц0 = 4я10-7 Гн/м — магнитная постоянная.
Уравнения состояния (2.7) задавались в приближении идеального газа с показателем адиабаты у = 1.67 и молярной массой ц = 40 г/моль. Транспортные (см. [4]) и радиационные (см. [5], [6]) свойства аргона вводились в расчетную программу в виде таблиц.
Для решения системы уравнений газовой динамики (2.1)—(2.7) использована явная схема Мак-Кормака с применением расщепления по пространственным координатам и FCT (Flux-Corrected Transport) метода коррекции потоков. Методика расщепления позволяет свести решение трехмерной задачи к последовательному решению наборов одномерных газодинамических задач, при этом согласованная последовательность применения одномерных операторов обеспечивает 2-й порядок аппроксимации, как по времени, так и по пространству. Уравнение электродинамики (2.8) решается численно методом установления с помощью расщепления по про-
TC электроды
странственным координатам, при этом алгоритм сводится к решению множества независимых одномерных задач. Каждая одномерная задача описывается системой линейных уравнений с трехдиагональной матрицей. Разрешение систем с такими матрицами осуществлялось с помощью метода циклической редукции (см. [7]). Более подробное описание вычислительного алгоритма системы уравнений (2.1)—(2.8) приведено в [8], [9].
Проведение расчетов на персональном компьютере требует больших временных затрат, поэтому для повышения скорости их выполнения актуальным является распараллеливание вычислительного алгоритма. Используемые в данной работе алгоритмы решения многомерных уравнений газодинамики и электродинамики основаны на методике расщепления, это позволяет естественным образом эффективно осуществлять процедуру распараллеливания вычислений. Усовременных персональных компьютеров видеокарты оснащены мощными процессорами GPU (Graphics Processing Unit), которые содержат сотни процессорных ядер, оптимизированных для интенсивной обработки потоков данных, и поэтому могут быть подключены к выполнению вычислений. Расчетный алгоритм данной задачи был реализован так, что параллельные вычисления по схемам расщепления осуществлялись на GPU, а CPU (Central Processing Unit) использовался для оставшихся последовательных вычислений. Для программной реализации был использован интерфейс и язык программирования OpenCL, который является открытым стандартом и обеспечивает реализацию параллельных вычислений на различных графических и центральных процессорах. Моделирование задачи проводилось на компьютере, оснащенном графическим процессором NVIDIA GeForce GTS 450 (чип GF106). Указанная система содержит 192 процессорных ядра, сгруппированных в 4 вычислительных модуля. На основе технологии SIMT (Single Instruction Multiple Thread) система обеспечивает параллельное выполнение 6144 вычислительных потоков. В сравнении с однопоточной реализацией на центральном процессоре использование техники параллельных вычислений на графическом процессоре позволило примерно на 2 порядка увеличить производительность используемых в данной работе расчетных алгоритмов.
3. РЕЗУЛЬТАТЫ ВЫЧИСЛИТЕЛЬНОГО МОДЕЛИРОВАНИЯ
При моделировании рассматривалось газодинамическое течение, которое в начальном состоянии имело температуру Т0 = 103 К и включало в себя цилиндрическую разрядную область с ко-локолообразным распределением температуры со значением в центре Ттах = 104 К. Диаметр начальной разрядной области задавался меньшим ширины канала а и равным 2а/3. Остальные параметры течения имели однородное распределение: давление р = 0.1 МПа, компоненты скорости и = и0 = 103 м/с, V = 0, w = 0. Величина индукции внешнего магнитного поля В принималась равной 2 Тл. Расчеты проводились для высоты канала к = 12.5 см и двух значений ширины канала а = 7.5 см и 15 см, при этом разностная сетка имела 400 х 30 х 50 (а = 7.5 см) и 400 х 60 х 50 (а = 15 см) интервалов по осям х, у и г, величины шагов кх = ку = кг = 2.5 мм. Система координат задавалась движущейся для отслеживания положения ТС по координате х, это позволяет минимизировать количество узлов в расчетной области при заданной длительности процесса.
Для варьирования интенсивности МГД-взаимодействия в качестве управляющего параметра использовано омическое сопротивление внешней нагрузки Параметр нагрузки К = Я1/(ЯР + + = (1 + определяет соотношение величин омических сопротивлений плазмы ЯР и
нагрузки Электродвижущая сила £,, генерируемая при движении электропроводного газа в магнитном поле, распределяется в виде падений напряжения на межэлектродном промежутке АУР = ^(1 — К) и внешней цепи = ^К. При прочих равных условиях интенсивность МГД-вза-имодействия возрастает с уменьшением сопротивления нагрузки, при этом увеличивается падение напряжения и напряженность электрического поля в межэлектродном промежутке, а это, в свою очередь, приводит к росту плотности тока, величины электромагнитной силы и мощности джоулевой диссипации.
Расчеты проведены для длительности процесса 1 мс в широком диапазоне значений соответствующих различной интенсивности МГД-взаимодействия в канале. Динамика формирования структуры ТС при а = 15 см исследовал
Для дальнейшего прочтения статьи необходимо приобрести полный текст. Статьи высылаются в формате PDF на указанную при оплате почту. Время доставки составляет менее 10 минут. Стоимость одной статьи — 150 рублей.