научная статья по теме KP1-ИТЕРАЦИОННЫЙ МЕТОД РЕШЕНИЯ УРАВНЕНИЯ ПЕРЕНОСА В ТРЕХМЕРНЫХ ОБЛАСТЯХ НА НЕСТРУКТУРИРОВАННЫХ СЕТКАХ Математика

Текст научной статьи на тему «KP1-ИТЕРАЦИОННЫЙ МЕТОД РЕШЕНИЯ УРАВНЕНИЯ ПЕРЕНОСА В ТРЕХМЕРНЫХ ОБЛАСТЯХ НА НЕСТРУКТУРИРОВАННЫХ СЕТКАХ»

ЖУРНАЛ ВЫЧИСЛИТЕЛЬНОЙ МАТЕМАТИКИ И МАТЕМАТИЧЕСКОЙ ФИЗИКИ, 2015, том 55, № 10, с. 1727-1740

УДК 519.634

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

© 2015 г. Н. И. Коконков, О. В. Николаева

(125047Москва, Миусская пл., 4, ИПМ РАН) e-mail: nika@kiam.ru Поступила в редакцию 27.01.2015 г.

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

Ключевые слова: уравнение переноса, трехмерные неструктурированные сетки, итерационные методы KPi и DSA.

DOI: 10.7868/S0044466915100154

ВВЕДЕНИЕ

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

Изначально для решения этих систем использовался метод итераций источника (см. [1]), в котором распределение, полученное на очередном итерационном шаге, описывает частицы, испытавшие еще одно рассеяние на атомах среды. Если среда является протяженной и/или слабо поглощающей, среднее число рассеяний, испытываемых частицей, оказывается велико, поэтому большим оказывается и число необходимых итераций источника.

В начале 1960-х гг. были получены первые двухшаговые итерационные методы, в которых на первом шаге выполняется итерация по источнику рассеяния, а на втором шаге, на основании полученного на первом шаге распределения, находятся ускоряющие величины. Эти величины получаются как решение грубого по угловым и/или пространственным переменным приближения уравнения переноса. Двухшаговые итерационные методы зачастую существенно уменьшают число итераций по сравнению с обычным методом итераций источника, поэтому второй их шаг называют ускорением итераций. Наиболее полная информация о развитии двухшаговых итерационных методов для 8м-уравнений из [1] — уравнений переноса с SN аппроксимацией (метод дискретных ординат) — содержится в обзорной статье [2].

На сегодняшний день наиболее популярными методами ускорения итераций являются методы, опирающиеся на огрубление уравнения переноса по угловым переменным, когда ускоряющие величины определяются при помощи решения уравнения типа диффузии. Часто используются метод квазидиффузии Гольдина из [3] — [6], метод средних потоков из [7], [8], методы предобусловливания — KP-метод из [9], [10] и синтетический метод из [11]—[20].

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

1727

как в классическом синтетическом методе (см. [11]), либо как линейная функция угловых переменных (КР1-метод).

В [12] было показано, что в одномерных областях синтетический метод быстро сходится для всех оптических толщин пространственных ячеек, если пространственная аппроксимация 8м-урав-нений для итерации источника и диффузионного уравнения для ускоряющей поправки строится с помощью одной и той же алмазной схемы. Полученный двухшаговый итерационный алгоритм был назван Diffusion synthetic acceleration (DSA). В [13] было показано, что для улучшения сходимости DSA-метода в задачах с сильной анизотропией рассеяния можно использовать не только не зависящую от угловых переменных, но и линейную по ним поправку к решению SN-урав-нений, полученную из решения уравнения типа диффузии. В [14]—[16] была предложена процедура получения уравнений синтетического метода для аддитивной поправки, согласованных с любой пространственной аппроксимацией SN-уравнений. В [10] была модифицирована "four step''-процедура из [16], которая применяется для получения согласованных с взвешенной алмазной схемой Р1-уравнений КР1-метода для регулярной сетки в двумерной области. В [17] был применен DSA-метод к нодальной схеме, в [18] — к билинейному разрывному методу конечных элементов. В [19] и [20] модифицирован DSA-метод для решения SN-уравнений методами длинных характеристик.

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

Предложенный КР1-метод используется для решения системы сеточных уравнений, отвечающих сеточным схемам из [21], аппроксимирующим уравнение переноса на неструктурированной тетраэдрической сетке. Для решения системы уравнений для ускоряющей поправки используется открытая реализация (см. [22]) обобщенного метода минимизации невязки в подпространстве Крылова из [23] с предобусловливанием (см. [24]). Тестирование эффективности предложенного КР1-метода выполняется на известных модельных задачах из [17], [25], [26]. В тестировании оценивается уменьшение как числа итераций источника, так и времени счета благодаря использованию ускоряющей поправки. Также оценивается доля общего времени счета, затрачиваемого на вычисление поправки.

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

1. КРАЕВАЯ ЗАДАЧА ДЛЯ УРАВНЕНИЯ ПЕРЕНОСА НЕЗАРЯЖЕННЫХ ЧАСТИЦ

Рассматривается стационарное уравнение переноса незаряженных частиц

Q) = Q) + q(r,а), г е V, fi еП, (1.1)

где дифференциальный оператор переноса имеет вид

L¥(r, а) = QV ¥(г, а) + at (г)¥(г, а). (1.2)

Интегральный оператор рассеяния определяется выражением

SY(r, О) = — IW,О■О' Щг, О' ДО. (1.3)

4п J

о

Искомая функция ¥(г, Q) определяет распределение незаряженных частиц в точке г выпуклой области V, движущихся в направлении единичного вектора fi еП, D. — единичная сфера. Здесь

at(г) — полное сечение, а(г, Q • ft') — сечение рассеяния, 4п — площадь единичной сферы Q, д(г, Ц) — плотность внутреннего источника.

В выпуклой области G введена декартова система координат, радиус-вектор г имеет координаты (x, y, z). Единичный вектор направления переноса ft в сферической системе координат задается двумя углами — полярным 9 и азимутальным ф, см. фиг. 1. В декартовой системе координат он определяется значениями проекций на оси x, y, z

Qx = sin 0 cos ф, Qy = sin0 sinф, Qz = c os0. (1.4)

Краевое условие для уравнения переноса (1.1) задается для тех точек г границы 5 V области V, в которых вектор направления переноса ft образует тупой угол с внешней нормалью n(r) поверхности в этой точке:

^ «) n.n(№Sv = «>l а5)

где ЭД¥(г, ft)] — оператор отражения.

Краевое условие (1.5) определяет распределение входящих в область V незаряженных частиц. Рассматриваются два типа краевого условия (1.5).

1. Вакуумное условие (частицы не входят в область извне):

^ ,П<0,геЖ = (16)

2. Зеркальное отражение

¥(г,")|,n<0,rev = ArW(T,ft'), íl' = я - 2(ft • n)n, (1.7)

где A e [0,1] — альбедо отражения. Можно видеть, что условие (1.6) может рассматриваться как частный случай условия зеркального отражения (1.7) при альбедо A = 0.

Сечение рассеяния a(r, ft • ft') в точке г выпуклой области Vзависит от скалярного произведения векторов ft' и ft, описывающих направление движения частицы до и после взаимодействия со средой. Зависимость сечения рассеяния от скалярного произведения ft • ft' задается разложением по полиномам Лежандра — ортогональной системе функций на отрезке [—1, 1]:

¿(г)

а(г, a • a') = ^(2/ + 1)ст; (г)р(а • a'),

/=0

где ¿(г) — порядок приближения, ст/(г) — коэффициенты разложения, P¡ (ft • ft') — полином Лежандра степени l.

В приближении первого порядка имеем соотношение

а(г,ft • ft') = а0(г) + 3а1(г)(ПxQ'x + QyQ'y + QzQ'z). (1.8)

2. СЕТОЧНЫЕ АППРОКСИМАЦИИ ДЛЯ УРАВНЕНИЯ ПЕРЕНОСА

При решении краевой задачи (1.1)—(1.5) вводятся сетки по направлениям переноса й (или угловые сетки) и по пространственным переменным х, у, z.

Сетки по направлениям переноса й вводятся набором узлов на единичной сфере

й, {ПхЛ, пуЛ, п ^}, , = 1,..., И.

Каждому узлу отвечает вес . Сумма всех весов полагается равной величине 4 п — площади единичной сферы. Для каждого узла вдоль направления Ц, определяется значение распределения частиц

¥ < (г) = ¥(г, Ц)

и плотности внутреннего источника

Цй (г) = Ц(г, Ц,).

Интегральный оператор (1.3) вдоль направления Ц, имеет вид

^(г, а) = — Го(г, ц ■ ащг, й^й'. (2.1)

п

Уравнение переноса (1.1) в этом случае записывается в форме

< (г) + а, (г)¥ < (г) = Тл (г), (2.2) где правая часть Ей (г) (см. (1.1)) содержит как внутр

Для дальнейшего прочтения статьи необходимо приобрести полный текст. Статьи высылаются в формате PDF на указанную при оплате почту. Время доставки составляет менее 10 минут. Стоимость одной статьи — 150 рублей.

Показать целиком