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

доктора физико-математических наук
Гинкин, Владимир Павлович
город
Обнинск
год
2004
специальность ВАК РФ
05.13.18
Диссертация по информатике, вычислительной технике и управлению на тему «Алгоритмы и комплексы программ для решения задач математической физики с использованием метода неполной факторизации»

Автореферат диссертации по теме "Алгоритмы и комплексы программ для решения задач математической физики с использованием метода неполной факторизации"

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

ГИНКИН Владимир Павлович

АЛГОРИТМЫ И КОМПЛЕКСЫ ПРОГРАММ ДЛЯ РЕШЕНИЯ ЗАДАЧ МАТЕМАТИЧЕСКОЙ ФИЗИКИ С ИСПОЛЬЗОВАНИЕМ МЕТОДА НЕПОЛНОЙ ФАКТОРИЗАЦИИ

05.13.18 — математическое моделирование, численные методы и комплексы программ

АВТОРЕФЕРАТ диссертации на соискание ученой степени доктора физико-математических наук

Обнинск - 2004 г.

Работа выполнена в Государственном научном центре Российской Федерации - Физико-энерг ическом институте им А.И.Лейпунского

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

член-корреспондент РАН,

доктор физико-математических наук, профессор

ЧЕТВЕРУШКИН Борис Николаевич

доктор физико-математических наук, доцент

НЕЧЕПУРЕНКО Юрий Михайлович

доктор физико-математических наук, профессор

ПОЛЕЖАЕВ Вадим Иванович

Ведущая организация - Российский научный центр «Курчатовский институт», г.Москва

Защита диссертации состоится « 29 » октября 2004 года в 9:00 часов на заседании диссертационного совета Д 201.003.01 при ГНЦ РФ-ФЭИ по адресу: 249033, г.Обнинск Калужской области, пл.Бондарснко, 1.

С диссертацией можно ознакомиться в библиотеке ГНЦ РФ-ФЭИ. Автореферат разослан «_»_2004 года.

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

Ю.А.Прохоров

<?¥ Л?

2005-4

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

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

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

Реальные задачи, как правило, являются нестационарными, нелинейными и несамосопряженными. В этих трех «не» заключена вся трудность решения пространственно-временных з

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

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

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

> ¿-'и*." I f '»«хЦ» «

^__fM ** <v. *

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

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

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

л.

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

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

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

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

2. Написаны специализированные комплексы программ: WIMS-ВОЛГА для стационарного расчета реакторов на тепловых нейтронах и W1MS-BOJ1HA для расчета кампаний реакторов типа ВВЭР, в которых впервые организован прямой расчет констант и выгорания по прецизионной программе WIMSD4 в каждом пространственном узле расчетной сетки в каждый момент времени. В рамках комплекса программ WIMS-ВОЛНА разработан алгоритм учета влияния формоизменения TBC на нейтронно-физические характеристики реактора.

3. Разработан комплекс программ VOLNA для расчета быстрых переходных процессов в реакторе в квазистатическом приближении в {hex, z) геометрии. На основе этого комплекса, а также комплексов программ теплогидравлического расчета GR1P-SM и расчета констант АРАМАКО, создан комплекс программ GVA для совместного трехмерного нейтронно-физического, теплогидравлического и термомеханического расчета аварий реакторов на быстрых нейтронах. С помощью этого комплекса программ выполнены трехмерные нестационарные расчеты аварий реактора типа БН-800, вызванных остановкой главных циркуляционных насосов (ГЦН) и несанкционированным движением группы компенсирующих стержней органов регулирования. Разработаны метод и программа расчета деформаций TBC и алгоритм учета теплового расширения активной зоны реактора на его нейтронно-физические характеристики.

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

программа вЮАТЧ для трехмерного нестационарного расчета тепловой конвекции и выполнены расчеты бенчмарков по конвекции воздуха в кубической полости, подогреваемой с одной из боковых сторон, на сетках большой размерности (до миллиона узлов). Показано, что для чисел Релея 106 и выше необходимо использовать еще большее количество узлов или применять специальное сгущение сеток для обеспечения приемлемой точности расчетных значений чисел Нуссельта.

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

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

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

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

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

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

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

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

Автор выносит на защиту:

• Варианты метода неполной факторизации: метод ¿-факторизации, метод параболических прогонок, вариант метода сопряженных градиентов с предобуславливателем по схеме неполной факторизации с периферийной компенсацией итерируемых членов СвРШ, комбинированные схемы неполной факторизации НРРР, НРВЗ, ШВ4, С1Р. Анализ сходимости блочных итерационных методов неполной факторизации, доказательство сходимости метода Л-факторизации, метода параболических прогонок, явной схемы Булеева с диагональной компенсацией итерируемых членов.

• Алгоритмы и программы стационарного и квазистационарного трехмерного расчета реакторов ВОЛГА и ВОЛНА. Комплекс программ WIMS.BOЛГA для расчета реакторов в (х, у, г) и {г, 9, г) геометриях. Комплекс программ 'ММБ-ВОЛНА для расчета кампаний реакторов типа ВВЭР в

{hex, Ч) геометрии с учетом влияния формоизменения TBC на нейтронно-физические характеристики реактора вследствие термомеханическнх деформаций.

• Алгоритм и комплекс программ VOLNA для расчета быстрых переходных процессов в реакторе в квазистатическом приближении в (hex, z) геометрии. Алгоритм и комплекс программ GVA для совместного нейтронно-физического, тепло1гидравлического и термомеханического расчета аварий реакторов на быстрых нейтронах типа БН-800 с учетом формоизменения TBC вследствие термомеханических деформаций. Результаты расчетов аварий реакторов типа БН-800, вызванных остановкой пйвных циркуляционных насосов и несанкционированным движением группы компенсирующих стержней органов регулирования.

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

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

Апробация работы. Основные положения и результаты диссертационной работы докладывались и обсуждались на ведущих тематических отечественных и международных конференциях: на IV Всесоюзном совещании по вычислительным методам линейной алгебры (Новосибирск, 1977), на Vli, IX, X, XI школах по моделям механики сплошной среды (Батуми, 1983; Якутск, 1987; Казань, 1989; Владивосток, 1991), на 13

Международном конгрессе по вычислительной и прикладной математике (IMACS'91, Дублин, 1991), на 7 международной конференции по численному анализу полупроводниковых приборов и интегральных схем (NASE CODE VII, Колорадо, США, 1991), наЗ, 4,5 международных семинарах по моделированию устройств и технологий (Обнинск, 1994; Южная Африка, 1996,1998), на 1, 2, 3, Российских национальных конференциях по теплообмену (РНКТ, Москва, 1994, 1998, 2002), на 9 международном совещании по безопасности ядерных реакторов (Москва, 1995), на III Минском международном форуме по тепло- и массообмену (Минск, 1996), на 1, 2 3 международных симпозиумах по достижениям в вычислительном тепломассопереносе (Турция, 1997; Австралия, 2001, Норвегия, 2004), на международной алгебраической конференции, посвященной памяти Д.К.Фадеева (С.-Петербург, 1997), на семинарах «Нейтроника-98», «Нейтроника-99» (Обнинск, 1998, 1999), на 3, 4, 5 международных конференциях «Рост монокристаллов и тепломассоперенос» (Обнинск, 1999, 2001, 2003), на международной конференции по вычислительному тепломассопереносу (Кипр, 1999), на IX и X Национальных конференциях по росту кристаллов (НКРК, Москва, 2000, 2002), на 4 Международной конференции по супервычислениям в ядерных приложениях (SNA 2000, Япония, 2000), на 4 Сибирском конгрессе по прикладной и индустриальной математике (ИНПРИМ-2000, Новосибирск, 2000), на международной конференции по твердым кристаллам (Польша, 2000), на 13 международной конференции по росту кристаллов (ICCG-13, Япония, 2001), на IV Международной конференции по неравновесным процессам в соплах и струях (NPNJ-2002, С.-Петербург, 2002), на 1 международной конференции по приложениям пористых сред (Тунис, 2002), на международной конференции по проблемам в тепловой конвекции (Пермь, 2003) и др. (всего на 67 международных и национальных совещаниях и конференциях).

Структура и объем диссертации. Диссертация состоит из введения, 5 глав и заключения (общие выводы по работе). Общий объем диссертации - 258

страниц, в том числе 123 рисунка и 46 таблиц, список литературы из 240 наименований на 25 страницах.

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

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

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

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

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

членов. Ниже приведено краткое описание метода сопряженных градиентов с предобуславливателем по этой схеме.

Пусть требуется решить уравнение Ay=f, где А - конечно-разностный оператор, аппроксимирующий уравнение эллиптического типа. Этот оператор действителен и симметричен, т.е. А = А1.

Метод сопряженных градиентов для решения этой системы уравнений имеет вид

Н . Н , )-1 /-1 У /-1 ,1 J-I 1 J-i ! 1

г = Ау -/, z =К г , g =Z ~ь g , у = у -ag',

Ь' =

J = l ( j-i i \ к".,'-1 ... -g]

и н\ J>1 W,g1]

W >S

Для решения пространственных конечно-разностных задач в качестве предобуславливающей матрицы К выберем матрицу из метода неполной факторизации. Метод неполной факторизации для этого случая имеет следующий вид: (А + в)у' =/+ Вун, j = 1,2,... где

(А + В) = MN, N = Мr, М г, транспонированная матрица,

\м' у' = г.

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

Ау = -ау^ -byk l -gy( l-гу^мГа^ум-Ъмум-gMyM "V*-.^.,*-, + РУ,

/ = 1.....я,; * = 1 = \,..лу

Здесь и далее тройной индекс /,*,/, для простоты записи, опускаем и пишем только ту часть индекса, которая отличается от i,k,l.

Выберем оператор М в виде нижней пятидиагонапьной матрицы Ш - рг-аг(ч -уг,_и+1-> ГДе Р> Р, У, Я - неопределенные

пока коэффициенты. Тогда оператор Л/2 будет верхней пятидиагональной матрицей

Мту = ру-аму^ ~РШУШ ~ГМУМ -К^У,^.,-

Эти выражения задают явную схему неполной факторизации. Найдем вид факторизованного оператора ММ7 : ММГу = р(ру-амум-0мум -ГмУм-*-М*-хУм*-^-

~ а¡+1,1-^1+1.1-1 ~ Рк+и-\Ук+и-\ ~УУ ~ ) ~~

~~ Х(Р,->М1У^.к+\ ~ак+^Ук+\ ~ Д-1,¿+2-^1-1,4+2 ~

По определению, это выражение равно (А + В)у. Отсюда находим явный вид итерируемого оператора В: Ву = (ММТ - А)у = (р2 +а2 + р1 +у2 - р)у- {арм - а)у^ -

~^Ркл ~а\л ~ЪУкл ~(РРМ ~Лак+, ~Ьк+,)Ук+, ~ -<>Р„ ~8)У,-> ~(РГы ~ем)Ум -(лР,-кк+1 ~г)У,-и+1 ~

~(РЛ.+,,к-, ~Ра,+кк-х ~гм,к->)У/+.,*-. +аГ,-и+1У,-ин +

+ РГк-и+>Ук-и+1 +гРк+и.1Ук,и.1+Рл^.гУММг+Щ-1М2У1-1М2 +

+ 1-1.4+1,/+!-У1-и+и+1 + '

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

Выберем компенсирующие члены так, чтобы итерируемое выражение для Ву приняло вид:

ВУ = - ) + (У - У,->) +

+ 7°мм - Фм) + - у,-,) +

-1./+1

(У* -и* I

+, - )+0> - )+ +- )+(У - У>-1М1)+ + Щ-хмг (У.-ХМ2-°Ука) + (У ~ У*-, > +

+ (У,-цч и+. - и)+Щ-,Умм1 (У-У >-,)+

+УК и

где в - итерационный параметр. При в = О компенсация отсутствует; при в = 1 компенсация полная, (то есть, Ву = 0, если у1 ^ - константа); при 0 < в < 1

компенсация частичная.

Перегруппировав члены в этом выражении для Вр, запишем его в виде: Ву = Ща, ,У(, +Ц1лХк ^Р^Л-^К/.-гм^У-^Мч -

-и+1 У 1-М+1

Потребуем теперь, чтобы это выражение тождественно совпало с

т

полученным выше выражением для Ву = (ММ - А)у. Вычитая их друг из друга получим:

О = (р2 + а2 + рг + Гг + Л2 - р - 2$(а1{у(1 + + + У,-^ "

-(РГМ -в(<*У,-1т+РУк-м+1 + -

-(яЛнл , -'•/и,*-. >■

Приравнивая нулю пять первых членов этого выражения, получим

. 'i

систему уравнений

ср2 +a2 + ^J +r2 + Л2-р-Ща^г» + /W*-. +Дм

(М- -4-,=0' 0-Л-, - « "Ж«,.,^+ А-/*-,+ =°>

<ЛЛ-и+, Ч-Л-иА-и+. =0>

решая которую, находим рекуррентные формулы для вычисления неопределенных коэффициентов, задающих матрицу М:

р = ^р-а2 - р2 -у2-Л2 +Ща, + + /U+A-, + Л,ч гЬ1>Ь1) ; у = (g + + + ))/рм ;

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

На численных примерах показано, что зависимость числа итераций J от числа узлов по одному направлению и для схемы CGPÍF пропорциональна 4п. Оптимальное значение итерационного параметра близко к единице.

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

свойством быстро гасить гладкие составляющие ошибки. Приведем описание этой схемы.

Пусть требуется решить уравнение А<р=/ , являющееся разностным аналогом двумерной задачи эллиптического типа:

А(Р = + Р,,^'

/ = 1,2.../я; к = 1,2..л;

>°; +а,у

Запишем метод неполной факторизации в следующем виде: Мг = / + В<рн

. У —

Ы<р = г

Структура матриц и В задается, й элементы этих матриц

определяются из тождества = А + В. Для схемы Ь-факторизации эти матрицы имеют вид:

= У,^,* "А,*^-. = V

В(Р,,к -«у 14-м+ 9» + + -«'„и

где

а,л = "а /(Гм.* ~ 4-и "4-иУ,а = Л. " "а <*.-и " 4-и " 4-и 4л=*и+аи4-и; 4*=<,и+а.Л-и;

Из выражения для В<р^ следует, что если не зависит хотя бы от

одного из индексов / или к, то Вр^ = 0, и схема становится безытерационной.

На основе сформулированной и доказанной автором в первом разделе главы теоремы 12 доказана сходимость схемы 1ьфакторизации и показано, что скорость сходимости ее в асимптотике пропорциональна т.

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

параболических прогонок и показано, что эта схема обладает свойством быстро уничтожать весь спектр ошибок, кроме гладкой составляющей. Поэтому автором была введена и численно исследована комбинированная схема неполной факторизации, в которой схема параболических прогонок чередуется со схемой Ь-факторизации (схема НФПП). Численно показано, что при этом имеет место значительное ускорение сходимости итерационного процесса по сравнению с каждой из входящих в нее схем в отдельности, так как схема параболических прогонок быстро сглаживает начальную ошибку, а схема Ь-факторизации быстро уничтожает гладкую составляющую ошибки. Этот эффект аналогичен действию многосеточных методов, в которых мелкая сетка используется для сглаживания ошибки, а крупная - для быстрого уничтожения гладких составляющих ошибки. Численно показано, что такой же эффект наблюдается при чередовании схемы Ь-факторизации с третьей и четвертой схемами Булеева (схемы НФВЗ и НФВ4).

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

Комбинированные схемы наиболее эффективны при решении задач с несимметричными матрицами коэффициентов, для которых метод сопряженных градиентов не работает. Но даже для случая симметричных матриц коэффициентов, как показал численный эксперимент, скорость сходимости схемы С1Р пропорциональна V», также как скорость сходимости схемы ССР1Р. Для несимметричных же матриц коэффициентов неявные схемы неполной факторизации (а также комбинированные схемы) сходятся

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

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

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

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

-1 (я)

-ЧЬ^чР+Е^чР-^^-^—а, g = 1, 2,..., N0, " К*

Ш (г) (*)

где (2= 2-1Ф ~ источник нейтронов деления, ¿^-эффективный

¡н

коэффициент размножения. В каждой геометрической области константы

Г^1, постоянны для данной группы

В квазистационарном расчете реактор на каждом временном шаге выводится в критическое состояние (соответствующее К = 1 ), путем

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

Главным результатом данной главы является описание созданных

автором комплексов программ WIMS-ВОЛГА и WIMS-ВОЛНА. Первый из них предназначен для трехмерного стационарного расчета реакторов в (x,y,z) и (г, в,г) геометриях, второй - для расчета кампаний реакторов типа ВВЭР в (hex.z) геометрии. Большое внимание уделяется описанию организации вычислений в комплексах программ, объединяющих программы подготовки констант, программы нейтроно-физического расчета, расчета теплогидравлики, выгорания топлива, перегрузок топлива, а также исследованию погрешности аппроксимации, влиянию изменений толщины водяных зазоров между TBC на нейтронно-физические характеристики реактора, и т.д.

В программах ВОЛГА и ВОЛНА на внутренних итерациях для решения групповых уравнений диффузии нейтронов используются эффективные схемы неполной факторизации. Оба комплекса в качестве константной базы используют прецизионный комплекс программ WIMSD4.

Отличительной чертой комплекса программ WIMS-ВОЛНА является прямой расчет выгорания в каждой пространственной точке реактора в каждый момент времени по программе WIMSD4. Это позволяет наиболее полно и корректно учесть изменения состава топлива и всех технологических параметров во взаимосвязи в течение кампании для реакторов типа ВВЭР с различными видами TBC и различными топливными загрузками, что важно при проведении поисковых исследований.

В комплексе программ WIMS-ВОЛНА для расчета кампаний реакторов типа ВВЭР реализован следующий алгоритм пространственно-временного расчета реактора:

-рассчитывается начальный набор констант SO, начальное распределение плотности нейтронного потока и начальное распределение числа делений FO(r);

- по известным значениям F0(r) рассчитывается выгорание по WIMSD4 и новый набор констант SI;

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

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

- по полученным значениям Fl(r) рассчитывается выгорание по WIMSD4 до следующего момента времени.

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

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

- энерговыделение по радиусу топливного сердечника считается постоянным;

- растечкой тепла вдоль оси Z, тепломассообменом между ячейками соседних ТВЭЛ и TBC, изменением газового зазора между топливом и оболочкой, влиянием примесей борной кислоты на теплогидравлические свойства теплоносителя пренебрегается;

-результаты теплогидравлического расчета ячейки ТВЭЛ распространяются без изменений на все ячейки данной TBC;

В рамках сделанных предположений распределение температуры

теплоносителя по высоте канала имеет вид:

где qt - тепловой поток в расчёте на единицу длины ТВЭЛа;

Gt - расход теплоносителя через /-ую TBC;

Ср - теплоемкость теплоносителя.

Линейный тепловой поток q{z) связан с объёмной плотностью теплового потока в топливном сердечнике, который, в свою очередь, зависит от удельной плотности числа делений и мощности реактора.

Описанная модель реализована в отдельной программе. Исходными данными для неё являются: поле удельных плотностей чисел делений, расходы теплоносителя по каждой TBC, температуры на входах TBC, среднее давление в теплоносителе и число ТВЭЛов в активной зоне. Результатом расчёта являются: трёхмерные поля температур топлива, оболочки и средних температур оболочек ТВЭЛов по кассетам.

Результаты тестирования комплекса WIMS-ВОЛНА на расчете кампании реактора ВВЭР-1000 со штатной загрузкой и трехгодичным топливным циклом показали хорошее совпадение с результатами по аттестованным программам.

Анализ результатов расчетных исследований основных источников погрешности в нейтронно-физическом расчете реактора по комплексу программ WIMS-ВОЛНА позволил отработать оптимальную методику расчета констант по программе WIMSD4 на всех стадиях топливного цикла. При расчете распределений энерговыделения по активной зоне реактора ВВЭР-1000 по программе ВОЛНА рекомендуется использовать семь точек на кассету. При этом погрешность получаемого решения оказывается не более 8% в распределении энерговыделения по кассетам активной зоны реактора.

Приведенные примеры расчетов по комплексу программ WIMS-ВОЛНА кампаний реакторов ВВЭР-1000 различных блоков Балаковской, Запорожской и Калининской АЭС при различных загрузках топлива и схемах его перегрузок, в которых испытывались новые усовершенствованные TBC разных типов и разных глубин выгорания показали, что коэффициенты неравномерности энерговыделения удовлетворяют технологическим требованиям. Результаты этих расчетов были использованы в термомеханических расчетах поведения

активных зон реакторов с новыми TBC в течение кампании.

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

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

Решается нестационарное диффузионное уравнение реактора

vw dt ** i-i i-i ' ' ' g = \,...,NG,

где Q = V vZ^'V''* - источник деления, w

Эта система NG уравнений дополняется уравнениями для концентраций ядер - предшественников запаздывающих нейтронов:

I = 1,..., NC;NC- число групп запаздывающих нейтронов (обычно NC=6).

Представим искомую функцию <pig)(r,t) в виде произведения двух функций (рЫ)(г,г) - P(t)y/lR)(r,t), где P(t) -амплитудный фактор, a y/®(r,i) групповая форм-функция. Тогда исходная система уравнений представляется в виде двух систем уравнений, эквивалентных исходному уравнению: -системы уравнений типа точечной кинетики для P(t)

-и системы пространственно-временных уравнений для групповых форм-(в)

функций Iff (r,t)

dw(t\ 1 dP{t)_M)

Л P(t) dt

- VD(g)Vy{g) +2TjVte> =

■f.

м

гдее = /§у

-„/«I L tel J

- скорость распада предшественников запаздывающих нейтронов в момент временила р(/)> Ж0> -^(О. Р, (0 вычисляются по формулам

1 ж; " «с , ,„ лк> ,„ ., ,.

К 8-1

/-«+1

,<*>I •Чг)

Л-

Г г я*

д4Е хТр&ГМ m=±m

r У ' 1-1 J

r rg=IV у g=i

интефирование ведется по всему объему реактора, - групповая функция ценности нейтронов, соответствующая невозмущенному состоянию реактора.

При этом (it) по физическому смыслу близко к значению реактивности, интерпретируемому как относительное изменение «динамического» коэффициента размножения p = (kd-1 )jkd, вызванное изменениями в

макросечениях; Л(г) - время генерации нейтронов; /? (/) - эффективная доля запаздывающих нейтронов /-той группы; fl(t) - общая эффективная доля запаздывающих нейтронов. Если считать, что P(i) представляет собой мощность реактора, то С является числом предшественников, умноженным на

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

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

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

широких пределах.

Описанная выше методика реализована в комплексе программ УОЬЫА. На основе этого комплекса был создан комплекс программ йУА (ОЮТ-БМ, УОША, АЯЛМАКО), позволяющий рассчитывать аварии в реакторах на быстрых нейтронах. В этом комплексе на каждом временном шаге осуществляется теплогидравлический расчет реактора по программе вШЕ-БМ, передача файлов с распределением температур и плотностей натрия из программы вЯ^М в комплекс программ расчета групповых нейтронно-физических макроконстант АРАМАКО, передача констант в комплекс нейтронно-физического расчета реактора УОЬЫА, и обратная передача из программы УОЬЫА распределений энерговыделений в программу СШР-вМ.

Комплекс СУ А состоит из двух частей:

1. блок стационарного расчета критического состояния реактора и функций ценности;

I

2. блок нестационарного расчета реактора.

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

1. задаются средние значения энерговыделений по активной зоне реактора;

2. по заданным значениям энерговыделений по программе СЮР-ЭМ рассчитываются распределения температуры и плотности натрия по активной зоне реактора;

3. по заданным значениям температур и плотностей натрия по программе АРАМАКО рассчитываются нейтронно-физические константы для каждого элементарного объема пространственной гексагональной сетки, отличающегося хотя бы одним из значений температур и ядерных плотностей материалов;

4. решается условно-критическая задача по программе УОЬЫА, находятся Кф распределения плотностей потоков нейтронов и энерговыделений;

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

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

На стадии нестационарного расчета реализуется неявная схема Эйлера. При этом необходимо, в зависимости от сценария развития аварии, определить критерии для выбора шага крупномасштабной временной сетки при пересчете пространственной форм-функции. Выбор шага крупномасштабной временной сетки tk для пересчета форм-функции осуществляется на основе опыта и диктуется необходимостью пересчета температурных полей по критериям, принятым в программе теплогидравлического расчета GRIF-SM. Как правило, этот шаг находится в диапазоне 0.001-0.1 сек.

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

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

При переходе к следующему временному шагу крупномасштабной сетки происходит усреднение всех рассчитанных зависимостей p(t), ß,{t), A(t), P(t) по предыдущему временному шагу & и экстраполяция полученных значений плотностей потоков нейтронов на середину следующего шага.

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

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

С помощью комплекса GVA на примере реактора БН-800 были рассчитаны и проанализированы аварии с остановкой главных циркуляционных насосов (ГЦН) с одновременным отказом органов аварийной защиты, а также аварии, обусловленные несанкционированным движением группы компенсирующих стержней. Приведены результаты параметрических исследований этих видов аварий с учетом радиальных расширений активной зоны. Показано, что величина радиального расширения активной зоны реактора существенно влияет на протекание переходного процесса. Приведены результаты параметрических расчетов аварий, вызванных самоходом группы компенсирующих стержней.

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

где Рх1 и Ру, - изгибные силы, действующие на сборку по направлениям х и у на слое i (i=l,2); Qm, - контактная сила, действующая на m-тую грань шестигранного чехла на уровне /. Кроме того, должно выполняться равенство контактных напряжений, действующих на смежные грани двух соседних сборок А и В в местах контактов и условие совместности деформаций.

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

быстро изменяющегося поля температур характеризуется изменением во времени размера TBC под ключ за счет температурного расширения, изменением координат осей TBC за счет свободного температурного изгиба и за счет изгибных сил в местах контактов Кроме того, изменяются прогибы граней чехлов TBC, но они не влияют на координаты осей TBC. Используя эти условия, записывается система уравнений для вычисления контактных напряжений Q. AlQi=Fl для/с!Р, ß, = 0 дпя/<=(Я/¥^

где Q - вся совокупность возможных мест контактов TBC между собой, Ч* - совокупность реальных мест контактов.

Решение полученной системы уравнений однозначно определяет множество Ч*: /с У при Q\ > 0. Так как множество Ч* заранее не известно, то решение системы уравнений ищется методом установления. Вначале по заданному приближению множества Ч'* решается система уравнений и находитсся Q?, затем, используя условие Q\ > 0, по найденному значению Qi вычисляется новое приближение для множества Если 4**" * Ч*\ то процесс вычислений ß/ повторяется. Описанный алгоритм термомеханического расчета реактора реализован в рамках комплекса программ GVA. Приведен результат демонстрационного расчета деформаций TBC в активной зоне реактора типа БН-800 при аварии, вызванной остановкой главных циркуляционных насосов с одновременным отказом органов аварийной защиты.

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

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

Система уравнений магнитной гидродинамики в приближении Буссинеска имеет следующий вид

рЦ-+р{\У)У = -««*( р+1+-¿-(нУ)Н+аыЧ+{ргт+рсс%,

ы

Нг\ 1 ^

8яг I 4л

рс,

дТ дI

+ (уу)Г

дх 4 '

где Д,. =

А

я.г

Г .

( Р{дс)р.г

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

р1-о! =

еж7 15/ ажу>

дг

~е г г дг ' дг

У.дУУ* 15 ,д*Г' г -е

V---е -

' дг дг 1 дг

г д<р г д<Р 9 д<Р + / = 1,..., 8.

Таблица. Коэффициенты и правые части уравнения.

] Ж7 г <• г' о!

1 V, V т V * У г V 1 г2 + г1/ к У г ¿1/ 1 дР р дг

2 V, V г V р V г К V 1 — + — + — г г* А ✓Л * А рг д<р

3 V, V г V г У 2 1 М у У ' Л! 1дР р дг

4 Нг V жг V т 9 У -1 -IV) ¿г »-а- ' н й р" +-1- А 0

$ Н, V я V •V Ут 1 У \дуа ___г___*__9 А г г д<р и & * /Ц 0

6 Н* V V т. У тг \ V дУ ___Г___2_ А г дг и Й * & 0

7 т а г а * а г У* % 0

8 с в, О» а У* с/ /м 0

где

' г г2 дг 4яр

т г дф Аяр

ДК 1а/р Га/, ют/\,ягя.'

\дг г д<р) \ & г др) I-

4яр

„ (дН, дНг\ „ (1 дН, дНЛ

Г."

г"

Применяя экспоненциальное преобразование, приведем это уравнение к диффузионному виду (индекс] для простоты опускаем):

15. dW l д 8W 1 д Ш m r. n

--r —----——tus —---—Tj£ -+ yW = F-Q

rXdr r dr r со dp * dq> rjdz • dz и

Используя интегро-интерполяционный метод и разнесенные сетки

построим конечно-разностные аналоги этих уравнений вида:

-crjQV^ -Ю-смцл{WM -W) + yDW = FD, I = 1,2,...л ; к = 1,2,...и ; / = 1,2,...«

с" ' ' <р' ' г

г . е D е D е D

Н Vi П-1 г,а

Гдеа = -^—7~; Ь= г* ' , ' с= л \ D =

r'Aqt ,Aq>k Л,

г K~t " 2

С V V

-rift. .'-isibd. ..Hiti

е. 2 «.i е 2

Vi ',i

A = е ' ; со = е 1 ; »/ = е г

Г ---Т*-------------т - -

ее ее ее

г.-> «V ».ч г,

Матрицы этих систем уравнений в общем случае не симметричны (если коэффициенты Л, со, т] отличны от единицы). Симметризуем полученное уравнение, прибавляя к обеим частям этого уравнения следующее выражение а(А-!)»;_, + ам(Л-1 -Щ^+Ь^-Щ^ +

В результате в левой части полученного уравнения исчезнут коэффициенты А, а, т\, и мы получим симметричное уравнение вида:

, -а Ж , -ЬИГ , -Ъ Ж , -<?Г , -с, Ж,, = /

(-1 /+1 /+| 4-1 4+1 4+1 /-1 /+1 /+1 ^ ■>'

/ = 1,...«,; к = 1,...пг; / = 1,...и3;

+ 1)^ , + (А"1 -1)+ Ь(со-1)^ , + +<-/+|('7 ' -О»',.,. Аналогично поступим с уравнениями для компонент скоростей.

Уравнение для давления получается из конечно-разностного представления уравнения неразрывности. Преобразованные уравнения линеаризуются с использованием неявного метода установления.

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

а > О, А > 0, с > 0, а > 0, д^а + Ь + с + а +Ь ,+с. ,.

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

Метод реализован в программе ОЮА1М, предназначенной для расчета тепловой конвекции в (х, у, г) геометрии. Приведены результаты расчетов по программе GIGAN численных бенчмарков о конвекции воздуха в кубической полости с подогреваемой боковой стенкой, для различных чисел Релея на сетках до миллиона узлов. Показано, что для получения хорошей точности решения, для чисел Релея более, чем 106, требуется либо использовать еще большее количество узлов, либо применять специальные сгущения в областях больших градиентов искомых решений.

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

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

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

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

а) задача теплопроводности;

б) задача Стефана.

Поверхности зон могут быть излучающими.

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

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

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

- условие теплового излучения на излучающих границах;

- условия 1 -го, 2-го и 3-го рода.

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

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

Шаг 1. По заданным значениям температур на излучающих поверхностях вычислить потоки тепла на излучающих границах. На неизлучающих границах выразить потоки тепла через температуры.

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

: »»ОС НАЦИОНАЛЬНА®

библиотека

СПе—

09

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

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

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

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

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

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

Система уравнений конвективного тепломассопереноса в этой модели имеет вид:

а ра й

е,—= )с

* д! д1 к

где и ^ означают доли жидкой и твердой фазы, соответственно: при А ¿А,,

=

А-А

-при А, <А<А1 = 1 -^

1 при А > А + А,-максимальное значение энтальпии твердой фазы, £ - скрытая теплота кристаллизации, К = коэффициент пористости, аг - эмпирическая

«константа, К^ - коэффициент сегрегации примеси. Температура выражается

через энтальпию следующим образом:

. (А-А,) Т + ' при А < А,

7 = •{Г* при А, <А<А, + £,

. А-А,

Т +-^ при А>А,+£,

СДА)

где С^ (А) - теплоемкость твердой фазы С^(А)-теплоемкость жидкой фазы Т -температура кристаллизации.

Усредненные свойства вещества, находящегося в промежуточном состоянии ^ < А < ^ + £, определяются следующим образом: X = Хчва + X [е1,

где параметр X принимает значения р, Л, С£>.

Оценка значения эмпирической константы а была произведена из условия совпадения результатов расчетов распределения примеси с данными предыдущей модели, применявшейся для объяснения аномальных распределений примеси в экспериментах по кристаллизации йе, легированного ва в космосе, где было принято возрастание значения вязкости вблизи границы фазового перехода в области температур ( Т +10К) в 60 раз. Такое сравнение приводит к значению а=1-1017.

Приведены результаты расчетов экспериментов по кристаллизации расплавов в космосе с использованием описанных моделей.

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

В заключении сформулированы основные результаты диссертации, полученные впервые. Приведем эти формулировки в кратком изложении.

1. Показано, что для решения конечно-разностных задач с симметричными

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

сопряженных градиентов с предобуславливателем по схеме неполной

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

предложенной автором. Скорость сходимости этого метода пропорциональна ¡1

п , где и - число узлов по одному направлению.

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

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

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

4. Разработан комплекс программ WIMS.BOJirA для трехмерного расчета реакторов в стационарном многогрупповом диффузионном приближении в (х, у, 2) и (г, в, z) геометриях. Он прошел необходимую верификацию, включен в отраслевой фонд алгоритмов и программ расчета ядерных реакторов ОФАП ЯР, и широко использовался для исследований критичности гомогенных и гетерогенных размножающих систем.

5. Разработан комплекс программ WIMS-ВОЛНА для расчета кампаний реакторов типа ВВЭР в квазистационарном многогрупповом диффузионном приближении в (hex, z) геометрии. Отличительной чертой комплекса является использование для расчета выгорания и констант в каждой пространственной точке реактора в каждый момент времени в процессе его эксплуатации с помощью прецизионной программы WIMSD4. При этом наиболее корректно учитываются реальные значения всех технологических параметров реактора. Разработан алгоритм учета влияния неравномерных водяных зазоров между TBC на нейтронно-физические характеристики реактора. С помощью этого комплекса выполнены расчеты кампаний реакторов ВВЭР-1000 различных

блоков Калининской, Запорожской и Балаковской АЭС при различных загрузках топлива и схемах его перегрузок с участием усовершенствованных TBC новой конструкции, в обосновании работоспособности которых использовался комплекс программ WIMS-BOJIHA.

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

7. На основе программ нейгронно-физического расчета VOLNA, теплогидравлического расчета GRIF-SM, и расчета констант АРАМАКО создан трехмерный комплекс программ расчета аварий реакторов на быстрых нейтронах GVA. В рамках этого комплекса разработаны алгоритм и программа расчета термомеханического поведения активной зоны реактора и алгоритм учета влияния теплового расширения активной зоны на нейтронно-физические характеристики реактора. Разработанный комплекс программ GVA позволяет моделировать сложные аварийные процессы в реакторах на быстрых нейтронах, приводящих к значительным изменениям концентраций компонентов реактора в переходном процессе, а именно таких, как кипение натрия или перемещение расплавов стали и топлива после предполагаемого расплавления активной зоны. С помощью этого комплекса выполнены тестовые расчеты аварий реакторов типа БН-800, вызванных остановкой главных циркуляционных насосов с одновременным отказом органов регулирования, и извлечением группы шести компенсирующих стержней, которые продемонстрировали работоспособность и эффективность комплекса. Они показали, в частности, что локальные возмущения энерговыделения в активной зоне могут в несколько раз превышать интегральное изменение мощности

i i I

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

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

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

«

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

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

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

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

Основные результаты диссертационной работы опубликованы в 8 статьях в реферируемых журналах, 29 докладах в сборниках трудов Международных и Российских конференций, 2 статьях в зарубежных журналах, 20 препринтах, трудах международных и отраслевых семинаров, сборниках трудов организаций и других работах. Разработанные автором методы параболических прогонок и /{-факторизации включены в виде отдельных разделов в монографию В.П.Ильина «Методы неполной факторизации для

решения алгебраических систем» Москва, Наука, Физматлит, 1995, стр.128-139, изданную также за рубежом в переводе на английский язык.

Статьи в рецензируемых журналах

1. Гинкин В.П. О влиянии релаксации на сходимость схемы ¿-факторизации при решении двумерных уравнений типа диффузии // ВАНТ, Физика и техника ядерных реакторов. - 1980. - Вып.4 (13). - С.196-202.

2. Гинкин В.П., Гадияк Г.В., Шварц H.JI. и др. Программа расчета стационарных характеристик МДП-транзистора MOS-1 // Автометрия. -1987. - №1. - С.70-75.

3. Гинкин В.П., Крылова J1.M., Кухтин В.П. Комплекс программ WIMS.ВОЛ ГА // ВАНТ, Физика и техника ядерных реакторов. - 1988. -Вып.З. - С.56-57.

4. Гинкин В.П., Васильев Ю.Ю., Гурин В.Н. и др. Расчеты критических параметров гомогенных и гетерогенных систем различных спектров на базе малогрупповой библиотеки микроконстант // ВАНТ, Физика и техника ядерных реакторов. - 1988. - Вып.4. - С.46-51.

5. Гинкин В.П. Эффективный метод решения одного класса задач эллиптического типа при отсутствии диагонального преобладания // Автометрия. -1988. - №5. - С.64-67.

6. Гинкин В.П., Кузеванов B.C., Иваненко И.Ю. и др. Численное исследование профилей параметров неравновесного двухфазного потока на основе двумерной модели II Известия ВУЗов, Ядерная энергетика. - 1994. - №4-5. -С.80-88; №6. - С.56-62.

7. Гинкин В.П., Ганина С.М. Новый метод расчета трехмерной конвекции на сетках большой размерности // Математическое моделирование. - 2003. -Т.15. - №6. - С.53-58.

8. Гинкин В.П., Картавых A.B., Мильвидский М.Г. и др. Численное моделирование процесса тепломассопереноса с позиции кластерной модели строения расплава // Поверхность. - 2004. - №6. - С.93-100.

Труды конференций

9. Гинкин В.П., Гадияк Г.В., Шварц Н.Л. и др. Моделирование МДП-транзисторов конечно-разностными методами // 4-я Международная конференция по численному моделированию полупроводниковых приборов и интегральных схем. - Дублин, Ирландия, 1985. - С.25-35.

Ю.Гинкин В.П., Жиганова И.Г., Гадияк Г.В. Исследование сходимости метода параболических прогонок в областях непрямоугольной формы //

Международная конференция по математическому моделированию приборов микроэлектроники. - Дублин, Ирландия, 1987. - Юс.

П.Гинкин В.П., ФирсоваЭ.В., ЕфановА.Д. и др. Влияние неравномерности распределения входной температуры на эффективность жидко-металлического теплообменника // VIII Всесоюзная конференция. - Минск, 1988.-9с.

12.Гинкин В.П., ХудаскоВ.В., ЗининаГ.А. и др. Методика трехмерного теплогидравлического расчета парогенераторов АЭС // VIII Всесоюзная конференция. - Минск, 1988. - Юс.

П.Гинкин В.П., ТрояноваН.М. Вариант метода неполной факторизации для решения трехмерных девятиточечных уравнений эллиптического типа // 13-й Международный конгресс по вычислительной и прикладной математике. - Дублин, Ирландия, 1991. - 6с.

И.Гинкин В.П., Гадияк Г.В., Шиманский A.A. и др. МОПИТ: открытая система сквозного моделирования технологии и характеристик полупроводниковых приборов // 7-я международная конференция по численному анализу полупроводниковых приборов и интегральных схем NASE CODE VII. -Колорадо, США, 1991. - 9с.

15.Гинкин В.П., Артемьев В.К., Бердников B.C. Численное исследование конвективного теплообмена на модели метода Чохрапьского // 1-я Российская национальная конференция по теплообмену (РНКТ1). - Москва, 1994. - Т.2. - С.26-30.

ló.Ginkin V.P., Artemyev V.K., Gusev N.V. et al. Numerical Simulation of GaAs Crystal Growth of a Big Diameter by Bridgman Method (Численное моделирование выращивания кристалла GaAs большого диаметра методом Бриджмена) // The 5-th Int. Conf. on Simulation of Devices and Technologies (ICSDT'96), May 13-17 1996. - Obninsk, 1996. - P.126-127.

17.Ginkin V.P. A Numerically-Analytical Method to Solve 3D Diffusion Multigroup Equation of Reactor for Domains with Step-Wise Coefficients (Численно-аналитический метод для решения многогруппового диффузионного уравнения реактора для областей с кусочно-постоянными коэффициентами) // The 5-th Int. Conf. on Simulation of Devices and Technologies (ICSDT'96), May 13-17 1996. - Obninsk, 1996. - P.159-160.

18. Гинкин В.П., Артемьев B.K. Численное моделирование трехмерной естественной конвекции // 2-я Российская национальная конференция по теплообмену (РНКТ2). 26-30 октября 1998. - Москва, 1998. - Т.З. - С.38-41.

19.Гинкин В.П., Артемьев В.К., Свириденко И.П. и др. Экспериментальные и расчетные исследования температурных режимов в макете ростовой установки при выращивании монокристаллов методом плавающей зоны // 2-й Российский симпозиум «Процессы тепломассопереноса и рост

монокристаллов и тонкопленочных структур» (НТ&СС97), 22-24 сентября 1997. - Обнинск, 1998. - С.17-29.

20.Гинкин В.П., Галкин В.А., Драков А.А. и др. Моделирование и оптимизация теплофизических процессов роста монокристаллов в условиях микрогравитации // 2-й Российский симпозиум «Процессы тепломассопереноса и рост монокристаллов и тонкопленочных структур НТ&СС97», 22-24 сентября 1997. - Обнинск, 1998. - С.140-147.

21.Гинкин В.П. Проблемы численного моделирования процесса роста монокристаллов // 2-й Российский симпозиум «Процессы тепломассопереноса и рост монокристаллов и тонкопленочных структур» (HT&CG'97), 22-24 сентября 1997. - Обнинск, 1998. - С.148-157.

22.Ginkin V.P. Problems of Numerical Simulation of the Single Crystal Growth Process (Проблемы численного моделирования процесса роста монокристаллов) // The 5-th International Conference on Simulation of Devices and Technologies, October 11-20 1998. - Cape-Tow, South Africa, 1998. -P.228-232.

23.Ginkin V.P. Problems of numerical simulation in crystal growth process (Проблемы численного моделирования в процессе роста кристаллов) // The International Conference on Computational Heat and Mass Transfer, April 26-29 1999. - Famagusta, Northern Cyprus, 1999. - P.32-37.

24.Ginkin V.P., Artemyev V.K., FolomeevV.I. et. al. Numerical benchmark for calculation of convective heat transfer during crystal growth by the floating zone method under the null gravity conditions (Численный бенчмарк для расчета конвективного теплопереноса в процессе выращивания кристалла методом плавающей зоны в условиях полного отсутствия гравитации) // The International Conference on Computational Heat and Mass Transfer, April 26-29 1999. - Famagusta, Northern Cyprus, 1999. - P. 195-200.

25.Ginkin V.P., Artemyev V.K., Prostomolotov A.l. et. al. The numerical and experimental study of the Bridgman crystal growth on a model device (Численные и экспериментальные исследования процесса выращивания кристаллов методом Бриджмена на модельной установке) // The International Conference on Computational Heat and Mass Transfer, April 26-29 1999. -Famagusta, Northern Cyprus, 1999. - P.201 -206.

26.Ginkin V.P., GorevN.A., Senchenkov A.S. et. al. Numerical and Experimental Study of Heat Transfer Processes in Modular Furnace "Polizon" (Численное и экспериментальное исследование процесса теплопереноса в модульной установке ПОЛИЗОН) // The 3-d International Conference "Single Crystal Growth, Strength Problems, and Heat Mass Transfer" (ICSC-99), September 21-24 1999. - Obninsk, 2000. - Vol.1. - P.15-25.

27.Ginkin V.P., Milvidsky M.G., RakovV.V. et. al. Numerical Simulation of Anomalous Doping Effect in Ge Single Crystal Grown by FZ-Technique Aboard the Space Crafts (Численное моделирование эффекта аномального

распределения примеси при выращивании кристаллов Ge методом плавающей зоны на борту космического корабля) // The 3-d International Conference "Single Crystal Growth, Strength Problems, and Heat Mass Transfer" (1CSC-99), September 21-24 1999. - Obninsk, 2000. - Vol.1. - P.26-36.

28.Ginkin V.P., Sviridenko I.P., Shishulin V.P. et. al. Numerical and Experimental Study of High Temperature Heat Zones in Growth Devices Based on the Bridgman and Zones Methods (Численные и экспериментальные исследования вцсоко температурных зон в ростовых установках методами Бриджмена и зонной плавки) // The 3-d International Conference "Single Crystal Growth, Strength Problems, and Heat Mass Transfer" (ICSC-99), September 21-24 1999. - Obninsk, 2000. - Vol.2. - P.477-491.

29.Ginkin V.P., Kulik A.V., Naumenko O.M. An efficient preconditioning procedure in the conjugate gradient method for 3D (hexj) geometry (Эффективный предобуславливатель в методе сопряженных градиентов для трехмерной (hex.z) геометрии) // The 4-th International Conference on Supercomputing in Nuclear Applications (SNA 2000), September 4-7 2000. - Tokyo, Japan, 2000. -P.PS2-G3-15.

30.Ginkin V.P. Methods of solution of convective heat mass transfer at single crystal growth problem (Методы решения задачи конвективного тепломассопереноса в проблеме роста монокристаллов) // The 2-nd International Symposium on Advances in Computational Heat Transfer (CHT'01), May 20-25 2001. - Palm Cove, Queensland, Australia, 2001. - Vol.2. - P.l 161-1168.

31.Ginkin V.P. Algorithm of solution of 3-D magnetic hydrodynamic equations for crystal growth problem (Алгоритм решения трехмерных уравнений магнитной гидродинамики для проблемы роста кристаллов) // The 4-th International Conference "Single Crystal Growth and Heat & Mass Transfer" (ICSC-2001), September 24-28 2001. - Obninsk, 2001. - Vol.3. - P.792-807.

32.Ginkin V.P., Makeev Kh.I., Shotaev A.N. et. al. Simulation of convective transfer of oxygen at growth of crystals of silicon by Czochralski method (Моделирование конвективного переноса кислорода при выращивании кристаллов кремния методом Чохральского) // The 4-th International Conference "Single Crystal Growth and Heat & Mass Transfer" (ICSC-2001), September 24-28 2001. - Obninsk, 2001. - Vol.4. - P.913-920.

33.Ginkin V.P., Maronchuk 1.1., Kartavykh A.V. et. al. Methodical research of thermal fields during crystallization in the "POLYZONE" facility by Bridgman method (Методические исследования тепловых полей в процессе выращивания кристаллов методом Бриджмена на установке ПОЛИЗОН) // The 4-th International Conference "Single Crystal Growth and Heat & Mass Transfer" (ICSC-2001), September 24-28 2001. - Obninsk, 2001. - Vol.4. -P.1003-1008. „1/!(,

34.Ginkin V.P. Porous solid model to describe heat-mass transfer near phase transition interface in crystal growth from melt simulations (Модель пористого

тела для описания тепломассопереноса вблизи границы фазового перехода при моделировании роста кристаллов из расплава) // The 1-st International Conference on Applications of Porous Media. June 2-8 2002. - Jerba, Tunisia, 2002. - P.I55-160.

35.Ginkin V.P., Zabud'ko M.A., Kartavykch A.V. et. al. Numerical modeling of heat-mass transfer process from the point of view of cluster model of a melt constitution (Численное моделирование процесса тепломассопереноса с точки зрения кластерной модели структуры расплава) // The 5-th International Conference "Single Crystal Growth and Heat & Mass Transfer" (ICSC-2003), September 22-26 2003. - Obninsk, 2003. - Vol.1. - P.31-42.

36.Ginkin V.P., Artemyev V.K., Makeev Kh.I. et. al. Hydrodynamic and growing parameters influence on oxygen distribution at Czochralski silicon crystal growth process (Влияние гидродинамики и ростовых параметров на распределение кислорода в процессе выращивания кристаллов кремния методом Чохральского) // The 5-th International Conference "Single Crystal Growth and Heat & Mass Transfer" (ICSC-2003), September 22-26 2003. - Obninsk, 2003. -Vol.2. - P.522-536.

37.Ginkin V.P. Efficient method to solve 3D-elliptical equations with symmetrical and nonsymmetrical coefficient matrices (Эффективный метод для решения трехмерных уравнений эллиптического типа с симметричными и несимметричными матрицами коэффициентов) // The 5-th International Conference "Single Crystal Growth and Heat & Mass Transfer" (ICSC-2003), September 22-26 2003. - Obninsk, 2003. - Vol.2. - P.631-636.

Прочие

Зв.Гинкин В.П. Метод А-факторизации для решения двумерных уравнений эллиптического типа // IV Всесоюзное совещание по вычислительным методам линейной алгебры / Вычислительные методы линейной алгебры. Новосибирск: ВЦ СОАН СССР. - 1977. - С.123-132.

39.Гинкин В.П., Корниенко Ю.Н., Кузеванов B.C. Расчетное определение локальных характеристик двухфазного потока на участке установления профиля газосодержания // Физическая механика неоднородных сред. Новосибирск: ИТПМ СОАН СССР. - 1984. - С.54-59.

40.Ginkin V.P., Artemyev V.K., GusevN.V. et. al. Numerical Simulation of the Process of rystal Growth from the Melt (Численное моделирование процесса роста кристалла из расплава) // The 4-th Int. Seminar on Simulation of Devices and Technologies (ISSDT'95), November 15-17 1995. - Berg-en-Dal, Kruger National Park, South Africa, 1995. - P.228-232.

41.Гинкин В.П., ГанинаС.М., Кузнецов И.А. и др. Расчетные исследования нестационарных процессов в быстром реакторе, обусловленных несанкционированным извлечением из активной зоны компенсирующих

стержней. // Отраслевой семинар «Нейтроника-97», 28-30 октября 1997. -Обнинск, 1997. -24с.

42.Гинкин В.П., Кулик А.В. Метод сопряженных градиентов в сочетании с неполной факторизацией и компенсацией итерируемых членов для трехмерной (hex,z) геометрии // Отраслевой семинар «Нейтроника-99», 26-28 октября, 1999. - Обнинск, 2000. - С. 160-172.

43.Ginkin V.P., Ganina S.M. Method and Code for Calculation of Three-Dimensional Convection on Grids of Big Dimension (Метод и программа расчета трехмерных уравнений конвекции на сетках большой размерности) // The 10-th Int. Meeting of the International Association for Hydraulic Engineering and Research (IAHR) "Thermal Hydraulics for Fast Reactors with Different Coolants", July 17-19,2001.-Obninsk, 2003.-P.167-183.

44.Ginkin V.P., Artemyev V.K., Kartavykh A.V. et. al. The mechanism of Marangoni convection influence on dopant distribution in Ge space-grown single crystals (Механизм влияния конвекции Марайгони на распределение примеси в монокристаллах германия, выращенных в космосе) // J. Crystal Growth. - 2001. - No.223. - Р.29-37.

45.Ginkin V.P., FolomeevVJ., Maronchukl.I. et. al. Method research and developing a method to control directed semiconductor crystallization in space (Методические исследования и развитие метода контроля направленной кристаллизации полупроводников в космосе) // J. Crystal Growth. - 2002. -No.236.-P.551-556.

Подписано к печати.29.06.2004 г. Формат60х84 1/16. Усл.п.л. 1,3. Уч.-иэд.л. 1,9. Тираж 100 экз. Заказ № Отпечатано в ОНТИ методом прямого репродуцирования с оригинала авторов. 24903 3, Обнинск Калужской обл., ФЭИ.

»18444

РНБ Русский фонд

2005-4 13686

4'

/

V

Оглавление автор диссертации — доктора физико-математических наук Гинкин, Владимир Павлович

Введение.

Глава 1. Метод неполной факторизации.

1.1 Основные теоремы.

1.2 Явные схемы.

1.3 Схемы для симметричных матриц.

1.4 Доказательство сходимости явной схемы Булеева для случая симметричных матриц.

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

1.6 Неявные схемы.

1.7 Неявные схемы для пятидиагональных матриц.

Вторая схема неполной факторизации Булеева.

Схема А-факторизации (НФ).

1.8 Схема параболических прогонок (ПП).

Схема параболических прогонок (ПП) с tj, = 0.

Схема параболических прогонок (ПП) с 0 = 0.

С Асимптотические оценки для элементов матриц в схеме ПП.

Спектральный анализ сходимости схемы ПП.

Комбинированная схема НФПП.

1.9 Неявные схемы для семидиагональных матриц.

1.10 Численные примеры.

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

Вторая схема Булеева.

Третья схема Булеева.

Четвертая схема.

Схема А-факторизации.

Схема параболических прогонок.

Комбинированные схемы.

Сравнение различных схем.

Неоднородная задача Дирихле.

Ф Однородная задача Неймана-Дирихле.

Несимметричные задачи.

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

Численные примеры.

Выводы к главе 1.

Глава 2. Алгоритмы решения стационарных и квазистационарных нейтронно-физических задач.

2.1 Трехмерный стационарный комплекс программ WIMS.ВОЛГА.

2.2 Трехмерный квазистационарный комплекс программ WIMS-BOJ1HA.

2.3 Алгоритм пространственно-временного расчета реактора.

2.4 Исследование погрешности аппроксимации при расчете активной зоны реактора ВВЭР-1000.

2.5 Расчет реальной зоны реактора ВВЭР-1000.

2.6 Расчеты кампаний реакторов типа ВВЭР.

Блок номер 3 Запорожской АЭС в 10 топливной кампании.

Блок номер 1 Балаковской АЭС в 8 топливной кампании.

Блок номер 2 Калининской АЭС в 10 топливной кампании.

2.7 Схема расчета активной зоны реактора с измененной формой ТВС.

Дискретизация области.

Методика учета неравномерного водяного зазора между ТВС.

Аппроксимация уравнений.

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

ВВЭР-1000 со смещенными ТВС.

2.8 Организация вычислений в комплексе программ WIMS-BOJ1HA.

Выводы к главе 2.

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

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

3.2 Условно критическая и сопряженная задачи.

Конечно-разностный аналог условно-критической и сопряженной задач.

Уточнение аппроксимации.

3.3 Нестационарная задача.

Решение уравнений кинетики.

Форм-функция. Уравнение, аппроксимация.

Алгоритм решения уравнения для форм-функции.

Константное обеспечение.

3.4 Краткое описание комплекса VOLNA, блок-схема.

3.5 Комплекс программ GVA.

3.6 Описание моделируемой аварии с остановкой ГЦН.

3.7 Расчет переходного процесса.

3.8 Исследование точности аппроксимации нестационарных уравнений реактора по времени.

3.9 Учет расширения активной зоны реактора в (hex.z) геометрии.

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

3.11 Аварии, обусловленных несанкционированным извлечением стержней.

3.12 Описание моделируемой аварии с движением стержней.

3.13 Результаты расчетов аварий, связанных с движением стержней.

3.14 Пространственно-временное деформирование активной зоны быстрого реактора.

3.15 Трехмерная нестационарная модель деформаций ТВС в активной зоне реактора.

3.16 Расчет деформирования активной зоны реактора типа БНв результате аварии.

Выводы к главе 3.

Глава 4. Алгоритм решения трехмерных уравнений конвективного тепломассопереноса и магнитной гидродинамики.

4.1 Нестационарные уравнения магнитной гидродинамики в приближении Буссинеска.

4.2 Преобразование уравнений.

4.3 Разностная аппроксимация.

4.4 Уравнение для давления.

4.5 Алгоритм решения.

4.6 Случай переменной вязкости и аксиальной симметрии.

4.7 Численные примеры.

Выводы к главе 4.

Глава 5. Математическое моделирование процессов кристаллизации расплавов в условиях микрогравитации.

5.1 Внешняя задача.

5.2 Задача теплопроводности.

5.3 Уравнения конвекции.

5.4 Уравнения радиационного теплообмена.

5.5 Методика решения внешней задачи.

5.6 Реализация алгоритма.

5.7 Внутренняя задача.

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

5.9 Аномальные эффекты распределения примеси в монокристаллах Ge, выращенных в космосе.

5.10 Кластерная модель процесса кристаллизации.

Вычислительный алгоритм.

Тестирование модели.

Выводы к главе 5.

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

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

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

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

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

Из множества разработанных автором приложений, в которых используются различные варианты метода неполной факторизации, в диссертацию, из-за ограничения объема, включены лишь некоторые. Это, прежде всего, комплекс программ WIMS.BOJirA для расчета ректоров в (х, у, z) и (г, в, z) геометриях, и комплекс программ WIMS-BOJ1HA для расчета кампаний реакторов на тепловых нейтронах типа ВВЭР в (hex, z) геометрии с прямым расчетом констант по прецизионной программе WIMSD4 в каждой пространственной точке реактора в каждый момент времени. С помощью этих комплексов программ было выполнено большое количество верификационных расчетов, и выполнены расчетные исследования топливных кампаний различных блоков реакторов на тепловых нейтронах применительно к проблеме деформаций тепловыделяющих сборок (ТВ С) в активных зонах реакторов типа ВВЭР-1000 в процессе эксплуатации.

Другой работой, включенной в диссертацию, является комплекс программ VOLNA для расчета быстрых переходных процессов в реакторах в квазистатическом приближении в (hex, z) геометрии, и созданный на его основе комплекс программ GVA (GRIF-SM, VOLNA, ARAMAKO) для расчета аварий в реакторах на быстрых нейтронах. Отличительной чертой комплекса GVA является совместный нейтронно-физический, теплогидравлический и термомеханический расчет реакторов на быстрых нейтронах. С помощью этого комплекса были просчитаны аварии в реакторе типа БН-800, вызванные остановкой главных циркуляционных насосов и несанкционированным движением группы компенсирующих стержней органов регулирования.

Из большого числа приложений в области тепловых расчетов автор включил в диссертацию лишь разработанный им новый метод решения уравнений магнитной гидродинамики, так как в нем отражены, с методической стороны, все основные особенности такого рода алгоритмов. На основе этого метода была написана программа GIGAN для расчета тепловой конвекции в трехмерной (х, у, z) геометрии и рассчитаны международные бенчмарки о конвекции воздуха в кубической полости, подогреваемой с одной из боковых сторон. Использование в программе GIGAN процедуры симметризации разностных уравнений и метода сопряженных градиентов с предобуславливателем по схеме PIF для решения линеаризованных систем с симметричными матрицами коэффициентов позволило решать задачи на сетках большой размерности (до миллиона узлов) на обычных PC.

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

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

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

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

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

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

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

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

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

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

Научная новизна результатов диссертации.

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

Ф блочных итерационных методов, сделан анализ сходимости схемы А-факторизации, схемы параболических прогонок и комбинированных схем неполной факторизации.

• Написаны специализированные комплексы программ: WIMS.BOJirA для стационарного расчета реакторов на тепловых нейтронах и WIMS-BOJ1HA для расчета кампаний реакторов типа ВВЭР, в которых впервые организован прямой расчет констант и выгорания по прецизионной программе WIMSD4 в каждом пространственном узле расчетной сетки в каждый момент времени. В рамках комплекса программ WIMS-BOJIHA разработан алгоритм учета влияния формоизменения ТВС на нейтронно-физические характеристики реактора.

• Разработан комплекс программ VOLNA для расчета быстрых переходных процессов в реакторе в квазистатическом приближении в (hex, z) геометрии. На основе этого комплекса, а также комплексов программ теплогидравлического расчета GRIF-SM и расчета констант АРАМАКО, создан комплекс программ GVA для совместного трехмерного нейтронно-физического, теплогидравлического и термомеханического расчета аварий реакторов на быстрых нейтронах. С помощью этого комплекса программ выполнены трехмерные f нестационарные расчеты аварий реактора типа БН-800, вызванных остановкой главных циркуляционных насосов (ГЦН) и несанкционированным движением группы компенсирующих стержней органов регулирования. Разработаны метод и программа расчета деформаций ТВС и алгоритм учета теплового расширения активной зоны реактора на его нейтронно-физические характеристики.

• Разработан новый метод решения уравнений магнитной гидродинамики в естественных переменных. На основе этого метода разработаны алгоритм и программа GIGAN для трехмерного нестационарного расчета тепловой конвекции и выполнены расчеты бенчмарков по конвекции воздуха в кубической полости, подогреваемой с одной из боковых сторон, на сетках большой размерности (до миллиона узлов). Показано, что для чисел Релея 106 и выше необходимо использовать еще большее количество узлов или применять специальное сгущение сеток для обеспечения приемлемой точности расчетных значений чисел Нуссельта.

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

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

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

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

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

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

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

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

Автор выносит на защиту:

• Варианты метода неполной факторизации: метод А-факторизации, метод параболических прогонок, вариант метода сопряженных градиентов с предобуславливателем по схеме неполной факторизации с периферийной компенсацией итерируемых членов CGPIF, комбинированные схемы неполной факторизации HFPP, HFB3, HFB4, CIF. Анализ сходимости блочных итерационных методов неполной факторизации, доказательство сходимости метода /г-факторизации, метода параболических прогонок, явной схемы Булеева с диагональной компенсацией итерируемых членов.

• Алгоритмы и программы стационарного и квазистационарного трехмерного расчета реакторов ВОЛГА и ВОЛНА. Комплекс программ расчета реакторов WIMS.ВОЛГА для расчета реакторов в (х, у, т) и (г, в, z) геометриях. Комплекс программ WIMS-ВОЛНА расчета кампании реакторов типа ВВЭР в (hex, z) геометрии с учетом влияния формоизменения ТВС на нейтронно-физические характеристики реактора вследствие термомеханических деформаций.

• Алгоритм и комплекс программ VOLNA для расчета быстрых переходных процессов в реакторе в квазистатическом приближении в (hex, z) геометрии. Алгоритм и комплекс программ GVA для совместного нейтронно-физического, теплогидравлического и термомеханического расчета аварий реакторов на быстрых нейтронах типа БН-800 с учетом формоизменения ТВС вследствие термомеханических деформаций. Результаты расчетов аварий реакторов типа БН-800, вызванных остановкой главных циркуляционных насосов и несанкционированным движением группы компенсирующих стержней органов регулирования.

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

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

Апробация работы. Основные положения и результаты диссертационной работы докладывались на IV Всесоюзном совещании по вычислительным методам линейной алгебры (Новосибирск, 1977), на встречах со специалистами ГДР по проблеме транспортирования отработавшего топлива (Берлин, 1980, 1982, 1984, 1986), на совещании стран-членов СЭВ (Ленинград, 1981), на VII, IX, X, XI школах по моделям механики сплошной среды (Батуми, 1983; Якутск, 1987; Казань, 1989; Владивосток, 1991), на 1 отраслевом совещании по методам и программам расчета (Обнинск, 1983), на межотраслевом семинаре «Сопряженные задачи теплообмена» (Москва, 1985), на IV международной конференции по численному моделированию полупроводниковых приборов и интегральных схем (Ирландия, 1985), на всесоюзном совещании по математическому моделированию приборов микроэлектроники (Новосибирск, 1987), на VIII Всесоюзной конференции по тепло- и массообмену, (Минск, 1988), на Всесоюзном совещании по математическому моделированию приборов микроэлектроники (Горно-Алтайск, 1989), на международном семинаре по математическому моделированию в микроэлектронике (Новосибирск, 1990), на Всесоюзном совещании по математическому моделированию в микроэлектронике (Ярославль, 1990), на Всесоюзном совещании по динамике реакторов (Гатчина, 1990), на международном семинаре «Теплофизика-90» (Обнинск, 1990), на 13 Международном конгрессе по вычислительной и прикладной математике IMACS'91 (Ирландия, 1991), на 7 международной конференции по численному анализу полупроводниковых приборов и интегральных схем NASE CODE VII (США, 1991), на XV Менделеевском съезде по общей и прикладной химии (Обнинск, 1993), на конференции по миграции радионуклидов в водных системах (Обнинск, 1993), на VI Российской конференции по радиационной защите ядерных установок (Обнинск, 1994), на 3, 4, 5 международных семинарах по моделированию устройств и технологий (Обнинск, 1994; Южная Африка, 1996, 1998), на 1, 2, 3, Российских национальных конференциях по теплообмену РНКТ (Москва, 1994, 1998, 2002), на 3 международном семинаре по математическому моделированию в микроэлектронике (Обнинск, 1994), на 9 международном совещании по безопасности ядерных реакторов (Москва, 1995), на III Минском международном форуме по тепло- и массообмену (Минск, 1996), на 5 международной конференции по моделированию устройств и технологий ICSDT'96 (Обнинск, 1996), на 2 международном симпозиуме ученых и исследователей России и США, выполняющих исследования по программе «НАУКА-НАСА» (Королев, 1996), на 1 и 2 международных симпозиумах по достижениям в вычислительном тепломассопереносе (Турция, 1997; Австралия 2001), на 2 Российском симпозиуме «Процессы тепломассопереноса и рост монокристаллов и тонкопленочных структур» HT&CG'97 (Обнинск, 1997), на международной алгебраической конференции памяти Д.К.Фадеева (С.-Петербург, 1997), на школе-семинаре МИФИ секции «Интегрированные математические модели и программные комплексы в ядерной энергетике» (Москва, 1998), на 5 Всероссийской научно-технической конференции с международным участием «Актуальные проблемы твердотельной электроники и микроэлектроники» (Дивноморское, Таганрог, 1998), на семинаре «Консультативная встреча специалистов IAEA по проблеме быстрых реакторов» (Обнинск, 1998), на семинарах «Нейтроника-98», «Нейтроника-99» (Обнинск, 1998, 1998), на XXXIII научных чтениях, посвященных разработке творческого наследия К.Э.Циолковского (Калуга, 1998), на 3,4, 5 международных конференциях «Рост монокристаллов и тепломассоперенос» (Обнинск, 1999, 2001, 2003), на международной конференции по вычислительному тепломассопереносу (Кипр, 1999), на рабочей группе по фазовым переходам с конвекцией, моделированию и валидации (Польша, 1999), на 1 и 2 Российских конференциях по космическому материаловедению (Калуга, 1999, 2003), на IX и X Национальных конференциях по росту кристаллов НКРК (Москва, 2000, 2002), на 4 Международной конференции по супервычислениям в ядерных приложениях SNA-2000 (Япония, 2000), на 4 Сибирском конгрессе по прикладной и индустриальной математике ИНПРИМ-2000 (Новосибирск, 2000), на международной конференции по твердым кристаллам (Польша, 2000), на конференции, посвященной памяти академика А.Н.Тихонова (Обнинск, 2000), на 13 международной конференции по росту кристаллов ICCG-13 (Япония, 2001), на XVII научном совещании «Высокочистые материалы с особыми физическими свойствами» (Суздаль, 2001), на международной конференции «Математические идеи П.Л.Чебышева и их приложение к современным проблемам естествознания» (Обнинск, 2002), на I Всероссийской конференции «Физико-химические процессы в конденсированном состоянии и на межфазных границах» ФАГРАН-2002 (Воронеж, 2002), на IV Международной конференции по неравновесным процессам в соплах и струях NPNJ-2002 (С.-Петербург, 2002), на 1 международной конференции по приложениям пористых сред (Тунис, 2002), на 1 научной школе-конференции «Актуальные вопросы теплофизики и физической гидрогазодинамики» (Украина, 2003), на международной конференции по проблемам в тепловой конвекции (Пермь, 2003).

Основное содержание диссертации изложено в 8 статьях в реферируемых журналах, 29 докладах в сборниках трудов Международных и Российских конференций, 2 статьях в зарубежных журналах, 20 препринтах.

Структура и объем работы. Диссертация состоит из введения, 5 глав и заключения (общие выводы по работе). Общий объем диссертации — 258 страниц, в том числе 123 рисунка и 46 таблиц, список литературы из 240 наименований на 25 страницах.

Заключение диссертация на тему "Алгоритмы и комплексы программ для решения задач математической физики с использованием метода неполной факторизации"

Основные результаты и выводы

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

Приведено описание различных вариантов метода неполной факторизации на сеточном и матричном языках. Дано доказательство сходимости явной схемы Булеева с диагональной компенсацией итерируемых членов и получены оценки оптимального значения итерационного параметра и скорости сходимости этой схемы. Развита теория блочных итерационных методов для анализа сходимости неявных схем метода неполной факторизации. Даны доказательства сходимости и оценки скоростей сходимости схем /z-факторизации и параболических прогонок.

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

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

Построена и численно исследована комбинированная схема неполной факторизации CIF для решения трехмерных конечно-разностных задач с симметричными и несимметричными матрицами коэффициентов. Показано, что для симметричных матриц по этой схеме зависимость числа итераций j от числа узлов п по одному направлению имеет вид j~0(4n), а для несимметричных матриц она существенно возрастает с ростом абсолютных значений конвективных членов в исходных уравнениях.

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

Разработан комплекс программ WIMS.ВОЛГА для трехмерного расчета реакторов в стационарном многогрупповом диффузионном приближении в (х, у, z) и (г, в, z) геометриях. Он прошел необходимую верификацию, включен в отраслевой фонд алгоритмов и программ расчета ядерных реакторов ОФАП ЯР, и широко использовался для исследований критичности гомогенных и гетерогенных размножающих систем.

Разработан комплекс программ WIMS-BOJIHA для расчета кампаний реакторов типа ВВЭР в квазистационарном многогрупповом диффузионном приближении в (hex, z) геометрии. Отличительной чертой комплекса является использование для расчета выгорания и констант в каждой пространственной точке реактора в каждый момент времени в процессе его эксплуатации с помощью прецизионной программы WIMSD4. При этом наиболее корректно учитываются реальные значения всех технологических параметров реактора. Разработан алгоритм учета влияния неравномерных водяных зазоров между ТВС на нейтронно-физические характеристики реактора. С помощью этого комплекса выполнены расчеты кампаний реакторов ВВЭР-1000 различных блоков Калининской, Запорожской и Балаковской АЭС при различных загрузках топлива и схемах его перегрузок с участием усовершенствованных ТВС новой конструкции, в обосновании работоспособности которых использовался комплекс программ WIMS-BOJIHA.

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

На основе программ нейтронно-физического расчета VOLNA, теплогидравлического расчета GRIF-SM, и расчета констант АРАМАКО создан трехмерный комплекс Программ расчета аварий реакторов на быстрых нейтронах GVA. В рамках этого комплекса разработаны алгоритм и программа расчета термомеханического поведения активной зоны реактора и алгоритм учета влияния теплового расширения активной зоны на нейтронно-физические характеристики реактора. Разработанный комплекс программ GVA позволяет моделировать сложные аварийные процессы в реакторах на быстрых нейтронах, приводящих к значительным изменениям концентраций компонентов реактора в переходном процессе, а именно таких, как кипение натрия или перемещение расплавов стали и топлива после предполагаемого расплавления активной зоны. С помощью этого комплекса выполнены тестовые расчеты аварий реакторов типа БН-800, вызванных остановкой главных циркуляционных насосов с одновременным отказом органов регулирования, и извлечением группы шести компенсирующих стержней, которые продемонстрировали работоспособность и эффективность комплекса. Они показали, в частности, что локальные возмущения энерговыделения в активной зоне могут в несколько раз превышать интегральное изменение мощности реактора, что, в свою очередь, может привести к расплавлению ТВС в местах этих возмущений.

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

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

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

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

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

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

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

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

Библиография Гинкин, Владимир Павлович, диссертация по теме Математическое моделирование, численные методы и комплексы программ

1. БулеевН.И. Численный метод решения двумерных уравнений диффузии 1. Атомная энергия. - 1959. - Т.б. - №3. - С.338-340.

2. OliphantT.A. An implicit numerical method for solving two-dimensional time-dependent diffusion problems = Неявный численный метод для решения двумерных, зависящих от времени задач диффузии) // Quart. Appl. Math. -1961. Vol.19. - P.221-229.

3. Булеев Н.И. Метод неполной факторизации для решения двумерных и трехмерных разностных уравнений типа диффузии II Вычислительная математика и математическая физика. 1970. - Т. 10. - №4. - С. 1042-1044.

4. Ильин В.П. Методы неполной факторизации для решения алгебраических систем. М.: Наука, 1995.-288с.

5. Гинкин В.П. Метод А-факторизации для решения двумерных уравнений эллиптического типа // Доклад на IV Всесоюзном совещании по вычислительным методам линейной алгебры / Вычислительные методы линейной алгебры. Новосибирск: ВЦ СОАН СССР, 1977. - С.123-132.

6. Гинкин В.П. О численном решении двумерных и трехмерных уравнений эллиптического типа. Препринт ФЭИ-767. - 1977. - 12с.

7. Гинкин В.П. О влиянии релаксации на сходимость схемы А-факторизации при решении двумерных уравнений типа диффузии // Физика и техника ядерных реакторов. ВАНТ, 1980. - Вып.4 (13). - С. 196-202.

8. Гинкин В.П. Метод параболических прогонок для решения двумерных уравнений эллиптического типа. Препринт ФЭИ-1153. - 1981. - 13с.

9. Гинкин В.П. Об одном варианте метода параболических прогонок для решения двумерных уравнений эллиптического типа. Препринт ФЭИ-1306. - 1982. - 14с.

10. Гинкин В.П. Решение эллиптических уравнений со смешанными производными методом параболических прогонок. Препринт ФЭИ-1365. - 1982. - 14с.

11. Гинкин В.П. Решение двумерных эллиптических уравнений методом параболических прогонок // Доклад на межотраслевом семинаре «Сопряженные задачи теплообмена» / ФЭИ, 1985. -Зс.

12. Гинкин В.П., Кончиц А.П. Элементы теории блочных итерационных методов. Препринт ФЭИ-1738. - 1985. -22с.

13. Гинкин В.П., Жиганова И.Г., ГадиякГ.В. Исследование сходимости метода параболических прогонок в областях непрямоугольной формы и гексагональной геометрии. Препринт ФЭИ-1841. - 1987. - Юс.

14. Гинкин В.П. Эффективный метод решения одного класса задач эллиптического типа при отсутствии диагонального преобладания // Доклад на всесоюзном совещании по математическому моделированию приборов микроэлектроники / Новосибирск, 1987. 7с.

15. Гинкин В.П. Эффективный метод решения одного класса задач эллиптического типа при отсутствии диагонального преобладания // Автометрия. 1988. - №5. - С.64-67.

16. Гинкин В.П. Метод неполной факторизации для решения трехмерных уравнений эллиптического типа // Доклад на международном семинаре по математическому моделированию в микроэлектронике / Новосибирск, 1990. 6с.

17. Гинкин В.П., Степанов И.Н. Распараллеливание трехточечной прогонки // Доклад на Всесоюзном совещании по математическому моделированию в микроэлектронике / Ярославль, 1990. 8с.

18. Алгоритмы решения сеточных уравнений для многопроцессорных ЭВМ : отчет о НИР / ФЭИ. Обнинск, 1991. - 19с. - Исполн.: Гинкин В.П., Ледовской В.М., Степанов И.А. -Инв. № 7869.

19. Гинкин В.П., Троянова Н.М. Использование метода неполной факторизации в трехмерной задаче нейтронно-физического расчета реакторов типа ВВЭР. Препринт ФЭИ-2104.- 1990.- 15с.

20. Гинкин В.П., Троянова Н.М. Метод неполной факторизации для решения трехмерных девятиточечных разностных уравнений эллиптического типа. Препринт ФЭИ-2181. -1991.- 17с.

21. Гинкин В.П., Троянова Н.М. Вариант метода неполной факторизации для решения трехмерных девятиточечных уравнений эллиптического типа // Доклад на 13-ом Международном конгрессе по вычислительной и прикладной математике / Дублин, Ирландия, 1991. 6с.

22. Гинкин В.П., Гущин Е.В. Блочный метод неполной факторизации для решения систем двумерных уравнений эллиптического типа. Препринт ФЭИ-2268. - 1992. - 8с.

23. Гинкин В.П., Ваньков К.А., Троянова Н.М. Метод сопряженных градиентов в сочетании с неполной факторизацией и компенсацией итерируемых членов. Препринт ФЭИ-2324. -1993.- 13с.

24. Гинкин В.П., Кулик А.В. Оптимальный предобуславливатель в методе сопряженных градиентов для трехмерной (hex.z) геометрии. Препринт ФЭИ-2786. - 1999. - 20с.

25. Transfer" (ICSC-2003), September 22-26 2003 / Proceedings. Obninsk: SSC RFIPPE, 2003. -Vol.1.-P.631-636.

26. Маркус M., Минк X. Обзор по теории матриц и матричных неравенств. М.: Наука, 1972. - 232с.

27. Schneider G.E., Zedan M. A modified strongly implicit procedure for the numerical solution of field problems = Модифицированная строго неявная процедура для численного решения пространственных задач // Numerical Heat Transfer. -1981. Vol.4. - P.1-19.

28. Дьяконов Е.Г. Разностные методы решения краевых задач. Издание МГУ, 1971. - Вып.1.

29. Самарский А.А. Теория разностных схем. М.: Наука, 1977.

30. Артемьев В.К., Булеев Н.И. О сходимости явной схемы неполной факторизации при решении двумерных уравнений диффузии. // Вопросы атомной науки и техники. Сер. "Физика и техника ядерных реакторов". 1983. - Вып.5 (34). - С.19-25.

31. Лебедев В.И., Финогенов С.А. О порядке выбора итерационных параметров в чебышевском циклическом итерационном процессе // ЖВМ и МФ. 1971. - Т.П. - №2. -С.425-438.

32. Булеев Н.И. Методы неполной факторизации для решения двумерных уравнений эллиптического типа // Вопросы атомной науки и техники. Сер. "Физика и техника ядерных реакторов". 1980. - Вып.З (13). - С.3-14.

33. Ильин В.П. Сверхнеявные итерационные методы для решения систем уравнений с блочно-трехдиагональными матрицами // Сибирский математический журнал. 1985. -T.XXVI. - №1. - С.83-89.

34. VargaR.S. Matrix Iterative Analysis = Матричный итерационный анализ. New Jersey, 1962.

35. Birkhof G., VargaR., Young D. Alternating direction implicit method = Неявный метод переменных направлений // Advances in Comput. 1962. - V.4. - P.190-274.

36. Ильин В.П. Численные методы решения задач электрооптики. Сибирское отделение: Наука. - Новосибирск, 1974. - 202с.

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

38. Марчук Г.И. Методы расчета ядерных реакторов. М.: Госатомиздат. - 1961. - 667с.

39. Марчук Г.И. Методы вычислительной математики. М.: Наука. - 1977. - 455с.

40. Гинкин В.П., Крылова JI.M., Кухтин В.П. Комплекс программ WIMS.BOJirA // Аннотация в ВАНТ. Серия "Физика и техника ядерных реакторов". ИАЭ, 1987. - 6с.

41. Апробирование методов расчета гетерогенных уран-водных сборок : отчет о НИР / ФЭИ. Обнинск, 1981. - 22с. - Исполн.: Гинкин В.П. и др. -Инв. № 6505.

42. Гинкин В.П., Гурин В.Н. Фактор-ЗК комплекс двумерного расчета гомогенных реакторов на базе трехгрупповой библиотеки микроконстант. - Препринт ФЭИ-1152. -1982.- Юс.

43. MIX двумерный комплекс программ для расчета гомогенных и гетерогенных реакторов в малогрупповом диффузионном приближении : отчет о НИР / ФЭИ. - Обнинск, 1983. -22с. - Исполн.: Гинкин В.П. и др. - Инв. № 6839.

44. Гинкин В.П., Будникова Г.Г. ВОЛГА программа трехмерного расчета реактора в малогрупповом диффузионном приближении. - Препринт ФЭИ-1583. - 1984. - 22с.

45. Комплекс программ трехмерного расчета гомогенных и гетерогенных реакторов на основе программ WIMS-D4 и ВОЛГА : отчет о НИР / ФЭИ. Обнинск, 1986. - 40с. -Исполн.: Гинкин В.П., Крылова Л.М. - Инв. № 7192.

46. Расчет зависимости К3ф водных решеток твэлов и ТВС реакторов ВВЭР от шага решетки, плотности воды и выгорания с помощью программ WIMS-D4 и MIX : отчет о НИР / ФЭИ. Обнинск, 1986. - 51с. -Исполн.: Гинкин В.П., Налимова М.А. - Инв. № 7269.

47. Гинкин В.П. Предложения по содержанию каталога-справочника точных критических экспериментов // Доклад на встрече специалистов СССР и ГДР по проблеме транспортирования отработавшего топлива / Обнинск. ФЭИ, 1988. - 14с.

48. Гинкин В.П., Троянова Н.М. Использование метода неполной факторизации в трехмерной задаче нейтронно-физического расчета реакторов типа ВВЭР. Препринт ФЭИ-2104. - 1990. - 15с.

49. Гинкин В.П., Троянова Н.М. Использование метода неполной факторизации в трехмерной задаче нейтронно-физического расчета реактора ВВЭР-1000 // Доклад на Всесоюзном совещании по динамике реакторов // Гатчина, 1990. 15с.

50. Стационарный расчет реактора БН-800М по программам ВОЛНА и ВОЛГА : отчет о НИР / ФЭИ. Обнинск, 1992. - 25с. - Исполн.: Гинкин В.П. и др.

51. Развитие программного комплекса ВОЛНА нестационарного пространственного расчета реактора в групповом диффузионном приближении : отчет о НИР / ФЭИ. Обнинск, 1994. - 27с. - Исполн.: Гинкин В.П. и др. - Инв. № 8545.

52. Нейтронно-физический расчет кампании реактора на тепловых нейтронах с помощью комплексов программ ВОЛНА, ВОЛГА и WIMSD4 : отчет о НИР / ФЭИ. Обнинск, 1994. - Исполн.: Гинкин В.П. и др.

53. Альбом-PC. Программа комплексного расчета нейтронно-физических характеристик реакторов ВВЭР. Общее описание и инструкция для пользователя. ВНИИАЭС НПО «Энергия» : отчет о НИР / Москва, 1991. Инв. № 03-2950/91.

54. Расчет нейтронно-физических характеристик реактора ВВЭР-1000 с помощью комплекса программ WIMS-ВОЛНА : отчет о НИР / ФЭИ. Обнинск, 1995. - 51с. - Исполн.: Гинкин В.П. и др. - Инв. № 9188.

55. Тестирование трехмерного комплекса программ WIMS.BOJIHA нейтронно-физического расчета кампании реакторов применительно к активным зонам серийных ВВЭР-1000 : отчет о НИР / ФЭИ. Обнинск, 1995. - 37с. - Исполн.: Гинкин В.П. и др. - Инв. № 9068.

56. Проведение нейтронно-физических расчетов реальных состояний активной зоны реактора ВВЭР-1000 : отчет о НИР / ФЭИ. Обнинск, 1996. - 33с. - Исполн.: Гинкин В.П. и др. - Инв. № 9258.

57. Разработка метода и программы расчета трехмерной кинетики с учетом изменений конфигурации ячеек и водяных зазоров во время эксплуатации ректора : отчет о НИР / ФЭИ. Обнинск, 1997. - 75с. - Исполн.: Гинкин В.П. и др.

58. Расчет прогибов ТВС в процессе 10-й топливной кампании на блоке №2 Калининской АЭС : отчет о НИР / ФЭИ. Обнинск, 1998. - Исполн.: Лихачев Ю.И., Гинкин В.П. и др. -Инв. №9878.

59. Расчет прогибов серийных и циркониевых ТВС в процессе 10-й топливной кампании на блоке №3 Запорожской АЭС : отчет о НИР / ФЭИ. Обнинск, 1998. - Исполн.: Лихачев Ю.И., Гинкин В.П. и др. - Инв. № 9879.

60. Расчетный прогноз деформации всех ТВС активной зоны блока №1 Балаковской АЭС в процессе перегрузки и 8-й топливной кампании : отчет о НИР / ФЭИ. Обнинск, 1998. -Исполн.: Лихачев Ю.И., Гинкин В.П. и др. - Инв. № 9880.

61. Halsall M.J. A summary of WIMSD4. Input options = Описание WIMSD4. Входные опции // Winfrith, Dorchester, Dorset, Jul. 1980. AEEW-M 1327.

62. Белл Д., Глестон С. Теория ядерных реакторов. М.: Атомиздат. - 1974.

63. Хохлов В.Ф., Николаев M.H., Савоськин M.M. Комплекс программ АРАМАКО для расчета групповых макро- и блокированных микросечений на основе 26-групповой системы констант в подгрупповом представлении М.: Атомиздат. - 1972. - Ядерные константы: вып.8.

64. Алгоритм решения трехмерного нестационарного уравнения реактора в многогрупповом диффузионном приближении : отчет о НИР / ФЭИ. Обнинск, 1992. - 39с. - Исполн.: Гинкин В.П., Троянова Н.М., Шиманская Т.М., Юрина Г.П., Быков В.И. - Инв. № 8122.

65. Развитие программного комплекса ВОЛНА нестационарного пространственного расчета реактора в групповом диффузионном приближении : отчет о НИР / ФЭИ. Обнинск, 1994. - 22с. - Исполн.: Гинкин В.П., Волков Ю.В., Троянова Н.М., Ваньков К.А. -Инв. № 8545.

66. Гинкин В.П., Троянова Н.М., Ваньков К.А. ВОЛНА программа трехмерного нестационарного расчета реактора в квазистатическом приближении. - Препринт ФЭИ-2360.- 1994.-23с.

67. Гинкин В.П., Безбородов А.А., Волков А.В., Ганина С.М., Кузнецов И.А., Троянова Н.М., Швецов Ю.Е. Spatial and time-dependent calculation of fast reactor transients =

68. Пространственно-временной расчет переходных процессов в быстром реакторе // Доклад на семинаре «Consulting Meeting of IAEA specialists on fast reactor problem», 4 июня 1998. Обнинск.

69. Безбородов A.A., Волков А.В., Ганина С.М., Гинкин В.П., Кузнецов И.А., Троянова Н.М., Швецов Ю.Е. Пространственно-временной расчет переходных процессов в быстрых реакторах. // Доклад на отраслевом семинаре «Нейтроника-98», 27-29 октября 1998. -Обнинск.

70. ЮО.Шишков JI.K. Методы решения диффузионных уравнений двухмерного ядерного реактора. М.: Атомиздат. - 1974.

71. Ринейский А.А. Методы и программы подготовки групповых констант для расчета ядерных реакторов и исследование погрешностей группового приближения : диссер. канд. физ.-мат. наук. ФЭИ. - Обнинск, 1994.

72. Шиманская Т.М., Зродников А.В. Эффективный алгоритм интегрирования уравнений кинетики реактора на основе численных методов Гира. Препринт ФЭИ-1478. - 1983.

73. Ю5.Плетчер Р.Г.Достижения в области исследования турбулентной вынужденной конвекции // Современное машиностроение. Серия А. 1989. - №6. - С. 1-31.

74. Conference on Computational Fluid Dynamics = Конференция по вычислительной гидроаэродинамике, April 1995. Oxford.

75. Advances in Computational Heat Transfer = Достижения в вычислительном теплопереносе // Symposium organized by International Center for Heat and Mass Transfer (CHT'97), 2630 May 1997 / Proceedings. Ce§me, Turkey, 1997.

76. Proceedings of the International Conference on Computational Heat and Mass Transfer (CHMT99) = Труды международной конференции по вычислительному тепломассопереносу // Eastern Mediterranean University, Famagusta, April 26-29 1999.

77. Ю9.Булеев Н.И. Пространственная модель турбулентного обмена. М.: Наука, 1989. - 344 с.

78. Ю.Артемьев В.К. Явный метод неполной факторизации с чебышевским адаптируемым ускорением сходимости. Препринт ФЭИ-2095. - 1990. - 18с.

79. Ш.Артемьев В.К. Неявный метод для решения уравнений Навье-Стокса в естественных переменных // Моделирование в механике. Новосибирск, 1992. - Т.6 (23). - №1. - С.17-22.

80. Артемьев В.К. Развитие численных методов решения задач динамики вязкой жидкости : диссер. канд. физ.-мат. наук. ФЭИ. - Обнинск, 1997.

81. З.Полежаев В.И., БунэА.В., Гончаров АЛ., ВерезубН.А. и др. Математическое моделирование конвективного тепломассообмена на основе уравнений Навье-Стокса. -М.: Наука, 1987.-271с.

82. Гончаров A.JL, Фрязинов И.В. Сеточный метод решения трехмерных уравнений Навье-Стокса в параллелепипеде // Дифференциальные уравнения. 1991. - Т.27. - №7. - С.1137-1144.

83. Мажорова O.C., Марченко М.П., Фрязинов И.В. Монотонизирующие регуляризаторы и матричный метод решения уравнений Навье-Стокса для несжимаемой жидкости. // Математическое моделирование. 1994. - Т.6. - №12. - С.97-116.

84. Гинкин В.П., Корниенко Ю.Н., Кузеванов B.C. Расчетное определение локальных характеристик двухфазного потока на участке установления профиля газосодержания // Доклад на VII школе по моделям механики сплошной среды / Обнинск. ФЭИ, 1983. - Зс.

85. Гинкин В.П., Корниенко Ю.Н., Кузеванов B.C. Расчетное определение локальных характеристик двухфазного потока на участке установления профиля газосодержания // Физическая механика неоднородных сред. Новосибирск: ИТПМ СОАН СССР, 1984. -С.54-59.

86. Гинкин В.П., Корниенко Ю.Н., Кузеванов B.C. Локальная модель двухфазного неравновесного потока в круглой трубе // Доклад на II Международной школе по моделям механики сплошной среды / Владивосток, 1991. 6с.

87. Гинкин В.П., Корниенко Ю.Н., Кузеванов B.C. Локальная модель двухфазного неравновесного потока в круглой трубе. Препринт ФЭИ-2337. - 1993. г 10с.

88. Гинкин В.П., Кузеванов B.C., Иваненко И.Ю., Корниенко Ю.Н. Численное исследование профилей параметров неравновесного двухфазного потока на основе двумерной модели // Известия ВУЗов. Сер. «Ядерная энергетика». 1994. - №4-5. - 6с.

89. Решение нестационарного уравнения теплопроводности с сильно меняющимися коэффициентами в конечном слоистом цилиндре : отчет о НИР / ФЭИ. Обнинск, 1981. -14с. - Гинкин В.П. и др. - Инв. № 3115.

90. Гинкин В.П., Анисимов В.В., ГрошевА.И., КащеевВ.М., ХудаскоВ.В. Локальная теплоотдача при турбулентном течении жидкости в кольцевом канале при постоянной температуре стенки. Препринт ФЭИ-1619. - 1984. - 20с.

91. Гинкин В.П., Канухина С.В., Иванова Т.В. Численное решение нестационарной сопряженной задачи теплообмена в режиме начала кипения в трубе. Препринт ФЭИ-1834.-1986.-14с.

92. Программа расчета температурного поля в 2-х мерной радиационной защите : отчет о НИР / ФЭИ. Обнинск, 1987. - 12с. - Исполн.: Гинкин В.П., Кураченко Ю.Б. Инв. №7341.

93. Гинкин В.П., ГрошевА.И., Казанцев А. А. РЕССА программа расчета трехмерной кондуктивно-конвективной задачи теплообмена. - Препринт ФЭИ-1875. - 1987. - 17с.

94. Трехмерная методика теплогидравлического расчета испарителя установки БН-350 : отчет о НИР / ФЭИ. Обнинск, 1987. - Исполн.: Гинкин В.П. и др. - Инв. № 5132.

95. Гинкин В.П., ХудаскоВ.В., ДорошенкоВ.А., ЗининаГ.А. Методика трехмерного теплогидравлического расчета парогенераторов АЭС // Доклад на VIII Всесоюзной конференции / Минск, 1988. 10с.

96. Гинкин В.П., Худаско В.В., Зинина Г.А., Замирайлова Л.Г. Турбулентный теплообмен на входном гидродинамическом и тепловом участке в круглых каналах // Доклад на X Всесоюзной школе по моделям механики сплошной среды / ФЭИ, 1989. Зс.

97. Гинкин В.П., Гориченко В.А., ХудаскоВ.В. Методика расчета теплообмена в сборке твэлов с учетом анизотропии коэффициента переноса тепла // Доклад на международном семинаре «Теплофизика-90» / Обнинск, 1990. Юс.

98. Артемьев В.К., Гинкин В.П. Численное моделирование трехмерной естественной конвекции // Доклад на второй Российской национальной конференции по теплообмену (РНКТ2), 26-30 октября 1998 / Сб. трудов. Москва: изд-во МЭИ, 1998. - Т.З. - С.38-41.

99. Гинкин В.П., Ганина C.M. Новый метод расчета трехмерной конвекции на сетках большой размерности // Математическое моделирование. 2003. - Т. 15. - №6. - С.53-58.

100. Булеев Н.И., Тимухин Г.И. О составлении разностных уравнений гидродинамики вязкой нестационарной жидкости // Численные методы механики сплошной среды. — СО РАН, 1972. Т.З. - №4.-С. 19-26.

101. Мюллер Г. Выращивание кристаллов из расплавов. М.: Мир, 1991. - 150с.

102. Келли А., Гровс Г. Кристаллография и дефекты в кристаллах. М.: Мир, 1974. - 496с.

103. Мильвидский М.Г., Освенский В.Б. Получение совершенных монокристаллов // Сб. статей «Проблемы современной кристаллографии». М.: Наука, 1975. - С.79-109.

104. Смирнов Б.И. Дислокационная структура и упрочение кристаллов. М.: Наука, 1981. -236с.

105. Новиков И.И., Розин К.М. Кристаллография и дефекты кристаллической решетки. М.: Металлургия, 1990.-336с.

106. Бартел И., БурингЭ., ХайнК., КухаржЛ. Кристаллизация из расплавов. М.: Металлургия, 1987. - 320с.

107. Математическое моделирование получения монокристаллов и полупроводниковых структур. М.: Наука, 1986. - 197с.

108. Вабищевич П.Н. Численные методы решения задач со свободной границей. Изд-во Моск. универ., 1987. - 164с.

109. Роуч П. Вычислительная гидродинамика. -М.: Мир, 1980. 616с.

110. Белоцерковский О.М. Численное моделирование в механике сплошных сред. М.: Наука, 1984.-520с.

111. Пасконов В.М., Полежаев В.И., Чудов JI.A. Численное моделирование процессов тепло-и массообмена. М.: Наука. - 288с.

112. Патанкар С. Численные методы решения задач теплообмена и динамики жидкости. — М.: Энергоатомиздат. 149с.

113. Артемьев В.К., Булеев Н.И. О решении уравнения Навье-Стокса в переменных "вихрь -функция тока". Препринт ФЭИ-1734. - 1985. - 28с.

114. ПейреР., Тейлор Т. Вычислительные методы в задачах механики жидкости. — Л.: Гидрометеоиздат, 1986. 352с.

115. Булеев Н.И. Пространственная модель турбулентного обмена. М.: Наука, 1989. - 344с.

116. Андерсон Д., Таннехил Дж., Плетчер Р. Вычислительная гидромеханика и теплообмен. — М.: Мир, 1990.-Т. 1-2.

117. Тарунин Е.Л. Вычислительный эксперимент в задачах свободной конвекции. Изд-во Иркут. универ., 1990. - 228с.

118. Sabhapathy P., Salcudean М.Е. Numerical analysis of heat in LEC growth of GaAr = Численный анализ теплообмена методом LEC при выращивании GaAr // J. Crystal Growth. 1989. - No.97. - P.125-135.

119. Бердников B.C., ПанченкоВ.И. Некоторые характеристики смешанной конвекции в лабораторной модели метода Чохральского // Сб. трудов «Теплофизика кристаллизации веществ и материалов». Новосибирск: ИТФ СО АН СССР, 1987. - С.5-15.

120. БердниковB.C., ПанченкоВ.И., СоловьевС.В. Тепловая гравитационно-капиллярная конвекция в методе Чохральского // Сб. трудов «Теплофизика кристаллизации высокотемпературной обработки материалов». Новосибирск: ИТФ СО АН СССР, 1990. -С. 162-199.

121. ПО.Гинкин В.П., Гориченко В.А., Глухов С.Б., Межерицкий Г.С., Худаско В.В. Методика расчета температурных полей в многослойной двумерной конструкции с лучистым теплообменом между слоями. -Препринт ФЭИ-1986. — 1989. — 17с.

122. Марченко М.П., Сенченков А.С., ФрязиновИ.В. Математическое моделирование процесса роста кристаллов из раствора-расплава методом движущегося нагревателя // Математическое моделирование. 1992. - Т.4. - №5. - С.67-79.

123. Марченко М.П., Фрязинов И.В. Комплекс программ КАРМА1 решения нестационарных задач выращивания монокристаллов в ампулах // Вычислительная математика и математическая физика. 1997. -Т.37. -№8. - С.988-998.

124. Авдеев А.Ю., Гончаров В.А. Решение задачи Стефана для нестационарного процесса роста кристаллов в условиях микрогравитации. // Изв. вузов. Сер. «Электроника». 1999. -№1. - С.9-15.

125. Гончаров В.А. Об одном методе решения двухфазной задачи Стефана с неплоской границей // Вычислительная математика и математическая физика. — 2000. — Т.40. №11. - С.1706-1715.

126. ArtemyevV.K., GinkinV.P., GusevN.V., ZininA.I., OzernykhI.L. Numerical Simulation of GaAs Crystal Growth of a Big Diameter by Bridgman Method = Численное моделирование выращивания кристалла GaAs большого диаметра методом Бриджмена // The 5-th Int.

127. Математическое моделирование приборов и технологий для системы МОПИТ : отчет о НИР / ФЭИ. Обнинск, 1992. -276с. - Исполн.: Гинкин В.П. и др. - Инв. № 8075.

128. Гинкин В.П., Гадияк Г.В., Шварц Н.Л., Обрехт М.С. Синица С.П. Программа расчета стационарных характеристик МДП-транзистора MOS-1 // Автометрия. 1987. - №1. -С.70-75.

129. Гинкин В.П., Артемьев В.К., Бердников B.C. Численное моделирование конвективного теплообмена на модели Чохральского // Доклад на 1-ой Российской национальной конференции по теплообмену (РНКТ), 21-25 ноября 1994 / Тез. докл. М.: МЭИ, 1994. -Юс.

130. Методика получения оценок времени нагревания и плавления образца кристалла и получающихся при этом температур и их градиентов : отчет по контракту №260 / ФЭИ. -Обнинск, 1994. 12с. - Исполн.: Гинкин В.П. и др.

131. Методика оценки риска полного расплавления затравочного стержня при подводимой мощности : отчет по контракту №260 / ФЭИ. Обнинск, 1994. - Юс. - Исполн.: Гинкин В.П. и др.

132. Гинкин В.П., Артемьев В.К., Денисов B.C., ЗининаГ.А., Корниенко Ю.Н. Численное моделирование процессов выращивания монокристаллов : обзор ФЭИ-0267 / М.: ЦНИИатоминформ, 1994. 13с.

133. Ginkin V.P. Problems of Numerical Simulation of the Single Crystal Growth Process // The Fifth International Conference on Simulation of Devices and Technologies, October 11-20 1998 / Proceedings. Cape-Town, South Africa, 1998.

134. August 2001 / Abstracts. Kyoto, Japan, 2001. - P.386.

135. Гинкин В.П., Забудько M.A., Картавых A.B., Мильвидский М.Г., Науменко О.М. Численное моделирование процесса тепломассопереноса с позиции кластерной модели строения расплава // Поверхность. 2004. - №6. - С.93-100.

136. Chen G., Roux В., Camel D., Tison P., Garandet G.P., Favier J.J., Senchenkov A.S., Moreau R. Numerical study of GEZON Experiment = Численные исследования эксперимента GEZON // Microgravity Sci. Technol. 1994. - Vol.7. -No.2.

137. Котов С.В., Лютиков А.Р., Хухрянский Ю.П. и др. О возможном механизме роста кристаллов из расплава в условиях невесомости // Письма в ЖТФ. 2002. - Т.28. - №14. -С.15-18.

138. Duffar Т., Serrano M.D., Moore C.D., CamasselJ., ContrerasS., DusserreP., Rivoallant A., Tanner B.K. Bridgman solidification of GaSb in space = Кристаллизация GaSb методом Бриджмена в космосе // J. Crystal Growth. 1998. - No.192. - P.63-72.