автореферат диссертации по информатике, вычислительной технике и управлению, 05.13.18, диссертация на тему:Методы оценки и повышения точности решения задач физики плазмы методом частиц в ячейках

На правах рукописи

< (

Мссяц Екатерина Александровна

МЕТОДЫ ОЦЕНКИ И ПОВЫШЕНИЯ ТОЧНОСТИ РЕШЕНИЯ ЗАДАЧ ФИЗИКИ ПЛАЗМЫ МЕТОДОМ ЧАСТИЦ В ЯЧЕЙКАХ

05.13.18 - математическое моделирование, численные методы и комплексы программ

АВТОРЕФЕРАТ диссертации на соискание ученой степени кандидата физико-математических наук

2 г МАЙ 2014

005548974

Новосибирск - 2014

005548974

Работа выполнена в Федеральном государственном бюджетном учреждении науки Институт вычислительной математики и математической геофизики Сибирского отделения Российской академии наук (ИВМиМГ СО РАН).

Научный руководитель: доктор физико-математических наук, профессор

Вшивков Виталий Андреевич

Официальные Медведев Сергей Борисович,

оппоненты: доктор физико-математических паук,

Федеральное государственное бюджетное учреждение науки Институт вычислительных технологий Сибирского отделения Российской академии наук (ИВТ СО РАН), старший научный сотрудник

Астрелин Виталий Тимофеевич,

кандидат технических наук,

Федеральное государственное бюджетное учреждение науки Институт ядерной физики им. Г.И. Буд-кера Сибирского отделения Российской академии наук (ИЯФ СО РАН), старший научный сотрудник

Ведущая организация: Федеральное государственное бюджетное учреждение науки Институт теоретической и прикладной механики Сибирского отделения Российской академии наук (ИТПМ СО РАН)

Защита состоится 30 июня 2014г. в 16.30 на заседании диссертационного совета Д 003.061.02 на базе Федерального государственного бюджетного учреждения науки Институт вычислительной математики и математической геофизики Сибирского отделения Российской академии наук, расположенном по адресу: 630090, г. Новосибирск, проспект Академика. Лаврентьева, 6.

С диссертацией можно ознакомиться в библиотеке и на сайте ИВМиМГ СО РАН, http://www.sscc.ru.

Автореферат разослан 28 апреля 2014 г.

Ученый секретарь

диссертационного совета Д 003.061.02 доктор физико-математических наук

Сорокин Сергей Борисович

ОБЩАЯ ХАРАКТЕРИСТИКА РАБОТЫ

Актуальность работы. Метод частиц в ячейках, возникший еще в середине прошлого столетия, на сегодняшний день широко применяется при моделировании поведения плазмы. Его также используют и при расчете динамических процессов в других средах (жидкости, газы, сплошные среды и др.), но именно в решении задач физики плазмы он получил наиболее широкое распространение. Как самостоятельный раздел вычислительная физика плазмы сформировалась одновременно с развитием вычислительной техники во второй половине XX века. Это привело к выделению математического моделирования в отдельное направление исследовательской работы. Ввиду больших возможностей диагностики моделируемых в численных экспериментах кинетических процессов и относительно невысоких экономических затрат на их проведение (по сравнению со строительством установок для проведения реальных экспериментов), численные методы играют все более важную роль в решении задач физики плазмы.

Существует ряд математических моделей описания плазмы: кинетическое описание, магнитогидродинамическое приближение, модели термоядерной плазмы и другие. В 1938 году Власов предложил свою концепцию описания широкого круга процессов в плазме. В основе модели Власова лежит кинетическое уравнение без столкновительного члена для функции распределения частиц но скоростям /(£,х, V). Электрическое и магнитное поля, входящие в уравнение через силу Лоренца, определяются из уравнений Максвелла. Система уравнений замыкается формулами для плотности зарядов и токов, которые выражаются через функции распределения частиц. Система уравнений Власова-Максвелла послужила основой для большого числа работ по теории колебаний, устойчивости, коллективных процессов в плазме.

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

Метод частиц в ячейках принадлежит группе вычислительных алгоритмов, объединенных способом дискретизации, которая носит название «методы частиц». В отличие от конечно-разностных методов, объединенных аппроксимацией дифференциальных или интегральных операторов исходных уравнений, методы частиц объединяет концепция представления среды в виде множества модельных частиц, которые являются носителями набора ха-

рактеристик этой среды (масса т. заряд д, импульс р и т.п.). Множество точек, которое составляют координаты частиц, принято называть «лагран-жевой» сеткой, в противовес «эйлеровой» сетке, привязанной к области. В этой терминологии методы частиц можно разделить на лагранжевы методы и смешанные эйлерово-лаграпжевы методы. Метод частиц в ячейках относится к группе смешанных алгоритмов, вычисления в которых проводятся в два этапа: на эйлеровой и на лаграпжевой сетке попеременно.

Такие преимущества, как отсутствие аппроксимацпонной вязкости и возможность осуществления сквозного счета при резких изменениях поведения решения в областях больших градиентов, поставили методы частиц для задач с большими объемными деформациями и пеустойчивостями на первое место. Вследствие этого методы частиц нашли широкое применение именно для моделирования таких сложных явлений, как неустойчивости в физике плазмы. Главным преимуществом метода частиц в ячейках по сравнению с эйлеровыми методами является простота реализации и экономичность, так как в методах частиц не вычисляется напрямую функция распределения. Простота реализации данного метода (в том числе и на параллельных ЭВМ), наглядность получаемых распределений частиц, большое число решенных на его основе задач составляют сильные стороны этого метода.

Основным недостатком метода частиц в ячейках всегда являлось наличие в решении нефизических эффектов, так называемых «численных шумов». Эта проблема очень сложна, так как причин возникновения шумов много. Влияние различных факторов на решение трудно разделить, особенно при моделировании неустойчивых сред. Шумовые гармоники взаимодействуют друг с другом, накладываются на гармоники неустойчивых решении, что может приводить к развитию численных неустойчивостей.

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

Цель диссертационной работы заключается в исследовании свойств метода частиц в ячейках и разработке общих алгоритмов уменьшения счет-

пых шумов данного метода. Для достижения этой цели было решено

• создать алгоритм подавления шума в методе частиц в ячейках,

• разработать методику выбора оптимальной формы ядра, которая обеспечивает по сравнению со стандартным Р1С-ядром уменьшение самосилы,

• создать комплекс программ для моделирования взаимодействия электронного пучка с плазмой в трехмерной геометрии и провести оценки точности получаемого решения.

Научная новизна работы заключается в следующем:

1. Разработан новый алгоритм уменьшения счетных эффектов метода частиц в ячейках, основанный на вычитании шумовой добавки, обеспечивающий подавление счетных флуктуаций при сохранении самого решения. Данный алгоритм реализован при решении задачи о распаде разрыва плотности ионов в неизотермической плазме.

2. Проведено исследование причины возникновения самосилы для метода частиц в ячейках на смещенных сетках. Для уменьшения этой ошибки разработан экономичный подход, основанный на использовании нового ядра модельной частицы с меньшей, по сравнению с РЮ-ядром, самосилой.

3. Показано, что использование нового ядра приводит к увеличению времени саморазогрева модельной плазмы и, как следствие, к лучшему сохранению импульса и полной энергии.

4. На основе метода частиц в ячейках создан комплекс программ для моделирования взаимодействия электронного пучка с плазмой в трехмерной постановке при параметрах, соответствующих условиям лабораторных экспериментов на установке ГОЛ-3, который позволяет проводить вычисления как в гидродинамическом, так и в кинетическом режимах развития пучковой неустойчивости.

5. Проведены оценки величины погрешности решения, получаемого данным комплексом программ. Для различных режимов приводятся оценки достаточного количества частиц.

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

дач физики плазмы методом частиц для уменьшения уровня нефизических флуктуаций.

Представленные в диссертации исследования проводились в рамках Интеграционных проектов СО РАН №113, №40, Научной школы Годунова -4292.2008.1 и по проектам, поддержанным Российским фондом фундаментальных исследований (№ 08-01-00615, 11-01-00249, 11-01-00178).

Достоверность результатов. Все численные алгоритмы проверялись по отдельности на тестовых расчетах. Программные комплексы также проходили тестирование, проводилось сравнение результатов моделирования с имеющимися аналитическими решениями или с решениями той же задачи конечно-разностными методами. Также проверялась сходимость численных методов решения при сгущении сетки и увеличении числа модельных частиц.

На защиту выносятся:

• алгоритм вычитания шумовой добавки.

• алгоритм поиска оптимальной формы ядра и новое С^СР-ядро, минимизирующее самовоздействие при вычислении полей на сдвинутых друг относительно друга сетках,

• созданный комплекс программ для моделирования взаимодействия электронного пучка с плазмой в трехмерной геометрии,

• оценки точности вычисления инкремента неустойчивости для созданного комплекса программ.

Апробация работы Основные положения диссертационной работы докладывались и обсуждались на семинарах «Математическое моделирование больших задач» под руководством д.ф.-м.н. Вшивкова В.А. (ИВМиМГ СО РАН), на семинаре «Математическое и архитектурное обеспечение параллельных вычислений» под руководством д.т.н. Малышкина В.И. (ИВМиМГ СО РАН, декабрь 2009), на семинаре «Задачи механики и математической физики» под руководством д.ф.-м.н. Медведева С.Б. (ИВТ СО РАН, июнь 2013), на объединенном семинаре ИВМиМГ СО РАН и кафедры вычислительной математики НГУ под руководством д.ф.-м.н. Ильина В.П. (ИВМиМГ СО РАН, октябрь 2013), на семинаре «Законы сохранения и инварианты» под руководством д.ф.-м.н. Медведева С.Б. (ИВТ СО РАН, февраль 2014), на Всероссийской конференции по математике и механике, посвященной 130-летию ТГУ (2008, Томск), на Международной конференции, посвященной 100-летию со дня рождения С.Л. Соболева (2008, Новосибирск), на Конференции молодых ученых ИВМиМГ СО РАН (2009, 2011, Новосибирск), на XIII Всероссийской молодежной конференции-школе «Современные проблемы математического моделирования. Математическое моделирование, чис-

ленные методы и комплексы программ» (2009, Дюрсо), на VI Всесибирском конгрессе женщин-математиков (2010, Красноярск), на XV Байкальской Всероссийской конференции «Информационные и математические технологии в науке и управлении» (2010, Иркутск-Байкал), на Всероссийской конференции по вычислительной математике КВМ-2011 (Новосибирск), на Российско-монгольской конференции молодых ученых по математическому моделированию, вычислительно-информационным технологиям и управлению (2011, Иркутск (Россия) - Ханх (Монголия)), на XIX Всероссийской конференции «Теоретические основы и конструирование численных алгоритмов и решение задач математической физики» (2012, Дюрсо), на XIII Всероссийской конференции молодых ученых по математическому моделированию и информационным технологиям (2012, Новосибирск), на Международной конференции «Mathematical modeling and computational physics» (2013, Дубна).

Основные результаты опубликованы в 15 работах, из которых 2 в журналах, рекомендованных ВАК. Список публикаций приведен в конце автореферата.

Личный вклад соискателя заключается в обсуждении постановок задач, разработке численных алгоритмов и методов решения, создании и тестировании программ, проведении расчетов и анализе полученных результатов. Все выносимые на защиту результаты принадлежат лично автору. Представление результатов совместных исследований и разработок согласовано с соавторами.

Структура и объем работы. Диссертация состоит из введения, 4 глав и заключения. Список использованной литературы содержит 133 наименования (включая 15 публикаций автора). Текст диссертации содержит 110 страниц машинописного текста, включая 36 рисунков и 4 таблицы.

КРАТКОЕ СОДЕРЖАНИЕ РАБОТЫ

Во введении обоснована актуальность темы диссертации, кратко приведены основные результаты, их новизна и практическая ценность.

В первой главе содержится обзор методов численного решения системы уравнений Власова-Пуассона (1), (4), (6) (В = О, Е = — Vy) и Власова-Максвелла (1), (2)-(6)

| + v.^ + A(E+[vxB]).^ = 0, (1)

at or та 3v

(ЭЕ

— =crotB-4Trj, (2)

<ЭВ

— = -crotE, (3)

сИу Е = 47гр, сИув = о.

I /а^У, = ! }асЫ.

(4)

(5)

(6)

Задание начальных условий

ут+\пт+х :ЭV

Э( Эх

дп 4

Э( дх

X

Корректировка положения частиц

Х"! = X"' +

Р"

V >е

Решение уравнений движения, Е=0

X

Здесь а - сорт частиц (ионы или электроны); /„(£, г, у)-функция распределения частиц сорта а\ тп,

га - заряд, масса, скорость, положение частицы; Е, В - напряженности электрического и магнитного полей.

Описано место метода частиц в ячейках среди других методов и приведен обзор литературы по проблеме численных шумов в методе частиц. В общем виде сформулирована схема метода частиц в ячейках, даны основные определения и понятия.

Вторая глава посвящена разработке алгоритма вычитания шумовой добавки. Предлагаемый алгоритм рассматривается на примере одномерной задачи распада разрыва плотности ионов в бесстолк-

новительной полностью ионизованной плазме, которая описывается системой уравнений Власова-Пуассона. В методе частиц в ячейках уравнение Власова заменяется уравнениями его характеристик, которые задают траектории движения модельных частиц

Корректировка скорости и 7 = и"1 -AV¡

Лагранжев этап

X

□

Эйлеров этап

Рис. 1. Блок-схема цикла модифицированого метода частиц в ячейках.

л

сЦ сИ

+ [и, X В;])

(7)

Здесь у - номер частицы, (х^, щ) - координата и скорость частицы.

Проведено сравнение решения, полученного по методу частиц в ячейках, с решением, полученным по конечно-разностной схеме (схема Лакса-Вендроффа) для той же задачи. Показано влияние счетных параметров модели на уровень численного шума.

Блок-схема метода частиц в ячейках с разработанным алгоритмом вычитания шумовой добавки показана на рис. 1. Суть его состоит в следующем.

На шаге т скорости и координаты частиц (и™, х™) корректируются таким образом, чтобы на те + 1 шаге, сделанном без действия электрического поля

(Е = 0), средняя скорость и плотность , п™+1), полученные по мето-

ду частиц в ячейках, совпали с решением уравнений для первых моментов уравнения Власова (У/™+\ (уравнение неразрывности и уравнение пе-

реноса). Д V-, - шумовая добавка в среднюю скорость в ячейке г, с?; - корректирующие сдвиги координат частиц в ячейке г.

Как показали результаты проведенных расчетов, алгоритм вычитания шумовой добавки позволяет получать методом частиц в ячейках более гладкое решение (рис. 2), не подавляя при этом развивающихся физических волн.

Средняя скорость

0,04 0,02 О

2,1

Плотность

- метод частиц

----- модифицированный метод частиц

" * "* \

0,5

1,5

■1,9

2 0 0,5

1,5

Рис. 2. График средней скорости и плотности для стандартного и модифицированного (с алгоритмом вычитания шумовой добавки) метода частиц в ячейках.

В третьей главе рассматривается одна из причин, приводящих к изменению движения частицы, и как следствие, к повышению столкновительности в плазме, а именно наличие пространственной сетки и возникновение так называемой «самосилы». Показано, что при вычислении полей на смещенной сетке (плотность и потенциал задаются в узлах сетки, а напряженность электрического поля вычисляется в центре ячейки (одномерный случай), в центре грани (двумерный случай) самосила не равна нулю. Случай со смещенной сеткой (при hx = hy = h) и рассматривается в данной главе.

Новая предложенная форма ядра частицы, названная QCP («quadruple-core-particle»), представляет собой следующую комбинацию

Rqср{х, у) = a0RPIC{x, у) + aiRQCpi(x, у) + a2RQCP2(x, у). (8) Здесь функции ядер имеют следующий вид 1. ядро PIC (классическое ядро метода частиц в ячейках)

(h - M),

\х\ < h,

Ы > h.

2. ядро QCP1

|®| < h, \у\ < h,

Rqcpi{x,y) = тт7 { (2h - |y|)|x|,

(2h-\x\)\y\, h < |x| < 2h, \y\ < h,

(2ft-|vlM \x\<h, h<\y\<2h, (10)

(2h - \x\)(2h - \y\), h<\x\< 2h, h < \y\ < 2h,

0

иначе.

3. ядро QCP2

' h(\x\ + \y\)-2\xy\, \x\<h, \y\ < h, ! (2h-\x\)(h-\y\), h < \x\ < 2h, \y\ < h,

RQCP2(x,y) = —

(А-М)(2Л-|у|), \x\ < h, h < \y\ < 2h,

(H)

иначе.

Получено, что при значениях параметров а о = 0, ai = 0.6912, а2 = 0.3088 самос.ила равна нулю. Но дальнейшие расчеты показали, что при моделировании методом частиц в ячейках с использованием та.кого ядра иногда возникают численные неустойчивости. Дополнительным требованием, необходимым для того, чтобы функция ядра частицы не приводила к развитию нефизических неустойчивостей, является положительность образа Фурье функции ядра модельной частицы (R(k) > 0). Поэтому в качестве оптимального было предложено ядро, имеющее минимальную самосилу при сохранении положительности образа Фурье функции ядра. Данные условия удовлетворяются при параметрах а0 = 0.25, ai = 0.25, ai = 0.5. Свойства QCP-ядра были проверены на решении задачи движения одной и двух частиц в отсутствии внешних полей. Все тесты показали лучший по сравнению с классическим PIC-ядром результат. Тест на эффект саморазогрева плазмы показал лучшее сохранение энергии при использовании QCP-ядра.

Раздел 3.6 посвящен дополнительной проверке свойств нового QCP-ядра на решении задали эволюции пылевой компоненты самогравитирующего про-топланетного диска. Численная модель динамики пылевой компоненты газопылевого диска или дисковой галактики сводится к гравитационной задаче N-тел. Она описывается системой, состоящей из уравнений Власова-Лиувилля и Пуассона

и

где f(t, г, v) - функция распределения частиц пыли по скоростям, а - ускорение частиц, которое определяется формулой а = —УФ, р - объемная плотность вещества (сумма плотности газа и частиц), Ф - гравитационный потенциал, 7 - гравитационная постоянная.

Так как в этой модели используется уравнение Пуассона, то вычисление силы происходит таким же образом, как и в электростатической задаче. Использование QCP ядра, при моделировании эволюции цротогшанетного диска на ряде тестов показало лучшее или аналогичное PIC сохранение полной энергии на временах порядка одного-двух оборотов диска (рис. 3). Резкий рост энергии для PIC-ядра, на рис. 3 а) объясняется недостаточно хорошим учетом энергии частиц, попавших в центр.

a)Ai~0.8 Прирост полной энергии б)Аи=4

т; ; ..........................'......—г 0.06 ;..........1—• » —

Рис. 3. Прирост полной энергии при начальной дисперсии скорости а) Аг — 0.8 и б) Ау = 4.

Четвертая глава посвящена описанию разработанного параллельного комплекса программ для моделирования взаимодействия теплого электронного пучка малой плотности с плазмой. Проводится оценка точности решения в зависимости от количества модельных частиц.

Данный комплекс программ ориентирован на исследование устойчивости и нагрева плазмы электронным пучком при параметрах, соответствующих экспериментам, проводимым на ГОЛ-3 ИЯФ СО РАН. Для изучения влияния пучковых нелинейностей на поведение неустойчивости в условиях развитой турбулентности необходимо численное моделирование, которое, с одной стороны, способно на больших временах отслеживать эволюцию возбуждаемой пучком турбулентности, а с другой, обеспечивает достаточно подробное описание кинетических эффектов, связанных с захватом пучка. Поэтому основной целью данной главы является изучение вопроса точности получаемого решения в зависимости от количества частиц. Формулируются критерии оценки точности решения и на примере решения задачи о взаимодействии электронного пучка с плазмой в разных режимах дается ответ на вопрос,

сколько частиц достаточно для качественно и количественно корректного результата.

Рассматривается задача в следующей постановке. В начальный момент времени в трехмерной области решения, имеющей форму прямоугольного параллелепипеда х 6 [О, Ь], у,г е [О, находится плазма, состоящая из электронов и ионов водорода, и пучок электронов. Задаются плотности пучка щ и электронов плазмы пе = 1 — щ, температура электронов плазмы Те и пучка Ть\ температура ионов считается нулевой 7; = 0. Начальное распределение частиц но скоростям максвелловское.

Для решения системы уравнений Максвелла используется метод Лэнгдона-Лазинского:

-от+1/2 _ тэт—1/2

— = -г<ЛлЕт, (13)

г

1 _

= Г01,Вт+1/2-Г+1/2.

Г

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

Для решения системы уравнений движения частиц используется схема с перешагиванием

»71+1/2 _ т—1/2

Рг,е Рг,е_____ _ ^ | ^г?

т+1/2 . т-1/2

V- + V-

* 1 Г> 1 Т } I

х В"

(14)

~г.е ~г,е т+1/2 /1Г\

-1—--= V,, е ■ (15)

Здесь индекс г, е - сорт частиц (ионы или электроны), г,р - координата и импульс частицы, р - у/\/1 — V2.

Рассматривается развитие отдельно взятой неустойчивой моды, поэтому длина области в направлении ж выбрана равной одной длине исследуемой волны. Граничные условия для всех функций по всем направлениям периодические. Таким образом из непрерывного спектра волновых чисел к вырезаются дискретные значения, из которых только основная мода с кх = 2-к/Ь, ку = кг = 0 попадает в область большого инкремента (электрическое поле ведет себя как Е ~ е7*). Требуется найти распределение ионов и электронов по энергиям и инкременты плазменных неустойчивостей для физических параметров, соответствующих этим трем режимам.

Основные результаты моделирования показаны на рис. 4, а также в табл. 1. Данный комплекс программ позволяет проводить моделирование эволюции

XXX

Рис. 4. Электроны пучка на фазовой плоскости (х. ух) для а - гидродинамического (п/, = 0.002, Ди = 0.007, 1р = 50), б - переходного (пь = 0.002, Ли = 0.028, 1р = 500) и е -кинетического (щ = 0.0002, Ау = 0.028,1р — 5000) режимов.

Таблица 1. Инкремент и его относительная погрешность (в %) в зависимости от количества частиц пучка в ячейке (1р). 7! - расчетное значение инкремента в гидродинамическом режиме, 72 - в переходном режиме, 73 - в кинетическом режиме. Аналитическое значение для тех же инкрементов 7? = 0.0722, ^ = 0.0232, 73° = 0.0027.

1р 71 е(%) 72 е(%) 7ч е (%)

25 0.066 9 - - - -

50 0.068 6 0.017 27 - -

250 0.069 4 0.02 14 - -

500 0.07 3 0.02 14 - -

2500 0.071 2 0.021 9 0.002 20

5000 0.071 2 0.022 5 0.002 20

10000 0.071 2 0.022 5 0.002 20

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

7° - аналитическое значение инкремента неустойчивой волны, 7^ - численное значение инкремента. Посчитано, что для воспроизведения не только качественной картины, но и количественных характеристик пучковой неустойчивости (величины инкремента неустойчивой моды с точностью не менее 20%) необходимо иметь не менее 100000 модельных частиц пучка в интервале скоростей, в котором происходит взаимодействие с неустойчивой волной.

В заключении диссертационной работы дано краткое изложение основных результатов.

ОСНОВНЫЕ РЕЗУЛЬТАТЫ РАБОТЫ

• Разработан алгоритм уменьшения счетных шумов метода частиц в ячейках. Данный алгоритм реализован при решении одномерной задачи распада разрыва плотности ионов. Он позволяет существенно уменьшить разброс в средней скорости и плотности в методе частиц в ячейках.

• Проведено исследование причины возникновения самосилы в методе частиц в ячейках. Предложен экономичный подход к уменьшению самосилы для метода частиц в ячейках на сдвинутых сетках в двумерном случае, основанный на использовании нового ядра модельной частицы. Подобраная форма ядра модельной частицы (С^СР-ядро) позволяет уменьшить величину самосилы в 4 раза по сравнению с РЮ-ядром.

• Показано, что использование нового <ЗСР-ядра приводит к увеличению времени саморазогрева модельной плазмы и к лучшему сохранению импульса и полной энергии. Использование С^СР-ядра при моделировании эволюции протопланетного диска на ряде тестов показало лучшее или аналогичное Р1С-ядру сохранение полной энергии на временах порядка одного-двух оборотов диска.

• На основе метода частиц в ячейках создан комплекс программ для моделирования взаимодействия электронного пучка с плазмой в трехмерной постановке при параметрах, соответствующих условиям лабораторных экспериментов на установке ГОЛ-3.

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

• Проведены оценки точности решения, полученного на основе созданного комплекса программ. Для различных режимов представлены оценки

минимального количества частиц, необходимого для получения величины инкремента неустойчивости с требуемой точностью. Показано, что для воспроизведения количественных характеристик пучковой неустойчивости необходимо иметь не менее 100000 модельных частиц пучка в интервале скоростей, в котором происходит взаимодействие с неустойчивой волной.

ПУБЛИКАЦИИ ПО ТЕМЕ ДИССЕРТАЦИИ

В журналах из списка ВАК:

1. Месяц, Е. А. Форма ядра частицы л проблема «самовоздействия» в методе частиц-в-ячейках / В.А. Вшивков, Е.А. Месяц // Научный вестник НГТУ. - № 1 (42). - 2011. - С. 47-56.

2. Месяц, Е.А. О выборе числа частиц в методе частиц-в-ячейках для моделирования задач физики плазмы / Е.А. Месяц, А.В. Снытников, К.В. Лотов // Вычислительные технологии. - Т. 18. - № 6. - 2013. - С. 83-97.

В других изданиях:

3. Mesyats, Е.А. A noise-reducing algorithm for particle-in-cell plasma simulation / E.A. Mesyats // Bulletin of the Novosibirsk Computing Center. Ser.: Numerical Analysis. - 2009. - № 14. -P. 21-30.

4. Месяц, Е.А. Разработка алгоритма, уменьшающего шум в методе частиц-в-ячейках / Е.А. Месяц //' Всероссийская конференция по математике и механике, посвященная 130-летию Томского государственного университета и 60-летию механико-математического факультета: сборник тезисов. - Томск: ТГУ. - 2008. - С. 138-139.

5. Месяц, Е.А. Исследование шумовых свойств метода частиц-в-ячейках / Е.А. Месяц // Дифференциальные уравнения. Функциональные пространства. Теория приближений. Международная конференция, посвященная 100-летию со дня рождения С.Л. Соболева. Тезисы докладов. - Новосибирск. -2008. - С. 529.

6. Месяц, Е.А. Об одном новом ядре для метода частиц-в-ячейках / Е.А. Месяц // Труды конференции молодых ученых ИВМиМГ СО РАН. - Новосибирск: ИВМиМГ СО РАН. - 2009. - С. 94-105.

7. Месяц, Е.А. Исследование численных свойств метода частиц-в-ячейках / В.А. Вшивков, Е.А. Месяц // XIII Всероссийская молодежная конференция-школа «Современные проблемы математического моделирования. Математическое моделирование, численные методы и комплексы программ». - Ростов-на-Дону: ЮФУ. - 2009. - С. 149-156.

8. Месяц, Е.А. Минимизация погрешности в методе частиц-в-ячейках /' Е.А. Месяц // VI Всесибирский конгресс жешцин-математиков (в день рождения Софьи Васильевны Ковалевской): Материалы Всероссийской конференции. - Красноярск: РИД СибГТУ. - 2010. - С. 284-287.

9. Месяц, Е.А. Форма ядра частицы и проблема «самовоздействия» в методе частиц-в-ячейках / Е.А. Месяц /'/ Труды XV Байкальской Всероссийской конференции «Информационные и математические технологии в науке и управлении»; ч. 1. - Иркутск: ИСЭМ СО РАН. - 2010. - С. 197-203.

10. Месяц, Е.А. Новое QCP-ядро с пониженной самосилой для метода частиц-в-ячейках и его использование при моделировании эволюции прото-планетного диска / Е.А. Месяц /,/ Труды конференции молодых ученых ИВ-МиМГ СО РАН. - Новосибирск: ИВМиМГ СО РАН. - 2011. - С. 62-74.

11. Месяц, Е.А. Новое ядро с пониженной самосилой в методе частиц-в-ячейках и его использование при моделировании эволюции протопланетного диска [Электронный ресурс] / Е.А. Месяц // Всероссийская конференция по вычислительной математике КВМ-2011. Тезисы докладов. - Режим доступа: http://www.sbras.ru/ws/show_abstract.dhtml?ru+2204-16170

12. Месяц, Е.А. Применение нового QCP-ядра с пониженным самовоздействием для метода частиц-в-ячейках при моделировании эволюции протопланетного диска / Е.А. Месяц // Российско-монгольская конференция молодых ученых по математическому моделированию, вычислительно-информационным технологиям и управлению, тезисы докладов. Иркутск (Россия) - Ханх (Монголия). - Иркутск: ИДСТУ СО РАН. - 2011. - С. 55.

13. Месяц, Е.А. Трехмерная численная модель насыщения двухиотоковой неустойчивости электронного нучка в плазме / Е.А. Месяц, A.B. Снытииков // Тезисы докладов XIX Всероссийской конференции «Теоретические основы и конструирование численных алгоритмов и решение задач математической физики», посвященной памяти К.И. Бабенко. Дюрсо. - Издательство ИПМ им. М.В. Келдыша. - 2012. - С. 73-74.

14. Месяц, Е.А. Трехмерная численная модель релаксации электронного пучка в плазме / Е.А. Месяц, A.B. Снытников // Материалы XIII Всероссийской конференции молодых ученых по математическому моделированию и информационным технологиям. - Новосибирск: ИВТ СО РАН. - 2012. -С. 28.

15. Mesyats, Е.А. Particle-in-cell simulation of kinetic instability of an electron beam in plasma / E.A. Mesyats, A.V. Snytnikov // Book of abstracts of the international conference «Mathematical modeling and computational physics». -Dubna. - 2013. - P. 129-130.

Подписано в печать 28.04.2014 Формат 60x84 1\16 Усл. печ. л. 1 Объем 16 стр. Тираж 100 экз. Заказ №78 Отпечатано Омега Принт 630090, г. Новосибирск, пр. Ак.Лаврентьева,6 email: omegap@yandex.ru

Федеральное государственное бюджетное учреждение науки Институт вычислительной математики и математической геофизики Сибирского отделения Российской академии наук

На правах рукописи 04201459272 ^¿¿ф

МЕСЯЦ Екатерина Александровна

МЕТОДЫ ОЦЕНКИ И ПОВЫШЕНИЯ ТОЧНОСТИ

РЕШЕНИЯ ЗАДАЧ ФИЗИКИ ПЛАЗМЫ МЕТОДОМ ЧАСТИЦ В ЯЧЕЙКАХ

05.13.18 - математическое моделирование, численные методы и комплексы программ

Диссертация на соискание ученой степени кандидата физико-математических наук

Научный руководитель профессор, д.ф.-м.н. В.А. Вшивков

Новосибирск - 2014

Оглавление

Введение......................................................................................5

1 Обзор методов численного решения системы уравнений Власова-Максвелла 10

1.1 Плазма, основные характеристики....................................................10

1.2 Кинетическое описание бесстолкновительной плазмы..............................12

1.3 Методы решения уравнения Власова, основанные на восстановлении функции распределения/(¿,х, V) ..............................................................13

1.4 Методы частиц..........................................................................15

1.5 Метод частиц в ячейках для числениого моделирования бесстолкновительиой плазмы ..................................................................................18

1.5.1 Уравнения движения модельных частиц, форма модельной частицы . . 18

1.5.2 Сеточные ядра..................................................................22

1.5.3 Ядра модельных частиц ......................................................23

1.5.4 Общая схема метода частиц в ячейках......................................25

1.5.5 Проблема численных шумов метода частиц в ячейках....................26

2 Алгоритм уменьшения счетных эффектов метода частиц в ячейках на примере моделирования распада разрыва плотности ионов в одномерной постановке.............................................31

2.1 Постановка задачи о распаде разрыва плотности ионов ..........................31

2.1.1 Исходные уравнения............................................................31

2.1.2 Начальные и граничные условия ............................................32

2.2 Схема метода частиц в ячейках......................................................33

2.3 Схема Лакса - Вендроффа для уравнения Власова................................36

2.3.1 Сравнение метода частиц и конечно-разностного метода..................37

2.4 Алгоритм уменьшения счетных эффектов метода частиц в ячейках (алгоритм вычитания шумовой добавки) ........................................................38

2.4.1 Алгоритм вычитания шума I..................................................39

2.4.2 Алгоритм вычитания шума II................................................39

2.4.3 Алгоритм вычитания шума III................................................40

2.5 Корректировка положения частиц....................................................41

2.6 Схемы, использованные в Алгоритме III............................................43

2.7 Результаты расчетов....................................................................44

2.7.1 Зависимость уровня шума от количества частиц..........................44

2.7.2 Сравнение схем для уравнений на и, п........................................47

2.8 Выводы..................................................................................49

3 Форма ядра частицы и проблема самовоздействия в методе частиц в ячейках .............................................50

3.1 Самосила и VSP-ядро в одномерном случае........................................51

3.2 Самосила в двумерном случае........................................................52

3.2.1 Самосила в двумерном случае, PIC-ядро....................................57

3.2.2 Самосила в двумерном случае, QCPl-ядро..................................57

3.2.3 Самосила в двумерном случае, С^СР2-ядро..................................58

3.2.4 Потенциал поля одиночного заряда..........................................58

3.2.5 Сравнение ядер PIC, QCP1, QCP2 ..........................................59

3.3 Новое ядро..............................................................................60

3.3.1 Новое ядро, тесты..............................................................61

3.4 Ядро QCP ..............................................................................64

3.4.1 Фурье-образ функции ядра частицы ........................................66

3.4.2 Выбор оптимальных параметров..............................................68

3.5 Саморазогрев модельной плазмы....................................................72

3.6 Моделирование эволюции протоплапетного диска с QCP-ядром..................74

3.7 Выводы..................................................................................77

4 Число модельных частиц и точность на примере задачи моделирования кинетической неустойчивости теплого электронного пучка малой плотности в плазме.......................................78

4.1 Основные уравнения....................................................................80

4.2 Методы и алгоритмы решения ......................................................81

4.3 Вычисление инкремента неустойчивости............................................84

4.4 Результаты расчетов..................................................................85

4.5 Фазовые плоскости ....................................................................91

4.6 Выводы..................................................................................93

Заключение Литература

Введение

Метод частиц в ячейках, возникший еще в середине прошлого столетия, на сегодняшний день широко применяется при моделировании поведения плазмы. Его также используют и при расчете динамических процессов в других средах (жидкости, газы, сплошные среды и др.). Но именно в решении задач физики плазмы он получил наиболее широкое распространение. Как самостоятельный раздел вычислительная физика плазмы сформировалась одновременно с развитием вычислительной техники во второй половине XX века. Это привело к выделению математического моделирования в отдельное направление исследовательской работы. Ввиду больших возможностей диагностики моделируемых в численных экспериментах кинетических процессов и относительно невысоких экономических затрат на их проведение (по сравнению со строительством установок для проведения реальных экспериментов), численные методы играют все более важную роль в решении задач физики плазмы [16].

В диссертации рассматривается приближение бесстолкновительиой плазмы и метод частиц в ячейках применительно к системам уравнений Власова-Пуассона и Власова-Максвелла. Для многих нелинейных процессов найти решение этой системы возможно только численно, а не аналитически. Основные методы численного решения кинетических уравнений бесстолкновительиой плазмы были разработаны к началу 70-х годов [17]. Это эйлеровы алгоритмы, метод водяного мешка, метод преобразований и метод частиц.

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

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

мепенне именно для моделирования таких сложных явлений, как неустойчивости в физике плазмы. Главным преимуществом метода частиц в ячейках по сравнению с эйлеровыми методами является простота реализации и экономичность, так как в методах частиц не вычисляется напрямую функция распределения /(¿,х, и). Развитая теория, простота реализации данного метода (в том числе и на параллельных ЭВМ), наглядность получаемых распределений частиц, большое число решенных на его основе задач в разнообразных областях физики плазмы составляют сильные стороны этого метода.

Актуальность работы. Основным недостатком метода частиц в ячейках всегда являлось наличие в решении нефизических эффектов, так называемых «численных шумов». Эта проблема очень сложна, так как причин возникновения шумов много [18]. Влияние различных факторов на решение трудно разделить, особенно при моделировании таких неустойчивых сред, как плазма. Шумовые гармоники взаимодействуют друг с другом, накладываются на гармоники неустойчивых решений, что может приводить к развитию численных иеустой-чивостей.

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

Цель диссертационной работы заключается в исследовании свойств метода частиц в ячейках и разработке общих алгоритмов уменьшения счетных шумов данного метода. Для достижения этой цели было решено

• создать алгоритм подавления шума в методе частиц в ячейках,

• разработать методику выбора оптимальной формы ядра, которая обеспечивает по сравнению со стандартным Р1С-ядром уменьшение самосилы,

• создать комплекс программ для моделирования взаимодействия электронного пучка с плазмой в трехмерной геометрии и провести оценки точности получаемого решения.

Научная новизна работы заключается в следующем:

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

2. Проведено исследование причины возникновения самосилы для метода частиц в ячейках на смещенных сетках. Для уменьшения этой ошибки разработан экономичный подход, основанный на использовании нового ядра модельной частицы с меньшей, по сравнению с Р1С-ядром, самосилой.

3. Показано, что использование нового ядра приводит к увеличению времени саморазогрева модельной плазмы и, как следствие, к лучшему сохранению импульса и полной энергии.

4. На основе метода частиц в ячейках создан комплекс программ для моделирования взаимодействия электронного пучка с плазмой в трехмерной постановке при параметрах, соответствующих условиям лабораторных экспериментов на установке ГОЛ-3, который позволяет проводить вычисления как в гидродинамическом, так и в кинетическом режимах развития пучковой неустойчивости.

5. Проведены оценки величины погрешности решения, получаемого данным комплексом программ. Для различных режимов приводятся оценки достаточного количества частиц.

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

Представленные в диссертации исследования проводились в рамках Интеграционных проектов СО РАН №113, №40, Научной школы Годунова - 4292.2008.1 и по проектам, поддержанным Российским фондом фундаментальных исследований (№ 08-01-00615, 11-01-00249, 11-01-00178).

Достоверность результатов. Все численные алгоритмы проверялись по отдельности на тестовых расчетах. Программные комплексы также проходили тестирование, проводилось сравнение результатов моделирования с имеющимися аналитическими решениями или

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

• алгоритм вычитания шумовой добавки,

• алгоритм поиска оптимальной формы ядра и новое QCP-ядро, минимизирующее самовоздействие при вычислении полей на сдвинутых друг относительно друга сетках,

• созданный комплекс программ для моделирования взаимодействия электронного пучка с плазмой в трехмерной геометрии,

• оценки точности вычисления инкремента неустойчивости для созданного комплекса программ.

Апробация работы. Основные положения диссертационной работы докладывались и обсуждались на семинарах «Математическое моделирование больших задач» под руководством д.ф.-м.н. Вшивкова В.А. (ИВМиМГ СО РАН), на семинаре «Математическое и архитектурное обеспечение параллельных вычислений» под руководством д.т.н. Малышкина В.И. (ИВМиМГ СО РАН, декабрь 2009), на семинаре «Задачи механики и математической физики» под руководством д.ф.-м.н. Медведева C.B. (ИВТ СО РАН, июнь 2013), на объединенном семинаре ИВМиМГ СО РАН и кафедры вычислительной математики НГУ под руководством д.ф.-м.н. Ильина В.П. (ИВМиМГ СО РАН, октябрь 2013), па семинаре «Законы сохранения и инварианты» под руководством д.ф.-м.н. Медведева C.B. (ИВТ СО РАН, февраль 2014), на Всероссийской конференции по математике и механике, посвященной 130-летшо ТГУ (2008, Томск), на Международной конференции, посвященной 100-летию со дня рождения C.JI. Соболева (2008, Новосибирск), па Конференции молодых ученых ИВМиМГ СО РАН (2009, 2011, Новосибирск), на XIII Всероссийской молодежной конференции-школе «Современные проблемы математического моделирования. Математическое моделирование, численные методы и комплексы программ» (2009, Дюрсо), на VI Всесибирском конгрессе женщин-математиков (2010, Красноярск), на XV Байкальской Всероссийской конференции «Информационные и математические технологии в науке и управлении» (2010, Иркутск-Байкал), на Всероссийской конференции по вычислительной математике КВМ-2011 (Новосибирск), на Российско-монгольской конференции молодых ученых по математическому моделированию, вычислительно-информационным технологиям и управлению (2011, Иркутск (Россия) - Ханх (Монголия)), на XIX Всероссийской конференции «Теоретические основы и конструирование численных алгоритмов и решение задач математической физики» (2012,

Дюрсо), на XIII Всероссийской конференции молодых ученых по математическому моделированию и информационным технологиям (2012, Новосибирск), на Международной конференции «Mathematical modeling and computational physics» (2013, Дубна).

Основные результаты опубликованы в 15 работах, из которых 2 в журналах, рекомендованных ВАК [1]- [15].

Личный вклад соискателя заключается в обсуждении постановок задач, разработке численных алгоритмов и методов решения, создании и тестировании программ, проведении расчетов и анализе полученных результатов. Все выносимые па защиту результаты принадлежат лично автору. Представление результатов совместных исследований и разработок согласовано с соавторами.

Структура диссертационной работы. Диссертация состоит из введения, четырех глав и заключения. В первой главе содержится обзор методов численного решения системы уравнений Власова-Пуассона и Власова-Максвелла, общая схема метода частиц в ячейках и обзор литературы по проблеме шумов в методе частиц. Вторая глава посвящена описанию алгоритма вычитания шумовой добавки. Разработка данного алгоритма проводится па примере одномерной задачи о распаде разрыва плотности ионов в бесстолкновительной полностью ионизованной плазме. Третья глава посвящена проблеме самовоздействия в методе частиц в ячейках. Разрабатывается новое QCP ядро, позволяющее минимизировать самовоздействие в двумерной постановке. В четвертой главе приводится описание созданного параллельного пакета программ для моделирования взаимодействия электрпного пучка с плазмой в трехмерной постановке. Проводится оценка точности получаемого решения в зависимости от количества модельных частиц. Показано, что данный комплекс программ охватывает широкий диапазон физических параметров и позволяет проводить моделирование пучковой неустойчивости как в гидродинамическом, так и в кинетическом режимах развития неустойчивости.

Глава 1

Обзор методов численного решения системы уравнений Власова-Максвелла

1.1 Плазма, основные характеристики

Подавляющая часть вещества во Вселенной находится в состоянии плазмы. Плазма представляет собой полностью или частично ионизованный квазинейтральный газ. В полностью ионизованной плазме взаимодействие между частицами является в основном электромагнитным, если плазма является еще и нерелятивистской, то взаимодействие можно считать электростатическим. Поля, действующие па плазму, делятся на внешние и внутренние (создаваемые зарядами и токами в сам