Содержание к диссертации
Введение
Глава 1. Используемый численный метод, схема и архитектура кодов NARAL 21
1.1. Методы конечных разностей для решения системы дифференциальных уравнений в частных производных 21
1.2. Конечно-элементный метод контрольных объемов 33
1.2.1. Применение метода конечных элементов к задачам гидродинамики и теплообмена 33
1.2.2. Способы аппроксимации скорости и давления с учетом особенностей метода контрольных объемов 36
1.2.3. Решение обобщенною конвективно-диффузионного переноса методом контрольных объемов... 38
1.2.4. Применение описанной модификации метода контрольного объема с одинаковым порядкомаппроксимации скорости и давления для решения поставленных задач гидродинамики и теплообмена 45
1.3. Базовые уравнения кода NARAL/FEM 53
Глава 2. Моделирование взаимодействия кориума с наполнителем и формирования инверсно стратифицированного бассейна расплава 61
2.1. Исходные данные для расчетов 64
2.2. Нагрев и плавление куска гематита, окруженного оксидными компонентами кориума 66
2.2.1 Математическая модель 66
2.2.2 Результаты расчетов 68
2.3. Одномерная модель теплового взаимодействия кориума с наполнителем 70
2.3.1 Математическая модель
2.3.2 Результаты расчетов 74
Глава 3. Численное исследование нестационарного теплового состояния кориума после инверсии стратификации расплава . 79
3.1. Одномерная модель задачи 81
3.2. Двумерная модель задачи 85
3.2.1 Расчеты с условием теплового излучения в серой среде 87
3.2.2 Расчеты с условием теплового излучения в прозрачной среде 93
3.2.3 Граничные условия и структура ванны расплава 94
3.2.4. Метод решения 97
Глава 4. Результаты расчетов теплового состояния кориума после инверсии стратификации расплава 101
4.1. Результаты расчетов по одномерной модели 101
4.2. Результаты расчетов по двумерной модели 111
Заключение 136
Список использованной литературы 138
- Применение метода конечных элементов к задачам гидродинамики и теплообмена
- Применение описанной модификации метода контрольного объема с одинаковым порядкомаппроксимации скорости и давления для решения поставленных задач гидродинамики и теплообмена
- Одномерная модель теплового взаимодействия кориума с наполнителем
- Расчеты с условием теплового излучения в серой среде
Введение к работе
Коды расчета процессов взаимодействия кориума с корпусом реактора
Мировая ядерная энергетика сегодня продолжает развиваться. АЭС сегодня вырабатывают 16,4% всей электроэнергии, производимой в мире. Ключевой проблемой, определяющей настоящее и ближайшее будущее развития ядерной энергетики является безопасность АЭС и, прежде всего, устойчивость АЭС к тяжелым авариям. Изучение тяжелых аварий АЭС проводят во многих странах , в том числе и в России. В настоящее время разработано много экспериментальных моделей [1,2].
Существуют следующие программные коды для расчета аварий АЭС по направлению развития их процессов (табл. В.1) [3].
Таблица В.1. Программные коды для расчета аварийных процессов на АЭС
Все крупные страны, развивающие атомную энергетику разрабатывают модели и вычислительные коды для описания этих процессов. При этом используется двухуровневный подход к созданию вычислительных кодов:
Уровень 1. Интегральные коды, предназначенные для расчетов всей совокупности аварийных процессов для каждого конкретного сценария в целом. В США Комиссия по регулированию ядерной деятельности (USNRC) разрабатывает и использует два кода: STCP (Source Term Code Package) и MELCOR. Первый интегральный код STCP охватывает 9 стадий развития аварийного процесса (кроме 10-ой стадии анализа последствия
за пределами АЭС). Второй интегральный код MELCOR предназначен для расчета всех 10 стадий развития аварийного процесса.
Уровень 2. Детальные коды, моделирующие те или иные аварийные процессы либо их отдельные стадии.
Из рассматриваемых стадий особенно важна 5-тая стадия анализа и расчета процесса разрушения корпуса реактора. Корпус реактора является одним из барьеров безопасности. При попадании на его днище обломков и расплава разрушенной активной зоны за счет остаточного тепловыделения происходит разогрев образовавшегося кориума и при тепловом и физико-химическом взаимодействии с корпусом реактора может произойти его разрушение.
Можно выделить две основных группы процессов, определяющих поведение расплава в корпусе реактора:
первичное нарушение целостности активной зоны, обусловленное взаимодействием различных компонентов (UO2 - Zr - Zr02 - В4С) при нагреве активной зоны, приводящее к перемещению материалов, образованию блокад и формированию предварительного состава смеси, попадающей на днище корпуса реактора.
тепловое и физико-химическое взаимодействие образовавшегося
кориума с внутрикорпусным конструкциями и корпусом реактора,
которое зависит от исходного состава первичных перемещаемых
материалов и сценария взаимодействия с корпусом.
В 1992 была инициирована российско-американская программа
«РАСПЛАВ» [2]. Программа включает исследование всего круга явлений,
имеющих место при взаимодействии расплава с корпусом реактора, а
именно:
Физико-химические проблемы взаимодействия расплава с корпусом реактора;
Теплофизические проблемы взаимодействия;
3. Механическое поведение корпуса с точки зрения его механических свойств с учетом приведенных выше механизмов взаимодействия с расплавом.
Предварительный анализ позволил выделить следующие диапазоны
температур:
Область относительно низких температур: нижняя граница-наименьшая температура образования первичных жидких эвтектик (« 1250 К) , верхняя граница - ниже температуры плавления корпусной стали (* 1800 К).
Область промежуточных температур - выше температуры плавления корпусной стали (« 1800К), но ниже температуры полного расплавления композиций UO2 - Zr02 - Zr(O) - нержавеющей стали (2000 - 2500 К в зависимости от концентрации кислорода и состава композиции).
Область наиболее высоких температур - область образования окончательного бассейна расплава в корпусе реактора (выше 2500К в зависимости от концентрации кислорода и состава расплава).
Главные цели программы «РАСПЛАВ» заключаются в следующем:
Оценить возможность удержания расплава активной зоны реактора в корпусе реактора на основании экспериментальных данных, полученных в интегральном эксперименте по взаимодействию расплавленных реальных материалов активной зоны корпуса реактора.
Разработать на основе анализа результатов интегрального эксперимента и ряда поддерживающих экспериментов методологию теоретического моделирования исследуемого явления.
Разработать соответствующий расчетный код.
Проводимые в рамках программы исследования позволили решить следующие проблемы:
оценить роль кинетики физико-химического взаимодействия корпуса с кориумом различного состава, имеющего эвтектический характер, в области исследуемых температур; исследовать тепломассообмен горячей части кориума с пристеночным слоем;
изучить физические характеристики расплавов (вязкости, излучательной способности, элементного и фазового составов, прочностных характеристик), образующихся на завершающей стали взаимодействия расплава с корпусом.
исследовать характеристики тепломассопереноса в различных температурных областях, для различного состава и фазового состояния кориума.
При анализе тяжелых аварий основным является вопрос проплавлений корпуса реактора. Необходимо достаточно точно представлять возможные последствия проплавлення корпуса с точки зрения нагрузок на контейнмент, а также основные временные характеристики процесса взаимодействия.
Рассматривается сценарий полного обрушения и расплавления активной зоны. Твердые или жидкие обломки активной зоны попадают на дно корпуса реактора и нагреваются до температуры плавления корпусной стали. Если остаточное тепловыделение велико, а теплоотвод от таблеток недостаточен, происходит разогрев и расплавление таблеток топлива. Остаточное тепловыделение распределяется по нескольким каналам потерь: на отвод тепла к корпусу и к верхней поверхности расплава, на нагрев топливных таблеток и расплава. После расплавления таблеток расплав разогревается до температур, превышающих температуру ликвидуса кориума и начинается стадия высокотемпературного взаимодействия.
Процессы, сопровождающие аварию с плавлением днища корпуса реактора, объединяют в себе тепло-и гидродинамику жидкого кориума с
образованием корок на поверхности расплава и гарнисажа на поверхности корпуса реактора , физико-химическое взаимодействие кориума со сталью корпуса с возникновением легкоплавких эвтектик и термомеханику корпуса, удерживающего массу расплавленного кориума. Одним из основных факторов, определяющих эти процессы, является уровень и профиль распределения тепловыделения. Он зависит от выхода осколков деления и перемещения активных составляющих кориума за счет температурной и концентрационной конвекции. Определяющими факторами для процесса взаимодействия кориума с корпусом реактора также являются внешние граничные условия, т.е. теплоотдача от днища реактора и поверхности кориума.
Особую задачу представляет собой исследование физико-химического состояния среды (кориума) в процессе развития аварии, эволюции компонентного и фазового составов, физико-химического взаимодействия расплава с конструкционными материалами, включая и корпус реактора. Известно, что в расплаве кориума неизбежно расслоение на две составляющие - оксидную и металлическую, причем наблюдается значительное отличие в концентрации урана в той или другой фазах.
К настоящему времени в России разрабатывается ряд компьютерных кодов и моделей для расчета процессов, сопровождающих взаимодействие кориума с корпусом реактора.
1. Пространственная модель теплопроводности (НИТИ)
Код расчета температур состояния корпуса реактора разработан на основе пакета программ «FEMTEM»[5]. Для расчета была дана модель, в которой рассматривается уравнение для изменяющейся во времени функции Т :
С(5Т/ д t) + div (КЛ gradT)+Q = О В области G с граничными условиями общего вида: КЛ ((ЭТУ 5п) + h(T-To) + q = 0 ,
где :
С, Q, h, Т, q - заданные пространственные функции
К - диагональная матрица коэффициентов
dTVdn- производная искомой функции в направлении нормали к
поверхности;
дТ/ dt-частная производная функции Т по времени t
Данное уравнение рассматривается в пакете «FEMTEM» в одномерной,
двухмерной (в декартовых или цилиндрических координатах) или
трехмерной геометрии.
В пакете «FEMTEM» входят модули:
- интерактивной подготовки текстового файла входных данных и
формирования пакета задачи;
интерполяции узловых значений температур на декартовую сетку и последующего построения изолиний;
графического представления геометрической модели и полей температур в пространственно-временных точках.
Основной язык программирования - ФОРТРАН -77.
Естественным недостатком кода является неопределенность
эффективного коэффициента теплопроводности кориума. Более того,
использование в расчетах постоянного значения этого коэффициента по
пространству (сложность аналитического или табличного представления
зависимости К(х,у), каким-либо образом имитирующего процесс
естественной конвекции) может привести к некорректным результатам
прогнозирования места проплавлення корпуса. Тогда оптимальной
пространственной моделью является математическая модель,
непосредственно описывающая естественную конвекцию кориума в полости днища корпуса реактора.
2. Двухмерная модель взаимодействия расплава кориума с корпусом реактора (код РАСПЛАВ) [6,7]
В программном коде РАСПЛАВ реализованы физические модели, описывающие теплогидродинамические процессы взаимодействия, включая проблемы охлаждения расплава.
Результатами численных исследований по программному коду РАСПЛАВ являются:
оценка времени установления стационарного режима естественной конвекции в расплаве топлива как в ламинарном, так и в турбулентном режимах;
картина температурных полей в корпусе, необходимая для расчетов термических напряжений и механических деформаций корпуса;
- получение расчетных корреляций зависимости потоков энергии к
верхней границе расплава и к контактной границе расплава с
корпусом в двух режимах конвекции — ламинарном и турбулентном,
как функции числа Рэлея;
- получение основных закономерностей и критериев подобия для
экстраполяции результатов и разработки быстрых упрощенных
моделей.
Уравнения кода «РАСПЛАВ». Система уравнений, самосогласованно описывающая процессы тепломассопереноса в расплаве топлива включает нестационарные уравнения теплопроводности и гидродинамики несжимаемой жидкости. Уравнение теплопроводности в двумерной геометрии с учетом цилиндрической симметрии записывается в виде:
\ot oz r J r or\ or) oz\ oz )
где ( r,z )- цилиндрические координаты ; t - время; T - температура ; р -плотность; с(Т)- теплоемкость; к(Т)-теплопроводность; к(Т)- объемная плотность тепловыделения ; и и v -соответствующие компоненты скорости.
Двухмерные нестационарные уравнения Навье-Стокса для несжимаемой жидкости с переменной вязкостью в цилиндрических координатах имеют вид
др д
ди ди ди
+ и— + v— = F(r)- — +
dt дг dz дг дг
,ди dv,
v(— + —) +-(2v—)
dz дг J r or
dp д
dv dv dv
— + u— + v— = F(z)- —+—
dt dr dz dz dz
2v— dz j
H
.du dv. v(— + —)
dz dr
+ — r
,du dv. v(— + —)
V dz dr
Уравнение неразрывности в принятых предложениях о несжимаемости записывается следующим образом:
d(ru) d(rv)
Здесь v - кинематическая вязкость; F(r), F(z) - компоненты внешних сил ( в этой задаче F(r)=0, F(z) = g(l+P(T-Tc)); g -ускорение свободного падения, Р - коэффициент объемного расширения )
Приведенная система данных уравнений описывает ламинарную конвекцию расплава. Уточнение модели для турбулентных конвекции может проводиться за счет изменения.теплофизических и гидродинамических характеристик среды (добавляется, например, турбулентная вязкость).
Модели турбулентности. Численное моделирование турбулентной конвекции в замкнутых объемах является в настоящее время трудной проблемой. Это связано прежде всего с необходимостью построения адекватных моделей турбулентности. Среди общего многообразия моделей турбулентности следует отметить прежде всего алгебраические модели, использующие те или иные соображения по дополнительной (турбулентной) вязкости и модели с дополнительными
дифференциальными уравнениями, наиболее известной из которых является к — є модель.
Граничные условия. На верхней поверхности, расплава, в соответствии со сценариями тяжелых аварий, возможен ряд вариантов теплообмена.
а. Радиационный теплообмен с окружающими конструкциями,
доминирующий при отсутствии теплоносителя в корпусе. При этом
выражение для теплового потока на границе имеет вид:
qu = a(T4-Ts4) ,
где Г— температура поверхности расплава; Т,— температура окружающей среды; а— постоянная Стефана-Больцмана.
б. Теплообмен на верхней границе при взаимодействии расплава
кориума с водой. Выражение для потока на границе имеет вид:
qu = h(T - Т w) ,
где Т — температура поверности расплава; Тж — температура воды; h — коэффициент теплопередачи.
При изотермической постановке задачи, на верхней и внешней границе (граница раздела расплав/металл) задается постоянное значение температуры Т = То, равное температуре Т пов.
Для уравнений гидродинамики задаются следующие условия: на границе топлива с корпусом — условие прилипания (равенство нулю скорости на границе), на оси симметрии и верхней границе - условие скольжения (равенство нулю нормальных компонент скоростей). В случае образования корки, на верхней границе задастся условие прилипания.
Задача о переносе тепла через корпус реактора. Для изучения проблемы теплопередачи через корпус необходимо рассматривать согласованную задачу тепломассопереноса в топливе и корпусе реактора.
Уравнение теплопроводности должно решатся с учетом фазовых переходов и переменности теплофизических свойств корпусной стали.
При решении согласованной задачи тепломассопереноса в топливе и корпусе реактора используется схема сквозного счета в которой фазовая граница особо не выделяется, а используется однородная схема.
Можно рассматривать задачу о распространениии тепла в корпусе, определяя граничные условия на внутренней поверхности корпуса, из решения гидродинамических уравнений для расплава топлива (с учетом возможной корки на границе взаимодействия). В упрощенной постановке задачи о проплавлений корпуса, позволяющей выявить основные закономерности взаимодействи, граничный поток вычисляется по известным корреляциям [8,9].
Код РАСПЛАВ можно считать более развитым по отношению к другим. Он решает самосопряженую задачу теплообмена кориума с корпусом реактора при его плавления (задача Стефана). Однако используемая в нем алгебраическая модель турбулентности давала только качественное но отнюдь не количественное согласие с опытом (результаты при ее верификации на эксперименте СОРО[10]). Кроме того, использование задачи Стефана при больших перепадах температуры солидуса и ликвидуса, которые характерны для проблемы взаимодействия кориума с корпусом реактора, может оказаться некорректным.
В течениии нескольких лет на кафедре АЭС МЭИ и ЭНИЦ разрабатывается расчетный программный комплекс (ПК) NARAL/FEM (являющийся дальнейшим развитием ПК NARAL/CAVITY [11-15]), для описания процессов, протекающих при свободной конвекции расплава A3 с учетом стенок окружающей конструкции.
Актуальность исследований с помощью этого ПК обусловлена следующими причинами:
необходимость изучения сценариев тяжелых аварий и решений по управлению ими,
необходимость оценки способности подреакторного устройства удержания расплава кориума.
Цель и постановка задачи настоящей работы. Настоящая работа посвящена созданию методики расчета и разработке программных модулей для расчета свободно-конвективного теплообмена в объеме расплава применительно к устройству локализации расплава.
В МЭИ, в научной группе под руководством проф. Н.Г.Рассохина
продолжается развитие расчетного программного кода NARAL/FEM [17-
28], который основан на концепции конечно-элементного метода
контрольного объема. Математическая модель включает в себя систему
осредненных по Рейнольдсу уравнений Навье-Стокса и уравнений
энергии. Для учета турбулентности используется или
двухпараметрическая к-є модель турбулентности или различные модификации алгебраических моделей. Учет плавления проводится в рамках модели эффективной энтальпии [58-64] как в бинарном расплаве или как в эффективном веществе с заданной температурой солидуса и ликвидуса.
Данная работа заключается в исследовании процесса и разработке метода расчета теплообмена в расплаве кориума с момента его поступления в подреакторное устройство (ловушку) после нарушения целостности корпуса реактора.
Подреакторная ловушка является новым оборудованием локализации тяжелой аварии и разрабатывается в ОКБ «ГИДРОПРЕСС» и «НИТИ» применительно к проекту ВВЭР в Китае. Она представляет собой бетонную шахту в подреакторном пространстве. Корзина в нижней части шахты, выполненная из чередующихся слоев стали и
диоксидциркониевого бетона, охлаждается водой. В корзину помещается крупнопористый наполнитель в виде кусков (блоков) гематита. Предполагается, что через некоторое время после поступления расплавленного кориума в ловушку куски гематита плавятся и расплав перемешивается с оксидными компонентами кориума. В результате плотность смеси оксидных компонент уменьшается и происходит инверсия стратификации расплава: если первоначально расплав металлов (стали и циркония) располагается на поверхности бассейна, то после инверсии расплавленные металлы оказываются внизу, под слоем оксидных компонент кориума.
В решаемой задаче естественным образом выделяются две стадии процесса:
плавление наполнителя и его растворение в оксидных компонентах кориума, после чего сразу же происходит инверсия расплава: в верхней части бассейна оказываются теперь уже более легкие оксидные компоненты, а вниз перемещается расплав металлов;
естественная конвекция в инверсно стратифицированном бассейне кориума при постепенно уменьшающемся остаточном тепловыделении и интенсивной передаче тепла стенкам шахты.
Детальные вариантные расчеты, выполненные ранее [17,20,23,26] применительно к задаче удержания в корпусе реактора, продемонстрировали перспективность совершенствования расчетов, учитывающих различные факторы, от которых зависит развитие ситуации в условиях с плавлением активной зоны и части корпуса реактора. В настоящей работе также производится детализация теплогидравлических расчетов и учета в них следующих основных факторов:
турбулентная свободная конвекция ;
отвод тепла в окружающую среду излучением;
- сопряженный конвективный теплообмен в ванне кориума и теплообмен
излучением над зеркалом расплава под сводом конструкции.
В работе на базе ПК NARAL/FEM применительно к ловушке проведена реализация расчетной методики в упрощенном коде, в котором учитывается теплообмен излучением, а для учета турбулентной конвекции используется достаточно простая, но проверенная на экспериментальных данных полуэмпирическая модель эффективной теплопроводности.
Основные результаты и их научная новизна:
Осуществлено численное исследование теплогидравлических процессов при взаимодействии расплава кориума с наполнителем и теплообменных процессов в ловушке после инверсии расплава . В частности:
Предложена модель теплового взаимодействия кориума с крупнопористым наполнителем с переменными распределенными стоками тепла для определения нестационарных профилей температуры компонент кориума и наполнителя. Выполненные расчеты позволили выявить несколько стадий теплового взаимодействия кориума и наполнителя вплоть до формирования инверсно стратифицированного бассейна расплава.
Проведена модификация ранее разработанной модели эффективной теплопроводности для расчета теплообмена после инверсии в бассейне кориума и определения конвективного теплового потока к стенкам подреакторного устройства (ловушки).
Впервые предложен численный алгоритм сопряженного «сквозного» расчета теплового состояния ловушки с учетом конвективного теплообмена в ванне расплава и переноса излучения в условиях поглощения и рассеяния излучения частицами аэрозоля.
- Впервые выполненная серия расчетов для вариантов с различным
составом кориума, соответствующих разным конструктивным
решениям, показала возможность длительного удержания расплава в
подреакторной ловушке.
Практическая значимость работы:
Проведена работа по доработке и верификации отечественного конечно-элементного кода NARAL/FEM, который с использованием реализованных в нем моделей позволяет исследовать теплогидравлические процессы, протекающие как при взаимодействии расплава кориума с наполнителем, так и при оценке состояния подреакторных устройств удержания расплавленного кориума после инверсии стратификации расплава.
Для облегчения работы с усложнившимися структурами данный программный код NARAL/FEM дополнен интерфейсными модулями предварительной подготовки расчетов и обработки полученных результатов.
Реализована и протестирована методика ускоренного сопряженного анализа теплового состояния в ванне расплава кориума. Используя существенное сокращение времени расчета, произошедшее благодаря применению модели эффективной теплопроводности, выполнены вариантные расчеты.
Для реальной геометрии подреакторной ловушки, применяемой на практике к проекту ВВЭР в Китае, на основе реалистичных оценок теплофизических свойств кориума и ловушки в нестационарной постановке впервые выполнены двумерные вариантные расчеты сопряженого «сквозного» теплообмена в ванне ловушки после инверсии расплава кориума. В расчетах выявлены сочетания варьируемых параметров (структура, высота ванны, условия среды над зеркалом
расплава кориума), при которых возможно длительное удержание кориума; выявлены ситуации с различной степенью тепловой нагрузки на стенки ловушки в зависимости от интенсивности наружного охлаждения и излучения с зеркала расплава. - Результаты данной работы могут быть использованы в термопрочностных расчетах для прогнозирования теплового и механического поведения корпуса подреакторной ловушки, заполненной расплавом кориума.
Структура диссертации. Диссертации состоит из введения, четырех глав и заключения.
Во введении проводится обзор различных отечественных и зарубежных программных кодов расчета тяжелых аварий на АЭС, определяется цель работы, представляется постановка задач диссертации.
В главе 1 приводится описание метода конечных элементов в концепции контрольного объема, численной схемы и архитектуры реализующего их кода NARAL/FEM, который используется для вычислений в данной работе.
В главе 2 рассматривается модель теплового взаимодействия кориума с крупнопористым наполнителем, учитывающая прогрев и плавление наполнителя. Приводится алгоритм решения одномерной задачи теплопроводности с переменными распределенными стоками тепла для определения нестационарных профилей температуры компонент кориума и наполнителя, формирующих инверсно стратифицированный бассейн расплава.
В главе 3 анализируется расчетное моделирование теплопереноса в инверсно стратифицированном бассейне расплава. В основу теоретической модели положено приближенное асимптотическое описание переноса тепла при турбулентной естественной конвекции в случае больших чисел
Рэлея, когда общая задача вырождается и возможно введение переменной эффективной теплопроводности среды.
В главе 4 представляются результаты расчетов теплообмена в подреакторном устройстве удержания расплава после инверсной стратификации.
В заключении сформулированы выводы по результатам работы.
Основное содержание проведенных с участием автора исследований отражено в работах [4,25,27,28,35,69] и далее представлено более подробно.
Применение метода конечных элементов к задачам гидродинамики и теплообмена
Методы конечных элементов все более широко используются для решения задач гидродинамики и теплообмена. Преимущества этих методов особенно проявляется при решении задач в областях сложной геометрии. Учебники Бэйкера [33] и Чэнга [34] посвящены исключительно предмету расчета течений жидкости с помощью метода конечных элементов, работы [36-42] являются примерами исследований в этой области.
Численные процедуры, в которых разностные уравнения получаются путем применения законов сохранения к контрольным объемам, окружающим узловую точку, имеют определенную привлекательность. Полная формулировка базируется на таких имеющих физический смысл величинах как потоки, источниковые члены, принципы сохранения и т.д. В то время, как конечно-разностные методы уже хорошо утвердились [43], разработка конечно-элементных методов на базе контрольных объемов происходит только в последнее время. Балига и Патанкар [30-32] впервые представили такие методы для решения двумерных задач с использованием треугольных элементов.
Как и многие другие конечно-элементные процедуры, метод Балиги и Патанкара [31 ] относится к типу с неодинаковым порядком аппроксимации скорости и давления. То есть давление рассчитывается на много меньшем числе узлов, чем скорость. Это делает метод достаточно сложным. Два различных набора элементов, называемых макро- и микроэлементами, предназначены для давления и скорости, соответственно. Уравнения движения и неразрывности дискретизируются с использованием двух различных семейств контрольных объемов и, следовательно, сохранение массы точно не выполняется по контрольным объемам для уравнений движения. Кроме того, из-за наличия двух наборов элементов и контрольных объемов, метод влечет за собой хранение большого объема геометрической (топологической) информации. Разностное уравнение для давления существенно отличается от разностного уравнения для любой другой переменной величины и нуждается в отдельном блоке решения. Наконец, поскольку метод неодинакового порядка аппроксимации скорости и давления использует более грубую сетку для давления, он может быть не очень точным для задач в которых давление сильно изменяется в пределах расчетной области.
В последующих работах [44, 45] Патанкаром и Пракашем предложен вариант конечно-элементного метода контрольного объема, использующий одну общую сетку конечных элементов. Предлагаемый метод равного порядка аппроксимации скорости и давления, основанный на [45], не имеет перечисленных выше трудностей. Скорость и давление вычисляются в одинаковых наборах узловых точек. Уравнения движения и неразрывности дискретизируются с использованием одного и того же набора контрольных объемов, и разностное уравнение для давления может быть решено с помощью той же процедуры решения, что используется для какой-либо другой переменной.
Конечно-элементный метод [30-32] использует треугольные элементы для решения двухмерных задач. Наибольшим вкладом, представленным в работах [30-32] является разработка соответствующих функций формы для моделирования конвективно-диффузионных процессов. Для каждого элемента предложенные функции формы изменяются экспоненциально в направлении вектора средней скорости и линейно в нормальном направлении. Такие функции формы правильно моделируют односторонний характер конвекции и двусторонний характер диффузии, и, поскольку функция формы распространяется в направлении локального вектора скорости, она в значительной степени уменьшает диффузию поперек потока.
Функции формы, разработанные Балигой и Патанкаром [30] очень похожи на противопоточную асимметричную процедуру дифференцирования Райтби [46]. Отрицательной стороной схемы является возможность получения отрицательных коэффициентов дискретизации, что может привести к превышению или занижению решения (то есть, когда средние значения расчетной величины оказываются больше или меньше его значения на границе даже в отсутствии каких-либо источников или стоков). Для практических задач, включающих турбулентность и горение, превышение или занижение решений неприемлемо, отрицательные значения могут возникать даже для таких величин как кинетическая энергия, скорость реакции и т. д., что не имеет физического смысла. Следовательно, для программ общего назначения необходима процедура дискретизации, гарантирующая положительные коэффициенты. Данная работа мотивирована именно этими соображениями. Схемы "против потока", используемые в данной работе, приводят к положительным коэффициентам в дискретном аналоге и устраняет возможность получения нефизичного решения. Однако, следствием процедуры, описанной в [45] является повышение диффузии в поперечном направлении [43]. Поэтому использование схемы "против потока" служит для оценки исследуемой задачи с точки зрения здравого смысла. Возможно, данная схема может быть полезной как альтернативная при использовании ее только для переменных, для которых последствия отрицательных коэффициентов дискретизации являются полностью неприемлемыми.
Шнейдер и Pay [50] предложили привлекательную схему, в которой используется концепции асимметричности для снижения поперечной диффузии. Но они ограничили степень асимметричности с тем, чтобы гарантировать положительность коэффициентов. Схема разработана в контексте метода конечных элементов в формулировке контрольных объемов с применением четырехугольных элементов. Метод оказался многообещающим даже несмотря на некоторые трудности реализации. Предлагаемая в данной работе скошенная схема "против потока", похожая на схему Райтби [46], соединяет все достоинства схемы "против потока" и экспоненциальных схем [30-32].
Применение описанной модификации метода контрольного объема с одинаковым порядкомаппроксимации скорости и давления для решения поставленных задач гидродинамики и теплообмена
Поле їна элементе P-L-К может быть выражено через узловые значения и, v, Du, Dv и р, значения Su и Sv на элементе и некоторые линейные функции координат г Когда таким образом определенное поле и используется в интегральном выражении (1.29), может быть получен разностный аналог для давления. Для типичной узловой точки он имеет вид Алгебраические детали вывода (1.33), аналогичные построению дискретного аналога для уравнения теплопроводности, изложены в работе [33]. Уравнение (1.33) применяется к внутренним узловым точкам. Для граничных точек, таких, как показанных на рис. 1.8, дискретный аналог имеет вид где тр- поток массы, покидающий область через границу g-P-a. Если в качестве граничного условия в точке Р задается давление, то уравнение (1.34) используется для определения тр. С другой стороны, если тр (или нормальная составляющая скорости на границе) задается, то по уравнению (1.34) определяется давление в точке Р. В принципе, процедура решения задач гидродинамики является полной. Уравнения (1.22), (1.23) и (1.33) должны решаться итерационно для определения скорости и давления во всех точках области. Однако, некоторые детали требуют дополнительного обсуждения. Прежде всего обсудим некоторые важные и тонкие особенности предлагаемого метода. Исключение неправдоподобных решений Как было упомянуто ранее основная трудность прямой процедуры равного порядка связана с появлением нереальных шахматообразных полей давления. В предлагаемом методе данная трудность устраняется путем установления зависимости поля скоростей и от разности давлений между смежными узловыми точками. Другими словами, хотя узловые скорости не воспринимают разницы между равномерным и шахматообразным полями давлений, это не справедливо по отношению к полю и . Движущей силой для поля и является разность давлений между смежными узловыми точками, и, следовательно, шахматное поле давлений неприемлемо для уравнения неразрывности. Полная система уравнений не допускает каких-либо нереальных колебаний давления в решении. Как было описано ранее в конечно-разностных методах, нереальные поля давления устраняются путем использования смещенных сеток. То, что сделано в предлагаемой процедуре по эффекту очень похоже на применение смещенных сеток. Поле її где-либо в элементе подобно смещенной скорости и уравнения (1.31), (1.32) являются эффективными уравнениями движения для этих смещенных скоростей. Для расчета й, v, Du и Dv в граничных узлах, исходя из массового баланса, используются выражения u =u, v =v, Du = Dv = 0. Рассмотрим граничный узел Р, показанный на рис. 1.10. Предположим, что точка Р является точкой входа потока, то есть масса втекает в область через сечение g-P-a. Если массовый расход очень велик, коэффициент а"ь для точки Р будет очень малым, так как большинство соседних точек располагается вниз по потоку. Это, как было обнаружено, ведет к расходимости итерационной процедуры. Для устранения этой проблемы член -rhp иР добавляется к обеим частям уравнения движения для узлов на втекающей границе. Аналогичный прием необходимо выполнить для получения коэффициент avnb для узлов на втекающей границе. Коррекция давления Уравнений (1.22), (1.23) и (1.33) достаточно для решения любой задачи гидродинамики. Однако, итерационная процедура с прямым использованием этих уравнений, как было найдено, сходится очень медленно. Сходимость была значительно ускорена путем подключенная процедуры коррекции давления, которую мы сейчас рассмотрим подробнее. Она разработана на базе схемы коррекции давления Патанкара [43]. Предположим р - текущее поле давления в итерационном процессе. Пусть поле скорости, полученное путем решения уравнений движения с полем давления р обозначено u , v . Пусть й и v - величины и и v, полученные с использованием поля u , v . Тогда, в общем случае, поле скорости итерационной и разработана на базе алгоритма SIMPLER Патанкара [43], также возможно применение процедур SIMPLE или SIMPLEC. То есть после того, как уравнения движения дискретизированы, дискретный аналог для давления (1.39) решается в первый раз. Используя это новое поле давлений, уравнения движения решаются вновь и получаются новые оценки узловых скоростей. Это сопровождается решением уравнений для поправок давления и расчетом нового приближения уравнения неразрывности, которому удовлетворяет поле и . В этом месте проверяется сходимость и полный цикл повторяется до тех пор, пока не будет достигнута сходимость всех переменных. При дискретизации любого другого конвективно-диффузионного уравнения именно удовлетворяющее уравнению неразрывности поле и используется в уравнении (1.5).
Одномерная модель теплового взаимодействия кориума с наполнителем
В настоящей главе приводится описание методики для расчета теплового взаимодействия кориума вплоть до формирования инверсно стратифицированной структуры ванны расплава в устройстве локализации тяжелой аварии.
Устройство локализации расплава (ловушка) кориума представляет собой бетонную шахту в подреакторном пространстве (см. рис.2.1). Корзина в нижней части шахты, выполненная из чередующихся слоев стали и диоксидциркониевого бетона, охлаждается водой. В корзину помещается крупнопористый наполнитель в виде кусков (блоков) гематита. Предполагается, что через некоторое время после поступления расплавленного кориума в ловушку куски гематита плавятся и расплав перемешивается с оксидными компонентами кориума. В результате плотность смеси оксидных компонент уменьшается и происходит инверсия стратификации расплава: если первоначально расплав металлов (стали и циркония) располагается на поверхности бассейна, то после инверсии расплавленные металлы оказываются внизу, под слоем оксидных компонент кориума. Ранее выполненные расчеты ОКБ «ГИДРОПРЕСС» и «НИТИ» показали, что при обычной стратификации (слой расплавленной стали - наверху) имеется серьезная опасность теплового разрушения стенок конструкции в зоне контакта с расплавленной сталью. Принятое разработчиками техническое решение, использующее крупнопористый наполнитель с относительно низкой плотностью, позволяет избежать неблагоприятного распределения теплового потока на боковых стенках.
В конструкции ловушки предусмотрена подача охлаждающей воды в объем над поверхностью бассейна кориума. Определение времени начала подачи воды относится к числу вопросов, которые должны быть решены расчетным путем с использованием разрабатываемых методик.
В данной работе предполагается следующая схема перемещения кориума в подреакторное пространство: расплавленный кориум в виде однородной смеси оксидных компонент (в основном, U02, Zr02 и Fe203) и не смешивающихся с ними расплавленных металлов (сталь и Zr) поступает в бетонную шахту под реактором и заполняет пространство между блоками крупнопористого наполнителя из гематита (Fe203) в нижней части ловушки. При этом сразу же происходит стратификация расплава, соответствующая плотности компонентов: более тяжелые оксиды располагаются в нижней части объема, а расплав металлов оказывается на поверхности. Получившаяся пористая двухфазная структура со значительным перепадом температуры между кориумом и наполнителем рассматривается как исходное состояние для последующего моделирования теплофизических процессов. В решаемой задаче естественным образом выделяются две стадии процесса: - плавление наполнителя и его растворение в оксидных компонентах кориума, после чего сразу же происходит инверсия расплава: в верхней части бассейна оказываются теперь уже более легкие оксидные компоненты, а вниз перемещается расплав металлов; - естественная конвекция в инверсно стратифицированном бассейне кориума при постепенно уменьшающемся остаточном тепловыделении и интенсивной передаче тепла стенкам шахты. Для расчета нестационарного теплового состояния кориума и определения тепловых потоков в стенки шахты на более теплонапряженной второй стадии процесса необходимо решать сопряженную задачу теплообмена, рассматривая одновременно как конвективный перенос тепла в бассейне кориума, так и теплообмен излучением в объеме над поверхностью бассейна. Названные физические особенности решаемой задачи определяют структуру работы, которая включает три относительно самостоятельных этапа: - моделирование теплового взаимодействия кориума с наполнителем вплоть до формирования инверсно стратифицированного бассейна расплава; - определение нестационарного теплового состояния кориума путем моделирования естественной конвекции в бассейне расплава; - разработка вычислительных моделей для расчета теплообмена излучением в объеме над поверхностью бассейна. Геометрическая схема решаемой задачи представлена на рис.2.1. Объемы, первоначально занимаемые наполнителем из гематита, оксидными компонентами кориума и расплавом металлов, рассчитаны по предварительным данным, предоставленным специалистами проектировщиками : масса гематита MJ=80-100T, его плотность рх = 5.25г/см , масса оксидных компонент кориума М2 «80 + 18 + 10 + 5 = 113т, плотность расплава - р2=7.7т/см , масса расплавленного металла М3 « 50 + 7 = 57 т, плотность его (если считать по стали) - р3 «6.5г/см . Характерный размер кусков гематита считался равным 0.2 м, приведенная интегральная излучательная способность поверхности расплава металлов є = 0.4, условная температура среды и стенок в верхней части подреакторного объема Те = 1500 К. При решении задачи на основании данных проектировщиков (ОКБ «ГИДРОПРЕСС» и «НИТИ») были приняты значения теплофизических параметров, приведенные в табл.2.1.
Расчеты с условием теплового излучения в серой среде
Как показывают выполненные расчеты, общая картина процесса выглядит следующим образом: расплав оксидных компонент кориума, заполняя пространство между кусками гематита, быстро (за две - три минуты) затвердевает и, несмотря на объемное тепловыделение, продолжает остывать еще в течение 20-2 5мин. В это же время происходит плавление гематита, которое завершается через 20мин. после начала процесса. В результате при 20 ? 60мин. расплавленный гематит окружен твердым кориумом, который расплавляется лишь к моменту времени t« 60 мин. - 90 мин.
На этом завершается рассматриваемая стадия процесса. Сразу после плавления оксидных компонент кориума они смешиваются с расплавленным гематитом, имеющим ту же температуру (см. рис. 2.4), но гораздо меньшую плотность. Получившийся однородный расплав имеет меньшую плотность, чем расплав металлов, который опускается на дно. Оставшиеся частично оплавленные куски гематита, находившиеся в слое металла, оказываются на поверхности бассейна и окончательно расплавляются в следующие 5-10 мин. В результате сравнительно быстро происходит переход к совершенно иной структуре среды: в верхней части бассейна расплава находятся расплавленные оксиды (включая и оксид железа - гематит), нижнюю часть шахты заполняет расплав металлов.
Заметим, что представленный одномерный расчет теплового взаимодействия кориума с засыпкой из гематита даже без использования соотношений теплового баланса для части шахты, заполненной оксидами, не требует большого времени вычислений. При 80 интервалах сетки по оси у и 20 интервалах по радиусу куска гематита время счета составляет немногим более одной минуты на PC Pentium-II-266.
Предложена одномерная модель теплового взаимодействия кориума с крупнопористым наполнителем, учитывающая его прогрев и плавление.
Приведен алгоритм решения задачи теплопроводности с переменными распределенными стоками тепла для определения нестационарных профилей температуры компонент кориума и наполнителя. Выполненные расчеты позволили выявить несколько стадий теплового взаимодействия кориума и наполнителя вплоть до формирования инверсно стратифицированного бассейна расплава.
Показано, что расплав оксидных компонент кориума, заполняя пространство между блоками наполнителя, быстро (за две-три минуты) затвердевает и, несмотря на остаточное тепловыделение, продолжает остывать еще в течение 20-25 мин. В это же время происходит плавление наполнителя, которое занимает примерно 20 мин. В результате расплавленный наполнитель оказывается окруженным твердым кориумом, который расплавляется лишь через 60 мин. после начала процесса. Сразу после плавления оксидных компонент кориума они смешиваются с расплавленным гематитом, и происходит инверсия расплава: в верхней части бассейна расплава находятся расплавленные оксиды, а нижнюю часть занимает расплав металлов.
Расчет теплообмена на этой, второй, стадии процесса удобно проводить с использованием отдельной специализированной двумерной модели, которая приведена в главе 3. Рассматриваемая задача, в принципе, вполне аналогична расчету теплового состояния кориума в корпусе реактора и для ее решения могут быть использованы ранее разработанные модели и алгоритмы. Вместе с тем, имеются и существенные отличия: инверсная стратификация кориума, когда слой расплавленного металла находится в нижней части бассейна; - более низкая температура плавления оксидов и возможное отсутствие твердых корок на границах оксидной части бассейна; значительно большая роль условий теплообмена на открытой поверхности бассейна. Названные обстоятельства делают необходимым анализ применимости разработанных ранее методик и их модификацию. В данной главе представлены результаты расчетов по моделированию теплопереноса в инверсно стратифицированном бассейне расплава. В основу теоретической модели положено приближенное асимптотическое описание переноса тепла при турбулентной естественной конвекции в случае больших чисел Рэлея, когда общая задача вырождается и возможно введение переменной эффективной теплопроводности среды. Аналогичная модель эффективной теплопроводности была использована авторами [26,29] ранее для расчета теплового состояния бассейна кориума в корпусе реактора. Тепловое состояние кориума описывается с помощью одномерной и двумерной моделей, которые дают изменяющееся во времени распределение температур инверсно стратифицированного кориума. При этом получаются также и значения плотности теплового потока к стенкам «корзины». Двумерная модель учитывает реальную форму конструкции и сопряженность конвективного теплообмена в ванне расплава с отводом тепла излучением с зеркала кориума в объеме под сводом конструкций. Поскольку граничное условие на поверхности бассейна достаточно "жесткое" и нелинейное, имеет смысл рассмотреть один из двух возможных способов организации решения сопряженной задачи: обычный расчет нестационарного поля температуры в кориуме с регулярным пересчетом граничного условия и специальными приемами, подавляющими пульсации решения вблизи поверхности бассейна; решение уравнения энергии сразу для всей области, включающей бассейн кориума и объем над бассейном.