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

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

Байдин Денис Федорович

ЭКОНОМИЧНАЯ ТРЕХМЕРНАЯ МЕТОДИКА РАСЧЕТА КРИТИЧЕСКИХ ПАРАМЕТРОВ АКТИВНОЙ ЗОНЫ БЫСТРОГО РЕАКТОРА С ЕСТЕСТВЕННОЙ БЕЗОПАСНОСТЬЮ

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

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

2 2 ЛЕК ?011

Москва-2011

005005962

Работа выполнена на кафедре вычислительной математики Московского физико-технического института (государственного университета)

Научный руководитель:

Доктор физико-математических наук АРИСТОВА Елена Николаевна

Официальные оппоненты:

Доктор физико-математических наук, профессор ТРОЩИЕВ Виталий Ефимович

Кандидат физико-математических наук БАСС Леонид Петрович

Ведущая организация:

РНЦ Курчатовский институт

« 4Ъ »

с£'р

. 2011 года в У час. на

Защита состоится

заседании диссертационного совета Д 212.156.05 при Московском физико-техническом институте (государственном университете) по адресу: 141700, г. Долгопрудный Московской обл., Институтский пер. д.9, ауд. 903 КПМ.

Автореферат разослан «.

Ученый секретарь диссертационного совета ( /

.2011 г.

А?

О.С. Федько

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

Актуальность темы

Задача создания реакторов нового типа, обеспечивающих высокий уровень безопасности, становится все более актуальной в связи с рядом аварий на АЭС и быстрым расходованием ограниченных запасов урана (и-235). Быстрые реакторы работают на смеси и-238 и Ри-239, воспроизводя в ходе своей работы плутоний. Реакторы такого типа, способные работать в саморегулируемых нейтронно-ядерных режимах, являются хорошими кандидатами для решения обеих этих проблем.

Построение сложных реакторных систем не обходится без их предварительного анализа посредством физического (натурного) и математического моделирования; второе из них является более вариативным и менее затратным.

Основой математической модели быстрых реакторов является многогрупповое уравнение переноса нейтронов. Решение уравнения переноса является одной из наиболее трудоемких частей моделирования реакторных задач. Это связано с большой размерностью задачи, т.е. большим количеством переменных, от которых зависит функция распределения (фазовые переменные и время). Исследование саморегулируемых нейтронно-ядерных процессов в активных зонах быстрых реакторов приводит к необходимости динамического моделирования процессов переноса нейтронов совместно с выгоранием, реакторной кинетикой и управлением. Поэтому необходимо разрабатывать высокоэкономичные методики решения уравнения переноса, чтобы при использовании современных суперкомпьютеров можно было проводить динамические расчеты за обозримое время. Большинство существующих на данный момент методов решения уравнения переноса имеют определенные недостатки: либо они вычислительно дороги, либо уравнение переноса заменяется более простой системой уравнений диффузии нейтронов. Создание эффективных численных методов для трехмерного моделирования переноса

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

Трудности решения уравнения переноса помимо большой размерности связаны с несколькими факторами:

-построение аппроксимации дифференциального оператора в левой части уравнения связано с дилеммой точность-монотонность;

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

-во многих физических приложениях уравнение переноса нужно решать совместно с другими уравнениями, такими как уравнения выгорания и реакторной кинетики в задачах переноса нейтронов; объединенная система уравнений может обладать сильной нелинейностью, хотя само уравнение переноса линейно относительно своих переменных;

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

Цель и задачи исследования

Настоящая работа посвящена разработке эффективных численных методов решения уравнения переноса с квазидиффузией. Метод квазидиффузии (В.Я. Гольдин. Квазидиффузионный метод решения кинетического уравнения // ЖВМ и МФ, 1964, т.4, № 6, с. 1078-1087) заключается в постепенном понижении размерности используемых уравнений. На первом этапе происходит усреднение уравнения переноса по угловым переменным, в результате которого получается многогрупповая система уравнений квазидиффузии. На втором - усреднение по

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

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

Методы исследования

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

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

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

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

На защиту выносятся следующие положения:

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

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

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

• На основе математической модели многогруппового уравнения переноса с квазидиффузией создана программа для расчета критических параметров сборки активной зоны реактора типа БН-800, способного работать в саморегулируемом нейтронно-ядерном режиме; с ее помощью получены численные результаты для параметров критической сборки.

Теоретическая и практическая ценность

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

Основные публикации

По теме диссертации опубликовано 7 работ, в том числе две [3,7] - в изданиях из списка, рекомендованного ВАК РФ.

Достоверность результатов

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

Апробация результатов диссертации

Результаты исследований, приведенных в диссертационной работе, были представлены и обсуждались на Всероссийских и Международных конференциях, семинарах ведущих институтов:

• 48-я научная конференция МФТИ, ноябрь 2005, Москва;

• XIII Международная конференция "Математика. Компьютер.

Образование", январь 2006, Дубна;

• 49-я научная конференция МФТИ, ноябрь 2006, Москва;

• The 20th International Conference on Transport Theory, July 2007, Obninsk, Russia;

• XV Международная конференция "Математика. Компьютер. Образование", январь 2008, Дубна;

• 52-я научная конференция МФТИ, ноябрь 2009, Москва;

• 53-я научная конференция МФТИ, ноябрь 2010, Москва;

• XVIII Международная конференция "Математика. Компьютер. Образование", январь 2011, Пущино;

• Семинары кафедры вычислительной математики Московского физико-технического института, научные семинары ФУПМ, 2009 - 2011, Москва.

Личный вклад соискателя

Ранее коллективом, в котором работает автор, была построена иерархия моделей более низкой размерности (ID, 1,5D, 2D); в создании последней из них автор принял активное участие. Частично методика расчета критических параметров была перенесена в ЗО-модель. Все научные результаты по методам расчета всей системы уравнений переноса с квазидиффузией в трехмерной геометрии, изложенные в диссертации, получены лично автором.

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

Связь с научными проектами

Работа выполнялась в рамках проектов Российского фонда фундаментальных исследований:

• фант 03-01000443а;

• грант 11-01-00389а.

Структура и объем диссертации

Диссертация состоит из введения, пяти глав, заключения и списка литературы. Общий объем диссертации составляет 110 страниц, список использованных источников содержит 62 наименования.

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

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

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

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

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

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

опосредованным через тензор квазидиффузии. Итерационный процесс решения системы уравнений с такой иерархической структурой быстро сходится. Таким образом, решение задач на собственные значения (задача нахождения критических параметров) для уравнения переноса сводится к ряду вычислительных задач: 1) решение уравнения переноса при известной правой части, усреднение полученного решения по углам и вычисление тензора квазидиффузии, 2) решение многогрупповой системы уравнений квазидиффузии и усреднение этой системы по энергии в одногрупповую систему, 3) решение полученной одногрупповой системы уравнений квазидиффузии, 4) нахождение параметров критической сборки реактора.

Поскольку саморегулируемый нейтронно-ядерный режим поддерживает в

реакторе квазистационарное распределение скалярного потока нейтронов \ур с

очень малой постоянной времени X, будем решать стационарное уравнение

переноса нейтронов, в котором члены рассеяния и деления преобразованы в

соответствии с методом квазидиффузии:

„ дФр „ дФр Л дФр _ А,.--

С1.--+ £1,-+(N0? +—)ФР =

1 & у ду г & ' у"

_1_ 4л

(1)

+ 5/2Wp)k\3'£QinjDfJ-l у V S

+ N&~^-cfvfmPV + Nafvf№Pd-f(P\V)4>]+QP-

Здесь угловые переменные

Qz = cos9, = sin9coscp, Q^, = sin9sin<p, dQ = sin0i/9ai(p,9е[0,л]. (2)

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

к

В соответствии с общепринятыми обозначениями:

• о, = о, + ст/ + ст„7 + с„2и + а„з„ - полное сечение всех процессов столкновения нейтронов в группер (р=\,...,Р), включая:

• а5 = + <3фп) - сечение рассеяния (упругое и неупругое);

• о„2П - сечение поглощения нейтрона с рождением двух нейтронов;

• _ суммарное сечение рассеяния и реакций п-*2п, п-*Ъп\

• (£1 • £1') - индикатриса рассеяния;

• О/ - сечение деления;

• а„7 - сечение поглощения нейтрона с излучением у кванта;

• р* - доля запаздывающих нейтронов в группе к;

• р/р - доля группы g в запаздывающих нейтронах;

• _ доля нейтронов деления группы к, попадающих в группу р;

• уР _ доля запаздывающих нейтронов группы g, попадающих в группу р\

• се- концентрация предшественников запаздывающих нейтронов;

• X = 1п2 / Тг, - период полураспада предшественников запаздывающих

нейтронов группы g.

В правой части (1) присутствуют скалярный и векторный потоки нейтронов, которые вычисляются не по функции распределения в соответствии с формулами

р р

V = | Ф"сЯ2, = = ,\У = X , (3)

й й р=1 Р=х

а из независимой системы уравнений квазидиффузии.

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

тКи0) = 1 /2*(0) + 3 /2и>(1)ц0 + 5 /2\У{2)(3\120 -1) /2, (4)

где = £!•£!' - косинус угла рассеяния.

Для нулевых коэффициентов разложения по полиномам Лежандра должны

быть выполнены законы сохранения числа частиц: р р р

2>ГоМеи=1> ЕИ'(0)п2,и =2' ЕИ(0)пЗп,* = 3" (5)

р=1 р=1

Для учета запаздывающих нейтронов нам нужна функция

(3

ХоЬ

V н У

Решение уравнения переноса позволяет найти групповые коэффициенты квазидиффузии

В? = | / (6)

й й

и коэффициенты граничного условия

с'= { Ф'Шй/ | Ф'<йг, (7)

(Пп)>0 / (Яп)>0

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

При решении задачи на собственные значения 1 = 0, /(р)(А.) = /(X) = 0. Однако, эти члены оставлены в уравнениях, поскольку они существенны для динамического моделирования, в котором методика решения уравнения переноса останется той же самой.

Нахождение критических параметров реактора сводится к решению системы уравнений для краевой задачи на определение собственного значения р, и, в дальнейшем, определение параметров активной зоны, обеспечивающих значение параметра реактивности р = 0 или 1/(1 -р) = 1.

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

д\Ухр д1Ур ки р р . р

дх ду & ^

ас Ф & ы\

ОХ Су & ¡Ы1

дх ду О! к=\

Граничные условия для системы (8)-(11):

\Уп = сп\\1 на внешней границе реактора, где п - внешняя нормаль к грани в плоскости ху.

IV; = с\\\>, 1Уг = с2г\\) на верхней и нижней границе реактора, при этом, если есть плоскость симметрии, то на ней сг2 = 0.

Все коэффициенты граничных условий с„, с\, с\ находятся по решению уравнения переноса в соответствии с (7).

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

Ш. Ж />г/ . X.

дх ду дг V (12)

= И{а- уРДА))^

дф^ ) оф \\1 ) дф ¥) -г-

° +-Г-+ (13)

дх ду дг

дф V) дф„ч ) дф у ) -=-

-г-+-т;-+-Т-+Ы(а,-а,угт)1У , (14)

дх ду дг

дфх:у ) , дфу^ ) | дф^ху )

дх ду дг

Здесь введена функция

Правые части dx , dy и dz включают в себя как просуммированные по всем группам правые части уравнений (9)-(11), так и добавки, связанные с возможно знакопеременным усреднением в левых частях уравнений.

В уравнениях (1), (8)-{11), (12)—(15) ищутся величины, стоящие в левой части уравнений, при этом правые части всех уравнений считаются известными с предыдущей итерации метода квазидиффузии.

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

Здесь коэффициенты аи£> считаются известными.

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

где 5 — проекция характеристики в плоскости, перпендикулярной оси 1. Эта замена соответствует выбору из всех плоскостей, проходящих через характеристику, той плоскости, которая параллельна оси г, и разложению движения на параллельное оси 2 и перпендикулярное ей.

Характеристический метод решения заключается в том, что по значению функции в точках 1, 2, 3 (рис. 1а) нужно определить значение в точке пересечения характеристики с границей расчетной ячейки Ф^, после чего

(16)

cos в— + sin0— + аФ =Q,

dz 6s

(17)

значение в точке 4 определяется интегрированием вдоль отрезка этой характеристики. Чтобы обеспечить второй порядок точности, нужно знать снесенное тем или иным способом по характеристике назад значение Ф^, построить квадратичную интерполяцию по значениям Ф|,Ф2,Ф3, и по ней определить промежуточное значение Ф^. Окончательно решение в узле 4 находится по формуле:

Ф4 = Ф; ехр (-ас//) + Q/ а(1- ехр(-окй)), (18)

где (II—длина отрезка характеристики внутри ячейки.

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

а) б)

Рис. 1. Характеристический метод (а) и его консервативный вариант (б).

Консервативная модификация характеристического метода решения заключается в восстановлении функции распределения на освещенных гранях (1-2, 2-3, рис. 16) по заданным значениям в узловых точках и интегралу от функции распределения по грани. Этим формируется правильная выпуклость аппроксимирующей функции. Значение Ф5 используется как в чисто

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

Рассмотрим случай, когда характеристика из точки 4 назад пересекается с нижней гранью ячейки (рис. 16). Входящий поток на грани 1-2, соответствующим образом преобразованный, целиком выходит через грань 1-4 (составляющая \>Лб) . Входящий поток на грани 2-3 разбивается на две части, одна из которых после учета поглощения и источника (\¥б4) в сумме с Wl6 дает полный поток через грань 1-4, а вторая часть формирует полный поток '\¥34. На отрезке 2-3 строится квадратичная или псевдоквадратичная интерполяция функции распределения, чтобы удовлетворить условиям:

Аналогично строится интерполяция F2^(¿;) на ребре 1-2.

Сравнение свойств исходного характеристического метода и его консервативной модификации проводилось на задачах, имеющих точное решение. В таблице 1 приведены численные результаты анализа сходимости на трехмерном расширении двумерной задачи. Эта задача аналогична реакторной, но на торцах реактора вместо нулевых падающих потоков задано точное решение для области бесконечно протяженной вдоль оси г, а параметры а и £> от 2 не зависят. Видно, что консервативный метод имеет второй порядок сходимости, а характеристический только первый. Аналогичные результаты по сходимости были ранее получены для двумерной геометрии.

(19)

о

Ошибка в норме С Метод 5x5x5 9x9x9 17x17x17 33x33x33

м характ. 7.02x10-2 4.38x10-2 2.47x10-2 1.32x10-2

консерв. 3.06x10-4 6.15x10-5 1.01x10-5 2.50x10-6

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

Для саморегулируемого нейтронно-ядерного режима не требуется тяжелая система компенсирующих стержней, поэтому активная зона может быть собрана с высокой степенью симметрии (рис. 2). Симметричность и гексагональная структура сборки позволяют строить угловую сетку одновременно практически равномерной и жестко связанной с пространственной (что по сути означает возможность использования свойств длинных характеристик). Это избавляет от необходимости использовать интерполяцию решения, что было характерно для двумерной цилиндрической геометрии. Благодаря переходу сборки самой в себя при повороте относительно центра на ж/Ъ (см. рис. 2), т.е. идентичности шести угловых сегментов (сборка, включая защиту, в горизонтальном сечении имеет форму шестиугольника), можно рассматривать 1/6 от общего числа характеристических лучей, что используется в программе расчета. Азимутальный угол в и угол места <р входят в уравнение (16) только как параметры, поэтому для различных в и <р решение может находиться независимо (в том числе, параллельно).

Рис. 2. Картограмма активной зоны быстрого реактора типа БН-800.

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

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

Многогрупповая система уравнений квазидиффузии (8)-(11) теперь может быть представлена в следующем виде:

д\У/ 81ГР ту; „ I

+ + (20)

оп от дг

. (22)

-л п I 2 г V '

СП ОТ 02

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

Построение разностной схемы для системы (20)-(22) проведем интегро-

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

остался только скалярный поток у/р. Проинтегрировав уравнение (20) по всему

объему ячейки (см. рис. За) получим:

и?: ( / +1,} +1 ,к) + IV/ а +1, у, к) + (/ + 1,у, к) \-5aS

-1,] -1 ,к)+у?:(.• -1,],к)+w:0,; -\,к) |

1.5ал/з

Для того, чтобы в уравнении (23) остался только скалярный поток, выразим все входящие в него нормальные векторные потоки и IV/ из оставшихся уравнений (21) и (22) - проинтегрируем эти уравнения по восьми парам полуячеек, прилегающих к граням данной ячейки (рис. 36,в). В полученных парах уравнений типа (24) и (25) для соответствующих полуячеек выразим векторные потоки через скалярные так, чтобы в объединенном уравнении, полученном из пары, исключить значение скалярного потока на общей грани.

руд+и+и*)++и,*)+д>"(|',у+1,*)

1.5ал/з

РУО'-1,7-1 Д) + (/ -1,7,к) + ,у -1Д) 1.5аТз

1.5а>/з

(Р:„ + >/зД£У'(|+ !,;,*) + (/% -Узр;гК(/,у + 1Д) | 1.5а>/з

, Д>'(/,Л* + 1)-Д>'(Ц,*-1) + = ^ л

Здесь а - сторона шестиугольника, Л - высота ячейки по оси г.

¡4-1.к+1

(25)

¡-1,1-1 а

¡-1,1 .к

¡Я,к*

• ¡О.к

•

• ¡.¡,к-1

1+1

¡-2,|-1,к

¡+1.)+2.к

в)

Рис. 3. Схемы для интегро-интерполяционного метода. 21

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

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

Для получения усредненной по энергиям одногрупповой системы уравнений квазидиффузии необходимо фактически просуммировать уравнения многогрупповой системы. При этом в конечной системе (12)-(15) фигурируют уже суммарные скалярные и векторные потоки. Получение решения этой системы ничем не отличается от получения решения многогрупповой системы квазидиффузии.

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

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

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

производился для стандартного 26-группового приближения по энергии нейтронов с учетом 44 нуклидов.

Рис. 4. Вертикальная конфигурация активной зоны реактора типа БН-800.

В табл. 2 приведено сравнение результатов расчета критических параметров в двумерной и трехмерной геометриях (результат данной работы).

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

Вариант расчета 20 ЗВ

Диффузионное приближение 15.17см 15.02 см

Квазидиффузия с переносом на каждой итерации источника 14.48 см 14.64 см

Квазидиффузия с переносом на второй итерации источника 14.49 см 14.65 см

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

X X

Рис. 5. Распределения скалярного потока нейтронов для седьмой группы.

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

ОСНОВНЫЕ РЕЗУЛЬТАТЫ ДИССЕРТАЦИИ

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

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

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

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

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

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

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

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

1. Е.Н.Аристова, Д.Ф. Байдин. Сравнение консервативного и характеристического методов решения уравнения переноса в r-z геометрии на основе перехода к переменным Владимирова в методе квазидиффузии // Труды XLVIII научной конференции МФТИ "Современные проблемы фундаментальных и прикладных наук". Часть VII. Управление и прикладная математика. / Моск. физ. - техн. ин-т. - М. - Долгопрудный, 2005. С. 162163.

2. E.H. Аристова, Д.Ф. Байдин, В.Я. Гольдин. Два варианта экономичного метода решения уравнения переноса в r-z геометрии на основе перехода к переменным Владимирова // Труды XIII международной конференции

"Математика. Компьютер. Образование", секция Вычислительные методы и математическое моделирование, 2006. Т. 2, С. 158-170.

3. E.H. Аристова, Д.Ф. Байдин, В.Я. Гольдин. Два варианта экономичного метода решения уравнения переноса в r-z геометрии на основе перехода к переменным Владимирова // Математическое моделирование, 2006, т. 18, № 7, с. 43-52.

4. E.H. Аристова, Д.Ф. Байдин. Два варианта экономичного метода решения уравнения переноса в r-z геометрии на основе перехода к переменным Владимирова с учетом логарифмической особенности // Труды 49-й научной конференции МФТИ "Современные проблемы фундаментальных и прикладных наук". Часть VII. Управление и прикладная математика. / Моск. физ. - техн. ин-т. - М. - Долгопрудный, 2006. С. 187-188.

5. D.F. Baydin, E.N. Aristova, V.Ya. Gol'din. Comparison of the efficiency of the transport equation calculation methods in characteristics variables // Transport Theory and Statistical Physics, 2008, vol.37, № 2, pp. 286-306.

6. E.H. Аристова, Д.Ф. Байдин. Экономичный метод решения уравнения переноса в 2D цилиндрической и 3D гексагональной геометриях для реализации метода квазидиффузии II Труды 53-й научной конференции МФТИ "Современные проблемы фундаментальных и прикладных наук". Часть VII. Управление и прикладная математика. / Моск. физ. - техн. ин-т. -М. - Долгопрудный, 2010. Т. 3, С. 91-92.

7. E.H. Аристова, Д.Ф. Байдин. Экономичный метод решения уравнения переноса в 2D цилиндрической и 3D гексагональной геометриях для метода квазидиффузии // Компьютерные исследования и моделирование, 2011, т. 3, № 3, с. 279-286.

Байдин Денис Федорович

ЭКОНОМИЧНАЯ ТРЕХМЕРНАЯ МЕТОДИКА РАСЧЕТА КРИТИЧЕСКИХ ПАРАМЕТРОВ АКТИВНОЙ ЗОНЫ БЫСТРОГО РЕАКТОРА С ЕСТЕСТВЕННОЙ БЕЗОПАСНОСТЬЮ

Автореферат

Подписано в печать 20.10.2011. Формат 60 х 84 '/,6. Усл. печ. л. 1,0. Тираж 100 экз. Заказ № 742. Федеральное государственное бюджетное образовательное учреждение высшего профессионального образования «Московский физико-технический институт (государственный университет)» Отдел оперативной полиграфии «Физтех-полиграф» 141700, Московская обл., г. Долгопрудный, Институтский пер., 9

Введение.

Глава 1. Квазидиффузионный подход к моделированию саморегулируемых нейтронно-ядерных режимов.

1.1. Физические аспекты саморегулируемых нейтронно-ядерных режимов (СНЯР).

1.1.1. СНЯР первого рода.

1.1.2. СНЯР второго рода.

1.2. Математическое моделирование саморегулируемых нейтронно-ядерных режимов.

1.2.1. Развитие математических моделей СНЯР.

1.2.2. Постановка задачи для многогруппового уравнения переноса.

1.2.3. Многогрупповые уравнения квазидиффузии.

1.2.4. Усреднение в одногрупповую систему уравнений квазидиффузии.

Глава 2. Решение уравнения переноса в собственных характеристических переменных.

2.1. Решение уравнения переноса в двумерной г-2 геометрии.

2.1.1. Постановка задачи.

2.1.2. Переход к переменным метода Владимирова.*!.

2.1.3. Угловая дискретизация.

2.1.4. Параболический естественный сплайн.

2.1.5. Характеристический и консервативно-характеристический методы решения.

2.1.6. Интегрирование по углам.

2.1.7. Анализ аппроксимации функции распределения на логарифмических разрывах.

2.1.8. Результаты численного сравнения.

2.2. Решение уравнения переноса в трехмерной х-у-г геометрии.

2.2.1. Постановка задачи.

2.2.2. Переход к переменным метода Владимирова.

2.2.3. Характеристический и консервативно-характеристический метод решения.

2.2.4. Интегрирование по углам.

2.2.5. Результаты численного исследования.

Глава 3. Решение уравнения переноса в трехмерной геометрии в рамках метода квазидиффузии для быстрого реактора.

3.1. Активная зона быстрого реактора.

3.2. Пространственная и угловая сетки.

3.3. Построение решения для уравнения переноса в трехмерной геометрии.

Глава 4. Решение системы уравнений квазидиффузии.

4.1. Решение многогрупповой системы уравнений квазидиффузии.

4.2. Метод последовательной верхней релаксации.

4.3. Результаты методических расчетов.

4.4. Эффективная одногрупповая система уравнений квазидиффузии.

Глава 5. Численное нахождение критических параметров активной зоны быстрого реактора.

5.1. Постановка задачи.

5.2. Описание итерационного процесса.

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

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

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

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

Широкий класс физических явлений и устройств может быть изучен с помощью моделей, основанных на приближении многокомпонентной сплошной среды с неравновесным переносом. Во многих задачах математической физики, таких как математическое моделирование процессов, протекающих в звездах, в задачах управляемого термоядерного синтеза, при моделировании распространения излучения в атмосфере и многих других возникает необходимость численного решения многомерного уравнения переноса излучения. При проектировании активных зон реакторов и при расчетах защиты реакторов, а также систем управления, встает задача нахождения решения уравнения переноса нейтронов, во многом родственного уравнению переноса излучения [1-3]. Уравнение переноса является линейным интегро-дифференциальным уравнением первого порядка относительно функции распределения частиц (фотонов или нейтронов). Эти два типа уравнений переноса отличаются структурой правой части, отвечающей за источники возникновения частиц. Соответственно, могут отличаться и постановки задач для переноса излучения и нейтронов: если для переноса излучения ставится начально-краевая задача, то в задачах переноса нейтронов помимо начально-краевой задачи возможна постановка задачи на нахождение собственных значений (например, задача расчета критических параметров реактора или защитных экранов). Однако многие проблемы решения для обеих разновидностей уравнения переноса являются общими.

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

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

Первые численные методы решения уравнения переноса были созданы в ходе работы над советским и американским атомными проектами и касались, главным образом, решения уравнения переноса в одномерной сферической геометрии. Практически это был первый опыт численного решения уравнений в частных производных. Первые предложенные методы интегрирования уравнения переноса можно разделить на два больших класса: это методы, которые в дальнейшем стали называться методами Карлсона (в советском атомном проекте аналогом этого метода был КН метод В.Я.Гольдина) [4-6] и характеристические методы, самым знаменитым из которых является метод Владимирова [7-9]. Если в £„ методах производится разностная аппроксимация непосредственно уравнения в частных производных, то характеристические методы базируются на сведении уравнения переноса к обыкновенному дифференциальному уравнению (ОДУ) вдоль некоторого набора характеристик или на представлении уравнения переноса в интегральной форме. Соответственно достоинства и недостатки у каждого из этих классов методов свои. К достоинствам метода нужно отнести консервативность (если разностной аппроксимация строилась для уравнения, записанного в дивергентной форме), возможность включения в расчет учета других физических процессов, легко достижимую аппроксимацию второго порядка точности. К его недостаткам следует отнести следующую из второго порядка аппроксимации по теореме Годунова

10] теоретическую и практическую немонотонность метода [11], а во многих практически важных случаях и неположительность схемы. В свою очередь, методы характеристик тоже можно разделить на два больших подкласса: метод длинных характеристик и метод коротких характеристик. Метод длинных характеристик может обеспечить практически любую точность решения при использовании, во-первых, соответствующего метода решения ОДУ (на текущий момент существует огромное множество методов высокого порядка аппроксимации) и, во-вторых, тщательно подобранного набора характеристических направлений [12]. Метод Владимирова является примером блестящего сочетания метода длинных характеристик и относительной экономичности метода, достижимой в одномерных сферической и цилиндрической геометриях [7-8]. Однако в более сложных геометриях метод длинных характеристик, для которого решение в заданной точке пространства получается численным интегрированием ОДУ по всей длине характеристики, начиная от границы, является слишком затратным, поскольку в общей постановке множества характеристик для различных узлов не пересекаются. Это означает, что для каждого узла сетки нужно интегрировать уравнение переноса вдоль характеристики через всю сложную расчетную область, начиная от границы. Такие затраты расчетного времени себя не оправдывают. Кроме того, метод характеристик без дополнительных условий не обеспечивает консервативности, а неудачный выбор характеристических направлений может приводить к так называемому «эффекту луча» [13], при котором среди выбранных направлений в данной точке может отсутствовать направление «на источник», обеспечивающий главную часть, например, плотности частиц. Угловая неоднородность решения в заданной точке является неотъемлемым свойством решения уравнения переноса при пространственной неоднородности распределения источников. «Эффект луча» оказывается особенно неприятным при расчете потоков в удаленных от источника неоднородностях, что приводит к сильной чувствительности методов характеристик к выбору угловой сетки. В методах коротких характеристик решение в заданном узле разностной сетки ищется интегрированием не от границы всей расчетной области, а интегрированием вдоль отрезка характеристики, приходящей в узел с освещенной грани расчетной ячейки, т.е. только вдоль части характеристики, лежащей внутри заданной ячейки. Этот метод сталкивается с необходимостью многократных интерполяций на освещенных гранях для получения значения в точке входа характеристики, что ухудшает качество получаемого численного решения по сравнению с методом длинных характеристик. Также может возникнуть необходимость в сложных алгоритмах обхода расчетной области.

Примерно одновременно с £„ методом и методом характеристик был предложен метод прямого интегрирования Рихтмайера [14]. Несколько позже для случая сферической геометрии был предложен метод характеристических трубок В.Е.Трощиева. Он объединяет достоинства S„ и характеристического методов, положителен и монотонен, имеет второй порядок точности и консервативен, однако тяжело распространяемый на случай учета других физических процессов и многомерных геометрий [15]. Тогда же был предложен дискретный S„ метод (DSn метод) [5]. Этот метод в отличие от метода характеристических трубок легко обобщается на многомерные геометрии, но немонотонен.

Простота реализации S„ метода в многомерных геометриях, второй порядок аппроксимации и консервативность сделали S„ метод весьма привлекательным в глазах многих поколений вычислителей [16-18].

Следующим классом предложенных схем повышенного порядка аппроксимации для решения уравнения переноса стали моментная .Diamond Difference (DD) схема, использующая только основное уравнение баланса, и родственные ей нодальные схемы, увеличивающие порядок аппроксимации с увеличением числа используемых уравнений баланса [19].

Варьируя форму дополнительных соотношений DD схемы,, удается улучшить качество сеточного решения, не отказываясь от второго порядка аппроксимации. При этом прибегают, например, к разбиению ячейки на части отрезками характеристик исходного уравнения [20] или отрезками, параллельными граням ячейки [21], к аппроксимации уравнения в интегральной форме [22]. Аналогичным образом удается улучшить качество нодальных схем с большим числом моментов. Однако эти улучшенные схемы не являются безусловно положительными и монотонными. Решение, полученное из неположительной схемы, может содержать отрицательные скалярные потоки, и сгущение сеток все равно может приводить к появлению у решения нефизических осцилляций большой амплитуды. Метод, соединяющий характеристический подход с сохранением консервативности, предложен в [23].

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

Чтобы соединить в рамках единой сеточной схемы высокий порядок аппроксимации и свойства положительности и монотонности, переходят к нелинейным схемам, т.е. взвешенным схемам с весовыми параметрами, зависящими от сеточного решения. Первоначально расчет ячейки выполняется со значениями параметров, отвечающих наибольшему из возможных порядков аппроксимации. Если полученное решение не удовлетворяет условию положительности и/или монотонности, проводится коррекция сеточного решения: ячейка пересчитывается со значениями параметров, гарантирующих полное или частичное выполнение рассматриваемых условий. При коррекции порядок схемы снижается. Наиболее известными из нелинейных схем являются взвешенная алмазная схема \\fDD [24], адаптивные А\¥00 [17,25], адаптивные нодальные схемы [26].

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

Еще один класс монотонных нелинейных методов высокого порядка точности носит название TVD схем (Total Variation Diminish) [27,28], основанных на введение ограничений на потоки. Другой подход к построению нелинейных схем высокого порядка аппроксимации с подсеточным разрешением разрывов предложен в [29]. Этот метод использует плавающий шаблон, что позволяет определить возможное положение разрыва решения внутри ячейки. Эти методы наиболее полно применяются для конструирования схем газовой динамики.

Все, что было сказано выше, относится к разностной аппроксимации дифференциального оператора в уравнении переноса. Исследование порядка аппроксимации уравнения переноса производится в предположении непрерывности и ограниченности частных производных функции распределения вплоть до некоторой степени я, отвечающей порядку главного члена погрешности аппроксимации. Однако решение задач для уравнения переноса, как правило, имеет особенности на внешних граничных поверхностях, в окрестности сосредоточенных источников, на характеристиках, касательных к поверхностям разрыва параметров среды. Это означает, что вблизи особенностей решение сингулярное, т.е. обладает большими по величине градиентами, или является недифференцируемым, или даже разрывным. Согласно Р.М.Шагалиеву [30], это приводит к тому, что в расчетах реальных гетерогенных задач порядок сходимости, оцениваемый по ошибке решения при сгущении сеток, в три раза меньше декларируемого порядка аппроксимации (в предположении гладкости функции и ее производных порядок сходимости и порядок аппроксимации для устойчивых разностных схем обязаны совпадать): для схем третьего порядка аппроксимации имеет место сходимость первого порядка, для схем второго порядка - сходимость порядка 0.7-0.6, а для схем первого - 0.3. При таком анализе становится ясно, что схемы первого порядка аппроксимации являются совершенно неудовлетворительными по порядку сходимости.

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

• положительность (неотрицательность решения при неотрицательных источниках и индикатрисе рассеяния);

• монотонность (сохранение в сеточном решении числа и расположения экстремумов точного решения);

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

Трудности решения уравнения переноса не исчерпываются только разностной аппроксимацией уравнения в частных производных. В общем случае функция распределения частиц зависит от семи переменных: трех пространственных, двух угловых, энергии и времени. Развитие современных суперкомпьютеров и вычислительных кластеров позволяет во многом разрешить проблему ограничения по требуемой памяти для задач со многими измерениями. Однако многомерные динамические расчеты в многогрупповом приближении, обеспечивающем необходимую точность, даже сейчас возможны в единичных случаях на многопроцессорных вычислительных системах. Эффективные параллельные алгоритмы решения задач переноса разрабатываются во ВНИИЭФ и в ИПМ РАН им. М.В.Келдыша [31-33]. Предлагаемые в этих работах алгоритмы основываются на распараллеливании задачи по пространственным областям и/или распараллеливании конвейерного типа, когда несколько сегментов вычислительной системы выполняют каждый свои действия в последовательной вычислительной структуре. Поскольку при любом из этих двух методов распараллеливания обмены данными между сегментами вычислительной системы достаточно существенны, то обычно используют архитектуру MPI.

Спектральное описание решения задач переноса нейтронов (и переноса излучения) также представляет значительные трудности. Связаны они с двумя факторами. Общим местом уже является сложная зависимость коэффициентов уравнения переноса от энергии частиц. Эта сложная зависимость коэффициентов поглощения и других, входящих в уравнение, включающих как непрерывный, так и линейчатый спектры, приводит к необходимости либо использовать чрезвычайно подробную сетку по энергии (не менее 10 точек на каждую линию), либо применять некоторые приближения для описания и расчета задачи. На этом этапе введение иерархии вычислительных моделей является наиболее оправданным. Например, в наиболее точных методах 'line-by-line' при учете всех линий поглощения атмосферными газами в задаче теплового баланса атмосферы Земли получается система уравнений переноса, содержащая порядка нескольких миллионов энергетических точек. Если учесть необходимую пространственную и угловую дискретизацию задачи, то такой расчет даже для одномерной пространственной геометрии становится возможен в единичных случаях, особенно при наличии рассеяния. Более простой моделью спектрального представления задачи является многогрупповое приближение, при введении которого используется эффективное усреднение по отрезкам частот, соответствующим некоторому разбиению энергетической шкалы. Строгое введение групповых коэффициентов поглощения возможно в ограниченном числе случаев: для оптически тонкого тела, для непрерывного спектра, для излучения, близкого к локальному термодинамическому равновесию. Тем не менее, практика расчетов показывает, что для большого числа задач может быть использовано многогрупповое приближение при соответствующем выборе весовой функции. Аналогичная проблема усреднения возникает при расчете групповых микроконстант, требующихся1 для реакторных задач, т.к. при различных взаимодействиях нейтронов с ядрами также возникают резонансные области с изменением величины микросечений этих реакций на несколько порядков вблизи резонанса [1]. Выстраивание правильной иерархии моделей и правильный выбор модели для конкретной задачи также является предметом математического моделирования.

Еще одна проблема возникает в задачах с рассеянием: при наличии сильной анизотропии рассеяния и, что для нас более существенно, в реакторных задачах, при нахождении собственного значения. Достаточно хороший обзор большинства существующих методов ускорения численного решения в области итераций по правой части уравнения переноса присутствует в [34]. Там подробно разобрано развитие методов, присутствующих, как в отечественной, так и в зарубежной литературе. Основой всех этих методов является метод итераций источника (81). Среди используемых в современной литературе методов ускорения итераций источника, применительно к реакторным задачам, следует отметить метод пространственного ребаланса [35], заключающийся в нормализации решения на каждой итерации путем введения мультипликативных поправок для нулевого углового момента решения, чтобы удовлетворить уравнению нейтронного баланса в некоторой области (ячейка сетки или объединение ячеек); КР-методы [36] или «синтетические методы» (SA), среди которых наиболее распространены DSA (diffusion) и P|SA, заключающиеся в разбиении итерации на два шага, причем на втором шаге считается поправка к решению первого шага из уравнения с оператором в левой части более простого вида (например, диффузионным); метод квазидиффузии [37,38], который выбран за основу в данной работе и будет описан ниже, и различные моментные методы; а также многосеточные методы.

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

В предлагаемой работе основой численного решения задач переноса нейтронов является метод квазидиффузии, предложенный В.Я.Гольдиным в 1964 году [37,38]. Он также относится к классу нелинейных методов решения уравнения переноса, но не к классу методов коррекции скалярного потока. Метод квазидиффузии заключается в постепенном понижении размерности задачи введением ряда дробно-линейных функционалов, слабо зависящих от решения. Понижение размерности задачи проходит в два этапа. На первом этапе происходит усреднение уравнения переноса по угловым переменным, в результате которого получается многогрупповая система уравнений квазидиффузии. На втором - усреднение по энергии, что приводит к эффективной одногрупповой системе уравнений квазидиффузии для скалярного и векторного потока (излучения или частиц), которая уже может быть объединена с уравнениями, описывающими другие физические процессы, происходящие в системе. При внешнем усложнении подхода метод квазидиффузии позволяет разрешить ряд трудностей. Во-первых, метод квазидиффузии нелинеен, так как вводит дробно-линейные функционалы для вычисления компонент тензора квазидиффузии, замыкающих систему уравнений меньшей размерности. Во-вторых, уравнения квазидиффузии выражают собой законы сохранения, поэтому консервативность получается автоматически. В-третьих, при умеренной анизотропии рассеяния метод квазидиффузии позволяет выразить главную часть интеграла- рассеяния внутри группы через групповые скалярный и векторный потоки, которые вычисляются самостоятельно из системы уравнений квазидиффузии, что обеспечивает быструю сходимость итераций по рассеянию. Также быструю сходимость итераций по рассеянию и делению обеспечивает переноса некоторых членов и правой части (источник) налево. В-четвертых, введение интегральных по энергии уравнений квазидиффузии позволяет эффективно объединять эту систему с уравнениями выгорания и реакторной кинетики. Для задачи на нахождение собственных, значений и/или критических параметров, активной зоны реактора значительно сокращается общее число итераций. Схематично структура метода квазидиффузии приведена на рис. 1.

Работы группы под руководством В.Я.Гольдина и их коллег [39-41] развивали метод квазидиффузии в одномерной и двумерной геометриях. Работа по исследованию устойчивости квазидиффузионного метода и по созданию методик решения двумерного уравнения переноса в методе квазидиффузии ведется параллельно в США [42,43]. Как в отечественных, так и в зарубежных работах было показано, что при слабой анизотропии рассеяния метод квазидиффузии быстро сходится. В англоязычной литературе, как многогрупповые уравнения квазидиффузии, так и эффективная одногрупповая» система уравнений квазидиффузии называются решения полученных разностных уравнений и способов ускорения сходимости итераций. Проводится сопоставление численных решений с точными решениями, там, где они существуют, исследуется сеточная сходимость, а также проводится сравнение результатов математического моделирования с результатами, полученными для задач меньшей размерности.

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

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

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

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

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

Основные результаты, представленные в диссертации, опубликованы в работах [44-50], были представлены и обсуждались на Всероссийских и Международных конференциях, семинарах института: • 48-я научная конференция МФТИ, ноябрь 2005, Москва;

• XIII Международная конференция "Математика. Компьютер. Образование", январь 2006, Дубна;

• 49-я научная конференция МФТИ, ноябрь 2006, Москва;

• The 20th International Conference on Transport Theory, July 2007, Obninsk, Russia;

• XV Международная конференция "Математика. Компьютер. Образование", январь 2008, Дубна;

• 52-я научная конференция МФТИ, ноябрь 2009, Москва;

• 53-я научная конференция МФТИ, ноябрь 2010, Москва;

• XVIII Международная конференция "Математика. Компьютер. Образование", январь 2011, Пущино;

• Семинары кафедры вычислительной математики Московского физико-технического института, научные семинары ФУПМ, 2009-2011, Москва.

Основные положения, выносимые на защиту:

1. Предложена консервативная модификация характеристического метода численного решения двумерного стационарного уравнения переноса в ячейке с точным перераспределением выходящих потоков по граням.

2. Разработан алгоритм решения стационарного уравнения переноса для трехмерной гексагональной геометрии на основе метода коротких характеристик и предложенного консервативно-характеристического метода решения уравнения переноса в ячейке.

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

4. На основе математической модели многогруппового уравнения переноса с квазидиффузией создана программа для расчета критических параметров сборки активной зоны реактора типа БН-800, способного работать в саморегулируемом нейтронно-ядерном режиме; с ее помощью получены численные результаты для параметров критической сборки.

1. Вейнберг А., Вигнер Е. Физическая теория ядерных реакторов. - М.: Изд-во Иностранной литературы, 1961. - 732 с.

2. Марчук Г.И. Численные методы расчета ядерных реакторов. М.: Атомиздат, 1958.-381 с.

3. Марчук Г.И., Лебедев В.И. Численные методы в теории уравнения переноса нейтронов. — М.: Атомиздат, 1981. 454 с.

4. Голъдин В.Я. Методы расчета переноса нейтронов и горения в термоядерном изделии (1948-1960гг.) Международный симпозиум, Дубна, 14-17 мая 1996г. В сб: «Наука и общество: история советского атомного проекта (49-50-е годы)», 1999. - Т.2. - С. 497-501.

5. Карлсон Б. Численное решение задач кинетической теории нейтронов. В сб.: «Теория ядерных реакторов». М.: Госатомиздат, 1963. - С. 243258.

6. Карлсон Б.Г., Латроп КД. Теория переноса. Метод дискретных ординат. В сб.: «Вычислительные методы в физике реакторов». / Под ред. X. Гринспена, К. Келбера, Д. Окрента. М.: Атомиздат, 1972. - С. 102157.

7. Владимиров B.C. Численное решение кинетического уравнения для сферы // Вычислительная математика. 1958. - Т.З - С. 3-33.

8. Владимиров B.C. Математические задачи односкоростной теории переноса частиц. Тр. Матем. ин-та им. Стеклова АН СССР. 1961. -158 с.

9. Голъдин В.Я., Данилова Г.В., Калиткин H.H. Численное интегрирование многомерного уравнения переноса // В сб. «Численные методы решения задач математической физики». М., 1966.-С. 190-193.

10. Годунов С.К, Рябенький B.C. Разностные схемы (введение в теорию). -М.: Наука, 1973.-440 с.распределенной памятью // Вопросы атомной науки и техники. Серия «Математическое моделирование физических процессов». 2000. -Вып. 1 - С. 21-26.

11. Adams M.L., Larsen Е. W. Fast iterative methods for discrete-ordinates particle transport calculations // Progress in Nuclear Energy. 2002. - V.40, Issue 1-P. 3-159.

12. Сычугова Е.П. Исследование устойчивости и эффективности метода пространственного ребаланса для ускорения сходимости итераций в задачах переноса частиц // Математическое моделирование. 2008. -Т.20, №9 - С. 75-93.

13. Волощенко A.M. КР|-схема ускорения внутренних итераций для уравнения переноса в трехмерной геометрии, согласованная со взвешенной алмазной схемой // Журнал вычислительной математики и математической физики. 2009. - Т.49, №2 - С. 344-372.

14. Гольдин В.Я. Квазидиффузионный метод решения кинетического уравнения // Ж. вычисл. матем. и матем. физ. 1964. - Т.4, №6 - С. 10781087.

15. Гольдин В.Я. О математическом моделировании задач сплошной среды с неравновесным переносом // В сб. «Современные проблемы матем. физ. и вычисл. матем.» М.: Наука, 1982. - 340 с.

16. Голъдин В.Я., Данилова Г.В., Калиткин Н.Н. Численное интегрирование многомерного уравнения переноса // В сб. «Численные методы решения задач математической физики» М., 1966. - С. 190-193.

17. Голъдин В.Я., Голъдина Д.А., Колпаков А.В. О решении двумерной стационарной задачи квазидиффузии: Препринт / ИПМ им. М.В.Келдыша. М., 1982. - №49 - 13 с.

18. Голъдин В.Я., Колпаков А.В. Нелинейный метод потоковой прогонки для решения многомерного диффузионного уравнения: Препринт / ИПМ им. М.В.Келдыша. М., 1982. - №22 - 13 с.

19. Constantinescu A., Anistratov D.Y. Stability Analysis of the Quasidiffusion Method for ID Periodic Heterogeneous Problems // Trans. Am. Nucl. Soc. -2006.-V.95-P. 565-567.

20. Wieselquist W.A., Anistratov D.Yu. The Quasidiffusion Method forTransport Problems in 2D Cartesian Geometry on Grids Composed of Arbitrary Quadrilaterals // Trans. Am. Nucl. Soc. 2007. - V.97 - P. 475-478.

21. Аристова Е.Н., Байдин Д.Ф., Голъдин В.Я. Два варианта экономичного метода решения уравнения переноса в r-z геометрии на основе переходак переменным Владимирова // Математическое моделирование. 2006. — Т.18, №7 - С. 43-52.

22. Baydin D.F., Aristova E.N., Gol'din V.Ya. Comparison of the efficiency of the transport équation calculation methods in characteristics variables // Transport Theory and Statistical Physics. 2008. - V.37, №2 - P. 286-306.

23. Аристова Е.Н., Байдин Д.Ф. Экономичный метод решения уравнения переноса в 2D цилиндрической и 3D гексагональной геометриях для метода квазидиффузии // Компьютерные исследования и моделирование. 2011. - Т.З, №3 - С. 279-286.

24. Феоктистов Л.П. Анализ одной концепции физически безопасного реактора: Препринт / ИАЭ. M., 1988. - №4605/4.

25. Феоктистов Л.П. Безопасность ключевой момент возрождения ядерной энергетики//Успехи физ. наук. - 1993. - Т. 163, №8 - С. 89-102.

26. Голъдин В.Я., Анистратов Д.Ю. Реактор на быстрых нейтронах в саморегулируемом нейтронно-ядерном режиме // Мат. моделирование -1995. -Т.7, №10 С. 12-32.

27. Гольдин В.Я., Пестрякова Г.А., Трощиев Ю.В. Усовершенствование математической модели саморегулируемого реактора // Мат. моделирование. 2002. - Т. 14, №12 - с. 39-47.

28. Гольдин В.Я., Трощиев Ю.В. Управление мощностью быстрого реактора в саморегулируемом режиме и его пуск // Атомная энергия. 2005. -Т.98, вып.1 — С. 18-24.

29. Гольдин В.Я., Пестрякова Г.А., Трощиев Ю.В., Аристова E.H. Быстрый реактор на оксидном уран-плутониевом топливе в саморегулируемом режиме // Атомная энергия. 2003. - Т.94, вып.З - С. 184-190.

30. Гольдин В.Я., Пестрякова Г.А., Трощиев Ю.В., Аристова E.H. Саморегулируемый нейтронно-ядерный режим в реакторе с жестким спектром и карбидным топливом // Математическое моделирование. -2002. Т. 14, №1 - С. 27-40.

31. Гольдин В.Я., Аристова E.H., Пестрякова Г.А., Стойнов М.И. Проект активной зоны реактора типа БН-800, работающего без запаса реактивности с минимальным управлением в течение длительного времени // Мат. моделирование. 2009. - Т.21, №10 - С. 76-84.

32. Николаев М.Н., Цибуля A.M., Цикунов А.Г. и др. Комплекс программ CONSYST/ABBN подготовка констант БНАБ к расчетам реакторов и защиты: Отчет ГНЦ РФ ФЭИ. - 1998. - №9865.

33. Четверушкин Б.Н. Построение тестов и некоторые вопросы численного решения уравнения переноса нейтронов // В сб.: Вычислительные методы в теории переноса под. ред. Г.И.Марчука. М.: Атомиздат, 1969. - С. 189-201.

34. Калиткин H.H., Кузьмина JI.B. Об естественных интерполяционных сплайнах // Мат. моделирование. 1994. - Т.6, №4 - С. 77-1 10.

35. Бакирова М.И., Карпов В.Я., Мухина М.И. Характеристико-интерполяционный метод решения уравнения переноса // Дифференциальные уравнения. 1986. -Т.22, №7 - С. 1141-1148.