Моделирование плавления нанокристаллов диоксида урана

В большинстве современных ядерных ректоров в качестве топлива используется диоксид урана. Из таблеток диоксида урана формируются топливные сердечники тепловыделяющих элементов (твэл) [1]. Во время эксплуатации таблетки подвергаются воздействию высоких температур и давления, облучению, механическим нагрузкам. При номинальной мощности реактора температура на оси твэла может превышать 2000 К, а давление газообразных продуктов деления урана внутри герметичного твэла может достигать 80-100 атм [2]. Поскольку диоксид урана обладает малой теплопроводностью, то в топливном сердечнике возникает высокий температурный градиент [3]. Всё это может привести к изменению микроструктуры топлива, его распуханию, рекристаллизации, спеканию зерен [2].

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

1. ОПИСАНИЕ МЕТОДА

Моделирование проводилось с использованием программного комплекса DL_POLY [4].

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

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

U eff (r ij , T) =

z i ( T ) z j ( T ) e ²

+

Г-ij

+ f ( T ) (b_i+b_j) exp

ai + a j

'll

( b+b, J

cicj

6 , (1) ^rij

где первое слагаемое соответствует кулоновскому взаимодействию, а второе и третье — потенциалу Борна-Майера [5]. Значения независимых от температуры параметров потенциала были взяты из работы [6].

Для восстановления параметров потенциала z_i(Т) и f(T) были использованы экспериментальные данные по тепловому расширению решётки UO₂ и изменению энтальпии [7]. Были получены зависимости дробного заряда иона кислорода и параметра f от температуры в виде:

z ( T ) = z o 1 +a— T ° , где a =^ L T o J I

^

f (T ) = f o 1 + ^—T , где P =

9,61 - lo 3 T < T _o

2,9 3 - lo 2 , T > T _o

| 9,3 - 10 ^{е 2}, T <  T o 2,69(1 , T >  T o ’

, z 0 = 1,193, T o = 2666 K

f , = 0,354, T o = 2627 K

Таблица 1. Моделируемы кристаллы

Размер, эл ячеек	4x4x4	5x5x5	6x6x6	7x7x7	8x8x8
Размер, нм	2,188	2,735	3,282	3,829	4,376
Количество атомов	768	1500	2592	4116	6114

Рис. 1. Зависимость параметра f и дробного заряда иона кислорода от температуры и кусочно-линейная аппроксимация с переходом вблизи точки 2670 K

На рис. 1 можно видеть, что найденные по этой методике значения параметров дробного заряда и f(T) имеют две характерные линейные области с переходом при температуре вблизи 2670 K – температуре перехода в суперионное состояние [7]. Более подробно вопросы обоснования температурной зависимости указанных параметров потенциала рассмотрены в наших предыдущих работах [8, 9].

• 2.    ОПРЕДЕЛЕНИЯ ТОЧКИ ПЛАВЛЕНИЯ

• 2.1    Радиальная функция распределения

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

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

Вид этой функции является графическим представлением упорядоченности расположения частиц в кристаллической решетке [10]. Математически РФР описывает вероятность нахождения частиц на расстоянии R от заданной точки [11]. Вычисляется РФР по формуле:

g ( R ) =

V N ( R , A R ) N 4 л К ² A R

где – N число частиц моделируемого кристалла N ( R, A R ) - число частиц находящихся в сферическом слое шириной A R •

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

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

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

На рис. 2. показаны РФР для нанокристалла 5х5х5 при различных температурах. Хорошо видно, что при температуре 1500 К второй пик становится едва различим, а при 1800 совсем исчезает, т.к. кристалла находится в жидком состоянии. Отсюда следует, что температура плавления может находиться между двумя этими температурами.

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

• 2.2    Метод структурного фактора

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

^

S _a ( k )

n

N ^ ^ A² X ^f co ^{s( k} • r n ) I , n = 1 )

(        ^ ^ A²

+ IX f z ^sin( k • r n ) I

V n = 1                 )

, (4)

где индекс a – соответствует атомам урана или кислорода, N – количество атомов в системе, f_a – атомный фактор рассеяния, вектор k – вектор обратной решетки, вектор r_n – вектор выбранного направления.

При 0 K кристалл представляет собой идеальную структуру, структурный фактор которой равен единице. Напротив, при высоких температурах для расплавленного состояния квадрат модуля структурного фактора стремиться к нулю, а его численные оценки колеблются в диапазоне 0,01-0,04. Величина структурного фактора определяет интенсивность отражения при рентгеновской дифракции на кристаллах для определенных направлений. При этом известно, что логарифмическая зависимость интенсивностей дифракционных максимумов представляет собой линейную функцию от температуры. Теоретически эта зависимость получается в приближении Эйнштейна, когда колебания атомов считаются независимыми осцилляторами при температурах выше комнатной [12]. Следовательно, в методе МД при неизменном характере колебаний атомов зависимость логарифма модуля структурного фактора от температуры также должна быть линейной.

На рис. 3 в логарифмическом масштабе показана зависимость квадрата модуля структурного фактора, усредненного по трем независимым направлениям (001), (010) и (100) для бесконечного кристалла. Усреднение проводилось с целью уменьшения амплитуды разброса точек и упрощения анализа результатов.

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

т, к

Рис. 3. Температурные зависимости квадрата модуля структурного фактора, усредненного по направлениям (001), (010) и (100), отдельно для решеток урана и кислорода с линейной аппроксимацией между точками фазовых переходов в области 2400 K и 3200K.

на зависимостей существенно меняется. Линейный характер зависимостей сохраняется до температуры плавления кристалла 3200K, дальнейшие изменения связаны с переходом в расплавленное состояние. При температуре 3600K решетка атомов урана разрушается, и структурный фактор скачком меняется от значений 0,5 до 0. Данный факт никак не может быть связан с процессом плавления, поскольку разрушение решетки происходит за несколько десятков шагов интегрирования. Аппроксимируя линейные участки методом наименьших квадратов и находя соответствующие точки пересечения, были получены значения для температур фазовых переходов 2471±400 К и

Рис. 4. Температурные зависимости квадрата модуля структурного фактора, усредненного по направлениям (001), (010) и (100) для кубического нанокристалла размером 5х5х5 элементарных ячеек

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

s * £ s	Температура плавление кислородной подрешетки, К		Температура плавления урановой подрешетки, К
	Метод структурного фактора	Метод энтальпии	Метод ст руктурного фактора	Метод энтальпии
4x4x4	700±25	700±25	1125±25	1100±25
5x5x5	930±260	1125±25	1500±300	1525±25
6x6x6	1250±25	1200±25	2000±25	2025±25
7x7x7	1500±50	1450±50	2200±50	2300±50
8x8x8	1500±50	1500±50	2400±50	2400±50

3262±280 К, которые согласуются с рекомендациями МАГАТЭ 2670 К и 3120 К соответственно [7]. Указанные погрешности расчетов вычислены на основе оценок точности аппроксимации по методу наименьших квадратов.

Используя описанный подход по оценке изменения структурного фактора с ростом температуры, можно судить о фазовых переходах и для нанокристаллов. Для примера рассмотрим нанокристалл размером 5х5х5. Результаты изменения структурного фактора представлены на рис.4. Полученная зависимость имеет характерные линейные (в логарифмических координатах) области. Температуры фазовых переходов определяли как точки пересечения аппроксимационных прямых. Полученные значения составили 930±260 K и 1500±300 K для подрешеток кислорода и урана соответственно.

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

• 2.3 . Зависимости энтальпии системы от температуры

На рис. 5 представлены зависимости приведённой энтальпии для кубических нанокристаллов различных размеров. Качественно эти зависимости имеют одинаковый характер для кристаллов различных размеров: в начале происходит плавный рост, затем наблюдаются два скачка энтальпии. Первый скачок на величину примерно 0,5 эВ/ат. связан с переходом кристалла в суперионное состояние, второй более резкий на величину примерно 1,5 эВ/ат. – с плавлением кристалла. Поскольку повышение температуры проводилось ступенчатым образом, а время моделирования на каждой температурной точке было ограниченно, то в качестве критических выбирались точки, которые лежат на середине отрезка резкого роста энтальпии. Рассчитанные таким образом температуры приведены в табл.2. Согласно полученным результатам обе эти температуры существенно снижаются с уменьшением размера кристалла. Так, температура плавления наименьшего из рассмотренных кристаллов (размером 2,2 нм) оказывается почти втрое ниже температуры плавления макрокристалла диоксида урана. С увеличением размера обе температуры быстро возрастают.

В литературе нет экспериментальных данных о температуре плавления наноразмерных кристаллов UO2. Однако для металлов известно, что с уменьшением размера кристалла до десятков нанометров его температура плавления резко уменьшается [13]. Таким образом, результаты моделирования показывают, что и для диоксида урана справедливо такое уменьшение.

ЗАКЛЮЧЕНИЕ

В настоящей работе путем молекулярно-динамического моделирования проведено исследование фазовых превращений диоксида урана. Для определения критических температуры было использовано несколько методов: 1) анализ вида функции радиального распределения; 2) анализ изменения структурного фактора рассеяния; 3) анализ зависимости энтальпии от температуры.

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

Рис. 5. Зависимость приведённой на атом энтальпии от температуры для кубических нанокристаллов различных размеров. Размеры указаны в элементарных ячейках

лов в диапазоне размеров 2,2 – 4,4 нм установлено снижение температур суперионного перехода и плавления с уменьшение размера.

Работа выполнена при поддержке Минобрнауки в рамках государственного задания на 20122014 гг., федеральной целевой программы “Научные и научно-педагогические кадры инновационной России на 2009-2013 годы” и “Исследования и разработки по приоритетным направлениям развития научно-технологического комплекса России на 2007-2012 годы”, а также при частичной поддержке гранта РФФИ: проект № “12-08-97076-р_поволжье_а”