Участник:S s serov/DSn-метод решения уравнения переноса

Автор текущей версии статьи: Серов Сергей, гр. 417, 2017-ый год.

1 Свойства и структура алгоритма

1.1 Описание алгоритма

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

Геометрическая модель решаемой задачи имеет следующий вид: перпендикулярно единственной участвующей в рассмотрении оси $x$ (будем считать, что она направлена слева направо) расположено несколько бесконечных по осям $y$ и $z$ плоских слоев некоторого вещества (например, урана) заданной толщины, прилегающих друг к другу.

$\mu \frac{dN(x,\mu)}{dx} + \sigma(x)N(x,\mu) = Q(x,\mu),$

где

$Q(x,\mu) = \frac{\sigma_S^0(x)}{2}\int_{-1}^1 N(x,\mu')d\mu' + \frac{1}{2}q(x).$

Используемые обозначения:

• $\mu$ --- косинус угла $\theta$ между направлением полета нейтрона и положительным направлением оси $x$;
• $N(x,\mu)$ --- поток нейтронов, летящих в заданном направлении $\mu$ в точке $x$;
• $SN(x) = \int_{-1}^1 N(x,\mu')d\mu'$ --- скалярный поток нейтронов в точке $x$;
• $\sigma(x)$ --- макроскопическое сечение взаимодействия нейтронов со средой в точке $x$, т.е. коэффициент обобщенного поглощения нейтронов этой точкой среды;
• $\sigma_S^0(x)$ --- макроскопическое сечение изотропного рассеяния нейтронов в точке $x$ среды, т.е. коэффициент обобщенного рассеяния нейтронов в этой точке среды;
• $q(x)$ --- независимый источник нейтронов в точке $x$.

Уравнение переноса является обыкновенным дифференциальным уравнением с интегральным членом в правой части. Целью решаемой задачи является нахождение потока нейтронов $N(x,\mu)$ или скалярного потока $SN(x)$, потому что из него путем подстановки в уравнение можно получить значение потока.

Для каждого слоя вещества необходимый набор параметров, задающих его, известен, а также известны краевые условия на левом и правом концах рассматриваемой системы. В рассматриваемой постановке задачи мы считаем, что в области $[x_{max}, +\infty)$ создан вакуум (то есть поток нейтронов равен нулю по всем направлениям): $N(x_{max},\mu) = 0, \forall \mu \in [-1, 1]$; а на левой границе нейтроны отражаются, изменяя направление своего полета с $\mu$ на $-\mu$, то есть краевое условие имеет вид: $N(x_{min},\mu) = N(x_{min}, -\mu), \forall \mu \in [0, 1]$. Также значения параметров $\sigma(x), \sigma_S^0(x)$ и $q(x)$ будем считать константными в каждом из слоев.

Для численного решения задачи вводится дискретная сетка по оси $x$ и по параметру $\mu$.

1.2 Математическое описание алгоритма

Входные данные:

• геометрическая модель задачи (количество слоев, их расположение на оси $x$ и толщина);

• $\sigma(x)$ для каждого слоя;

• $\sigma_S^0(x)$ для каждого слоя;

• $q(x)$ для каждого слоя.

Параметры:

• $IR$ --- количество ячеек дискретной сетки по оси $x$;

• $M_0$ --- количество ячеек дискретной сетки по параметру $\mu$.

Выходные данные:

• $SN(x)$ для каждой ячейки введенной сетки.

Введем дискретную сетку по $x$ и по $\mu$:

$\mathbb{X} = \{x_i, i \in \{1, \dots, IR + 1\} | x_{i+1} = x_i + \Delta_x, x_1 = x_{min}, \Delta_x = \frac{x_{max} - x_{min}}{IR}\},$

$\mathbb{M} = \{\mu_m, m \in \{1, \dots, M_0 + 1\} | \mu_{m+1} = \mu_m + \Delta_{\mu}, \mu_1 = -1, \Delta_{\mu} = \frac{1 - (-1)}{M_0}\}.$

Середины ячеек соответствующих сеток будем обозначать как $x_{i+1/2}$ и $\mu_{m+1/2}$.

Алгоритм решения задачи выглядит так:

1. Рассмотреть значения $\mu \lt 0$. Для каждого $m \in \{1, \dots, M_0 + 1\}:\mu_m \lt 0:$
   1. Использовать метод итерации источника. Учесть краевые условия в точке $x = x_{IR+1}$. Для каждого $i \in \{IR + 1, \dots, 2\}$ последовательно от $IR+1$ в порядке убывания$:$
      1. Рассчитать правую часть уравнения переноса по формуле$:$

$Q(x_i,\mu_m) = \frac{\sigma_S^0(x_i)}{2}SN(x_i) + \frac{1}{2}q(x_i), SN(x_i) = \int_{-1}^1 N(x_i,\mu_m)d\mu_m = \sum_{m=1}^{\frac{M_0}{2}} N(x_i,\mu_m) \cdot \Delta_{\mu} = \Delta_{\mu} \cdot \sum_{m=1}^{\frac{M_0}{2}} N(x_i,\mu_m).$

1. 1. 1. Решить обыкновенное дифференциальное уравнение с известной правой частью относительно $N(x,\mu_m)$ и вычислить $N(x_{i-1},\mu_m)$, $\bar{N}(x,\mu_m) = N(x_{i-1/2},\mu_m)$ и $SN_{Left}(x_{i-1/2})$ по выбранной разностной схеме. При необходимости использовать правило переключения схем.
2. Рассмотреть значения $\mu \geq 0$. Для каждого $m \in \{1, \dots, M_0 + 1\}:\mu_m \geq 0:$
   1. Использовать метод итерации источника. Учесть полученные на предыдущем шаге значения $N(x_1,\mu_m), \mu_m \lt 0,$ при вычислении краевых условий в точке $x = x_{1}$. Для каждого $i \in \{1, \dots, IR\}$ последовательно от $1$ в порядке возрастания$:$
      1. Рассчитать правую часть уравнения переноса по формуле$:$

$Q(x_i,\mu_m) = \frac{\sigma_S^0(x_i)}{2}SN(x_i) + \frac{1}{2}q(x_i),$ где

$SN(x_i) = \int_{-1}^1 N(x_i,\mu_m)d\mu_m = \sum_{m=1}^{\frac{M_0}{2}} N(x_i,\mu_m) \cdot \Delta_{\mu} = \Delta_{\mu} \cdot \sum_{m=1}^{\frac{M_0}{2}} N(x_i,\mu_m).$

1. 1. 1. Решить обыкновенное дифференциальное уравнение с известной правой частью относительно $N(x,\mu_m)$ и вычислить $N(x_{i+1},\mu_m)$, $\bar{N}(x,\mu_m) = N(x_{i+1/2},\mu_m)$ и $SN_{Right}(x_{i+1/2})$ по выбранной разностной схеме. При необходимости использовать правило переключения схем.
2. В качестве ответа вычислить $SN(x_{i+1/2}) = SN_{Left}(x_{i+1/2}) + SN_{Right}(x_{i+1/2})$ для каждого $i \in \{1, \dots, IR\}$.

Теперь опишем метод решения обыкновенного дифференциального уравнения с известной правой частью и используемые для этого разностные схемы. Также учитываем, что в указанной постановке задачи параметры $\sigma(x), \sigma_S^0(x)$ и $q(x)$ заданы константными для каждого слоя.

Уравнение имеет вид:

$\mu \frac{dN(x,\mu_m)}{dx} + \sigma(x)N(x,\mu_m) = Q(x_i,\mu_m).$

Проинтегрируем его по отрезку $[x_i, x_{i+1}]$:

$\mu_m(N_{i+1} - N_i) + \sigma \int_{x_i}^{x_{i+1}} N(x,\mu_m)dx = \int_{x_i}^{x_{i+1}} Q(x,\mu_m)dx.$

Далее по теореме о среднем:

$\int_{x_i}^{x_{i+1}} N(x,\mu_m)dx = \bar{N}\cdot \Delta_x, \int_{x_i}^{x_{i+1}} Q(x,\mu_m)dx = \bar{Q}\cdot \Delta_x, \Delta_x = x_{i+1} - x_i.:$

Тогда уравнение (называемое уравнением баланса) примет вид:

$\mu_m(N_{i+1} - N_i) + \sigma \bar{N} \Delta_x = \bar{Q} \Delta_x.$

Далее приведем формулы для шага 2 алгоритма (то есть для $\mu_m \geq 0$). Найдем решения этого разностного уравнения по двум разностным схемам. \\ St-схема (шаговая) ($\bar{N} = N_{i+1})$:

$\bar{N} = \frac{\bar{Q} \Delta_x + \mu_m N_i}{\sigma \Delta_x + \mu_m}, N_{i+1} = \frac{\bar{Q} \Delta_x + \mu_m N_i}{\sigma \Delta_x + \mu_m}.$

DD-схема (алмазная) ($\bar{N} = \frac{1}{2}(N_i + N_{i+1}), N_{i+1} = 2\bar{N} - N_i$):

$\bar{N} = \frac{\bar{Q} \Delta_x + 2\mu_m N_i}{\sigma \Delta_x + 2\mu_m}, N_{i+1} = \frac{2\bar{Q} \Delta_x + (2\mu_m - \sigma \Delta_x) N_i}{\sigma \Delta_x + 2\mu_m}.$

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

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

3 Литература

1. Басс Л.П., Волощенко А.М., Гермогенова Т.А. Методы дискретных ординат в задачах о переносе излучения, Препринт ИПМ АН СССР, М., 1986.
2. Карлсон Б.Г., Латроп К.Д. Теория переноса. Метод дискретных ординат. В сб.: Вычислительные методы в физике реакторов, М.: 1972, с. 102-157.
3. Reed W.H. New difference schemes for the neutron transport equation, Nucl. Sci. Eng., 1971, 42, №2, 309-314.
4. Гольдин В.Я., Калиткин Н.Н., Шишова Т.В. Нелинейные разностные схемы для гиперболических уравнений. ЖВМ и МФ, 1965, 5, №5, с. 938-944.
5. Rhoades W.A., Engle W.W. A new weighted-difference formulation for discrete ordinates calculations, Trans. Am. Nucl. Soc., 1977, 27, 776-777.