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

Введение

Глава I. Метод решения задач эволюции магнитного поля в рамках ЭМГ-модели.

1.1 Математическая модель. Постановка задачи.

1.2 Численный метод расчета уравнений эволюции магнитного поля.

1.2.1 Основные масштабы величин и обезразмеривание.

1.2.2 Расчет эволюции напряженности магнитного поля.

1.2.3 Разностная схема для второго этапа расчета

1.2.4 Вычисление изменения электронной температуры.

1.2.5 Разностная схема для третьего этапа расчета.

1.3 Результаты тестовых расчетов

1.3.1 Задача на установление распределения ДДг)

1.3.2 Задача на установление распределения B^z)

1.3.3 Задачи с использованием автомодельных решений.

1.4 Алгоритм перестроения подвижной сетки.

Глава II. Результаты расчетов задач движения плазмы в размыкателях тока.

2.1 Выбор граничных условий для напряженности магнитного поля.

2.2 Результаты расчетов по ЭМГ-модели.

2.2.1 Задача с фиксированными границами, отрицательная полярность импульса.

2.2.2 Задача с фиксированными границами, положительная полярность импульса.

2.2.3 Задача с граничными условиями II рода для напряженности магнитного поля на электродах

2.2.4 Влияние начальной температуры на динамику плазменного размыкателя тока.

Актуальность темы исследований.

С семидесятых годов XX века ведутся работы по изучению поведения вещества под воздействием высокоэнергетического излучения. Программы исследования управляемого термоядерного синтеза повлекли необходимость создания установок по преобразованию запасенной в накопителях энергии для достижения высоких температур. Одной из таких установок является генератор импульса напряжения (ГИН). Электрическая энергия, запасенная в конденсаторе, преобразуется в энергию электрического тока, проходящего через подсоединенный проводящий контур. При этом возникла проблема передачи и коммутации потока энергии большой плотности, что приводит к существенным потерям энергии во вспомогательных элементах установки. Было предложено применение плазменного размыкателя (прерывателя) тока (ППТ, plasma opening switch, POS) для коммутации электрического тока. В дальнейшем решалась задача получения большой мощности преобразователя, для чего его размеры уменьшались.

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

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

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

Цели настоящей работы.

Основными целями представленной работы являлись:

1. построение эффективного численного метода расчета динамики проникновения магнитного поля в материал плазмы в плазменном размыкателе тока в рамках ЭМГ-модели;

2. анализ существующих методов построения разностной сетки для повышения качества расчетов на подвижной лагранжевой сетке;

3. разработка программного комплекса для численного моделирования динамики плазмы в условиях, соответствующих экспериментальным установкам "Стенд-300" и РС20 (РНЦ КИ) для совместного использования с существующими программными комплексами;

4. проведение исследования динамики плазменного размыкателя тока при параметрах модели, соответствующих экспериментам в РНЦ К И на установке "Стенд-300".

Научная новизна.

1. Предложен метод эффективного численного расчета динамики проникновения магнитного поля в материал плазмы в плазменном размыкателе тока в рамках ЭМГ-м од ели;

2. Создан комплекс вычислительных программ для численного моделирования динамики плазмы в рамках ЭМГ-модели, расширивший область применения существующих программных комплексов; peaлизован метод перестроения сетки, который позволяет проводить расчеты на более длительных интервалах времени;

3. Проведены численные эксперименты по исследованию динамики плазменного размыкателя тока при параметрах модели, соответствующих экспериментам в РНЦ КИ на установке "Стенд-300".

Научная и практическая ценность.

Разработанный программный комплекс позволяет эффективно проводить численное моделирование динамики плазмы с учетом эффектов электронной магнитной гидродинамики в широком диапазоне параметров математической модели плазменного размыкателя тока. В сотрудничестве с Отделением прикладной физики Российского научного центра "Курчатовский институт" (ОПФ РНЦ КИ) было проведено более сотни вычислительных экспериментов по моделированию движения материала плазмы в установках "Стенд-300" и РС20.

Апробация работы.

Результаты, полученные в диссертации, докладывались на конференциях МФТИ (2001 и 2002 г.г.), на Девятой международной конференции "Математика. Компьютер. Образование." (Дубна, 2002 г.), 19th IAEA Fusion Energy Conference (Lyon, France, 2002 г.), на конференции BEAMS'2002 (Albuquerque, New Mexico USA, 2002 г.). Результаты также докладывались на семинарах Отделения прикладной физики Российского научного центра "Курчатовский институт" (2001-2002 г.г.) и научных семинарах кафедры вычислительной математики МФТИ (2002-2002 г.г.).

Публикации.

Основные результаты проведенной работы опубликованы в:

1. тезисах конференции Девятой международной конференции "Математика. Компьютер. Образование." (Дубна, 2002 г.) [1];

2. трудах конференции Девятой международной конференции "Математика. Компьютер. Образование." (Дубна, 2002 г.) [2];

3. 14th International Conference on High-Power Particle Beams and 5th International Conference on Dense Z-Pinches (Albuquerque, New-Mexico USA, 2002) [3j;

4. 19th IAEA Fusion Energy Conference (Lyon, France, 2002) [4].

Обзор существующих работ по теме исследований.

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

В обзоре А. С. Кингсепа и др. [5] подробно описывается модель электронной магнитной гидродинамики (ЭМГ) плазмы и рассматриваются присущие этой модели эффекты. По сравнению с уравнениями модели магнитной гидродинамики (МГД) в уравнениях ЭМГ учитывается более сильное влияние электронных течений на эволюцию магнитного поля, температуры и других параметров, описывающих состояние компонент плазмы. В границах применимости ЭМГ-приближения [6] математическая модель содержит уравнения движения по крайней мере двух компонент плазмы - электронов и ионов. При этом описание ионных течений оказывается близким к идеальной гидродинамике.

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

В обзоре [8] исследовались основные уравнения ЭМГ-приближения в условиях идеальной проводимости вещества и отсутствия инерции электронов. Эта модель носит название HMHD (Hall magnetohydrodynamics), по сравнению с МГД уравнение эволюции магнитного поля дополнено учетом эффекта Холла: дВ Jjxb\ = j.

Уравнение записано в безразмерном виде. Проводилось моделирование динамики плазмозаполненного диода с учетом аксиально-симметричного неоднородного начального распределения плотности электронов плазмы (Vne =/= 0, пе - концентрация электронов). Показано более быстрое проникновение магнитного поля от катода к аноду (в направлении В х Vne) при учете эффекта Холла (ЭМГ-модель) по сравнению со случаем идеальной МГД.

В работе [9] изучалась динамика магнитного поля на границе двух сред, различающихся по плотности электронов и проводимости. Проводилось моделирование на двумерной эйлеровой сетке в декартовых координатах. Применялось уравнение для изменения магнитного поля в рамках ЭМГ с учетом эффекта Холла и конечной проводимости при отсутствии градиента электронного давления. Показано влияние эффекта Холла на распространение магнитного поля вдоль границы.

Дальнейшее развитие это исследование получило при численном моделировании динамики плазмы Z-пинча [10]. Применялись двумерные уравнения двухжидкостной МГД с учетом эффекта Холла. Численное моделирование проводилось с использованием явного метода на эйлеровой сетке в цилиндрических координатах. В работе описаны эффекты более быстрого распространения магнитного поля вдоль катода и анода из-за влияния эффекта Холла, получены оценки скорости скольжения плазменной оболочки вдоль анода.

Изучение динамики плазменного размыкателя тока в режиме проводимости в рамках двухжидкостной МГД модели с учетом эффекта Холла проводилось в работе [11]. Были получены критерии режимов работы размыкателя с преобладанием эффектов МГД, Холла и смешанного типа. Численное моделирование осуществлялось на эйлеровой сетке с прямоугольными ячейками в цилиндрических координатах. Значение напряженности магнитного поля на границе плазмы со стороны генератора тока

2 J задавалось в виде = (в единицах системы СГС), где /о - ток генератора, г - расстояние от оси размыкателя; при этом полагалось, что ток генератора не зависит от времени, то есть к началу расчетов ток через размыкатель уже достиг этого значения. Граничные условия для магнитного поля на электродах и на противоположной генератору области задавались в виде равенства нулю тангенциальной компоненты электрического поля. В работе описано проникновение магнитного поля вдоль катода из-за эффекта Холла. Показано влияние комбинации эффекта Холла и МГД эффектов на формирование и эволюцию области малой плотности в при-электродном пространстве.

Необходимо также отметить работу [12], которая представляет один из существующих на кафедре Вычислительной математики МФТИ программных комплексов для численного моделирования динамики плазмы в рамках МГД-модели. В результате представленной работы этот комплекс был модернизирован с учетом поставленных целей исследований.

В работе [12] описывалась реализация эксперимента, в котором изучалась возможность применения метода совместного численного решения уравнений одножидкостной двухтемпературной МГД модели и кинетической модели, описывающей динамику состояний ионизации вещества плазмы. Реализация комбинированото метода проверялась при двумерном моделировании Z-пинча средней плотности на сетках от 60 х 134 до 100 х 206 лагранжевых ячеек. Показано хорошее соответствие результатам физических экспериментов на импульсных плазменных установках. Применение методов параллельных вычислений позволило существенно ускорить расчеты экспериментальных вариантов.

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

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

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

Реализация комбинированных подходов типа ALE (Arbitrary Lagrangian-Eulerian) [14],[15] предполагает проведение на каждом временном шаге лагранжевого этапа расчета, при котором узлы сетки несколько изменяют свое положение. Далее может быть произведен выбор новой сетки и пересчет на нее с учетом конвективных потоков при изменении границ ячеек. В случае, когда происходит возврат к исходной сетке, фактически реализуется эйлеровый подход. Лагранжевый метод расчетов реализуется без изменения сетки, полученой в результате лагранжевого этапа расчета. В общем случае пространственное расположение узлов новой сетки может не быть связано ни с движением среды, ни с исходной эйлеровой сеткой, но выбираться из других соображений, например, сохранения равномерности сетки или детализации определенной области. Методы типа ALE позволяют избежать катастрофических деформаций ячеек сетки, но применение пересчетов может приводить к потере информации о деталях моделируемых процессов, в частности о течении среды. Алгоритмы учета конвективного переноса могут быть реализованы с аппроксимацией первого, второго и более высоких порядков.

В работе [16] описана реализация подхода ALE при моделировании среды с многими контактными границами (MMALE, Multi-material ALE). Сравниваются различные алгоритмы пересчета на новые сетки. Проведен анализ применения в некоторых приложениях.

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

Для улучшения качества двумерной расчетной сетки наибольшей популярностью [18] пользуются эвристические стратегии на основе диффузионного сглаживания [19]. В простейшем случае рекурсивно производится перемещение узлов в усредненное по координатам положение относительно связанных с ними узлов. Хорошие результаты достигаются при отсутствии невыпуклых ячеек в сглаживаемой области. Но вблизи невыпуклой ячейки сетки алгоритм может приводить и к более сильному искажению или выворачиванию ячек. Положительным качеством этих методов является линейная зависимость объема вычислений от числа узлов при наличии простой топологии. В сложных случаях требуются дополнительные затраты на определение связей между узлами.

Более перспективными представляются методы, в которых вводится мера деформации ячеек и производится минимизация величины их искажения в заданной метрике. Метрика на основе взвешенного значения якобианов в узлах прямоугольной ячейки предложена в [20]. При этом используются различные методы оптимизации, например, метод градиентного спуска.

С учетом такого подхода в [18] предложена комбинация сглаживания на основе решения уравнения Лапласа для координат узлов сетки и локальной оптимизации положения отдельных узлов с использованием специально выбранных метрик для треугольных и четырехугольных ячеек. Комбинирование позволило избежать применения дорогостоящей процедуры оптимизации при малых в метрике искажениях геометрии ячеек. В результате в тестовых расчетах оптимизация применялась для одного процента от всех узлов.

Для сглаживания регулярной сетки из четырехугольных ячеек в [21] применен генетический алгоритм. Генетические алгоритмы [22],[23] являются алгоритмами поиска решения на основе механизмов естественного отбора. Чтобы найти оптимальное решение задачи, алгоритм начинает работу с заданным набором "хромосом" (предполагаемых решений) и развивает различные, но дающие наилучшие результаты, наборы "хромосом" (наборы решений) в последовательности "поколений" (итераций). В каждом "поколении" критерий меры пригодности (функция качества решения) определяет годные "хромосомы", и некоторые из них выбираются для воспроизведения в следующем поколении. Число воспроизведенных копий полагается прямо пропорциональным мере качества, таким образом происходит отбор наиболее пригодных решений и исключение остальных. В работе [21] "скрещивание" моделируется перестановками частей "хромосом", "мутации" реализованы методом случайного варьирования компонент "хромосом" в пространстве искомых решений. Приведено сравнение подходов на основе экспертных систем, генетических алгоритмов и нейронных сетей для решения поставленной задачи. Для построения функции качества решения использована мера деформации ячейки сетки на основе [20]. Отмечается хорошее качество получаемого решения для исходных сеток, содержащих невыпуклые ячейки. Недостатком метода является сильная зависимость качества результата и времени получения результата от выбора алгоритмов отбора, скрещивания и мутаций, а также большие затраты на промежуточное хранение данных о "поколениях" наборов решений.

В работе [24] предложен прямой метод сглаживания сетки. Система уравнений для координат каждого узла строится с учетом топологии и "жесткости" связей между узлами, и затем разрешается с учетом граничных условий в виде координат жестко закрепленных узлов на границе и внутри сетки. При этом новые координаты всех узлов сетки вычислиются одновременно, исключая итерационные процессы [18] с проходом по каждому узлу. Отмечается зависимость решения системы только от координат фиксированных узлов и заданных начальных связей между узлами сетки.

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

При соответствующем выборе меры деформации ячейки сетки появляется возможность устранить их выворачивание. В работе [25] предлагается метод введения "барьерной" компоненты в используемый при построении сетки функционал для исключения нежелательной деформации ячеек.

В работе [26] предложен метод построения сетки, основанный на следующей идее. Лагранжевая сетка, полученная после очередного шага расчета, отражает как физические перемещения среды, так и нефизичные искажения. Поэтому следует ограничить оптимизацию расположения узлов сетки в расчетной области, чтобы по возможности сохранить полезную информацию о движении среды. Если предположить, что нефизичные искажения ячеек имеют более короткую длину волны, то их можно отделить от физического перемещения методом усреднения по малому количеству соседних ячеек. Аналогичный подход применялся ранее в [27] при введении сглаживающего оператора на основе линейного уравнения теплопроводности для коррекции значений сеточной функции, .вычисляемых с помощью немонотонных разностных схем (метод Л. А. Чудова).

На первом шаге метода [26] для каждого узла сетки выбиралось перемещение с тем, чтобы минимизировать функционал, в случае четырехугольных ячеек сетки задаваемый следующим образом: где А обозначает площадь треугольника, и 1г) обозначают длины соответствующих сторон треугольника, индексы к и / нумеруют прилежащие к узлу (i, j) треугольники, получаемые разрезанием соответствующей ячейки сетки по диагонали. Этот функционал тесно связан с методом Уинслоу [28], в котором используется решение стационарного уравнения теплопереноса с соответствующими граничными условиями для задания гладкого отображения координат узлов двумерной прямоугольной равномерной сетки в логическом пространстве в координаты сетки в физическом пространстве. Значение функционала минимизировалось локально для каждого узла сетки при фиксированном исходном положении остальных узлов. Таким образом были получены образцовые матрицы якобианов, локально определяющие желаемые свойства сетки. Они не могут быть использованы для построения непрерывного однозначного отображения на всей сетке из-за возможнык конфликтов локальных целей оптимизации. Поэтому на втором шаге строилось отображение, наиболее близкое к совокупности образцовых отображений в смысле наименьших квадратов отклонения якобианов в каждом узле новой сетки от образцовых. Для этого производилась минимизация дискретного функционала где J с обозначает якобиан для узла п в ячейке с, а шаблон St (с) включает все узлы этой ячейки. Областью определения функционала является множество всех узлов сетки, включая граничные. При этом для отслеживания

0.1)

0.2)

16 границ между различными средами для граничных узлов возможно параметрическое задание направления перемещения. Функционалы 0.1 и 0.2 обладают свойством "барьера" [25] в том смысле, что при оптимизации исключается уменьшение до нуля или появление отрицательной площади каждого из треугольников, на которые разбивается сетка, при условии отсутствия таких треугольников в начальной сетке. В 0.2 роль "барьера" играет значение максимальной нормы якобиана получаемого отображения. Таким образом, исключается вырождение и "выворачивание" ячеек сетки.

Заключение

Сформулируем основные выводы работы.

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

2. Проанализированы алгоритмы перестройки расчетных лагранжевых сеток. Реализован алгоритм перестройки сетки, основанный на минимизации соответствующего функционала.

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

4. Показано, что проникновение магнитного поля в плазму за счет ЭМГ-эффектов может сопровождаться возникновением петель тока в приэлектродной области, что качественно отличается от характера проникновения поля в МГД-приближении.

1. А. С. Кингсеп, К. В. Чукбар, В. В. Яньков. Электронная магнитная гидродинамика. / В кн.: Вопросы теории плазмы. Под ред. Б. Б. Кадомцева. М.: Атомиздат, 1987. - Вып. 16. - с. 243-291.

2. А. С. Кингсеп. Введение в нелинейную физику плазмы. М.: Изд-во МФТИ, 1996. - 208 с.

3. Т. Б. Б. Кадомцев. Коллективные явления в плазме. М.: Наука, 1988.

4. В. В. Вихрев, О. 3. Забайдулин. Проникновение магнитного поля в плазму вдоль границы двух сред из-за эффекта Холла. // Физика плазмы. 1994 т. 20. - N 11. - с. 968-972.

5. В. В. Вихрев, О. 3. Забайдулин, А. Р. Терентьев. Моделирование динамики плазменной оболочки Z-пинча вблизи электродов. // Физика плазмы. 1995. -т. 21. - N 1. - с. 23-30.

6. О. Z. Zabaidullin. Numerical study of a plasma opening switch conduction phase. // Physics of Plasmas. 2000. - Y. 7. - N 4. - P. 1321-1330.

7. А. С. Кингсеп, В. E. Карпов, А. И. Лобанов, И. Марон, А. А. Старобинец, В. И. Фишер. Численное моделирование динамики медленного Z-пинча. // Физика плазмы. 2002. т. 28. - N 4. - с. 319-328.

8. JI. И. Седов. Механика сплошной среды. Т. 1. М.: Наука, 1983. - 528 с.

9. С. W. Hirt, A. A. Amsden, and J. L. Cook. An Arbitrary Lagrangian-Eulerian Computing Method for All Flow Speeds. // Journal of Computational Physics. -1974. Y. 14. -N 3. - P. 227-253.

10. С. Херт. Произвольный лагранжево-эйлеровый численный метод. / В кн.: Численные методы в механике жидкостей. М: Мир, 1973. - с. 156-164.

11. Freitag, Jones, and Plassmann. A Parallel Algorithm for Mesh Smoothing. // SIAM Journal on Scientific Computing. 1999. - V. 20.

12. D. A. Field. Laplacian Smoothing and Delaunay Triangulation. // Communications in Applied Numerical Methods. 1988. - V. 4. - P. 709-712.

13. M. Holder and J. Richardson. Genetic Algorithms, Another Tool for Quad Mesh Optimization? / In proceedings of the Seventh International Meshing Roundtable.- Sandia National Lab, 1998. P. 497-504.

14. J. H. Holland. Adaptation in Natural and Artificial Systems. 3rd printing, MIT Press, Cambridge, 1994 (1975).

15. D.E. Goldberg. Genetic Algorithms in Search, Optimization and Machine Learning. Addison-Wesley, Reading, MA, 1988.

16. Bala Balendran. A direct smoothing method for surface meshes. / In proceedings of 8th International Meshing Roundtable. South Lake Tahoe, CA, U.S.A., 1999.- P. 189-193.

17. A. Charakhch'yan and S. Ivanenko. A variational form of the Winslow grid generator. // Journal of Computational Physics. 1997. - N 136. - P. 385-398.

18. P. Knupp, L. Margolin, and M. Shashkov. Reference Jacobian Optimization-Based Rezone Strategies for Arbitrary Lagrangian Eulerian Methods. // Journal of Computational Physics. 2002. - N 176. - P. 93-128.

19. В. И. Полежаев, В. M. Пасконов, JI. А. Чудов. Численное моделирование процессов тепло- и массообмена. М: Наука, 1984 - 274 с.

20. A. Winslow. Numerical solution of the quasilinear Poisson equations in a nonuniform triangle mesh. // Journal of Computational Physics. 1966. N 1. - P. 149-172.

21. А. С. Кингсеп, В. И. Косарев, А. И. Лобанов, А. А. Севастьянов. Численное моделирование сжатия плазмы легким лайнером. / / Физика плазмы. 1997. - т. 23. - N 10. - с. 953-959.

22. Ю. Г. Калинин, А. С. Кингсеп, В. И. Косарев, А. И. Лобанов. К расчету осесимметричного сжатия легкого лайнера. // Матаматическое моделирование.- 2000. т. 12. - N 11. - с. 63-80.

23. А. А. Самарский, Ю. П. Попов. Разностные методы решения задач газовой динамики. М.: Наука, 1980. - 392 с.

24. В. В. Демченко, В. И. Косарев. Численные и аналитические решения двумерных задач МГД и высокотемпературной плазмы. / VI Всесоюзный съезд по теоретической и прикладной механике. Аннотации докладов. Ташкент, 1986. -с. 243.