Численные модели кода евклид/v2 для анализа аварий с разрушением активной зоны быстрых реакторов

Бутов А.А., Жданов В.С., Климонов И.А., Кудашов И.Г., Мосунова Н.А., Стрижов В.Ф., Усов Э.В., Чухно В.И. (ИБРАЭ РАН, Москва, Новосибирск, Россия)

Введение

Для обоснования безопасности перспективных установок с жидкометаллическим теплоносителем требуется система расчетных кодов, позволяющая моделировать тяжелые запроектные аварий с плавлением активной зоны. Для этих целей в ИБРАЭ РАН на базе ЕВКЛИД/V1 разработан интегральный код ЕВКЛИД/V2, который должен позволить расширить спектр режимов работы, включая аварии с разрушением активной зоны с оксидным и нитридным топливом. Текущая версия кода ЕВКЛИД/V2 позволяет моделировать поведение реакторных установок на быстрых нейтронах с жидкометаллическими теплоносителями в стационарных и переходных режимах работы, а также при различных проектных авариях, путем проведения связанных нейтронно-физических, термомеханических, теплогидравлических расчетов, расчетов переноса ПД по контуру, поведения аэрозолей в помещениях АЭС, расчетов термического разрушения твэлов, ТВС и активной зоны.

При моделировании тяжелой аварии в составе кода ЕВКЛИД/V2 взаимодействуют твэльный модуль БЕРКУТ, тяжелоаварийный модуль SAFR/V1 и теплогидравлический модуль HYDRA-IBRAE/LM, нейтронно-физический модуль DN3D, модуль переноса ПД по контуру AEROSOL/LM. Расчетные сетки каждого из модулей согласованы друг с другом, что позволяет реализовать механизм обмена между модулями, в рамках интерфейсных возможностей кода ЕВКЛИД/V2. Ниже представлено краткое описание основных расчетных моделей, которые могут быть использованы для расчета тяжелых аварий в реакторных установках на быстрых нейтронах.

1. Модели теплогидравлических процессов в ТВС и а.з. для моделирования тяжелых аварий

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

Для расчетов теплогидравлических процессов в составе кода ЕВКЛИД/V2 используется аттестованный в Ростехнадзоре код HYDRA-IBRAE/LM [1]. Моделирование течения жидкого металла в одно- и двухфазных режимах проводится с использованием двухжидкостной модели с замыкающими соотношениями, зависящими от режима течения. Базовая система уравнений представлена ниже.

Уравнения сохранения массы:

Подробное описание карты режимов и замыкающих отношений можно найти в работах [2, 3].

В процессе расчета HYDRA-IBRAE/LM передает информацию о значениях коэффициентов теплообмена и температуре теплоносителя твэльному модулю БЕРКУТ [5], от которого получает информацию о величине теплового потока в теплоноситель. Распределение температуры в твэле модуль БЕРКУТ определяет путем решения уравнения теплопроводности на основе данных по энерговыделению в активной зоне, получаемых из модуля DN3D [4].

2. Модели термического разрушения твэла

1.1. Базовые модели термического разрушения твэла. 

Расчет тяжелой аварии невозможен без моделирования плавления отдельного твэла. Поведение твэла до момента плавления рассчитывается с использованием твэльного модуля БЕРКУТ [5]. После достижения температуры плавления происходит передача расчета тяжелоаварийный модуль SAFR/V1 [6,7].

Тяжелоаварийный модуль SAFR/V1 решает задачу разрушения одиночного твэла вследствие его плавления. Модуль состоит из двух основных блоков: блока решения тепловой задачи при плавлении и блока решения задачи перемещения расплава. Первый блок предназначен для решения уравнения теплопроводности при фазовом переходе. Уравнение записывается в следующем виде:

Граничные условия могут быть трех типов: первого, второго и третьего рода. При работе SAFR/V1 отдельно от кода ЕВКЛИД/V2 граничное условие третьего рода определяется путем задания температуры окружающей среды и коэффициента теплоотдачи. При работе в составе ЕВКЛИД/V2 температура теплоносителя и коэффициент теплоотдачи передается из теплогидравлического кода HYDRAIBRAE/LM. Возможно задание условия по радиационному теплообмену.

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

Информация о количестве образовавшейся массы передается из теплового блока модуля SAFR/V1. В процессе расчета учитывается трение расплава с газовым потоком, с поверхностью твэла, наличие контактного угла на границе раздела расплав – твердое тело, теплообмен расплава со стенкой и с газовым потоком. Пример расчета эксперимента с плавлением твэла, выполненный на стенде НГТУ им. Алексеева [8, 9], представлен на рис. 1.

Рис. 1. Расчет эксперимента с плавлением твэла.

Более детальное описание подходов для расчета термического разрушения твэла, включающих описание численной реализации базовой системы уравнений и способов расчета интенсивности теплообмена и трения расплава с поверхностью твэла  и с теплоносителем, представлены в работах авторов [6,7].

1.2. Модели диссоциации нитридного топлива 

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

Известно, что при относительно низких давлениях азота над поверхностью нитрида урана, конгруэтного плавления последнего не наблюдается [10]. Происходит разложение на уран и азот. Конгруэтное плавление начинается при температуре 3123 К при давлении азота над поверхностью выше 2,5·105 Па [11].

Для расчета диссоциации СНУП-топлива рассматривается общая реакция вида:

Предполагается, что диссоциация происходит на граничной поверхности, для которой должно быть поставлено следующее условие (рис. 2):

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

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

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

Пример расчета эксперимента [12] по диссоциации таблетки из нитрида урана представлен на рис. 3

Рис. 3. Расчет эксперимента по диссоциации нитрида урана

3.Модели движения компонент разрушенных частей твэла, ТВС и всей а. з.

Разрушение твэла приводит к выбросу разрушенных частей твэла в межтвэльное пространство их оседание на поверхность чехла и других твэлов, что приводит к дальнейшему развитию аварийной ситуации. Для моделирования подобных процессов на каждом расчетном шаге модуль SAFR/V1 передает в код HYDRAIBRAE/LM данные о температуре поверхности твэла, данные о скорости и толщине пленки расплава, а также информацию о доле расплава, которая уносится с поверхности твэла в поток теплоносителя. На основе полученных данных HYDRAIBRAE/LM рассчитывает изменение проходного сечения гидравлического канала и изменение скорости теплоносителя. Далее модуль HYDRA-IBRAE/LM возвращает в SAFR/V1 данные о коэффициентах теплообмена и трения между стенкой твэла и потоком теплоносителя, скорости и температуре теплоносителя, а также данные о доле расплава в теплоносителе.

В результате уноса расплава с поверхности твэла, часть расплава остается в потоке теплоносителя. Скорость уноса расплава определяется соотношениями из работ [13,14]. Дальнейшее перемещение расплава в потоке теплоносителя и теплообмен расплава с теплоносителем моделируется с помощью уравнений многокомпонентной модели, реализованной в модуле HYDRA-IBRAE/LM. В результате движения расплава в потоке теплоносителя, часть расплава может осаждаться на твэлах, пэлах и чехлах ТВС, приводя к их дальнейшему разрушению.

Для расчета перемещения компонент стандартная двухжидкостная модель кода HYDRA-IBRAE/LM была дополнена уравнениями сохранения массы, энергии и импульса, а также уравнение для отдельных компонент:

Для объемных долей компонент справедливо равенство:

Также была расширена карта режимов, которая стала включать в себя отдельные компоненты: жидкую сталь и жидкое топливо [15].

В каждой расчетной ячейке кода HYDRA-IBRAE/LM происходит пересчет площади проходного сечения канала, площади теплообмена, смоченного периметра гидравлического диаметра, шага решетки твэлов на основе информации, полученной из модуля разрушения твэла SAFR/V1. По новой геометрии также пересчитываются замыкающие соотношения.

Так площадь проходного сечения канала на новом шаге 1 nA  рассчитывается по формуле:

4.Опция для расчета процессов разрушения а. з. в трехмерном приближении

Недостатком одномерных моделей является их относительно невысокая точность и высокая доля используемых эмпирических соотношений, используемых для моделирования гидродинамических и тепловых процессов. Для увеличения точности расчета и снижения эмпиризма общей тенденцией становится использование трехмерных моделей для расчета процессов в реакторной установке. Также развиваются подходы для анализа тяжелых аварий в быстрых реакторах в трехмерном приближении с использованием расчетных кодов. Подобная опция есть, например, в коде COREMELT [16] и коде SIMMER-IV [17].

Для того чтобы учесть подобные тенденции в расчетном коде ЕВКЛИД/V2 и в его теплогидравлическом модуле HYDRA-IBRAE/LM [1 была реализована трехмерная многокомпонентная модель для расчета переноса и теплообмена компонент разрушенного твэла и ТВС в активной зоне и верхней камере смешения в реакторной установке с жидкометаллическим теплоносителем. Ниже представлено описание созданной модели и примеры тестовых расчетов, выполненных на ее базе.

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

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

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

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

Полная масса стекающих материалов рассчитывается по уравнению:

Полная энтальпия стекающих материалов удовлетворяет уравнению:

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

где k b – средняя степень выгорания k-го материала.

Средняя температура и состояние материалов в бассейне рассчитываются по следующему алгоритму.

Рассчитывается полная энтальпия бассейна расплава: 

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

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

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

Отсюда масса летучего продукта деления i-го типа, выходящего в единицу времени:

Заключение

Представлено краткое описание моделей разрушения активной зоны интегрального кода ЕВКЛИД/V2, которые могут быть использованы для расчета тяжелых аварий в реакторных установках на быстрых нейтронах с жидкометаллическим охлаждением. Представлены результаты расчета некоторых экспериментов по разрушению одиночных твэлов, демонстрирующих работоспособность моделей.

Список литературы

1. Алипченков В.М., Анфимов А.М., Афремов Д.А. и др. Базовые положения, текущее состояние разработки и перспективы дальнейшего развития теплогидравлического расчетного кода нового поколения HYDRA-IBRAE/LM для моделирования реакторных установок на быстрых нейтронах // Теплоэнергетика. 2016. №
2. С. 54– 64. 2. Усов Э.В., Прибатурин Н.А., Кудашов И.Г и др. Один из этапов верификации теплогидравлического кода HYDRA-IBRAE/LM/V1 для расчета течения натриевого теплоносителя в сборках твэлов // Атомная энергия. 2015. т.118. № 6. С. 309–313.
3. Усов Э.В., А. А. Бутов, Г. А. Дугаров и др. Система замыкающих соотношений двухжидкостной модели кода HYDRA-IBRAE/LM/V1 для расчета процессов при кипении натрия в каналах энергетического оборудования // Теплоэнергетика. 2017. № 7. С. 48–55
4. Koltashev D.A., Stakhanova A.A. Neutronic calculation of fast reactors by the EUCLID/V1 integrated code // IOP Conf. Series: Journal of Physics: Conf. Series, 2017, v. 781, 012003.
5. Veprev D.P., Boldyrev A.V., Chernov S.Y., Mosunova N.A.. Development And Validation Of The Berkut Fuel Rod Module Of The Euclid/V1 Integrated Computer Code // Annals of Nuclear Energy, 2017, v. 113, p. 237–245.
6. Усов Э.В., Бутов А.А., Чухно В.И. и др. Моделирование процессов плавления твэла и затвердевания расплава, образующегося при термическом разрушении твэла быстрого реактора, с помощью модуля SAFR/V1 интегрального кода ЕВКЛИД/V2 // Атомная энергия. 2018. т.124. № 3. С. 123–126.
7. Усов Э.В., Бутов А.А., Чухно В.И. и др. Моделирование перемещения расплава по поверхности твэла быстрого реактора при тяжелой аварии с помощью модуля SAFR/V1 интегрального кода ЕВКЛИД/V2 // Атомная энергия. 2018. т.124. № 4. С. 197–200
8. Усынин Г.Б., Аношкин Ю.И., Семенычев М.А. Исследование плавления твэлов на имитаторах с топливными композициями // Атомная энергия. 1991. т. 70. вып. 2. С. 108–110.
9. Усынин Г.Б., Аношкин Ю.И., Власичев Г.Н. и др Модельное изучение процессов, возникающих при перегреве твэлов // Атомная энергия. 1986. т. 61. вып. 5.  С. 347–350.
10. Кириллов П.Л., Терентьева М.И., Денискина Н.Б. Теплофизические свойства материалов ядерной техники. М: ИздАт. 2007. 200 с.
11. Olson W.M., Mulford R.N.R. The Decomposition Pressure and Melting Point of Uranium Mononitride // The Journal of Physical Chemistry. 1963. V. 67. No. 4. P. 952−954.
12. Baranov V.G., Tenishev A.V., Kuzmin R.S., et. all. Thermal stability investigation technique for uranium nitride // Ann. Nucl. Energy. 2014. V. 87(2). P. 784.
13. Ishii M. and Mishima K. Droplet Entrainment Correlation in Annular Two-Phase Flow // International Journal of Heat and Mass Transfer. 1989. Vol. 32. P.1835–1846
14. I. Kataoka, M. Ishii. Entrainment and Deposition Rates of Droplets in Annular Two-phase Flow. Proceedings of the ASME/JSME Thermal Engineering Joint Conference. v. 1. 1983.
15. А. А. Бутов, В. С. Жданов, И. А. Климонов и др. Физические модели для расчета разрушений твэла и активной зоны реактора с жидкометаллическим теплоносителем, реализованные в коде ЕВКЛИД/V2 // Теплоэнергетика. 2019. № 5. с. 5–15
16. Ю.М. Ашурко, А.В. Волков, К.Ф. Раскач. Разработка программных модулей для расчета запроектных аварий в быстрых реакторах с учетом пространственновременной кинетики // Атомная энергия, 2013, т. 114, № 2. С. 63–67.
17. Yamano H., Tobita Y., Fujita S. A three-dimensional neutronics-thermohydraulics simulation of core disruptive accident in sodium-cooled fast reactor // Nucl. Eng. Des., 2009, v. 239, i. 9, p. 1673–1681.