АЛГОРИТМ РЕШЕНИЯ МНОГОГРУППОВОГО СТАЦИОНАРНОГО УРАВНЕНИЯ ПЕРЕНОСА НЕЙТРОНОВ И ГАММА-КВАНТОВ НА СЕТКАХ, СОГЛАСОВАННЫХ СО СТРУКТУРОЙ РАСЧЕТНОЙ ОБЛАСТИ

( THE ALGORITHM FOR SOLVING THE MULTIGROUP STEADY-STATE TRANSPORT EQUATION FOR NEUTRONS AND GAMMA-QUANTA ON GRIDS, AGREED WITH THE STRUCTURE OF CALCULATION DOMAIN
Preprint, Inst. Appl. Math., the Russian Academy of Science)

СОДЕРЖАНИЕ

1. Введение                                                                                          

2. Постановка задачи и основные формулы                                      

3. Конечно-разностная формулировка задачи и метод решения     

4. Численные результаты                                                                  

5. Заключение                                                                                    

6. Литература                                                                                    

1. Введение.

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

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

2. Постановка задачи и основные формулы.

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

,                                (1)

Правая часть уравнения (1) для нейтронов группы  () имеет следующий вид:

 (2a)

Правая часть уравнения (1) для гамма - квантов группы  имеет следующий вид ():

                                       (2б)

где

                                                                   (3)

В (1) - (3) использованы следующие обозначения:

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

,                                                    (4)

где  - внешняя нормаль к границе.

Система уравнений (1) - (4) описывает стационарное распределение нейтронов и гамма - квантов с учетом процесса деления, а также заданных внутренних и граничных источников.

Если в правой части (2б) член с делением отсутствует (), система уравнений (1) - (4) описывает стационарное распределение нейтронов и гамма - квантов в зависимости от заданных внутренних и граничных источников.

Задача на  (однородная задача) описывается системой уравнений переноса нейтронов с нулевыми внутренними и граничными источниками и множителем  перед членом с делением.

Предполагается, что макроскопическое сечение рассеяния ,  задается в виде ряда по полиномам Лежандра  степени  (в  приближении, ):

,                     (5)

с условиями нормировки

                                                     (6)

                                           (7)

Для упрощения выражения (5) используется теорема сложения полиномов Лежандра: Пусть  и  - два единичных вектора со сферическими координатами  и  соответственно, и пусть , , . Тогда справедлива следующая формула:

        (8)

где  - присоединенные функции Лежандра первого рода с нормировкой:

                                   (9)

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

                        (10)

Для них имеем:

,                    (11)

В случае аксиально-симметричной двумерной R-Z геометрии (Рис.1) решение уравнения переноса (1) обладает свойством симметрии ():

                                   (12)

Рис.1. Система координат для уравнения переноса в R-Z геометрии

Уравнение переноса (1) в этом случае имеет следующий вид ():

,    (13)

где

,,

Система уравнений (13) решается в пространственной области , которая заключена между внешней образующей тела вращения и осью симметрии Z и представляет собой сечение тела вращения плоскостью, проходящей через ось Z (Рис. 2).

Рис.2. Область D и четырехугольная пространственная сетка

Подставляя (5) с использованием (8) в интеграл рассеяния, входящий в правую часть уравнения (1) для частиц группы , получим выражение для источника рассеяния:

 (14)

В (14) используются нормированные функции Лежандра (10) и угловые моменты  решения уравнения переноса (13) в группе , определяемые по формуле:

                             (15)

Выражение (14) получено с использованием формулы:

и свойства четности функции  по углу , а именно:

Запишем правую часть уравнения (13) для нейтронов группы  в аксиально-симметричной двумерной R-Z геометрии с учетом (14) и (15):

(16)

Правая часть (1) для гамма - квантов группы  имеет вид:

(17)

Граничное условие на внешней поверхности задается для потока частиц, входящих внутрь области :

                        (18)

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

                                       (19)

Для  имеем граничное условие:

         (20)

3. Конечно-разностная формулировка задачи и метод решения.

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

Угловая сетка задается  квадратурой, которая разбивает полусферу направлений на сектора следующим образом. На интервале  задаются узлы , . Для каждого  интервал изменения угла  ,  делится на подинтервалы , , где , . Диаграмма расположения точек для угловой квадратуры  на октанте  представлена на Рис. 3.

Рис. 3. Диаграмма расположения точек угловой квадратуры  на октанте.

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

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

Интегрируя уравнение (13) с весом  по четырехугольнику с координатами вершин  в пределах разностной угловой ячейки    площадью , применяя формулу Гаусса - Остроградского и сокращая на , получим следующее балансное уравнение (индекс группы  опущен):

         (21)

где

.                          (22)

,  ,                                

,    ,         ,              

Рис. 4. Значения функции  в четырехугольной пространственной ячейке.

Параметры  для заданного значения  и для всех  вычисляются из дополнительных соотношений:

,                       (23)

Площадь  четырехугольной пространственной ячейки в плоскости R-Z вычисляется через координаты вершин по формуле:

                                                       (24)

Объем  области, образованной вращением этой ячейки вокруг оси симметрии Z, вычисляется по формуле:

                                                (25)

Для замыкания системы (21) - (25) вводятся дополнительные “взвешенные” соотношения. По угловой переменной  дополнительное соотношение имеет вид:

                                          (26)

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

Граничные условия (18) задаются на внешней границе для тех сторон четырехугольной ячейки, для которых , что эквивалентно условию . Чтобы аппроксимировать граничное условие (19) на примыкающих к оси Z сторонах четырехугольников и в центре координат, это условие заменяется соотношением:

                                 (27)

Для аппроксимации граничного условия (20) для  используется следующее балансное уравнение:

,(28)

где .

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

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

Для фиксированного значения  сначала решается уравнение (28) для направления  снаружи пространственной области  во внутрь () с учетом граничных условий (18) и дополнительных соотношений по пространственным переменным. Полученное во всех пространственных ячейках решение  используется затем в качестве граничного значения для направления . Сеточные уравнения системы (21) - (26) решаются последовательно для каждого значения , , используя граничные условия (18), а также условия (19) для направлений изнутри области  наружу ().

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

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

Тип 1.

Освещена одна сторона, например, AB (Рис. 5).

Рис. 5. Освещена одна сторона четырехугольной ячейки.

В этом случае значение  известно. Следующие значения не известны: , ,  и . Используется дополнительное соотношение:

                                         (29)

и два соотношения одного из двух вариантов:

1) ;                                     (30а)

2) ;                                            (30б)

Соотношения (30) понижают порядок аппроксимации схемы в целом. Подставляя (26), (29), (30а) в балансное уравнение (21), получим формулу для вычисления :

Для вычисления  из граничного условия (28) при  используются соотношения (29) и (30а) при :

Тип 2.

Освещены две стороны, например, AB и AD (Рис. 6).

Рис. 6. Освещены две стороны четырехугольной ячейки.

В этом случае значения  и  известны, а значения ,  и  не известны. Используются следующие дополнительные соотношения:

                                        (31)

                                       (32)

Подставляя (26), (31) и (32) в (21), получим формулу для вычисления :

Для вычисления  из граничного условия (28) при  используются соотношения (31) и (32) при :

Тип 3.

Освещены три стороны, например, AB, BC и AD (Рис. 7).

Рис. 7. Освещены три стороны четырехугольной ячейки.

В этом случае известны следующие значения: , , . Значения  и  не известны. Используется следующее дополнительное соотношение:

                                        (33)

Подставляя (26) и (33) в (21), получим формулу для вычисления :

Для вычисления  из граничного условия (28) при  используется соотношение (33) при :

Скалярный поток нейтронов вычисляется по формуле:

                                                (34)

4. Численные результаты.

Приведем результаты численного решения одногруппового стационарного уравнения переноса нейтронов. Рассматривается однородная сфера радиуса . Поток частиц извне отсутствует. Параметры задачи: макроскопическое полное сечение, макро-сечения рассеяния и размножения, а так же внутренний независимый источник нейтронов взяты из работы [1]:

,      ,    ,            (35)

Требуется вычислить распределение частиц, решая уравнение переноса (13) в сферической области.

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

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

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

В одномерном варианте этой задачи решалось сферически-симметричное уравнение переноса:

                    (36)

,       ,                               

с теми же параметрами (35).

В Таблице 1 приведено распределение скалярной плотности потока частиц (34) по некоторым секторам в двумерном расчете по программе NPTQ и значения скалярной плотности потока для одномерных расчетов с совпадающими сетками по радиусу и по углу ( квадратура Гаусса) по программам KIN1D и американской программе ANISN.

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

Из Таблицы 1 видно, что результаты, полученные в двумерной геометрии на сетке, состоящей из выпуклых четырехугольников, и в одномерной геометрии близки друг к другу. Вблизи точки (0,0) имеется разброс значений потока в пределах величины 0.5. Это объясняется тем, что численное решение предполагает разрыв в этой точке для сохранения баланса частиц в системе. Полное число нейтронов в системе равно 145.87 для двумерного расчета и 148.4 для одномерного.

Аналогичные результаты, полученные в работе [1], используя ту же сетку по пространственным переменным и квадратуру  [2] по угловым переменным, представлены в Таблице 2. Результаты расчетов, полученные в двумерной геометрии на сетке, состоящей из выпуклых четырехугольников, и в одномерной геометрии менее близки друг к другу. Вблизи точки (0,0) разброс значений потока составляет 1.23, т.е. существенно больше. Для сравнения полное число нейтронов в системе из работы [1]: 151.2 и 148.2 соответственно.

Таблица 2. Результаты расчетов [1] по двумерной программе и по одномерным программам.

5. Заключение

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

В дальнейшем предполагается обобщить данную схему на нестационарный случай.

6. Литература

1. В.Е. Трощиев, В.А. Шумилин. Разностная схема решения двумерного уравнения переноса на нерегулярных четырехугольных сетках. ЖВМиМФ, № 2, Т. 26, 1986.

2. Carlson B.G. The numerical theory of neutron transport. Methods in Computational physics. Vol. 1, N.Y.: Acad. Press, 1963.