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

ОБЪЕДИНЕННЫЙ ИНСТИТУТ ЯДЕРНЫХ ИССЛЕДОВАНИЙ

ОД 11-98-295

1* 1\ •-<> 1С93 На правах рукописи

УДК 517.968.23; 519.642; 537.612

возницки

Збигнев И.*

ИТЕРАЦИОННЫЕ РЕШЕНИЯ МНОГОГРУППОВЫХ УРАВНЕНИЙ ДИФФУЗИИ НЕЙТРОНОВ

Специальность: 05.13.16 — применение вычислительной техники, математического моделирования и математических методов для научных исследований

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

*Институт атомной энергии, 05-400 Отвоцк-Сверк, Польша

Дубна 1998

Работа выполнена в Институте атомной энергии Сверк, Польша, в Лаборатории вычислительной техники и автоматизации Объеди ненного института ядерных исследований, г.Дубна.

Официальные оппоненты доктор физико-математических наук

профессор Вячеслав Иванович Лебедев

доктор физико-математических наук профессор Александр Витальевич Крян

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

Ведущая организация Институт атомной энергетики, г.Обни

Защита диссертации состоится "____года

в на заседании Диссертационного совет

Д047.01.04 при Лаборатории вычислительной техники и автомат зации Объединенного института ядерных исследований по адрес 141980, г.Дубна Московской области.

С диссертацией можно ознакомиться в библиотеке ОИЯИ. Автореферат разослан 1998г.

Ученый секретарь

Диссертационного совета Д047.01.04

кандидат физико-математических наук [/^/Ос*., З.М.Иванче!

(лХи

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

1. Актуальность темы диссертации

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

Во многих физических задачах нейтроны могут рассматриваться как среда нового типа - нейтронный газ. При этом уравнение Больцмана становится систематическим "каталогом" всех возможных способов входа, выхода, образования и поглощения частиц этого газа, причём•в каждой точке простра! 1ства имеет место распределение скоростей нейтронов по всем направлениям.

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

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

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

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

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

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

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

2. Цель работы

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

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

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

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

3. Научная новизна работы

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

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

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

4. Практическая ценность работы

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

также полезным при классификации предфакторизационных методов, других авторов , представленной в таблице 1.1.

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

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

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

Алгоритмы метода AGA с двойной верхней релаксацией были реализованы для двухмерного и трёхмерного случаев в компьютерных программах HEXAGA-II и HEXAGA-III, решающих уравнения диффузии нейтронов на треугольной и гексагональной разностных сетках на основе стандартной стратегии внешних-внутренних итераций. Программы HEXAGA, используемые интенсивно для расчётов связанных с проектированием и эксплуатацией ядерных реакторов в Польше и других странах (в Бельгии, Болгарии, Германии), делают возможным получать результаты вычислений в несколько раз быстрее чем другие компьютерные программы этого типа, и их преимущество возрастает с ростом числа точек разностной сетки. Сравнение. результатов расчётов HEXAGA-II и английской программы SNAP для реакторной задачи ВВЭР-1000 данное в следующей таблице (4 сек CYBER 73 s ю сек ЕС-1040).

Число сет. Число Время на

точек на энерго- CPU одну энерго-

Программа одной гек- -простр. ^е Г f время простр. Компьютер

сагональной точек сек. точку

кассете милисек.

1 288 1.11715 6 0.0208

12 1024 1.11375 14 0.0136

HEXAGA-II 27 2208 1.11279 34 0.0154 CYBER 73

48 3844 1.11240 60 0.0156

108 8462 1.11211 237 0.0280

1 92 1.11387 26 0.283

SNAP 6 552 1.11068 150 0.272 ЕС-1040

24 2208 1.11124 834 0.405

96 8832 1.11175 6435 0.729

Эти результаты были представлены на конференции по физике реакторов в Knoxvi1le, США, в 1985г.

В Польше программа HEXAGA-II была использована в расчётах реактора ядерной электростанции "Жарновец", которые проводились в 1980 годах, а также при физических расчётах исследовательских реакторов EWA и MARIA, необходимых при разработке отчёта безопасности и текущей эксплуатации реактора. Проверка результатов вычислений, полученных при помоши программы HEXAGA-II на основании данных по эксплуатации реактора EWA показала великолепную согласованность. Надо отметить, что для реактора EWA и ядерной электростанции "Жарновец", имеющих гексагональную структуру сетки топливных элементов, имело место полное совпадение с геометрией треугольной сетки программы HEXAGA-II. В поперечном сечении реактора MARIA каналы твэлов размещены в квадратной сетке, а именно применение мелкой треугольной сетки программы HEXAGA-II позволило воссоздать круглый вид каналов твэлов и детали сложной геометрическо-мате-риальной конфигурации реактора MARIA, более точным образом, чем это было бы возможно в программе с прямоугольной ориентацией сетки. Программа HEXAGA-II в настоящее время используется в расчетной системе разработанной для нужд реактора MARIA.

Для того, чтобы представить подготовку входных данных программы HEXAGA-II для реакторных задач при мелкой сетке были разработаны две вспомогательные программы HEXI-22 и HEX1-23. Эти вспомогательные программы предназначены для производства первой части входных данных программы HEXAGA-II (связанных с описанием материально-геометрической конфигурации точечной сетки) для санной реакторной задачи, в которой шаг однородной треугольной сетки уменьшается в два раза в случае программы HEXI-22 и в три раза в случае программы HEXI-23. Эти программы используют такие же входные данные, как программа HEXAGA-II без ввода добавочной входной информации. Каждые выходные данные программ HEXI могут быть использованы как новые входные данные программ HEXI. Таким образом, если мы имеем входные данные программы HEXAGA-II, описывающие данную рактор-ную задачу с минимальным числом точек сетки, и используем любую комбинацию входных/выходных данных программы HEXI-22 и/или HEXI-23 мы можем производить

входные данные программы HEXAGA-II для расположения точек сетки для этой

„¡„к

задачи с шагом сетки уменьшенным согласно следующим множителям: 2 3 для любых целых чисел i,k г О.

В случае трехмерной программы HEXAGA-III были разработаны аналогичные вспомогательные программы HEXI-32 и HEXI-33.

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

тических вопросов физики реакторов ВВЭР. Кроме того, автор тесно сотрудничал с Институтом ядерных исследований и ядерной энергетики в Софии, где была разработана болгарская версия программы HEXAGA, известная под названием НЕХАВ, и применяемая при расчётах связанных с эксплуатацией энергетических ядерных реакторов ВВЭР работающих в Болгарии.

Актуальные версии программ HEXAGA работают на основе, предложенной в диссертации, стратегии глобальных-внешних-внутренних итераций.

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

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

1. NPY Seminar on Numerical Solution of Multi-Dimensional Neutron Diffusion Equations, Warsaw, Poland, March 1969;

2. Seminar on Reactor Physics Calculations, Budapest, Hungary, October 1969;

3. The INEA Third International Advanced Summer School in Reactor Physics, Hercegnovi,' Yugoslavia, September 1970;

4. Reaktortagung, Dusseldorf, West Germany, March 1976;

5. Conference on Numerical Methods of Mathematical Physics, Novosibirsk, Soviet Union October 1980;

6. ANS/ENS International Topical Meeting on Advances in Mathematical Methods for Nuclear Engineering Problems, Munich, Germany, April 1981;

7. International Meeting on Advances in Nuclear Engineering Computational Methods, Knoxville, USA, April 1985;

8. International Topical Meeting on Advances in Reactor Physics, Mathematics and Computations, Paris, France, April 1987;

9. International Conference on Linear Algebra and Applications, Valencia, Spain, September 1987;

10. International Conference on the Physics of Reactors: Operation, Design and Computations, Marseille, France, April 1990;

11. Copper Mountain Conference on Iterative Methods, Copper Mountain, Colorado, USA, April 1992;

12. The Second Conference of the International Linear Algebra Society, University of Lisbon, Lisbon, Portugal, August 1992;

13. '92 Shanghai International Conference of Numerical Algebra and its Applications, Fudan University, Shanghai, P.R.China, October 1992;

14- Two Matrix Theory Meetings in Technion, Haifa, Israel, June 1993;

15. Third SIAM Conference on Linear Algebra in Signals, Systems and Control, Seattle, USA, August 1993;

16- The International Conference on Reactor Physics and Reactor Computations, Tel Aviv, Israel, January 1994;

17- ICIAM'95 - The Third International Congress on Industrial and Applied Mathematics, Hamburg, Germany, July 1995;

18- WCNA'96 - The Second Congress of Nonlinear Analysts, Athens, Greece July 1996;

19. PHYSOR'96 - International Conference on the Physics of Reactors, Mito, Ibaraki, Japan, September 1996;

20. Czech-U.S. Workshop on Iterative Methods and Parallel Computing, Milovy, Czech Republic, June 1997;

21. 15th IMACS World Congress on Scientific Computation, Modelling and Applied Mathematics, Berlin, Germany, August 1997;

22. XII Conference on Applied Mathematics PRIM'97, Pali6, Yugoslavia, September 1997;

23. Sixth SIAM Conference on Applied Linear Algebra, Snowbird, Utah, USA, October 1997;

24. 7th Russian Conference on Grid Generation Problems and Numerical Methods in Mathematical Physics, Durso (Krasnodar region), Russia, September 1998.

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

Диссертация состоит из введения и шести глав, содержит 20 таблиц, 60 рисунков, список цитированной литературы и изложена на 170 страницах машинописного текста.

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

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

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

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

В шестой главе сформулированы основные выводы диссертации.

7. Содержание работы

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

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

Рассмотрим итерационное решение системы линейных уравнений

Ах = Ь, (1!

где А е Рпхп является неособенной натрицей и х,Ь е к".

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

А = К - N (2)

с неособенной матрицей М, приближённое решение хш в этом случае имеет вид

ИхаП) = КхШ + Ь. 1>-0, (3)

что эквивалентно

с(и11 = мЛ1ха) ♦ И"1!.. £ * О, (4)

(0)

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

Анализ сходимости указанного метода основан на спектральном радиусе итерационной матрицы р(М_1Ы) = тах (А,I, где Л. являются собственными значениями матрицы М_1Ы. Итерационный метод сходится к единственному решению

х = А_1Ь (5)

для каждого х<0), если, и только если, е(м"'ю < 1. Для больших значений 1 на каждом итерационном шаге величина погрешности решения уменьшается приблизительно на коэффициент р(М_1Ю. Чем меньше р СМ"1Ы), тем быстрее сходимость. Следовательно, оценка итерационного метода зависит от двух условий: М надо выбрать таким образом, чтобы можно было легко найти обратную матрицу, а р{и"'ы) должно быть настолько малым, насколько это возможно.

Определение 1.1 Разложение А = Н - N называют сходящимся расщеплением матрицы А, если А и М - неособенные матрицы и рШ^Л) < 1.

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

Лемма 1.1 Пусть А = М - N - расщепление матрицы А. Если А и М являются неособенными матрицами, то

М'^А"1 = А"1™"1, (6)

матрицы м"1Ы и А^Ы, а также ММ"1 и ЫА-1, коммутируют.

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

называют предфакторизационером матрицы А или предкондиционером в случае метода сопряжённых градиентов. Матрица М"1 является матрицей предобращения матрицы А и называется предобразителем матрицы А приближающим матрицу А"1. Следовательно, N = M - А - это остаточная матрица, полученная при помощи предполагаемой предфакторизации матрицы А. Если N ненулевая матрица, то осуществляется частичная (неполная) факторизация матрицы А и решение имеет итерационный вид. В случае, когда матрица N становится нулевой, осуществляется точная факторизация матрицы А и решение уравнения (1) получается при помощи прямого метода исключения Гаусса.

Идея решения систем линейных уравнений с разреженными матрицами при помощи предфакторизационных методов не нова. Первые такие подходы под названием "примитивные итерационные методы", ограниченные симметрическими матрицами, были рассмотрены Varga'ой в i960 году. Весьма сходная методика в то же время была независимо введена Булеевым и Oliphant'ом. На самом деле методы Булеева были пионерскими и являлись существенным вкладом в итерационное решение пятиточечных и семиточечных уравнений конечных разностей, которые приближают двух- и трехмерные краевые задачи математической физики.

Развитие методов предфакторизации AGA, известных также под названием "двухпроходные итерационные методы AGA", началось в 19В8 году в виде первой и наиболее простой версии, называемой методом EWA. С тех пор эти методы успешно применяются для решения систем линейных уравнений, возникающих из конечно-разностного приближения задач диффузии нейтронов. Они удобны на практике, и в случае, когда применяются для расчёта ядерных реакторов, дают весьма неплохие результаты.

Допустим, что матрица А определена как

А = К - L - и, (7)

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

M = [I - CL+H)D"1]D(I - D-1(U+Q) ] . (8)

где H и Q соответственно добавочная строго нижняя треугольная и добавочная строго верхняя треугольная матрицы. Предполагается, что D неособенная диагональная матрица неявно определённая при помощи соотношения

D = К - diag{ (L+H)D~JCU+Q) } . (9)

Легко проверить, что

N = offdlagi (L+HJD1 CU+QJ > - H - Q, (10)

где выражение diag{В} означает диагональную матрицу с диагональными элемен-

тами, тождественно равными диагональным элементам матрицы В, и

offdiagiB} = В - diag{B>.

Итерационный метод теперь можно записать следующим образом

(ttl) „ (t) „-i, t „ х +ИЬ, í г о

и

3 = [I - D 1 (U+Q) ]D '[I - (L+H)D"']N.

(12)

Поскольку (L+H1D и D (U+Q) соответственно строго нижняя треугольная и

-i

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

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

(13) представляют собой общий вид широкого класса так называемых: предфакто-ризационных или двухпроходных итерационных методов. Каждый из них однозначно определён с помощью выбора матриц Н и О. Матрицу 5 называют предфакториза-ционной или двухпроходной итерационной матрицей.

Классификация предфакторизационных итерационных методов по отношению к выбору матриц Н и 0 представлена в таблице 1.1.

Метод АСА определяет широкий класс алгоритмов, каждый из которых задается путём выбора матриц Нд и <2Д таким образом, чтобы Нд + <зд * О (исключая случаи, в которых Нд = -I. и <3Д = -и) при предположении, что размещение ненулевых элементов матрицы ЫА не совпадает с размещением ненулевых элементов матрицы Ь + Нд + и + <3Д. Предположение о том, что матрицы НА и 0Д соответственно строго нижняя треугольняя матрица и строго верхняя треугольняя матрица, всегда позволяет определить в явном виде значения ненулевых элементов матриц Нд и 0Д, для произвольно принятой модели размещения ненулевых элементов в этих матрицах, непосредственно следует из неявного вида матрицы (Ь+НдШ^(и+<2д). Другими словами, при постулированной модели распределения ненулевых элементов в матрицах Нд и <3А, все ненулевые элементы матриц Нд, (Зд а следовательно, и 0Д и Лд, могут быть вычислены непосредственно при помощи приравнивания их соответствующим элементам неявной формы произведения матриц

(t+i)

для уменьшающихся индексов в проходе обратной подстановки. Уравнения

х

(L+HA)D¡ (U+QA).

-1

Таблица 1.1 Классификация предфакторизационных методов (A=K-L-U=M-N)

н q M=[I-(L+U)D"1]D[I-D"1(U+Q)} D=K-diag{(l+h)d~1(U+Q)> n=offdiag{(L+H)D~1 (U+Q)} Метод

-l -и mj=k dj=k nj=l+u Jacobi

-l 0 м^ш-к^и] dc=k n0=l Gauss-Seidel (backward order)

о -и m0=k[i-k"1l] dc=k ng=u Gauss-Seidel (forward order)

0 0 m^ii-ld^idE-ii-d^u] De=k-diag{LD21U> NE=offdiag{LD^1U> EWA Woinicki (1968)

о l= о =ит ^[I-uX'lDvlI-D'/u] Dv=K-djag{i;TDy Hv=offdiag{UTDy1U> Varga (1960)

0 О MB=tI-LD^' ]DB[I-D^U] DB=( 1+e IK-diagiLDjj'u} NB=offdiag(LD^1U}+eK (e=0, MB=ME, Nb=He) Buleev (1960)

(g-l)L 0 M0=t I-gLDg1 ]D01I-D^U] D0=K-d iag{gLDg1U} H0=offdiagigLD^)1\])-(.g-l )L (g=0, M0=MC, N0=NC; g=l, H0 N0=Ne) Oliphant (1962)

aL aU Ms=[ I-iL+aDDg1 ]DS [ I-Dg1 (U+aU)] Ds=K-d i ag {(L+aL)Dg1(U+aU)} ns=of/diag{{L+aLJDg1(U+aU)} (a=0, hs=he, ns=ne) -aL- aU Stone (1968)

"д Qa ma=[i-(l+ha)d;1jda[i-d;1 DA=K-diag{(L+HA)D^1(U+Qa (U+Qa)] )} NA=offdiag{ (L+HA )D;1 (U+Qa ) > -H*- QA AGA WoZnicki (1970)

HD Qd Md=[i-(l+hd)dp1]dd[i-d^1 dd=k-diag{(l+hd)d;)1(u+qd (U+Qd)] )} ND=0 (offd i ag < ( L+Hd )Dp1 ( U+Qd ) } =H D+Qc ) Direct" method Woinicki (1973)

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

В случае алгоритмов АСА, уравнения (13) могут быть записаны следующим образом.

y,ul) = 01[(LA+HA)D;,y<ttI) ♦ Naxíu + с] - (uj-1 )уlt) x(ul> = ц п;,[(0Л)«,м1 + y't+1,l - (u2-l )x(tl t*>

(14)

где ш, и иг являются параметрами релаксации.

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

Рдж"+1) = вдх«" ♦ 5Ах"-и . ь. г>0. (15)

соответствующие двухкратному расщеплению матрицы

А = РЛ + + Ба, (16)

где

рА = гАг (1-и1(1-+нЖюА[1-о,г(и+<1А)}, (17)

НА = {(J^t - (<J2-1)[I-<J1(L+Ha)D¡1]Da - lu1-l)[I-íJ:,D'1lV+QA)n,

(18)

(и,-1)(иг-1)

= -¡Tu- D*' " Ul * * (1S)

12

Так как для (J¡ * О и ш2 * О, РА является неособенной матрицей, уравнение (14) может быть записано эквивалентно, как

х"*1' = р;4*(" + PlV1"1' * р;гь. t > О. (20)

Этот метод называется методом двойной верхней релаксации AGA или, для краткости, методом АСА двойной SOR (AGA double SOR method). Свойства сходимости этого петода, основанные на таком же подходе как в случае однократного расщепления матрицы А = М - N, подробно проанализированы в [12] и рассмотрены в главе 4.

Очевидно, что этот метод сводится к однократному методу верхней релаксации SOR - или к методу прямово рода (the forward type), когда ut = 1 и и2 * 1 или к методу обратново рода (the backward type), когда ut * 1 и u2 = 1 для которых SA = О.

Эффективность двойной верхней релаксации AGA зависит от выбора параметров релаксации и о2. Скорость сходимости имеет максимальное значение тогда, когда cot = ш2 = иор1.

Решение задачи оценки оптимального параметра релаксации в этих методах является - с теоретической точки зрения - существенной трудностью, потому что итерационные матрицы в методах EWA и AGA имеют также отрицательные и комплексные собственные значения. Предлагаемый в работе [12] подход к предварительной (a priori) оценке оптимального значения uopt основан на численном анализе задачи Дирихле для эллиптических уравнений. Численные эксперименты подтверждают высокую эффективность предлагаемой процедуры для оценки оптимального параметра релаксации в методах AGA для практических применений.

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

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

Определение 1.2 Для матрицы А разложение А = М - N называют:

(a) регулярным расщеплением матрицы А, если Н4 г О и N г О,

—i -i

(b) неотрицательным расщеплением матрицы А, если М г о, М Н г о и NM"1 г О,

(c) слабым неотрицательным расцеплением матрицы А, если М-1 г 0 и: или M-1N г О (первого рода) или NM"1 £ О (второго рода),

(d) слабым расщеплением матрицы А, если М является неособенной матрицей и: или M-1N £ о (первого рода) или NM-1 г о (второго рода). В частности данное слабое расщепление может быть и первого, и второго рода одновременно,

(e) сходящимся расщеплением матрицы А, если М является неособенной матрицей и p(M_1N) < 1.

Определение, принятое в пункте СЬ), эквивалентно определению слабого регулярного расщепления матрицы А, предложенного Ortega'ой и Rheinboldt'oM. Однако, И"1 г о и только M_1N ь О (без условия NM"1 а О) определяется другими авторами как слабое регулярное расщепление матрицы А, но, в этом случае необходимо ввести добавочные предположения в теоремах сравнения, когда соотношение Mi1 i Vl~2 = О используется в качестве гипотезы. Следует отметить, что применение терминологии Ortega'и и Rhetnboldt'a "слабое регулярное" в пункте (Ь) является причиной путаницы при использовании названия расщепления в пункте (с). Употребление термина "неотрицательное" позволяет различить эти два случая и избежать неоднозначности. Определение слабого расщепления патрицы А для случая первого рода было введено Marek'ом и Szyld'oM.

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

Теорема 1.1 Пусть А = М - N будет слабым неотрицательным расщеплением матрицы А. Если А"1 £ О, тогда

2. если M_1N г 0, тогда A > M_1N и если NM"1 г о, тогда NA"1 £ NM"1

. р (A_1N) p(NA-1)

3. р(М JN) = —-5-----— < 1.

1 + р(А N) 1 + р(А Н)

Обратно, если дСН^Ю < 1, тогда А"1 г о.

Зависимость в пункте 3, подобная результату полученной Varga'ой для регулярных расщеплений матрицы А, была доказана при помощи коммутативных свойств A~IN и M"1N (или NA"1 и NM"1) представленных при помощи леммы 1.1. Как следствие этой теоремы можно сформулировать следующий вывод.

Вывод 1.1 Каждое слабое неотрицательное расщепление матрицы А при А"1 г О является сходящимся расщеплением и наоборот, для каждого сходящегося слабого неотрицательного расщепления матрицы А, А"1 г О.

Много теорем сравнения для различных расщеплений монотонной матрицы А, (т.е. А"1 г О) автор доказал при простых, а также при более сложных, предположениях. Некоторые из этих результатов приведены ниже.

Теорема 1.2 Пусть А = Mj - Hj = Мг - N2 суть два слабые неотрицательные расщепления матрицы А разных родов, то есть, или Mj1N1 г О и N^a1 г О, или И^1 г о и M¡'n2 & О, где А"1 г О. Если М^1 г

тогда

е(м^Ч) ^ е1*Ч\). (21)

В частности, если А"1 > О и М^1 > М^1, тогда

е(м;\) < р^1^). (22)

Теорема 1.3 Пусть А = М! - ^ = М2 - N3 суть два слабые неотрицательные расщепления симметрической матрицы А, где А 1 г о. Если по крайней мере одна из матриц М1 и Мг является симметрической матрицей, и если М^1 г Мг1, тогда

рШ^Х) ^ р(М>2). В частности, если А"1 > О и М^1 > М21, тогда р(М"1М1) *

Отдельный класс естественных условий, важных в приложениях представлен

ниже.

1. А 2: М

Теорема 1.4 Пусть А = М, - Nt = М2 - Пг суть два слабые неотрицательные расщепления матрицы А, одного и того же рода, то есть, или Mí'N! 2 о и М^Лз г О, или NjMj1 г О и N^^1 i 0, где А"1 £ О. Если Nj - Nt, тогда

etM^Nt) i e(M¡'N2). .(23)

Теорема 1.5 Пусть А = Mj - N, = М2 - N2 суть два слабые неотрицательные расщепления матрицы А, одного и того же рода, то есть, или Mj'n, — О и M^n;, 2: О, или NjMj1 2: О и N2M2' í О, где А"1 г О.

Если М^1 г (Mj1)1, тогда

píH^Nj) — p(M2'N2). В частности, если А"1 > 0 и М^1 > (М^1)т, тогда £(»!>,) í e(M¡%).

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

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

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

Young'он. В случае медленно сходящихся задач алгоритм Sigma-SOR имеет очень сильные преимущества по сравнению с общеиспользуемым адаптационным методом, в котором оптимальный параметр релаксации uopt определяется "динамически" во время выполнения итерационного процесса.

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

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

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

с

-VDg(rWg(r) + £g(r)0g(r) - [ Eg.^(r)0g,(r) =

я'*д

с

= I ^ГдЛг)фчЛг), g=2,2,..G, (24)

q' = l

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

D,(r) dn + ag(r)<pg(r) =0, г е Г, (25)

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

В области Q строится сетка fih, аппроксимирующая область fi. Затем на этой сетке происходит аппроксимация производных, входящих в уравнение (24) и краевых условий (25). Наиболее часто используются стандартные конечно-разностные схемы низкого порядка, которые сводят задачу (24) и (25) к системе линейных алгебраических уравнений. Эта система может быть преобразована в следующую матричную задачу на собственные значения

Е <р = |xfV (26)

л

Порядок квадратной матрицы

Е = А - Sd - Su (27)

равен числу точек сетки пространства-энергии, OxN (G число групп энергии и N число точек пространственной сетки), где матрицы A, Sd и Su представляют

соответсвенно члены диффузии-захвата, рассеяния вниз и рассеяния вверх. Кошюнснти вектора <р представляют приближённые значения потока нейтронов ф в узловых точках сетки пространства-энергии. Спектр деления и член продукции представлены соответсвенно с иомоиыо матриц X и Г, где верхний индекс Т обозначает транспозицию матрицы.

Свойства уравнения (26) были исчерпывающе исследованы в литературе и оно имеет единственный положительный собственный вектор <р, и соответствующее одиночное положительное собственное значение X! (интерпретированное как эффективный коэффициент размножения кеГГ), большее абсолютного значения любого другого собственного значения. Матрица Е неособенная и монотонная, т.е. её обратная матрица неотрицательна, и задача на собственные значения (26) может быть записана в следующем виде

Ар = Вр, (28)

где

В = Е_,ХГТ г О. (29)

Широкий класс методов для определения собственного значения с наибольшим модулем А, и соответствующего собственного вектора р, основан на степенном методе, в котором последовательные приближения для Л, и ?>, определяются с помощью вычислительного процесса

<р (1*1) = ХШВ*,П') (30)

и

А (1*1) = А (1)--, (31)

где 1 является номером итерации и две нормы как максимальная |*| так и эвклидовая |«| применяются наиболее часто. Так как наибольшее (по модулю) собственное значение А, неотрицательной матрицы В является положительным и простым, то степенной метод является сходящимся процессоп для почти случайно выбранного неотрицательного начального вектора р(0), т.е.

\(1) —» А,, <р(1) —> р1 ПРИ 1 —> со.

Скорость сходимости в степс! гном методе зависит от поддонинангного отношения

б(В) = и^/А,, (32)

предполагая, что собственные значения А, матрицы В упорядочены таким образом, что

А! > |Х2| г )А3| ^ ... . (33)

Чем меньше поддопинантное отношение, тем сходимость более быстрая. Предпола-

гая, что для достаточно больших значений 1, I) приблизительно равно Л,, можно равнозначно рассмотреть задачу

<р(1+1) = В1^1)<р(1), (34)

где

В(Л,) = | Е_1ХГТ (35)

Л1

имеет собственные значения уБ)1 = А^Л! и согласно вышеуказанному порядку 1>вд = р[В(Л,)] = 1 > |ув>2| г: |1/в з| г.....

и

б(В) = б[В() ] = |1>в>21 . (36)

В рассуждениях приведенных выше, было принято предположение, что матрица Е-1 дана в явном виде. Однако в случае одномерных задач, когда нет процесса рассеяния нейтронов вверх, т.е. Б" = О, часто возможно вычислить непосредственно матрицу Е"1. В дву- и трёхмерных задачах прямое обращение Е оказалось сильно затруднено и, следовательно, должна быть применена специальная итерационная стратегия.

Уравнение (26) вместе с уравнением (27) можно записать как

Ар = Б<1Р + + (37)

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

( м ^ _

АО)

А <р(1+1) = 5*(р(1 + 1) + Би<р(1) + -1-^Т<р(1) (38)

или эквивалентно

<р(1+1) = А~Х15Л<р(1+1) + + ^уХРХ1р(1)]. (39)

Факт, что Б11 является точно нижней треугольной блочной матрицей позволяет использовать в вычислительном процессе компоненты р, полученные в 1+1 итерации. После выполнения расчётов для итерации 1+1, новое значение Х(1+1) оценивается с помощью формулы (31). Вышеуказанное уравнение может быть представлено в следующем виде

<р(1+1) - Т(А(1))р(1), (40)

где

тиш) = II - А'^ГЧ"1 [Би + дщИ-1]- (41)

Как можно заметить, значения элементов матрицы КАОЛ меняются при переходе от данной итерации к очередной итерации в зависимости от того, кейс меняется значение ЛП) и для достаточно больших значений индекса 1, \(1) становится близким к А] и

Т(А.Ш) -> ТСЦ) = [I - А"18,1]"1А"1(5и + 1хкт] г О. (42)

Предположим, что собственные значения матрицы КА,) упорядочены следующим образом

= ) ] = 1 > |утг1 £ |УТ1з1 г: ...,

тогда для достаточно больших значений индекса 1, скорость сходимости итерационного процесса (2.20) зависит от поддоминантного отношения

5 [Т( ) ] = \1>т>г\- (43)

Очевидно, что при отсутствии верхнего рассеяния, т.е., когда = 0,

Т(А,) = ВС А!).

Итерационный процесс представленный при помощи уравнения (2.20) известен как "итерация источника деления" или как "внешняя итерация". Однако, при проведении итераций необходимо знать вид обратных подматриц Ад, образующих блочно-диагональную структуру матрицы А. В одномерных задачах Ад, являются трёхдиагональными матрицами, которые легко обратимы с помощью хорошо известной процедуры прогонки, прямого исключения - обратного подставления. Для дву-и трехмерных задач у матриц Ад более трёх диагоналей, более того их порядок на много больше, чем в одномерном случае; так что обращение матриц Ад становится очень трудным или просто оно невозможно из-за ошибок округления. В этом случае приближённое обращение матриц Ад обычно проводится итерационным методом используя ряд внутренних итераций. Стратегии, основанные на различных уровнях итераций, представлены в следующих главах.

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

Неоднородные решения многогрупповых задач диффузии нейтронов связаны с итерационным решением системы линейных уравнений

Аф = с, (44)

которые могут быть представлены в следующем виде

Кфа*и = N фм + с, 15 0, (45)

где фа) означает очередную итерацию и

А = М - N (46)

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

ф = Л~1с (47)

для каждого ф'0> тогда и только тогда, когда М является неособенной матрицей

и соответствующая итерационная патрица

5 = М-'и (48)

имеет р(5) < 1. Уравнение (45) может быть записано в эквивалентном виде

= °§фш + М_1с, г г о (49)

. (О)

или с использованием ф

.ИМ) „1*1 .(0) £-1 , —

ф = <§ ф + М(1)с, I г 0, (50)

где

I

М^ = (I + § + ГЧ . .. + = ^э'м"1. (51)

1=0

Из уравнения (3.3) следует, что отсюда

А"1 = (I - Н Чп"1!»"1 = (I - S)"V\ (52)

И"1 = (I - g)A~'. (53)

Подставляя уравнение (53) в (51), получают

Sai = (I " S1*1)*"1- (54)

Так как при предположении р(S) < 1, S1 достигает нулевой матрицы при £ —> ю, и последовательно Й"^ —> А"1, и решение уравнения (50) стремится к единственному решению определённому с помощью уравниения (47) для любого ф<0). В анализе сходимости итерационных методов (асимптотическая) скорость сходимости

R(i?) = - lnp(S) (55)

является наиболее простой практической мерой скорости сходимости для матрицы & и особенно полезной для сравнения эффективности различных итерационных методов.

Матрицу Н~| можно рассматривать в некотором смысле как неполную обратную форму матрицы А, приближающую А"1 поспе t итераций. Матрица была названа автором предобратной матрицей t-степени (prelnvertloner of t-degree) матрицы А, где M"^ = A"1, a = M"1 является предобратной матрицей 0-степени матрицы А. Свойства сходимости обычно изучаются с помощью проверки вектора ошибок определяемого при помощи

е<1+1> = ф - ф1м>. (56)

В частности, для ф10) = О, получают

е = А с - M(t)c,

отсюда с помощью уравнения (54)

(t+i>

Таким образом, решение уравнения (50) можно сравнить с решением, полученным

К,-

для предобратной матрицы М"' , приближающей А"1 после t итераций. Когда

ф(0> = О, его вектор ошибок по отношению к единственному решению (47) определяется уравнением (57).

Принимается расщепление матрицы А, как дано в (46), Уравнение (37) может быть записано в следующем виде

Н(Р = Np + sd<p + S"<p + ¿XFTf> (58 ) Л

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

И <р(1+1) = N <р(1) + sVfi + í J + S "<p(l) + ^jjjXTTip(l). (59)

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

<p(l + l) = M'^sVCJ+í,) + Щ(1 ) + S а<р(1) + ^jjjXFT<p(l)l (60)

или эквивалентно

ра*1) = 4(Ml))tp(l), (61)

где

VM1» = [I - M'VW + M_1(SU + ^yXF1)] (62)

и итерационная матрица V = М_1М связана с расщеплением матрицы А. Для достаточно больших значений индекса 1, А(1) ~ At и

УШ1)) V(At) = [I - M~V]_1[V + «"'(S" + | )CFT)], (63)

Ai

и предполагая, что собственные значения матрицы VUj) упорядочены следующим образом

"v,i = PlV(Aj)] = 1 > \vva\ г |i>V3| г ...,

скорость сходимости в итерационном процессе (61), представляющем собой метод называемый автором стратегией глобальных итераций (1) для однократных расщеплений, зависит от поддоминантного отношения

StVU^J = |1>у21, (64)

где индекс i относится к глобальным итерациям, которые являются просто итерациями степенного метода, и М"1 представляет собой предобратную матрицу 0-степени матрицы А.

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

R(V) = - 1пб[V), (65)

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

Спектр собственных значений матрицы ) зависит от предположенного расщепления матрицы А. Если А = М - N представляет неотрицательное расщепление матрицы А (т.е., М~' г о, íT'n г О и NM"1 г 0) или в частности регулярное расщепление матрицы А (т.е., М-1 2 0 и М г о), то У(Лг) является неотрицательной матрицей, которая обеспечивает сходимость итерационного процесса (61) для.всех р(0).

Обычно 6лочно диагональная матрица А определяется при помощи следующего разложения

А = К - L - U, (66)

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

М = К - L и N = U. (67)

Вышеуказанные уравнения представляют регулярное расщепление матрицы А и соответствующая итерационная матрица имеет вид

J?! = M_1N = [I - K'^rVu 2 о. (68)

Метод последовательной верхней релаксации SOR, тесно связанный с методом Гаусса-Зейделя, представлен с помощью следующих матриц расщепления

ми = ¿k[i - (jk_1l] и nu = i[uu - (u-l)k] (69)

и связанная итерационная матрица может быть записана следующим образом

£и = MuNu = [I " (oK^Li'V'tuU - (ш-1)1], (70)

где и параметр релаксации. Очевидно, что для u = 1, вышеуказанные уравнения сводятся к уравнениям (67) и (68), представляющим метод Гаусса-Зейделя. Как хорошо известно, < 1 для всех 1 < о < 2 и, когда А является циклической

последовательно упорядоченной матрицей с индексом 2 (2-cyclic consistently ordered matrix), существует

2

"opt = -/ <71)

1 + / 1-р(£г)

доводящее до минимума значения р(£и). Эффективный метод оценки uopt представлен в [15].

Собственные значения матрицы V(Aj), определённой при помощи уравнения (63), удовлетворяют уравнению

[I - M'VrVlN + Su + Í XFT)x = Ux (72)

или эквивалентно

[N + us" + Su + y XFTJx = vHx (73)

и подставляя уравнение (67) в (73), получаем

K~lli>(L + Sd) + U ♦ Su + i ХГт]х = ух. (74)

Ai

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

V(j(A,)у = [I - M^Vf'K^CM,,, + Su + I XFT)ly = Т1У (75)

или эквивалентно

[Hu + T)Sd + Su + i XFT]y = TjHjjy (76)

и после подстановки уравнений (69) в уравнение (76), получаем

if'lTja + Sd) + U + Su + i XFT]y = у. (77)

В уравнении (77) можно заметить неявную зависимость r¡ от релаксационного параметра ы, где очевидно, n = v при <о = 1. В случае матрицы мы заинтересованы в нахождении о, дающего минимум значения доминантного собственного значения, где для случая циклической матрицы с индексом 2, его значение определено с помощью уравнения (71). В рассматриваемых задачах на собственные значения связанные матрицы определены таким образом, что их доминантные собственные значения равны единице, так, что в случае матрицы УЫ(Л,), rj, = 1 и для соответствующего собственного вектора yt соотношение

К-1 [ (L + S11) + U + Sü + i X^y^ yt (78)

выполнено для всех и) * О. Таким образом, в этом случае главная задача связана с нахождением значения и, которое минимизует второе собственное значение rj2, эквивалентное поддоминантному отношению ВIVU( Aj)). Как замечено в численных эксперименте«, существует «best, минимизующее fftVyiA.,)) и его значение обычно больше, чем а>ор1, минимизующее спектральный радиус итерационной матрицы fu в методе верхней релаксации SOR. Для 1 * и з ubest 5 IV^l At) ] уменьшается монотонно, для tOb„t < и s ucrlt наблюдают сильное возрастание 6(Vu(А,)) и при (Jcrlt появляется расхождение итерационного процесса, где обычно <JcrU < г. Этот эффект обсуждается и демонстрируется в реакторных расчётах, приведённых в главе 5.

Интересно упомянуть, что в реакторной задаче 2, рассмотренной в главе 5, только эта стратегия является наиболее эффективным методом решения по сравнению со всеми остальными методами.

Определяя член источника селения при погода

Ш) = хщ»Г9<1). (79) можно переписать уравнение (59) как

(и - = т * я")^!} + га.». шо) Для стационарных значения 1(1) цикл внеяннх итераций р=1.г....,Р может быть выполнен согласно следующей схепе

(И - = 1Я * + 1(1) (81)

р р

или эквивалентно

<ра*р = с9(1*^-) * [г - нГЧЛ-'н"1*«.), не)

где

в = [I - И ^Г'к"'(Ы + 5") = [I - И^'Б"]"1»' + Н'^"]. (83)

Легко заметить, что

для

р=1, р(1*1) = С<р(1) * [I - и"1^] "и^сш.

для р=2, = С<р(1*~) * и - и ^"¡''к^гси

= Сг<р(1) * II * СНГ - М"18<,Г1К"11П>

для р=р. ри-и) = бТрш + [1 + с + ____ + с?*"1 ш - иГ'з*1]1)

и эти выражения могут быть записаны в следую ¡¡.ем обпей виде

(р (1*1) = ЧШ1)}р(1), (84)

где

пит = ^ ♦ (85)

и

р

к"* = Тс5"1» - ■1'Л-Ье-1. (86)

(р-1) и

Р=1

Предполагая, что для достаточно больших значений 1. У(А(1)) - и собст-

венные значения патрицы упорядэтены следующий образом

= Л = 1 > £ 2----

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

erEViAj.il = |?т.21. (87)

где индекс р относится к внешнии итерация«.

Как можно заметить для сходящейся матрицы G,

р —1 -1 Gr -> О и И('_и -> Е (88)

при Р —» со, отсюда, V(X1) —> B(Aj). Таким образом Патрицу м"^ ^ можно рассматривать как предобратную матрицу СР-1)-степени матрицы Е в глобальных-внешних итерациях. Поскольку, при Р = 1, эта стратегия сводится к стратегии глобальных итераций, отсюда, матрицу М ^ можно рассматривать как предобратную матрицу нулевой-степени матрицы Е в стратегии глобальных итераций.

Эта стратегия используется в методе EQUIPOISE, для которого не было дано доказательство. Однако, как можно заметить из уравнения (85), требование что матрица G должна быть сходящейся матрицей, т.е. p(G) < 1, является необходимым и достаточным условием сходимости глобальных-внешних итераций. Предполагая, что неотрицательная матрица G, определённая в уравнении (83), имеет собственные значения г, удовлетворяющие следующему уравнению

Gx = тх (89)

или

(Н + xSd + s")x = тМх (90)

и, используя расщеплённые матрицы уравнения (67), получают

K_1[t(L + Sd) + U + Sulx = тх. (91)

Так как диагональные доминантные свойства матрицы Е = А + Sd + Su заключают в себе то, что неотрицательная матрица

®Е= К"1[L + Sd + U + Su] (92)

представляющая итерационную матрицу в методе Якоби неприводимая и её спектральный радиус меньше единицы, то спектральный радиус матрицы G, представляющий в некотором смысле итерационную матрицу в методе Гаусса-Зейделя, также меньше единицы и, кроме того, p(G) < р (ВЕ).

В случае расщеплённых матриц уравнения (69), определяющие метод верхней релаксации S0R, принимают вид

С<оУ = [I " «¿Vr^CNu + Su)y = (93)

что даёт

(Nw + £Sd + s")y = е«(оУ

или эквивалентно

K"l[£(L + sd) + и + su]y = (94)

Вышеуказанное уравнение показывает неявную зависимость ? от параметра релаксации о. Как хорошо известно, для всех 1 < и < 2, |{|1 < 1 , т.е. Gu является сходящейся матрицей; кроме того, существует и минимиэующее доминантное собственное значение нв случае, когда матрица Е является циклической матрицей с индексом 2

ú = --(95)

1 + V ,-Tt

где Tt = p(G) = рг(ЖЕ). В общем случае обе матрицы Sd and Su неотрицательны и матрица 2Е не имеет никаких свойств цикличности с индексом 2, а также нет никакой точной формулы для £>. Однако оценка й при помощи вышеуказанной формулы является достаточным приближением для приложений. В реакторных задачах без рассеяния нейтронов вверх, для которых Su = О, ВЕ является блочной треугольной матрицей и её собственные значения связаны только с элементами диагональных подматриц 5Bg = K^tLg + Ug), которые являются циклическими матрицами с индексом 2.

Однако, минимизация g(G) не включает минимизации 6(V(Aj)] и обычно видно на практике, что значение й минимизуюшее меньше значения ubest миними-

зующего б(VÍA!)].

Матрично-блочная структура уравнения (58) позволяет, при фиксированном индексе 1, ввести следующий цикл внутренних итераций, t = 1,2,...,Тд, для каждой группы g

Mg¥>gU+f) = N (i+iil) + CgU), (96)

g 'д

где

g-i g с

cg(i) = + + l SÍ.gíPkCi) + (97}

k=l k=g+l k=l

включает члены суммарного рассеяния и деления для данной группы g. Уравнение (97) записывается в эквивалентном виде

(PgU+f) = + М~ cg(1) = v <р u+iiL) + м^сgU) (98)

1 /Ч * п

■f) =

'g 'g

и мы имеем

для t=l, í>gU+f) = Vg(i) + ШЧШ g

для t=2, pgU+f) = VgipgU+±) + M^CgCl) . g 9

= l/'í>g(J) + (I + VJM^CgU)

ДЛЯ t=Tg, <p,U + l) = v£%U) + (I + V, + v| + .... f vV'iM^CgCJ), или эквивалентно

(PgÍJ + l) = Vg(J) + M^CgÍJ), (99)

где

Kg н vjg = (M¿\)Te (100)

Tg

t-i -i ч

является предобратной матрицей (Т9-1^-степени матрицы Ад. Определяя

TV,

V

Vr.

Ü1

и c(l) =

CjO) с г(1)

(102)

групповые уравнения (99) можно выразить в сжатом виде при помощи одного уравнения

tp(l + l) = Vip(l) + Me (1)

<p(i+i) = vip(i) + m'^s^cj+j; + s"<p(i) + ^JjjXfVi.)] ,

и отсюда

где

ip(i +1) = Vl\<l))ip(l),

Ч(.Х(1)) = (I - M"1Sd]"1(P + H_1(SU + ^i-jXF1')].

(103)

(104)

(105)

(106)

Предполагая, что для достаточно больших значений 1, V(A(1)) к V(Aj), и собственные значения матрицы V{Aj) упорядочены в следующем виде

¿4,1 = £>[ VC At ) J = 1 > l?vi2l 2 |í>Vi3l £----

скорость сходимости в итерационном процессе (105), представляющем собой метод называемый автором стратегией глобальных-внутренних итераций (1, t) для однократного расщепления, зависит от поддоминантного отношения

eivUj)] = l?Vi2i, (107)

где индекс t относится к внутренним итерациям, и Tq может иметь разные значения для каждого g.

Как можно заметить, для каждой сходящейся матрицы Vg

V —> 0 и Й"1 —> А"1 (108)

при Г

<я для каждого g=1,2,... ,G. Отсюда

V(A,) -> Т(Aj), (109)

где T(Aj) определено при помощи уравнений (42). Таким образом, если значения Tv возрастают для всех g=l,2, . . . ,а, то спектр собственных значений матрицы VU (1)) стремится к спектру собственных значений матрицы Т( А ПЛ.

При 7g = 1 для каждого g = 1,2,...,а, эта стратегия сводится к стратегии глобальных итераций, м"1 = М"1 и может рассматриваться как предобратная матрица нулевой степени матрицы Л в стратегии глобальных итераций.

Эта стратегия глобальных-внутренних итераций эквивалентна стратегии, известной под названием стратегии внешних-внутренних итераций (где индекс 1 относится к внешним итерациям, а индекс í к внутренним итерациям), использующейся в большинстве' компьютерных программ.

Обычно метод SOR, определённый при помощи расщеплённых матриц (69), используют для ускорения сходимости внутренних итераций и поскольку Ад являются циклическими матрицами с индексом 2, оптимальный параметр релаксации uopt, минимизуюишй спектральный радиус итерационной матрицы V в £ы может быть определён при помощи формулы (71). Однако, как можно заметить на практике использование u = <oopt в итерационном процессе не приводит к минимизации eiVtAj)), которая происходит при использовании и = ubesl > toopt во всех группах.

Матрица \(\(1)) уравнений (106) связана с циклом внутренних итераций выполненных во всех группах, g = 1,2,... ,G, и может быть представлена следующим образом

vua» = G + [I - M"Vr!M~J^jjyXFT, (110)

где

G = [I - M-1Sd]_1 [í? + M'Vj. (Ill)

При вышеуказанном виде матрицы \(.\(1)) схема глобальных итераций может быть записана как

<р(1+1) = G<р(1) + [I - H"ISd]"1M"1^I-jXFTipCJ).

Теперь вводя в эту схему внешние итерации р = 1,..,Р, и используя связь (79), получаем

<p(l+Z) = g<р(1+^) + [i - m"1sd]"lh"1fa) =

= cpf(v + íi + g + g2+... + gp_1][i - m'vj'^fo; и при p = P

<p(l+l) = ЧШ1))<р(1), (112)

где

У(АШ) = GP + jj2jXF1 (113)

И

p

r1 = Y G^U - rVrlM"1. (114)

(p-i) l

p=1

Предполагая, что для достаточно больших значений I, VCXi 1 >) » VO^), и собственные значения матрицы упорядочены следующим образом

>Лм = ptVUj)] = 1 > |?Vi2l £ |?V3| £ ....

скорость сходимости в итерационном процессе (112), представляющем общий вид метода, называемого автором стратегией глобальных-внешних-внутренних итераций (1,р, t) для однократного расщепления зависит от поддоминантного отношения

elVUj)] = |S/V>2I (115)

и матрица является предобратной матрицей (Р-1)-степени матрицы Е, с

процессом внутренних итераций, включённым в матрицу М"1.

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

Очевидно, что описанные стратегии, являются частными случаями стратегии глобальных-внешних-внутренних итераций, полученной при предположении Р = 1 и/или Тд = 1 для всех g = 1,2, ...,G. Конечно, если Р —» и, то

Gp -> О, М"1 , -> Е"1 И VU,) BU,). (116)

(р-1) * 1

Как можно заметить структура матрицы V(At), выраженной явно при помощи уравнения (113), строго связана как с выбором матриц расщепления И и М, так и значений Р и Тд, определяющих степень так называемых "уместных" предобратных матриц; каждое из них влияет на поддоминантное отношение 5[v(x,)). Однако их явное взаимодействие создаёт серёзные трудности в теоретических исследованиях поведения поддоминантного отношения в зависимости от предполагаемого расщепления матрицы А = Н - N и значений параметров Р и Тд. Отсюда становится ясным почему на практике используют эмпирический подход для оценки итерационных параметров в разных стратегиях. Обычно предполагают, что параметры итерационных стратегий зависят от типа реактора.

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

В системе линейных уравнений

Аф = с (117)

матрица А может быть выражена в форме

А = Р - R + S (118)

называемой двухкратным расщеплением матрицы А и если Р является неособенной матрицей, то это расщепление ведёт к следующей трёхслойной схеме

t > О

Рф(м) = Rфа' - Зфа-и + с.

(119)

или эквивалентно

= P"VU - Р-^11"11 + Р_1с, t > 0.

(120)

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

г,(t+1) ф P-1R, -P_1S" У» ■ p-V

I, 0 Ф + 0

Обозначая

Ф

ФШ

Гл(и 1 9 P-V

,(t-l) У ■ v = 0

(121)

(122)

получают

W =

rP_1R, -P^S"1

=

+ v,

t > 0.

(123)

(124)

Таким образом, необходимое и достаточное условие, которое обеспечивает сходимость итерационного метода (124) к единственному вектору решения ф = а"*с для всех векторов ф'0> и ф'1', это условие того, что спектральный радиус итерационной матрицы Ж, дШ), должен быть меньше единицы.

Определяем

IM =

тогда

Р, о о, Р

Д = IM - IN

IN =

r, -s P, 0

P-R, S -P, P

(123)

(126)

представляет однократное расщепление матрицы Д, выведенной из двухкратного расщепления матрицы А, данного с помощью уравнения (118), где

W = М N.

(127)

и

По аналогии с результатами главы 3 очевидно, что

М"1 = (Л + Ж + ж2* . .. + ТсМм"1 = У V1 м"1 (128)

(1) и»

О

или эквивалентно

= (0 - Ж1*1)Я"1 (129)

является предобратной матрицей ¡-степени матрицы Д, и при р(Ж) < 1 и ( -) »

й"1. (130)

Очевидно, при 5=0, двухкратное расщепление матрицы А сводится к однократному расщеплению матрицы А.

Конструкция итерационных стратегий, основанных на двухкратном расщеплении матрицы А, подробно обсуждается в подглавах 4.1, 4.2, 4.3 и 4.4.

Очевидно, что все итерационные решения основанные на использовании двухкратного расщепления матрицы А не связаны с задачей (26) а со следующей эквивалентной задачей двухкратного порядка на собственные значения

ЕФ = (131)

где

Е = Д - Б* - 5й. (132)

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

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

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

Н( V) = - 1п5[V].

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

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

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

аи = аи + ьи х ciJ' (133)

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

Все вычисления были выполнены при помощи следующих программ HEXAGA, HEXSOR и HEXSLOR решающих двухмерное многогрупповое уравнение диффузии нейтронов (24) путём разностного приближения основанного на 7-точечной сеточно-граничной (mesh-edged) разностной формуле в однородной треугольной сетке, наложенной на 120-градусную параллелограммную область. В случае натурального упорядочения точек сетки, подматрицы Ад в матрице коэффициентов А, имеют только семь ненулевых диагоналей, принимающих трехдиагональную блочную структуру, подходящую для осуществления 1-линейных алгоритмов. Для решения неоднородных задач с семью диагональными матрицами в пределах внутренних итераций, предфактори-зационный метод AGA, точечный метод Гаусса-Зейделя и 1-линейный метод Гаусса--Зейделя используются как алгоритмы матричного расщепления соответственно в программах HEXAGA, HEXS0R и HEXSLOR. Сходимость глобальных итераций во всех программах ускоряется только при помощи процессов последовательных верхних релаксаций во внутренних итерациях. В случае программы HEXAGA процесс двухкратной последовательной верхней релаксации, представляющий итерационные стратегии двухкратных расщеплений обсуждаемые во главе 4, осуществлён как

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

Вычисления для каждой задачи продолжались до тех пор, пока максимальное, по всех точках сетки и всех группах, относительное изменение потоков нейтронов между глобальными итерациями, Ср, а также относительные изменения значений kerf, e,t, были меньше чем предписанные числа, где keff определён следующим образом

Krrll + l) s ЛП+JJ = \>l)[l+l)dv

и интегрирование источников деления ф{1+1) = ) приближается при помощи

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

е,„ = гпах

pn(J + 2,p=P, í=J) - <pnU,p=P,t=T)

ке„и+1) - kerf(i)

для всех n = i.2,...HxG (134)

Еь = шах

rU + l)

(135)

где 1, р, 4 итерационные индексы соответственно в глобальных, внешних и внутренних итерациях. Предполагается, что значения Р и Т (где Тд = Т для всех 1 2> ц з с) постоянны во всём итерационном процессе для данной задачи. Индекс 1 играет роль итерационного индекса в степенном методе, использующимся для вычисления доминантного собственного значения данной матрицы, структура которой зависит от предполагаемых значений Р, Г и применяемого (однократного или двухкратного) расщепления матрицы А.

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

ер * Ю"5 и ск = 10~6, (13В)

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

Т = 1, 2, 3, 5, а, 12 и 20.

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

у

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

Были выполнены компьютерные вычисления для ряда моделей связанных с шестью реакторными задачами (два быстрые и четыре термические реакторные примеры) взятыми, главным образом, из литературы. Первые четыре задачи были вычислены с арифметикой однократной точности (single precision arithmetic) на персональном компьютере, а две последние задачи на компьютере Convex С32Ю с арифметикой однократной точности для программы HEXAGA и с арифметикой двухкратной точности (double precision arithmetic) для программы HEXSLOR.

Во всех анализированных стратегиях предполагается одно и то же самое число внутренних итераций и одинаковое значение параметра релаксации ы, в каждой группе энергии для целого итерационного процесса. Таким образом, каждый численный эксперимент рассматривается как отдельная задача определённая при помощи предположенного числа внутренних итераций, приходящихся на одну внешнюю итерацию и предположенного значения параметра релаксации (J. В компьютерных вычислениях для каждой задачи используется нулевой начальный вектор р10> и в случае двухкратного расщепления матрицы А, выполненного по алгоритму HEXAGA, предполагается р(1> = 0. Однако, для вычисления стартового вектора источника деления 0111 = FTpl0) все компоненты вектора р10) были положены равными единице. Таким образом, стартовое распределение источников деления было предположено кусочно плоским в расщепляемых подобластях.

Наблюдаемые числа глобальных итераций и соответствующая компьютерная работа (выраженная при помощи числа флопов, приходящихся на одну точку сетки) в зависимости от числа внутренних итераций, приходящихся на одну внешнюю итерацию, предположенная как параметр, представлены графически в зависимости от коэффициента релаксации w дпя каждой итерационной стратегии, используемой для решения данной реакторной задачи. Значения оптимальных параметров релаксации wopt, минимализирующих спектральный радиус блочно-диагональных итерационных матриц V (или If в случае программы HEXAGA) и удовлетворяющих следующему неравенству

, . hexaga , hexslor „ hexlor

1 < GJopt < (Jopt < (Oopt < 2,

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

последовательно для программ

HEXAGA, HEXSOR и HEXSLOR.

Значение uopt определено при помощи формулы, данной с помощью уравнения (71) для программ HEXSOR и HEXSLOR, Однако, надо отметить, что итерационные матрицы V в HEXSOR не являются циклическими последовательно упорядоченными матрицами индекса 2, вычисление uopt при помоши уравнения (95) даёт достаточно хорошую оценку значения tJ0Pt. необходимую для практических применений. В случае использования HEXSLOR итерационные патрицы V суть циклические последовательно упорядоченные матрицы индекса 2 и оценка uopl может быть сделана при помощи алгоритма sigma-SOR, описанного в работе (15). Значение (Jopt, минимализирующее спектральный радиус итерационной матрицы W, появляющейся в процессе двухкратной верхней релаксации, использованное в программе HEXAGA, было определено при помощи процедуры описанной в работе [13).

Из численных результатов полученных для рассмотренных задач, можно заключить, что существует некоторый диапазон числа внутренних итераций, приходящихся на одну внешнюю итерацию Т £ 1 * т, при Т г j, для которого решение может быть получено при сравнимом итерационном усилии (представленном при помощи общего числа внутренних итераций) и приближённо пропорциональном величине арифметической работы (выраженной при помощи числа флопов на одну точку сетки) в задачах при простых моделях рассеяния. Как можно заметить, для результатов полученных для задачи Зад.6, этот диапазон равен 3 s г ^ 20 для HEXAGA и 2 — Т — 30 для HEXSLOR, но минимальное число флопов получено при Г = 1 для обеих программ, где HEXAGA даёт решение при числе флопов на одну точку сетки в два раза меньше.

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

Из результатов, полученных при помощи HEXSLOR для задач, рассмотренных в этой работе, скорее всего трудно найти формулу, приближающую достаточно хорошо значение t^est- в задаче Зад.6, а также в других задачах, «be!St изменяет своё значение в зависимости от Г, следовательно, может оказатся трудным отличить (Jbest от 4>pt как, например это видно в трёх вариантах задачи Зад.1 или в Зад.5, при 2 внутренних итерациях на одну внешнюю итерацию.

Однако, в случае алгоритма HEXAGA можно заметить несколько регулярных отклонений Wt,est от uopt в рассмотренных задачах и, кроме того, u)opl скорее всего нечувствительно к использованному числу внутренних итераций на одну внешнюю итерацию. Следующая формула

1.2 - <J0Dt

<"W = 4>Pt + pi»!)-437)

даёт очень хорошее приближение к верному значению (Jb<.st для большого класса реакторных задач, где p(Wi) есть спектральный радиус матрицы Ж при ы = 1, обозначенный как toopt оптимальный параметр релаксации, минимализирущий спектральный радиус матрицы Jfu в процессе двойной верхней релаксации, использованной в HEXAGA. Оценка uopt в HEXAGA может быть получена при помощи эффективной техники описанной в [131.

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

Анализ большого числа численных результатов полученных для рассмотренных задач, представляющих разные типы быстрых и термических реакторов, позволяет сделать вывод, что решения при наименьших затратах на вычисления, могут быть получены при помощи итерационной стратегии, где используется только 1 внутренняя итерация на одну глобальную итерацию в программе HEXAGA, алгоритм • которой основан на одной из простейших моделей предфакторизационных методов AGA (описанных в главе 1) при использовании процесса двойной верхней релаксации как эффективной техники ускорения сходимости. Это означает, что выбор матриц расщепления имеет доминирующее влияние на эффективность метода, использованного для решения задач теории диффузии нейтронов и наблюдаемое поведение сходимости решений при помощи алгоритма HEXAGA указывает на очевидный значительный потенциал двухкратных расщеплений, происходящих от предфакторизационных методов AGA.

Для всех рассматриванных задач, HEXAGA даёт решения при от 2 до 5 раз меньшем числе флопов по сравнению с результатами HEXSLOR.

В случае программ VALE и DXY, основанных на многочленном ускорении Чебышева, результаты доступны только для задач Зад.1, Зад.З и Зад.6, но их сравнение с результатами, полученными при помощи HEXAGA, показывает, что программа HEXAGA требует в два раза меньше вычисления. Результаты DXY для

задачи Зад.6 доступны только при критерии остановки итерационного процесса с<р ± 5х10~3 и было бы очень интересно сравнить результаты, полученные для этой задачи при с^ = icf5.

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

В заключение стоит упомянуть, что KEXAGA даёт решения при использовании арифметики с однократной точностью для больших реакторных задач, решение которых при помощи HEXSLOR может быть получено только при использовании арифметики с двухкратной точностью, как было показано в задаче Зад.6. Задача для априорной оценки ubCEt в HEXAGA достаточно хорошо решена для приложений. Кроме того, проверка симметрии решения в задачах Зад.1 и Зад.2, для которых HEXAGA даёт решения при максимальной относительной ошибке 60-градусной симметрии esym, меньше на десять порядков по сравнению с решениями программы HEXLOR, что показувает большую погрешность решении программы HEXAGA. Но это позволяет предположить, что в случае решений несимметрических задач при помощи HEXAGA, с одной стороны характеризующихся повышенной скоростью сходимости, с другой стороны они могут быть на много более достоверными, чем решения, полученные при помощи расщепления SLOR использованного в HEXSLOR. Это свойство решений, полученных при помощи HEXAGA является особенно важным в сопряжённых вычислениях потока нейтронов, использованных для получения ответа на малые возмущения.

8. Основные результаты, выносимые на защиту

На защиту выносится предложенный автором метод решения многомерных задач теории диффузии нейтронов на основе метода трехуровневых (глобальных-внешних-внутренних) итераций, который является обобщением хорошо известного метода внешних-внутренних итераций. Введённая автором матричная формулировка, позволяет точно и наглядно представить взаимозависимость между внутренними и внешними итерациями на уровне глобальных итераций. Эта матричная формулировка позволила также доказать сходимость метода EQUIPOISE, основанного на стратегии глобальных-внешних итераций. Расщепление матриц при помощи двойной верхней релаксации, использованной в предфакторизационном методе AGA, позволяет существенно увеличить скорость сходимости решаемых задач по сравнению с другими существующими методами. Вычислительные эксперименты доказали высокую достоверность результатов программы HEXAGA.

На основании проведенных исследований на защиту выносятся следующие результаты:

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

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

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

4. Доказана теорема сходимости этого процесса.

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

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

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

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

1. Z.I.Woinicki, EWA-II Two-Dimensional, Two-Group Diffusion Fast Code, Kernenergie 14, 325-329, 1971.

2. Z.I.Woánicki, Two-Sweep Iterative Methods for Solving Large Linear Systems and their Application to the Numerical Solution of Multi-Group, Multi-Dimensional Neutron Diffusion Equation, Doctoral Dissertation, Institute of Nuclear Research, áwierk-Otvock, Report No.1447/CYFRONET /РМ/А, 1973.

3. Т.Апостолов и З.Возницки, Диффузионная двумерная программа HEXAGA-II для многогрупповых расчетов гексагональных решеток, Атомная Энергия, т. 38, вып. 6, 372-374, 1975.

4. Zbigniew Woánicki, Dwuprzebiegowe nietody iteracyjne AGA rozwi^zywania duzych ukladóu równaA liniowych, Rocznik Polskiego Towarzystwa Matema-tycznego, Seria III: Matematyka Stosowana VI, 5-16, 1976.

5. Z.I.Woánicki, The АСА Two-Sweep Iterative Methods and their Application to the HEXAGA-II Programme for Solving the Two- Dimensional, Multi-Group Neutron Diffusion Equation for an Uniform Triangular Mesh, Proc. Reaktortagung, Düsseldorf, Germany, March-April, 1976, 83-84.

6. Т.Г.Апостолов, П.Т.Петков и З.Возницки, Двухмерная многогрупповая диффузионная программа ЕВА-II-30, Ядерна Энергия 3, 8-15, София, 1976.

7. Z.I.Woznicki, AGA Two-Sweep Iterative Methods and their Applications in Critical Reactor Calculations, Nukleonika 23, 941-968, 1978.

8. Z.1.Woánicki, HEXAGA-I1-120,-60,-30 two-dimensional multigroup neutron diffusion programmes for a uniform triangular mesh with arbitrary group scattering. Rep. KfK 2789, Kernforschungscentrum Karlsruhe, Germany, 1979.

9. Z.I.WoinLckl, Efficient methods of accelerating reactor diffusion codes, Proc. International Topical Meeting on Advances in Mathematical Methods for the Solution of Nuclear Enginnering Problems, Hilton International

München, April 27-29, 1981, Vol.1, 385-399.

10. Z.I.Woinicki, HEXAGA-III-120,-30 three-dimensional muJtigroup neutron diffusion programmes for a uniform triangular mesh with arbitrary group scattering. Rep. KfK 3572, Kernforschungscentrum Karlsruhe, Germany, 1983.

11. Z.I.Woinicki, Two- and three-dimensional benchmark calculations for triangular geometry by means of HEXAGA programmes, Proc. Internat. Mee'ting on Advances In Nuclear Engineering Computational Methods, Knoxville, Tennessee, April 27-28, 1985, American Nuclear Society, 147-156.

12. Z.I.Woinicki, AGA Two-Sweep Iterative Methods and their Applications for the Solution of Linear Equation Systems, Proc. International Conference on Linear Algebra and Applications, Valencia, Spain, September 28-30, 1987 (published in Linear Algebra and its Applications 121, 702-710, 1989).

13. Z.I.Woinicki, Estimation of the opt imam relaxation factors in partial factorization iterative methods, SIAM J. Matrix Anal. Appl., 14, 59-73, 1993.

14. Z. I .Woinicki, On numerical analysis of conjugate gradient method, Japan J. Indust. Appl. Math., 12, 487-519, 1993.

15. Z.I.Woinicki, The sigma-SOR algorithm and the optimal strategy for the utilization of the SOR iterative method, Mathematics of Computations, 206, 619-644, 1994.

16. Z.I.Woinicki, Nonnegative splitt ing theory, Japan J. Indust. Appl. Math., 11, 289-342, 1994.

17. Z.I.Woinicki, The Numerical Analysis of Eigenvalue Problem Solutions in the MuJtigroup Neutron Diffusion Theory, Institute of Atomic Energy, otwock-ävierk, Report IEA-6/A (125p.), June 1995.

18. Z.I.Woinicki, The Numerical Analysis of Eigenvalue Problem Solutions in the Hultigroup Neutron Diffusion Theory (the abbreviated version of the work 1371), Proc. International Conference on Reactor Physics and Reactor Computations, Tel Aviv, Israel, January 23-26, 1994, 641-653.

19. Z.I.Woinicki, On reliability of solution symmetry in neutron diffusion theory problems, Proc. International Conference on the Physics of Reactors PHYSOR'96, Mito, Ibaraki, Japan, September 16-20, 1996, A21-A30.

20. Z.I.Woinicki, Comparison Theorems for Regular Splittings on Block Partitions, Linear Algebra and its Applications 253, 199-207, 1997.

21. Z.I.Woinicki, Comparison Theorems for Splittings of Monotone Matrices, The invited talk in WCNA'96 - The Second Congress of Nonlinear Analysts, Athens, Greece, July 10-17, 1996. Published in Nonlinear Analysis, Theory, Methods & Applications, 30(2), 1251-1262, 1997 Elsevier Science LTD.

22. Z.I. Woinicki, Conditions for Convergence and Comparison, Proc. 15th IMACS World Congress on Scientific Computation, Modelling and Applied Mathematics, Berlin, August 1997, Vol.2, Numerical Mathematics, 291-296, Edited by A.Sydow, Wissenschaft & Technik Verlag.

23. Z.I. Woinicki, Matrix Splitting Principles, The invited talk in XII Conference on Applied Mathematics, Paliö, Yugoslavia, September 8-12, 1997. To be published in Proceedings in 1998. The final version of the paper submitted to Lin. Alg. Appl.

24. Z.I.Woinicki, H.J^drzejec, How Can Ve Accelerate the Convergence in the SLOR Method, The invited talk in 7th Russian Conference on Grid Generation Problems and Numerical Methods in Mathematical Physics, Durso (Krasnodar region), Russia, September 3-15, 1998. To be published.

Рукопись иосгутла в издательский отдел 16 октября 1998 гола.

ВЕДЕНИЕ лава 1. ИТЕРАЦИОННЫЕ МЕТОДЫ РЕШЕНИЯ СИСТЕМ ЛИНЕЙНЫХ УРАВНЕНИЙ

1.1. ПРЕДФАКТ0РИЗАЦИ0ННЫЕ МЕТОДЫ AGA

1.2. ТЕОРИЯ РЕГУЛЯРНОГО РАСЩЕПЛЕНИЯ

1.3. ТЕОРИЯ НЕОТРИЦАТЕЛЬНОГО РАСЩЕПЛЕНИЯ

1.4. ПРОЦЕССЫ ПОСЛЕДОВАТЕЛЬНОЙ ВЕРХНЕЙ РЕЛАКСАЦИИ

1.5. ЧИСЛЕННЫЙ АНАЛИЗ МЕТОДА СОПРЯЖЁННЫХ ГРАДИЕНТОВ

1.6. АЛГОРИТМ Sigma-SOR лава 2. ФОРМУЛИРОВКА ЗАДАЧИ

2.1. УПРОЩЕННАЯ ЗАДАЧА

2.2. ОСНОВНАЯ ЗАДАЧА лава 3. ИТЕРАЦИОННЫЕ МЕТОДЫ ОДНОКРАТНОГО РАСЩЕПЛЕНИЯ

3.1. ГЛОБАЛЬНЫЕ ИТЕРАЦИИ ОДНОКРАТНОГО РАСЩЕПЛЕНИЯ

3.2. ГЛОБАЛЬНЫЕ-ВНЕШНИЕ ИТЕРАЦИИ ОДНОКРАТНОГО РАСЩЕПЛЕНИЯ

3.3. ГЛОБАЛЬНЫЕ-ВНУТРЕННИЕ ИТЕРАЦИИ ОДНОКРАТНОГО РАСЩЕПЛЕНИЯ

3.4. ГЛОБАЛЬНЫЕ-ВНЕШНИЕ-ВНУТРЕННИЕ ИТЕРАЦИИ ОДНОКРАТНОГО РАСЩЕПЛЕНИЯ лава 4. ИТЕРАЦИОННЫЕ МЕТОДЫ ДВУХКРАТНОГО РАСЩЕПЛЕНИЯ

4.1. ГЛОБАЛЬНЫЕ ИТЕРАЦИИ ДВУХКРАТНОГО РАСЩЕПЛЕНИЯ

4.2. ГЛОБАЛЬНЫЕ-ВНЕШНИЕ ИТЕРАЦИИ ДВУХКРАТНОГО РАСЩЕПЛЕНИЯ

4.3. ГЛОБАЛЬНЫЕ-ВНУТРЕННИЕ ИТЕРАЦИИ ДВУХКРАТНОГО РАСЩЕПЛЕНИЯ

4.4. ГЛОБАЛЬНЫЕ-ВНЕШНИЕ-ВНУТРЕННИЕ ИТЕРАЦИИ ДВУХКРАТНОГО РАСЩЕПЛЕНИЯ лава 5. ЧИСЛЕННЫЕ ПРИМЕРЫ

5.1. ВЫЧИСЛИТЕЛЬНЫЕ АСПЕКТЫ

5.2. РЕАКТОРНЫЕ ЗАДАЧИ

5.3. ОБСУЖДЕНИЕ ЧИСЛЕННЫХ РЕЗУЛЬТАТОВ

5.4. ДОСТОВЕРНОСТЬ СИММЕТРИЧЕСКИХ РЕШЕНИЙ

5.5. ЗАКЛЮЧИТЕЛЬНЫЕ ЗАМЕЧАНИЯ лава 6. ВЫВОДЫ

ЛАГОДАРНОСТИ

ИТЕРАТУРА

ПИСОК КОНФЕРЕНЦИЙ

1РИЛ0ЖЕНИЕ - ТАБЛИЦЫ И РИСУНКИ К ГЛАВЕ

Бурное развитие атомной энергетики стимулировало активный поиск эффективных методов решения задач переноса нейтронов как основы моделирования физических процессов, протекающих в ядерных реакторах. Теория переноса излучения - одна из основных проблем современной науки - быстро развивается на основе достижений теоретической физики и, как каждая новая важная область прикладной физики, стимулирует развитие вычислительной математики. Весьма общая математическая формулировка задач теории переноса задаётся с помощью линеаризованного уравнения Больцмана [1,23,24,57], относящегося к классу интегро-дифференциальных уравнений математической физики, являющихся математической моделью для описания переноса нейтронов. В основе такого подхода лежит фундаментальная теория микроскопических процессов, однако подобно теории упругости и классической гидродинамике, эта теория может быть достаточно корректно сформулирована и на макроскопическом уровне.

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

Прогресс атомной науки и техники стимулировал развитие различных приближённых методов решения уравнения Больцмана. В настоящее время наиболее продвинутым методом в теоретическом и алгоритмическом аспектах является диффузионное приближение, которое выводится с помощью метода сферических гармоник. Границы применимости диффузионного приближения указаны, например в работах [5,23,24,57]. Диффузионное приближение является наиболее широко используемым методом анализа критичности ядерных реакторов. Рассмотрение критичности, вообще говоря, сводится к задаче о собственных значениях для многогрупповых уравнений диффузии нейтронов, решение которой даёт собственное значение реактивности - эффективный коэффициент размножения, и вектор собственных значений - поток нейтронов порождающих распределения мощности в реакторах большой мощности; эти величины связаны и с неравномерной загрузкой топлива и с его выгоранием.

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

Решение уравнений одномерной диффузии является сравнительно простой задачей, позволяющей использовать прямой численный метод прогонок [2]. Этот метод явился важным вкладом в теорию дифференциальных уравнений и привлёк к себе большое внимание, особенно после дискуссий в 1953-1954 гг. на известных семинарах Гельфанда при Московском университете [25,26,27].

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

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

Сеточный или конечно-разностный метод аппроксимации дифференциального уравнения системой алгебраических уравнениий описан, например, в работах [3, 54,58-65].

В вычислительной математике предложен целый ряд итерационных алгоритмов для решения многомерных уравнений диффузии, с ними можно ознакомиться в монографиях [57,61,54,69]. Обзору сеточных методов решения уравнения диффузии посвящена статья [71]. В последнее время развиваются итерационные методы, в которых сочетается подход спектральных и вариационных оптимизаций. В.И.Лебедев [29] сформулировал условия на операторы задач, для которых итерационный процесс имеет неулучшаемую оценку числа арифметических операций. Развивается еще один метод выбора оптимальных параметров итерации, основанный на вероятностном подходе. Ряд интерестных результатов в этой области получен Ю.В.Воробьевым [72]. Итерационные подходы рассматрываются также в книгах по методам расчёта реакторов [24,57], где обзор и систематизация итерационных методов даны в книге Г.И.Марчука и В.И.Лебедева [24]. Итерационные методы решения несовместных систем были предложены в работе Ю.А.Кузнецова [68]. Итерационным методам решения эллиптических уравнений на примере диффузионных реакторных уравнений посвящена содержательная монография [2]. До сих пор не утратил своего большого значения, ставший уже классическим и вышеупомятуный, метод верхней релаксации Янга-Франкела [3]. Исследования этого метода обобщены в монографиях [2,3,24,28,57], а также в статье автора [8], где предложена эффективная процедура для оценки наилудшего параметра релаксации. Новый подход для ускорения сходимости линейной верхней релаксации предложен в работе [73].

В последние сорок лет в области численного решения многогрупповых многомерных уравнений диффузии нейтронов много усилий было направлено на развитие эффективных итерационных методов и внедрение различных компьютерных программ (см., например, ссылки в [4,7]). Стандартный метод решения основан на использовании внешних-внутренних итераций [2,71], останавливаемых, когда их число достигает некоторого значения, установленного для каждой группы, либо когда удовлетворяется критерий сходимости. Выполнение цикла внутренних итераций во всех группах соответствует одной внешней итерации, в которой член деления нейтронов пересчитывается. Для увеличения скорости сходимости внешних итераций обычно используются методы ускорения полиномами Чебышева [2,3,24,29].

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

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

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

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

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

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

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

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

Теоремы теории неотрицательного расщепления [9], как результат многолетних исследований автора вопросов монотонности, существенно важны для анализа сходимости итерационных методов решения систем линейных уравнений с монотонными матрицами.

Алгоритмы метода AGA с двойной верхней релаксацией были реализованы для двухмерного и трёхмерного случаев в компьютерных программах HEXAGA-II [11,32] и HEXAGA-III [12], решающих уравнения диффузии нейтронов на треугольной и гексагональной разностных сетках, на основе стандартной стратегии внешних--внутренных итераций. Программы HEXAGA, используемые интенсивно для расчётов связанных с проектированием и эксплуатацией ядерных реакторов в Польше и других странах (в Бельгии, Болгарии, Германии), делают возможным получать результаты вычислений в несколько раз быстрее чем другие компьютерные программы этого типа, а их преимущество возрастает с ростом числа точек разностной сетки. Сравнение результатов расчётов HEXAGA-II и английской программы SNAP для реакторной задачи ВВЭР-1000 данное в следующей таблице (4 сек CYBER 73 = 10 сек ЕС-1040).

Число сет. Число Время на точек на энерго- CPU одну энерго-

Программа одной гек- -простр. kef f время простр. Компьютер сагональной точек сек. точку кассете милисек.

1 288 1. 11715 6 0. 0208

12 1024 1. ,11375 14 0. 0136

HEXAGA-II 27 2208 1. ,11279 34 0. 0154 CYBER 73

48 3844 1. ,11240 60 0. 0156

108 8462 1. ,11211 237 0. 0280

1 92 1. ,11387 26 0. 283

SNAP 6 552 1. ,11068 150 0. 272 ЕС-1040

24 2208 1. ,11124 894 0. 405

96 8832 1. .11175 6435 0. 729

Эти результаты были представлены на конференции по физике реакторов в Knoxville, США, в 1985г [13], где результаты для программы SNAP даны в [56].

В Польше программа HEXAGA-11 была использована в расчётах реактора ядерной электростанции "Жарновец", которые проводились в 1980 годах, а также при физических расчётах исследовательских реакторов EWA и MARIA, необходимых при разработке отчёта безопасности и текущей эксплуатации реактора. Проверка результатов вычислений, полученных при помощи программы HEXAGA-II на основании данных по эксплуатации реактора EWA показала великолемную согласованность [74]. Надо отметить, что для реактора EWA и ядерной электростанции "Жарновец", имеющих гексагональную структуру сетки топливных элементов, имело место полное совпадение с геометрией треугольной сетки программы HEXAGA-11. В поперечном сечении реактора MARIA каналы твэлов размещены в квадратной сетке, а именно применение мелкой треугольной сетки программы HEXAGA-11 позволило восстановить круглый вид каналов твэлов и детали сложной геометрическо-мате-риальной конфигурации реактора MARIA, более точным образом, чем это было-бы возможно в программе с прямоугольной ориентацией сетки. Программа HEXAGA-11 в настоящее время используется в расчётной системе разработанной для нужд реактора MARIA.

Для того, чтобы представить подготовку входных данных программы HEXAGA-II для реакторных задач при мелкой сетке были разработаны две вспомогательные программы HEX1-22 и HEXI-23 [11]. Эти вспомогательные программы предназначены для производства первой части входных данных программы HEXAGA-11 (связанных с описанием материально-геометрической конфигурации точечной сетки) для данной реакторной задачи, в которой шаг однородной треугольной сетки уменьшается два раза в случае программы HEX 1-22 и три раза в случае программы HEXI-23. Эти программы используют такие же входные данные, как программа HEXAGA-II без ввода добавочной входной информации. Каждые выходные данные программ HEXI могут быть использованы как новые входные данные программ HEXI. Таким образом, если мы имеем входные данные программы HEXAGA-11 описывающие данную реакторную задачу с минимальным числом точек сетки и используем любую комбинацию входных/выходных данных программы HEXI-22 и/или HEXI-23 мы можем производить входные данные программы HEXAGA-11 для расположения точек сетки для этой задачи с шагом сетки уменьшенным согласно i к следующим множителям: 2 3 для любых целых чисел 1.к > О.

В случае трёхмерной программы HEXAGA-III были разработаны аналогичные вспомогательные программы HEXI-32 и HEXI-33 [12].

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

Софии [50], где была разработана болгарская версия программы HEXAGA, известная под названием НЕХАВ, и применяемая при расчётах связанных с эксплуатацией энергетических ядерных реакторов ВВЭР работающих в Болгарии.

Актуальные версии программ HEXAGA работают на основе, предложенной в диссертации, стратегии глобальных-внешних-внутренних итераций.

Диссертационная работа состоит из введения, шести глав и приложения.

6. выводы

Опыт автора по численным методам линейной алгебры связан в основном с эшением многомерных и многогрупповых уравнений диффузии нейтронов, играющих /щественную роль в математическом описании физических явлений, имеющих место ядерных реакторах. Решение многомерных уравнений диффузии нейтронов, пред-гавляющих класс эллиптических дифференциальных уравнений в частных производ-ых второго порядка, сводится к решению эквивалентных разностных уравнений ипа Аф = с; где матрицы А имеют разреженную структуру и монотонны т.е., довлетворяют условию, что обратная матрица существует и она неотрицательна А."1 > 0).

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

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

Матричный формализм введённый автором для описания методов EWA и AGA казался также полезным при классификации предфакторизационных методов других второв. В таблице 1.1 из работы [41], приведённой на стр. 8, отдельные етоды определяются с помощью выбора матриц Н и Q, где методы Якбби и Гаусса-Зейделя могут рассматриваться как некоторые дегенеративные виды предфактори-ационных методов. Классификация предфакторизационных методов, представленная втором на международной конференции в Париже в 1987г., получила первый приз конкурсе на лучшую работу.

Более ранние теоретические результаты автора, полученные для регулярных асщеплений, приведённые в подглаве 1.2, были предметом его кандидатской дис-ертации, защищенной в 1973г. [38]. Varga [30], как первый в литературе пред-ета, обратил внимание на значение результатов кандидатской диссертации втора и с тех пор, многие авторы ссылаются на эти результаты, используя их в воих работах. Итак, они были отмечены основоположником методов расщепления ссылки 8,9,10,11,15,16 в работе [9]).

В подглаве 1.3 описаны более важные результаты теории неотрицательного асщепления [9] полученные, после защиты кандидатской диссертации, во время ноголетних исследований вопросов монотонности, которые ведутся автором осо-енно интенсивно в последние годы [34,35,52,53]. Теоремы этой теории сущест-енно важны для анализа сходимости итерационных методов решения систем линей-ых уравнений с монотонными матрицами.

В подглаве 1.4 описано применение верхней релаксации в предфакториза-июнных алгоритмах AGA [10], соответствующей двухкратному расщеплению матри-ы решаемой задачи, которые оказались исключительно эффективными для решения искретных многомерных эллиптических задач появляющихся в физике ядерных ректоров. Алгоритмы метода AGA с двойной верхней релаксацией были внедрены в вумерных и трехмерных компьютерных программах HEXAGA-II [11] и HEXAGA-III 12], решающих уравнения диффузии нейтронов в треугольной и гексагональной еометрии разностной сетки.

Анализ критериев сходимости, использующихся для остановки итерационного роцесса метода сопряжённых градиентов с точки зрения корректности полученных езультатов, представленный в работе [14], описан в подглаве 1.5. На основе нализа полученных численных результатов можно сделать вывод о том, что метод опряжённых градиентов не во всех задачах имеет преимущество по отношению к ругим методам для классических решений, а выбор соответствующих критериев ходимости недооценивается многими авторами. В этой работе показано, что ритерии, использованные некоторыми авторами, не гарантируют получения езультатов с ожидаемой пользователем точностью. Исследования эффективности етода сопряжённых градиентов привели к открытию компенсационного эффекта для екоторого класса стартовых векторов, употребление которых даёт очень сильное величение скорости сходимости итерационного процесса. Представленный в рабо-е анализ вектора ошибок позволил выяснить причины компенсационного эффекта.

В диссертации изложены многие теоретические результаты, которые были олучены неожиданно, как побочный продукт исследований посвящённых другой ематике. Примером этому может служить алгоритм Sigma-SOR [8,] дающий заранее очную оценку оптимального параметра релаксации и описанный в подглаве 1.6. аиболее важным теоретическим результатом этой работы является показание ого, что минимальные значения спектрального радиуса и поддоминантного отно-ения в методах SOR выражаются той самой формулой, полученной первоначально oung'oM [3]. В случае медленно сходящихся задач алгоритм Sigma-SOR имеет чень сильные преимущества по сравнению с общеиспользуемым адаптационным етодом [3], в котором оптимальный параметр релаксации a)opt определяется динамически" во время выполнения итерационного процесса.

Некоторые исследования разностных схем в треугольной геометрии, выпол-енные автором, привели к модификации 7-точечной разностной схемы, позволя-)щей сократить в три раза число точек разностной сетки, а ошибка аппроксимации получаемых решений близка ко второму порядку. Модификация разностной :хемы и полученные результаты были представлены на конференции по реакторной >изике в Мюнхен в 1981г. [33].

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

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

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

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

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

2. Доказательство сходимости метода EQUIPOISE, основанного на стратегии лобальных-внешних итераций.

В 1977 Nakamura в своей книге [54] написал на стр.133, что: "Доказательство сходимости EQUIPOISE является интересной задачей, хотя не было проведено ни одно удовлетворительное исследование." (" The proof for convergence of EQUIPOISE scheme is an interesting problem, although not any satisfactory study has been done.").

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

Свойство сходимости итерационных стратегий, рассмотренных в главах 3 и 4 , поясняется на большом числе численных экспериментов для различных типов реакторных задач. Все вычисления были выполнены при помощи двухмерных программ HEXAGA, HEXSOR и HEXSLOR, основанных соответственно на алгоритмах метода AGA с двойной верхней релаксацией, точечном SOR и 1-линейном SOR (SLOR), решающих по 7-точечной разностной формуле в однородной треугольной сетке.

Анализ большого числа численных результатов полученных для рассмотренных задач, представляющих разные типы быстрых и термических реакторов, позволяет сделать вывод, что решения при наименьшей затрате на компьютерные вычисления могут быть получены при помощи итерационной стратегии, использующей 1 внутреннюю итерацию на одну глобальную итерацию в программе HEXAGA, алгоритм которой основан на одной из простейших моделей предфакторизационных методов AGA (описанных в главе 1) при использовании процесса двойной верхней релаксации как эффективной техники ускорения сходимости. Это приводит к выводу, что выбор матриц расщепления имеет определяющее влияние на эффективность метода, используемого для решения задач теории диффузии нейтронов и наблюдаемое поведение сходимости решений при помощи HEXAGA указывает на значительный потенциал двухкратных расщеплений происходящих от префакториза-ционных методов AGA.

Для всех рассматриванных реакторных задач, программа HEXAGA даёт решения от 2 до 5 раз более быстрые по сравнению с результатами программы HEXSLOR.

4. Априорная оценка значения параметра релаксации wbest, максимизирующего корость сходимости метода AGA, который используется в алгоритме программы EXAGA.

В случае алгоритма HEXAGA можно заметить несколько регулярных отклонений wbest от wopt в рассмотренных задачах и, кроме того, w t скорее всего нечувствительно к использованному числу внутренних итераций на одну внешнюю итерацию. Следующая формула best = "opt + --даёт очень хорошее приближение к верному значению cobest для большого класса реакторных задач, где g(W1) есть спектральный радиус итера

-1 ционной матрицы W при и = 1, обозначенный как "W^ = Мд NA; и caopt оптимальный спектральный радиус матрицы Ww в процессе двухкратной верхней релаксации, использованной в программе HEXAGA. Оценка cjopt в программе HEXAGA может быть получена при помощи эффективной техники описанной в работе [10].

5. Графическая интерпретация собственных значений матриц, характеризу-|щих рассмотренные итерационные стратегии.

6. Доказательство высокой достоверности результатов программы НЕХАСА.

Как можно заметить во всех рассмотренных случаях, НЕХАСА даёт решения при максимальной относительной ошибке 60-градусной симметрии сзут меньшей на десять порядков по сравнению с результатами, полученными при помощи НЕХБЬОК. В решениях, полученных по критерию с^ ^ 10~5, потоки в точках сетки, расположенных симметрично по отношению к более короткой диагонали шаблона, имеют те же самые значения при 15 значащих цифрах для НЕХАСА и 5 значащих цифрах для НЕХЗУЭЯ. В случае НЕХБШН значения езут сравнимы по величине с у с использованными для остановки итерационного процесса. Последние результаты, использованы в исследовании достоверности задачи термического реактора, вытекающей из задачи Зад.5.

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

В заключение стоит упомянуть, что НЕХАСА даёт решения при использовании арифметики с однократной точностью для больших реакторных задач, для которых решение при помощи НЕХЗЬОИ может быть получено только при использовании арифметики с двухкратной точностью, как было показано в задаче Зад.6.

БЛАГОДАРНОСТИ

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

Автор обязан значительной существенно-редакционной помощью своему другу р. Г. Енджеецу, за что горячо его благодарит.

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

Автор благодарит также Проф.Е.П.Жидкова и Др.3.М.Иванченко за ценные амечания, полученные им во время пребывания в ОИЯИ в Дубне.

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

1. J.Duderstadt and L.S.Hamilton, Nuclear reactor analysis, Wiley, New York, 1976.

2. E.L.Wachspress, Iterative solution of elliptic systems and applications to the neutron diffusion equation of reactor physics, Prentice-Hall, Englewood Cliffs, N.J, 1966.

3. L.A.Hageman and D.Young, Applied iterative methods, Academic Press, New York, 1981.

4. D.R.Ferguson and K.L.Derstine, Optimized iteration strategies and data management considerations for fast reactor finite difference diffusion theory codes, Nucl. Sei. Eng., 64, 593-604, 1977.

5. R.Fröhlich, Positivity theorems for the discrete form of the multigroup diffu- sion equations, Nucl. Sei. Eng., 34, 57-66, 1968.

6. M.L.Tobias and T.B.Fowler, The EQUIPOISE method A simple procedure for group-diffusion calculations in two and three dimensions, Nucl. Sei. Eng., 12, 513-518, 1962.

7. C.H.Adams, Current trends in methods for neutron diffusion calculations, Nucl. Sei. Eng., 64, 552-562, 1977.

8. Z.I.Woznicki, The sigma-SOR algorithm and the optimal strategy for the utilization of the SOR iterative method, Mathematics of Computations, 206, 619-644, 1994.

9. Z.I.Woznicki, Nonnegative splitting theory, Japan J. Indust. Appl. Math., 11, 289-342, 1994.

10. Z.I.Woznicki, Estimation of the optimum relaxation factors in partial factorization iterative methods, SIAM J. Matrix Anal. Appl., 14, 59-73, 1993.

11. Z.I.Woznicki, HEXAGA-11-120, -60, -30 two-dimensional multigroup neutron diffusion programmes for a uniform triangular mesh with arbitrary group scattering, Rep. KFK 2789, Kernforschungscentrum Karlsruhe, Germany, 1979.

12. Z.I.Woznicki, HEXAGA-111-120, -30 three-dimensional multigroup neutron diffusion programmes for a uniform triangular mesh with arbitrary group scattering, Rep. KfK 3572, Kernforschungscentrum Karlsruhe, Germany, 1983.

13. Z.I.Woznicki, On numerical analysis of conjugate gradient method, Japan J. Indust. Appl. Math., 12, 487-519, 1993.

14. G.H.Golub and C.F.Van Loan, Matrix computations, Johns Hopkins Press, Baltimore, 1983.

15. Benchmark Problem Book, Rep.ANL-7416, Suppl.2, 228-237, 1977.

16. Benchmark Problem Book, Rep.ANL-7416, Suppl.3, 1985.

17. Э. J.M.Barry and J.P.Pollard, Performance of a nonstationary implicit scheme for solving large systems of linear equations, Nucl. Sci. Eng., 92, 2733, 1986.

18. Э. I.К.Abu-Shumays, Vectorization of transport and Diffusion Computations on the CDC Cyber 205, Nucl. Sci. Eng., 92, 4-19, 1986.1. Rep.ORNL-5792, 1981.

19. Z.I.Woznicki, On reliability of solution symmetry in neutron diffusion theory problems, Proc. International Conference on the Physics of Reactors PHYSOR'96, Mito, Ibaraki, Japan, September 16-20, 1996, A21-A30.

20. G.I.Bell and S.Glasstone, Nuclear Reactor Theory, Van Nostrant Reinhold Company, 1970.

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

22. С.К.Годунов и B.C.Рабенький, Разностные схемы. Введение в теорию, М., Наука, 1973.

23. И.С.Березин и И.П.Жидков, Методы вычислений, Т.2. М., Физматгиз, 1962.

24. I.M.Gelfand and S.V.Fomin, Calculus of variations, Prentice-Hall, Engle-wood Cliffs, 1963.

25. R.S.Varga, Matrix iterative analysis, Prentice-Hall, Englewood Cliffs, 1962.

26. В.И.Лебедев, Итерационный метод с чебышевскими параметрами для определения наибольшего собственного значения и соответствующей собственной функции, ЖВМ и МФ, т. 17, вып.1, 1977.

27. G.Csordas and R.S.Varga, Comparison of regular splittings of matrices, Numer. Math., 44, 23-35, 1984.

28. Z.I.Woznicki, EWA-II Two-Dimensional, Two-Group Diffusion Fast Code, Kernenergie 14, 325-329, 1971.

29. Т.Апостолов и З.Возницки, Диффузионная двумерная программа HEXAGA II для многогрупповых расчетов гексагональных решеток, Атомная Энергия, т. 38, ВЫП. 6, 372-374, 1975.

30. Z.I.Woznicki, Comparison Theorems for Regular Splittings on Block Partitions, Linear Algebra and its Applications 253, 199-207, 1997.

31. Z.I.Woznicki, The Numerical Analysis of Eigenvalue Problem Solutions in the Multigroup Neutron Diffusion Theory, Institute of Atomic Energy, Otwock-Swierk, Report IEA-6/A (125p.), June 1995.

32. Z.I.Woznicki, AGA Two-Sweep Iterative Methods and their Applications in Critical Reactor Calculations, Nukleonika 23, 941-968, 1978.

33. Н.И.Булеев, Численный метод решения двумерных и трехмерных уравнений диффузии, Математический Сборник, т.51, No.2, 227-238, 1960.

34. Н.И.Булеев, 0 численном решении двумерныцх уравнений элиптического типа, Численные Методы Механики Сплошной Среды, т.6, No.2, 18-28, 1975.

35. Н.И.Булеев, Новый вариант метода неполной факторизации для решения двумерных разностных уравнений диффузии, Численные Методы Механики Сплошной Среды, т.9, No.l, 5-19, 1978.

36. R.S.Varga, Factorization and normalized iterative methods, Boundary Problems in Differential Equations, edited by R.E.Langer, University of Wisconsin Press, Madison, 121-142, 1960.

37. T.A.Oliphant, An extrapolation procedure for solving linear systems, Quart. Appl. Math. 20, 257-267, 1962.

38. H.Stone, Iterative solution of implicit approximations of multi-dimensional partial differential equations, SIAM J. Numer. Anal. 5, 530-558, 1968.

39. J.M.Ortega and W.Rheinboldt, Monotone iterations for nonlinear equations with applications to Gauss-Seidel methods, SIAM J. Numer. Anal. 4, 171190, 1967.

40. I.Marek and D.Szyld, Comparison theorems for weak splittings of bounded operators, Numer. Math. 58, 389-397, 1990.

41. T.Г.Апостолов, П.Т.Петков и З.Возницки, Двухмерная многогрупповая диффузионная программа EBA-II-30, Ядерна Энергия 3, 8-15, София, 1976.

42. Zbigniew Woznicki, Dwuprzebiegowe metody iteracyjne AGA rozwi^zywania duzych ukladow röwnari liniowych, Rocznik Polskiego Towarzystwa Matema-tycznego, Seria III: Matematyka Stosowana VI, 5-16, 1976.

43. М.Н.Зизин, Л.К.Шишков, Л.H.Ярославцева, Тестовые нейтронно-физические расчеты ядерных реакторов, М., Атомиздат, 1980.

44. В. M.Makai, Response Matrix of Symmetric Nodes, Nucl.Sci.Eng., 86, 302-314, 1984.

45. Г.И.Марчук, Методы расчета ядерных реакторов, М. , Атомиздат, 1961.

46. Г.И.Марчук, Методы вычислительной математики, М., Наука, 1989.

47. А.А.Самарский. Е.С.Николаев, Методы решения сеточных уравнений, М.,Наука, 1978.

48. Э. А.А.Самарский, Теория разностных схем, М., Наука, 1983.

49. А.А.Самарский, Введение в численные методы, М., Наука, 1987.

50. А.А,Самарский, А.В.Гулин, Численные методы, М., Наука, 1989.

51. Н.Н.Яненко, Метод дробных шагов решения многомерных задач математической физики, Новосибирск, Наука, 1967.

52. Н.С.Бахвалов, Основы вычислительной математики: Курс лекций, М. , Изд-во МГУ, 1970.

53. В.В.Воеводин, Численные методы алгебры. Теория и алгоритмы, М., Наука, 1966.

54. В.В.Воеводин, Ю.А.Кузнецов, Матрицы и вычисления, М., Наука, 1984.

55. Ф.Р.Гантмахер, Теория матриц, М., Наука, 1967.

56. Ю.А.Кузнецов, Итерационные методы решения несовместных систем линейных уравнений. Некоторые проблемы вычислительной и прикладной математики, Новосибирск, Наука, 1975.

57. Д.К.Фаддеев, В.Н.Фаддеева, Вычислительные методы линейной алгебры, Л., Наука, 1975.

58. V.P.Il'in, Iterative Incomplete Factorization Methods, World Scientific, 1992.

59. А.Хассит, Теория диффузии в двух- и трёхмерной геометрии. В кн. Вычислительные методы в физике реакторов, М., Атомиздат, 1972, с.50.

60. Ю.В.Воробьев, Случайный итерационный процесс, ЖВМ и МФ, Т. 5, No.6, 1964.

61. K.Pytel, W.Mieleszczenko, Z.Woznicki, Физические расчёты реактора EWA с топливом ВВР-М2, Внутренний отчёт ИАЭ номер: 145/R-V/90, Институт атомной энергии, Отвоцк-Сверк, 1990 (по-польски).

62. NPY Seminar on Numerical Solution of Multi-Dimensional Neutron Diffusion Equations, Warsaw, Poland, March 1969;

63. Seminar on Reactor Physics Calculations, Budapest, Hungary, October 1969;

64. The INEA Third International Advanced Summer School in Reactor Physics, Hercegnovi, Yugoslavia, September 1970;

65. Reaktortagung, Düsseldorf, West Germany, March 1976;

66. Conference on Numerical Methods of Mathematical Physics, Novosibirsk, Soviet Union October 1980;

67. ANS/ENS International Topical Meeting on Advances in Mathematical Methods for Nuclear Engineering Problems, Munich, Germany, April 1981;

68. International Meeting on Advances in Nuclear Engineering Computational Methods, Knoxville, USA, April 1985;

69. International Topical Meeting on Advances in Reactor Physics, Mathematics and Computations, Paris, France, April 1987;

70. International Conference on Linear Algebra and Applications, Valencia, Spain, September 1987;

71. International Conference on the Physics of Reactors: Operation, Design and Computations, Marseille, France, April 1990;

72. Copper Mountain Conference on Iterative Methods, Copper Mountain, Colorado, USA, April 1992;

73. Two Matrix Theory Meetings in Technion, Haifa, Israel, June 1993;

74. Third SI AM Conference on Linear Algebra in Signals, Systems and Control, Seattle, USA, August 1993;

75. The International Conference on Reactor Physics and Reactor Computations, Tel Aviv, Israel, January 1994;

76. ICIAM'95 The Third International Congress on Industrial and Applied Mathematics, Hamburg, Germany, July 1995;

77. WCNA'96 The Second Congress of Nonlinear Analysts, Athens, Greece July 1996;

78. PHYSOR'96 International Conference on the Physics of Reactors, Mito, Ibaraki, Japan, September 1996;

79. Czech-U.S. Workshop on Iterative Methods and Parallel Computing, Milovy, Czech Republic, June 1997;1. 15th IMACS World Congress on Scientific Computation, Modelling and Applied Mathematics, Berlin, Germany, August 1997;

80. XII Conference on Applied Mathematics PRIM'97, Palic, Yugoslavia September, 1997;

81. Sixth SIAM Conference on Applied Linear Algebra, Snowbird, Utah, USA, October 1997;