ФИЗИКА ПЛАЗМЫ, 2009, том 35, № 4, с. 291-302
= ТОКАМАКИ
УДК 533.951.7
ТУРБУЛЕНТНАЯ ДИНАМИКА ВИНТОВЫХ ВОЗМУЩЕНИЙ В ПРИСТЕНОЧНОЙ ПЛАЗМЕ ТОКАМАКА Т-10
© 2009 г. Р. В. Шурыгин, А. В. Мельников
РНЦ "Курчатовский институт", Москва, Россия Поступила в редакцию 23.06.2008 г.
Проведен численный расчет турбулентной динамики пристеночного слоя плазмы токамака Т-10 на основе решения нелинейных МГД-уравнений в рамках 4-х полевой |ф, п, ре, р} редуцированной двухжидкостной гидродинамики Брагинского. Показано, что при переходе от омического к электронно-циклотронному нагреву плазмы амплитуды турбулентных флуктуаций падают вследствие увеличения продольной диссипации, связанной с ростом электронной температуры. Однако фазовые соотношения флуктуаций потенциала различных мод изменяются таким образом, что турбулентная сила Рейнольдса увеличивается. Рост этой силы приводит к увеличению полоидальной скорости в направлении ионного диамагнитного дрейфа. Так как в уравнение для радиального баланса сил для ионов полоидальная скорость и скорость диамагнитного дрейфа ионов входят с разными знаками, то величина радиального электрического поля уменьшается. Данные расчетов качественно согласуются с результатами эксперимента на Т-10. Кроме того, получена зависимость радиального электрического поля от плотности плазмы, температуры ионов и концентрации нейтралов.
РАСЯ: 52.35.Ra, 52.55.Fa
1. ВВЕДЕНИЕ
Расчет турбулентных транспортных потоков в магнитных ловушках можно отнести к одной из важных, но пока еще окончательно не решенных задач в проблеме УТС. Она проистекает как из потребностей теории, так и из необходимости разрешения инженерно-технологических вопросов, возникающих при проектировании установок с магнитным удержанием. В последнее десятилетие пришло понимание того, что для решения указанной задачи необходимо использовать прямые численные методы, способные реалистично и адекватно моделировать нелинейные физические процессы, протекающие в плазме токамака. Во многих лабораториях созданы эффективные численные алгоритмы, из которых наибольшее развитие получили коды, описывающие поведение плазмы на основе решения редуцированной системы двухжидкостных МГД-уравне-ний [1—4].
В данной работе проводится численное моделирование турбулентной динамики пристеночной плазмы токамака Т-10 на основе решения нелинейных МГД-уравнений в рамках редуцированной двухжидкостной гидродинамики Брагинского. При создании численного кода были использованы результаты работы [5], где проведена строгая процедура редуцирования полных МГД-уравнений с учетом конечной температуры ионов (Т ~ Те) и неоднородности плотности. При этом появляется возможность корректного учета
как скорости электрического дрейфа, так и скорости диамагнитного дрейфа ионов. В результате в уравнении для обобщенного вихря, наряду с основной нелинейностью, связанной с конвективным переносом вихря, возникает ряд новых нелинейных членов, связанных с флуктуациями плотности, электростатического потенциала и ионного давления. В нашей работе мы показываем, что дополнительные нелинейности играют важную роль при получении эволюционного уравнения для усредненного (по угловой переменной) полоидального момента импульса. Показано, что полученный в результате усреднения тензор турбулентных напряжений Рейнольдса включает новые члены, связанные с флуктуация-ми ионного давления, а также нулевой гармоникой электрического потенциала, которые приводят к дополнительной генерации полоидального момента импульса.
В современных экспериментах на токамаках широко используется метод дополнительного СВЧ-нагрева плазмы на электронно-циклотронном резонансе (ЭЦР). В данной работе проводится численный расчет поведения турбулентных полей плазмы в плоском пристеночном слое токамака Т-10 при омическом (ОН) и ЭЦР-нагреве. Переход от ОН к ЭЦР-нагреву моделировался увеличением электронной температуры ТЬе на внутренней границе слоя, примыкающей к основной плазме. Измерения профиля электростатического потенциала на токамаке Т-10 показы-
вают его снижение в направлении от стенки к центру (Ег = —д ф 0 /дг < 0) [6]. При переходе от ОН к ЭЦР-нагреву это снижение ослабляется, то есть модуль радиального электрического поля также уменьшается \Ег \ЕСК < \Ег |он, при этом поле остается отрицательным.
В разд. 2 выписана основная система уравнений. В разд. 3 приводится расчет радиального профиля электростатического потенциала в турбулентном режиме пристеночной плазмы тока-мака Т-10 для двух указанных режимов нагрева, что является основной целью настоящей. В разд. 4 изучаются зависимости радиального электрического поля от плотности плазмы, температуры ионов на внутренней границе расчетного слоя и концентрации нейтралов. В разд. 5 обсуждается природа структурной турбулентности на основе анализа двумерных полей флуктуаций электростатического потенциала и плотности.
2. ОСНОВНЫЕ УРАВНЕНИЯ
В рамках редуцированной двухжидкостной гидродинамики Брагинского проводилось прямое численное моделирование нелинейных МГД-уравнений. Предполагается, что указанные уравнения справедливы, так как выполнено условие сильной столкновительности X е/дЯ ^ 1 (где
X е = УТе1 V е1 — длина свободного пробега электронов), то есть рассматриваемая плазма достаточно холодная и плотная. Расчеты проводились для узкого плоского пристеночного слоя токамака в приближении плоской декартовой геометрии. Введем систему координат (х, у, г), где ось х направлена вдоль большого радиуса тора, у — в по-лоидальном направлении, г — вдоль тора. В этом случае магнитное поле токамака можно представить в виде
В = Во(х)е, + Воу(х)е
В0, - В01 1 -
Х_
Я0
= В0&
В,
0у
д(г)
где и г — большой и малый радиусы токамака, ^ — коэффициент запаса устойчивости, е = г/В0.
Введем винтовую переменную У = у — аг (параметр а указан ниже) и далее будем рассматривать случай винтовой симметрии /(х, у, г, 0 = =/(х, У, ?). Тогда безразмерная система МГД-уравнений будет [7, 8]
^ + = N0Уц/ ^ + + V„А(1) Dt Ь дУ
Оп = у/ + 4в д(^Ф- Ре) + ^ + Л„А пП, (2)
Dt 11 Ь дУ п п
3 ^ = 5 . + + 2 Л 2 " 2 Ь е
+ Х||е^2Ре + XleVipе + ^
ре>
3 = _т.0У||/ + 5_^ВК{ + 2 Л 2 " 2 Ь '
+Х||(^2р1 +х1_1У21р1 + s¡
Р"
П</ =-У|(ф- Рe/Nо),
„ = ^ V ±(лУ ±ф + У л Ь
Р = Ре + Р, Г - р дф 1Т дРе , т2 дп
(3)
(4)
(5)
(6)
дУ
дУ
К _ Р дФ + Т д(Р1 - Ре) т2 дп к I _ р0--+ т ю--т ю-■
' дУ дУ дУ
дУ
■- т0
Здесь ^о, Ре0, Рю, Те0, Т10 — нулевые, не зависящие от У, гармоники плотности, давлений и температур,
А, В} = е,УА хУВ,
Л = д +1 Dt дt Ь
У/ = ь-У/ = Ьу- + Ь,-, Ь = В, ду д, В
2й 0.51 vei 8в =—, П0 =-,
Я Шее
х = 3.16 Шее х _3.9Ю„-
л||е > Л||(
Vei Vil
Бш, Бп, Бре, 8р1 — члены, описывающие источники, выражения для которых приведены ниже. Приведение системы (1)—(6) к безразмерному виду проделано с помощью преобразования переменных
t — t|t*, (х,у,У) — (х/й,у/й,У/й), ф-> еф Т*, В ^ Ь = В В*,
п — п1п*, Ре,1 — Ре,^Р* , К — kyd,
Р = Р */й, Ш0 = Р2ю„-,
где t* = 1/®0, Р* = К*/Юе,-,К* = ^Т*/т, ©ее, I = = еВ/(те,е), Р* = п*Т*.
Решения системы (1)—(6) рассматриваются в узком плоском слое пристеночной плазмы шириной d = а — г0, где а и г0 — радиусы внешней и внутренней границ слоя. Значения нормировочных величин для плотности, магнитного поля и элек-
тронной температуры принимались равными п* = = 1013 см-3 , В* = 2 Тл и Т* = 200 эВ соответственно. В расчетах используется безразмерная радиальная переменная х = (г — г0)/(а — г0), 0 < х < 1.
Система нелинейных уравнений (1)—(6) описывает электростатическую турбулентность, возникающую на первоначальном этапе расчета как результат совокупного действия диссипативных дрейфовых неустойчивостей, связанных с градиентами плотности, давления ионов и электронов, а также тороидального магнитного поля — мода). Однако в режиме развитой турбулентности при описании динамики нелинейных дрейфовых волн, формирующих долгоживущие когерентные структуры, традиционное понятие линейной неустойчивости, основанное на анализе линеаризованных уравнений, во многом теряет физический смысл. Численные расчеты показывают, что турбулентный транспорт проявляет себя как сложная комбинация нелинейной конвекции и диффузии, создавая принципиальную невозможность адаптации результатов линейной теории к вычислению коэффициентов турбулентного переноса.
Важно отметить, что в предложенной модели
наряду с электрическим дрейфом УЕ = сВ_1Ъ х Уф, корректно учитывается диамагнитный дрейф
ионов Ур1 = с(внВ)-1Ь х Ур, а также неоднородность плотности (V п) в выражении для вихря. Это приводит к более сложному определению вихря (обобщенный вихрь и в (6)) и появлению в уравнении (1) дополнительных нелинейных членов Иж вида
Нж =Д23 [VI {ф, р}- (ф,у!р } 2Ь у '
Р,,Щ + n}}.
Для построения численной схемы в данной работе используется квазиспектральный подход, основанный на методе Галеркина. Все зависимости f = {n, w, ф, pe, p} выбирались в виде суммы волн с одной спиральностью, то есть
f(x, y, z, t) = Z fk,K exp[i(kyy - kzz)] =
= X fky exP[ iky(y -a z)\
Выбор волн с одной спиральностью подразумевает, что параметр а является постоянной величиной, а = kjky = const (для мод в тороидальной геометрии c = n/R0, ky = m/r получим
а = е/^ге8, е = гге5/К0), что позволило нам перейти к винтовой переменной У = у — а1.
Важно отметить, что переход к винтовой переменной сохраняет физическое содержание рассматриваемой МГД-задачи, в которой принципиально наличие трехмерных возмущений, вытянутых вдоль силовых линий магнитного поля. При этом в математическом смысле задача становится двумерной, что позволяет при численной реализации существенно экономить время счета.
Будем предполагать, что для любой функции / = {п, и, ф, ре, р} из системы (1)—(6) существует периодичность по винтовой переменной У
/(у) = /(У + 2 п).
Тогда на интервале 0 < У < ЬУ = 2пг0/ё справедливо разложение в ряд Фурье
f(x,Y, t) = f0(x, t) +
+
Z
L = 1
fsL Sin
2Кшл
Ly
Y 1 + fcL cos
2 nmY
Lv Y
m = m(L), k0y = dlrres, Ly = 2п/k0r
Расчетная пристеночная область может быть представлена прямоугольником r0 < r < a, 0 < y < < 2nr0, или — в безразмерных единицах — 0 < х < 1,
Для дальнейшего прочтения статьи необходимо приобрести полный текст. Статьи высылаются в формате PDF на указанную при оплате почту. Время доставки составляет менее 10 минут. Стоимость одной статьи — 150 рублей.