Теплофизика высоких температур, 2019, T. 57, № 6, стр. 870-881

ВВЕДЕНИЕ

В настоящее время расплавы системы Li–Pb привлекают интерес благодаря их использованию в ядерной энергетике. Свойства этих расплавов изучены довольно подробно [1–3]. Исследованы плотность этих расплавов [4–7], термодинамические свойства [4, 5, 7–12], структура [13–17], диффузионная подвижность [17–19] и широкий спектр физических свойств. Считается, что из-за большой разницы электроотрицательностей в расплавах происходит перенос электронного заряда с атомов лития на атомы свинца, так что связь становится частично ионной. Наиболее сильно этот эффект наблюдается при составе Li₄Pb [16, 20–23].

Опубликовано значительное число работ по моделированию расплавов Li–Pb, в том числе методом ab initio [12, 22]. В первых работах применялись преимущественно методы теории жидкостей, использующие концепцию парного взаимодействия, а также методы молекулярной динамики и Монте-Карло с эффективными парными потенциалами [24], в том числе с учетом кулоновского взаимодействия [16, 25]. Этот подход практически непригоден для описания термодинамических свойств растворов. Более близкие к реальным результаты получены с применением потенциалов погруженного атома (Embedded atom model, ЕАМ) [26–28].

Однако все еще остается много неясного в поведении расплавов Li–Pb. Методами компьютерного моделирования недостаточно исследованы особенности, связанные с наличием вклада ионной связи. Концепция трансферабельности потенциала ЕАМ переносит все особенности взаимодействия компонентов на парный потенциал для пар 1–2 (1 – Li, 2 – Pb). Такая схема может хорошо работать в системах, близких к идеальности, но для систем с сильным взаимодействием компонентов (кулоновским или ковалентным) она может быть недостаточной.

В настоящей работе рассмотрены два варианта описания расплавов литий–свинец: 1) “трансферабельный” вариант и 2) вариант с явным включением кулоновского взаимодействия в схему модели погруженного атома. Для проверки трансферабельного варианта выполнено моделирование растворов Li–Pb с иным выбором потенциалов ЕАМ. Использовались (без изменения) потенциалы ЕАМ для лития [29–31] и свинца [31–33], разработанные специально для описания жидких металлов.

ТРАНСФЕРАБЕЛЬНЫЙ ВАРИАНТ

В этом варианте моделирования методом молекулярной динамики (МД) парные вклады в потенциал (для пар 1–1 и 2–2) и потенциалы погружения ЕАМ взяты такими же, как в случае соответствующего чистого компонента (так называемое свойство транферабельности), и дополнительно вводится только парный потенциал взаимодействия 1–2. Форма и параметры потенциала пар 1–2 подбираются с учетом термодинамических свойств двойной системы – обычно теплоты образования (смешения) и данных по плотности растворов. Концепция транферабельности представляется адекватной в случае систем, близких к идеальным, когда состояние частиц в растворе мало отличается от их состояния в однокомпонентной жидкости. Однако в случае систем с сильным взаимодействием компонентов (например, в растворах металлов IV–V групп с щелочными металлами и в растворах Al в 3d-металлах, когда связь перестает быть чисто металлической) такая трансферабельность может нарушаться.

Растворы Li–Pb являются примером такой системы. Теплота смешения ΔH_f расплавов Li–Pb измерена в [8, 9] (рис. 1). Довольно глубокий минимум изотермы ΔH₁₀₀₀ (около –27.5 кДж/моль) находится вблизи состава Li₄Pb. “Трансферабельный” вариант ЕАМ был применен ранее для моделирования расплавов системы Li–Pb в работах [26–28]. В [27] удалось получить неплохое согласие значений ΔH_f МД-моделей с экспериментальными данными по теплоте смешения [8, 9], но со сдвигом минимума от 80 до 60 ат. % лития (рис. 1). Неплохое согласие получено в [27] также для изотермы мольного объема (рис. 2), причем здесь минимум находится на правильном месте (80% лития [6, 7]). Соответствует данным [6, 7] и коэффициент теплового расширения в области граничных растворов. В МД-работе [28] согласуется с экспериментом теплоемкость раствора Li₁₇Pb₈₃, но минимум изотермы теплоты смешения ΔH_f расплавов сдвинут к 50% лития. Что касается структуры расплавов, то неплохое согласие с дифракционными данными [13, 14] по суммарному структурному фактору расплава Li₁₇Pb₈₃ получено в МД-работе [27] и умеренное согласие с экспериментом [15] в [28].

Для пар 1–2 (т.е. Li–Pb) в настоящей работе применяется потенциал 8–4:

Параметры потенциала ε и r₀ подобраны по экспериментальным данным для теплоты образования растворов Li–Pb при 1000 К [8, 9] (рис. 1). Оптимальные значения параметров для пар 1–2 составляют: r₀ = 2.553 Å и ε = 0.0880 эВ. Совокупность этих потенциалов для системы Li–Pb далее обозначена как ЕАМ-1. Близкие результаты получаются с потенциалом Морзе:

при значениях параметров ε = 0.134 эВ, r₀ = 3.102 Å и α = 4.58. Использование для пар 1–2 потенциала Леннард-Джонса или сплайновых функций типа [27] дают либо завышенные значения модуля всестороннего сжатия и скорости звука в моделях Li–Pb, либо заниженные межчастичные расстояния Li–Pb.

В трансферабельном варианте модели растворов Li–Pb имели размер 2000 или 16 000 атомов в основном кубе и релаксировались с потенциалом ЕАМ-1. Применялся алгоритм Верле в режимах NpT или NVT. Радиус обрыва взаимодействия равен 7.50 Å для пар 1–1 и 9.01 Å для пар 1–2 и 2–2. Радиус сферы Верле равен 10.5 или 17.5 Å.

Важное значение в методе молекулярной динамики имеет выбор шага по времени. В компьютерной системе расчетов используются обычно три основных единицы: длины L = 1 Å, энергии E = 1 эВ и массы m (масса атома, в данном случае – лития). Единица времени является здесь производной величиной и рассчитывается по выражению t₀ = (mL²/E)^1/2. При расчете t₀ по средней массе атома лития (6.94) получается t₀ = 2.682 × 10^–14 с. В случае растворов Li–Pb массы атомов компонентов отличаются более чем в 29 раз, что усложняет расчеты, поскольку обычная длина шага в алгоритме Верле (Δt ≈ 0.01t₀) может оказаться неподходящей для одного из компонентов. Поэтому здесь длина шага выбрана из условия, по которому при прогоне длиной 1000 шагов в режиме NVE (т.е. в режиме изолированной системы) энергия модели менялась не более, чем в пятом–шестом знаке (в пределах 0.004% за 1 пс). При расчете t₀ по средней массе атома лития оптимальная величина Δt для растворов Li–Pb равна 0.05t₀ = 1.341 фс.

Результаты расчетов энергии и теплоты образования ΔU_MD моделей растворов Li–Pb с потенциалом ЕАМ-1 показаны на рис. 1 и в табл. 1. В таблице приведены экспериментальные значения теплоты образования ΔH_f [9] (вторая колонка) и мольного объема растворов Li–Pb [6] (третья колонка). В четвертой колонке приведены значения энергии моделей U_MD с потенциалом ЕАМ-1 при 1000 К и нулевом давлении, в пятой – теплота образования ΔH₁₀₀₀ = ΔU₁₀₀₀ моделей раствора (разность между энергией раствора и аддитивной величиной, рассчитанной по энергии чистых компонентов). При 1000 К электронные поправки к энергии и давлению (см. ниже) невелики и ими можно пренебречь. На рис. 2 показаны экспериментальные и расчетные мольные объемы сплавов Li–Pb при 1000 К. Настоящие данные с потенциалом ЕАМ-1 (линия 3) в целом близки к результатам МД-расчетов [27] и согласуются с экспериментальными данными [6] при концентрациях X₁ ≤ 0.6. При X₁ > 0.6 расхождения с экспериментом увеличиваются. Расчетные объемы чистых компонентов практически совпадают с экспериментальными данными (14.96 см³/моль для Li [34] и 20.39 для Pb [35]).

Расхождение с экспериментом при больших концентрациях лития обусловлено недостаточной точностью приближения трансферабельности. Ни один из упомянутых выше потенциалов ЕАМ не дает реально наблюдаемого минимума изотермы теплоты смешения при составе Li₄Pb [8, 9].

На рис. 3 показана зависимость плотности d расплава Li₁₇Pb₈₃ от температуры по разным данным. Расхождения результатов расчета плотности модели с потенциалом ЕАМ-1 и экспериментальных данных [5] невелики и составляют 0.6–2.6%. Модели расплава Li₁₇Pb₈₃ устойчивы при давлении p ≈ 0 до температуры немного выше 2000 К. Плотность модели с потенциалом работы [28] занижена.

При высокотемпературных расчетах следует учитывать электронный вклад E_eT в энергию и давление, рассчитываемый по методу [30, 33]. Для расплава Li₁₇Pb₈₃ он увеличивается от нуля при 298 К до 3.129 кДж/моль при 2000 К. На рис. 4 показана зависимость энергии этой модели от температуры, которая хорошо описывается выражением

Данные [28] располагаются немного ниже и почти параллельно кривой, полученной по (1). Дифференцированием (1) получается теплоемкость C_p = 3.4194 × 10^–7T ² + 1.1801 × 10^–3T + 27.274. Она немного возрастает с температурой (от 27.95 Дж/(моль К) при 500 К до 31.00 при 2000 К). Эксперимент показывает, что теплоемкость расплава Li₁₇Pb₈₃ слабо убывает от 32.9 Дж/(моль К) при 513 К до 32.8 при 623 К [1, 2]. Эти цифры выше примерно на 2–4 Дж/(моль К), чем полученные МД-данные. В МД-работе [27] теплоемкость монотонно убывает при нагревании от 37 Дж/(моль К) при 600 К до 30 при 1200 К (электронный вклад не учитывался). В [28] значения теплоемкости монотонно убывают от 34.1 Дж/(моль К) при 600 К до 30.2 при 1000 К.

В целом данные МД-моделирования с потенциалами ЕАМ [27, 28] и ЕАМ-1 показывают, что трансферабельный вариант приводит к значительным расхождениям с экспериментом для расплавов, богатых литием, в частности к смещению минимума изотермы теплоты смешения от состава Li₈₀Pb₂₀ к эквиатомному составу. В обычном потенциале ЕАМ не заложено инструментов, позволяющих варьировать положение минимума теплоты смешения двойных систем.

ЕАМ С УЧЕТОМ КУЛОНОВСКОГО ВЗАИМОДЕЙСТВИЯ

В работах [36, 37] структура растворов литий–свинец рассмотрена в различных вариантах теории жидкостей в приближении парного взаимодействия с межчастичными экранированными кулоновскими потенциалами (потенциалами Юкавы) вида

где a_ij – подгоночные коэффициенты, λ – параметр экранирования. В работах [16, 25] применены среднесферическое и гиперсетевое приближения для расчета свойств расплава Li₄Pb. Ионные оболочки описаны в [16] жесткосферным потенциалом с диаметром 2.65 Å для каждого компонента. Наилучшее согласие с дифракционными данными по структуре получено при эффективных зарядах ионов лития и свинца 0.5 и –2.0 соответственно и λ = 1.10. МД-моделирование расплава Li₁₇Pb₈₃ с этими парными потенциалами выполнено в [24]. Включение кулоновского взаимодействия в рассмотрение позволяет улучшить в целом согласие между расчетными структурными характеристиками растворов и дифракционными данными. По всей вероятности, учет кулоновского взаимодействия был бы полезным и при использовании потенциалов ЕАМ. Однако включения кулоновского взаимодействия при моделировании в схеме ЕАМ до сих пор не проводилось.

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

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

2) Избыточные заряды ионов компонентов q₁ и q₂ удовлетворяют условию электронейтральности: q₁X₁ + q₂X₂ = 0 (X_i – атомные доли компонентов).

3) Экранированный электрический потенциал i-го иона имеет вид (q_i/r) exp(–λr), является короткодействующим, и поэтому при расчетах кулоновского взаимодействия не требуется применять процедуру Эвальда.

4) При переходе части электронов данного атома в состояние экранированного заряда q_i соответствующая доля электронов изымается из эффективной электронной плотности схемы ЕАМ ψ_i(r) = p_1iexp(–p_2ir) или добавляется к ней, т.е. p_1i = $p_{{{\text{1}}i}}^{0}$(1 – ${{{{q}_{i}}} \mathord{\left/ {\vphantom {{{{q}_{i}}} {q_{i}^{0}}}} \right. \kern-0em} {q_{i}^{0}}}$), где $q_{i}^{0}$ = 1 для лития и 4 для свинца. Впрочем, этот эффект очень слаб и его можно не учитывать.

5) Равновесные значения зарядов ионов определяются конкуренцией между кулоновской энергией раствора E_Coul (E_Coul < 0) и энергией перезарядки ионов E_Trans (E_Trans > 0). Изменение энергии при перезарядке i-го иона ΔE_i может быть задано функцией, которая при q_i = 1 принимает значение первого потенциала ионизации I_i, а при q_i = –1 равна энергии сродства к электрону со знаком минус (–Q_i). Кроме того, ΔE_i = 0 при q_i = 0. Можно, например, выбрать функцию ΔE_i в виде

где k, l – натуральные числа. Первые потенциалы ионизации Li и Pb равны 5.39 и 7.42 эВ, а энергии сродства соответственно 0.168 и 0.364 эВ.

6) Эффективные заряды ионов q_i рассчитываются как величины, минимизирующие полную энергию модели при заданных давлении и температуре. Оптимальные заряды $q_{i}^{*}$ получаются путем коррекции минимизирующих зарядов для лучшего согласия плотности модели с экспериментальными данными.

7) Экранированный кулоновский потенциал и кулоновские межчастичные силы добавляются к потенциалам и силам, рассчитываемым обычным образом на основе ЕАМ-1. Добавка E_Trans на межчастичные силы не влияет.

Далее такой гибридный потенциал обозначен как ЕАМ-2.

Расчет энергии перезарядки E_Trans выходит за рамки метода МД и является предметом теории металлов. Чтобы оставаться в пределах метода компьютерного моделирования, следует параметризовать вклад E_Trans, выразив его через значения зарядов ионов. Оптимальное выражение для E_Trans получено (с учетом п. 2) сравнением энергии моделей и реальных растворов в виде

где n₁ – число ионов лития. Это выражение пригодно при 0 < q₁ < 0.3 и –1 < q₂ < 0.

В варианте с включением кулоновского взаимодействия модели растворов Li–Pb имели размер 2000 или 16000 атомов в основном кубе и релаксировались с гибридным потенциалом ЕАМ-2, равным сумме потенциала ЕАМ-1 и экранированного кулоновского потенциала (2), в котором a_ij = 14.399$q_{i}^{*}q_{j}^{*}$ (14.399 – коэффициент пересчета из компьютерной системы единиц в СГС). Оптимальный параметр экранирования оказался равным 1.10, т.е. таким же, как значение, найденное в [16]. В режимах NpT или NVT применялся алгоритм Верле. Радиус обрыва взаимодействия равен 7.50 Å для пар 1–1 и 9.01 Å для пар 1–2 и 2–2. Радиус сферы Верле составляет 17.5 Å.

Для каждой исследованной концентрации растворов Li–Pb выполнена релаксация моделей с различными допустимыми эффективными зарядами 0 < $q_{{\text{1}}}^{*}$ < 1 и –1 < $q_{{\text{2}}}^{*}$ < 0 (с учетом п. 2). Минимизирующие заряды найдены из условия минимума энергии модели при температуре 1000 К и нулевом давлении, а оптимальные заряды определены с учетом реального значения плотности раствора [6]. Минимизирующие и оптимальные значения зарядов приведены в табл. 2. Они значительно меньше найденных в [16] (0.5 и –2.0), однако значения в [16] приводят к неприемлемым расхождениям с экспериментом по плотности расплавов. При концентрациях лития менее 50% учет кулоновского взаимодействия не требуется, так как минимизирующие заряды близки к нулю.

Релаксация моделей проведена с включением зарядов $q_{i}^{*}.$ Полученные значения расчетной теплоты смешения показаны на рис. 5 и в табл. 1. Включение кулоновского взаимодействия позволяет значительно улучшить согласие расчетной теплоты смешения с экспериментальными данными (колонки 2 и 11). На рис. 2 и в табл. 1 (колонки 3 и 13) показаны значения мольного объема растворов при 1000 К. Отсюда также видно, что включение кулоновского взаимодействия значительно улучшает согласие с экспериментом [6].

Несколько хуже выглядят результаты для производных свойств. С использованием потенциала ЕАМ-2 рассчитаны теплоемкости C_V и C_p, а также изотермический модуль всестороннего сжатия K_T моделей при 1000 К. Результаты приведены в табл. 3. Прямые измерения теплоемкостей и модуля, видимо, не проводились. Расчетные значения C_p в [27] (метод МД) располагаются выше настоящих данных примерно на 5–6 Дж/(моль К). Полученная в [28] методом МД величина 30.3 Дж/(моль К) для раствора Li₁₇Pb₈₃ при 1000 К заметно выше полученного здесь значения 25.8.

В табл. 3 приведены расчетные значения скорости звука в моделях при 1000 К в сравнении с экспериментальными данными [6, 38]. МД-данные выше экспериментальных при 30–60% Li в среднем на ∼10%, но убывают с приближением к чистым компонентам. Расчеты с потенциалом Морзе дают примерно то же самое, а потенциал Леннард-Джонса приводит к еще более высоким скоростям звука. Данные [6], возможно, занижены. Неплохо согласуются с экспериментом [6] значения коэффициента теплового расширения α. Производная давления (∂p/∂T)_V зависит от концентрации с положительными отклонениями от аддитивности, а коэффициент Грюнайзена γ = = (V/C_V)(∂p/∂T)_V меняется почти аддитивно.

Структура. Структура расплавов исследована в [13, 14] методом рассеяния нейтронов в интервале температур 995–1225 К с использованием изотопа ⁷Li. На рис. 6 показаны парциальные парные корреляционные функции (ППКФ) модели расплава Li₁₇Pb₈₃ при 600 К. Координаты пиков ППКФ равны 5.05, 2.90 и 3.23 Å для пар 1–1, 1–2 и 2–2 соответственно. Первый пик для пар 1–2 с высотой 2.63 ниже пика для пар 2–2, но расположен левее. При температуре 1000 К координаты пиков ППКФ равны 4.92, 2.87 и 3.19 Å для пар 1–1, 1–2 и 2–2 при таких же соотношениях высот пиков (1.26, 2.24 и 2.59). Это указывает на слабую тенденцию к преимущественному окружению атомов частицами противоположного сорта.

Расплавы Li–Pb моделировались также в [27, 28]. В [28] высота пика ППКФ для пар 1–2 модели Li₁₇Pb₈₃ при 600 К больше, чем в настоящей работе: 3.38 вместо 2.63, хотя координаты пиков в основном совпадают (2.93 Å) (рис. 6). Высоты пиков ППКФ моделей Li–Pb в [27] почти не зависят от состава и равны 2.2–2.4 при 1000 К, т.е. близки к полученным в настоящей работе. На рис. 7 показаны ППКФ для пар 1–2 в модели расплава Li₈₀Pb₂₀ при 1000 К. Учет кулоновского вклада в ЕАМ-2 приводит к небольшому росту высоты первого пика ППКФ g₁₂(r) при 1000 К: от 2.68 при q₁ = 0 до 2.75 при q₁ = 0.13. Высота пика ППКФ пар 1–2 при 1000 К в модели Li₈₀Pb₂₀, полученная в [27], равна 2.18 (рис. 7). Эти различия обусловлены различной жесткостью потенциала пар 1–2: в [28] он более жесткий, а в [27] менее жесткий, чем в ЕАМ-2.

На рис. 8 показаны зависимости от концентрации координат r₁ и высот g₁(r) первых пиков ППКФ для пар 1–2 в рассматриваемых моделях. Координата растет с увеличением концентрации свинца, а высота пика проходит через максимум при 60–80 ат. % Li, т.е. в интервале концентраций с наибольшим кулоновским вкладом.

В работах [13, 14] рассчитана функция

полученная путем анализа дифракции нейтронов, которая является фурье-образом структурного фактора Бхатиа–Торнтона S_CC(r)/(X₁X₂).

На рис. 9 приведена функция 4πr²n₀ρ_CC(r) (n₀ – число частиц в единице объема) для модели расплава Li₄Pb при 1073 К в сравнении с данными [13, 14]. Согласие можно считать удовлетворительным только для первого минимума. Осциллирующая форма зависимости означает гетерокоординацию частиц, слоистость ближнего порядка вокруг частиц расплава и может являться признаком наличия вклада ионной связи, при которой чередуются слои частиц с зарядами разного знака. В МД-модели данной работы гетерокоординация быстро затухает с расстоянием из-за сильного экранирования кулоновского взаимодействия. В [27] наличие гетерокоординации вообще не отмечено.

На рис. 10 показана функция q[S(q) – 1], где q и S(q) – вектор рассеяния и структурный фактор рассеяния нейтронов на расплаве ⁷Li₁₇Pb₈₃ при 773 К в сравнении с дифракционными данными [13]. В настоящих расчетах амплитуды рассеяния на частицах ⁷Li и Pb приняты равными –1.90 и 7.60. Согласие с данными [13] здесь очень хорошее, за исключением небольшого смещения 2–4 пиков вправо.

Диффузия. Коэффициенты самодиффузии $D_{i}^{0}$ в равновесных компьютерных моделях с потенциалом ЕАМ-2 рассчитывались по уравнению 〈Δr²〉 = 6Dt, где 〈Δr²〉 – средний квадрат смещения частиц за время t. В табл. 4 приведены коэффициенты самодиффузии компонентов 10⁵D в см²/с при температурах до 2000 К. С увеличением концентрации лития скорость диффузии возрастает. На рис. 11, 12 показаны зависимости D(T) в логарифмических координатах для лития и свинца. Коэффициенты для прямых линий lg D = a₁lg T + a₂ имеют следующие значения:

Спрямление этих зависимостей означает, что температурная зависимость коэффициентов самодиффузии хорошо описывается выражением D = сT ^b. Такая форма зависимости выполняется вдоль изобар многих металлов [31]. Коэффициенты самодиффузии лития в расплаве Li₁₇Pb₈₃, полученные в [27], в 2–4 раза ниже приведенных в табл. 4. Коэффициенты, рассчитанные в [28], близки к полученным результатам при температурах выше 800 К, а при меньших температурах занижены.

Отметим, что графики расчетных коэффициентов самодиффузии из табл. 4 не спрямляются в координатах ln D–1/T и имеют вид парабол с положительной кривизной. Для приближенных оценок аппроксимация D = D₀exp(–E_act/RT) дает следующие величины:

Свойства расплавов Li–Pb в условиях ударного сжатия. Потенциал ЕАМ-2 позволяет рассчитывать свойства расплавов при высоких температурах и давлениях. Для чистых компонентов Li и Pb такие расчеты выполнены в [29–33] при температурах до 14 000 К (Li) и 25 000 К (Pb) и давлениях до 100 ГПа (Li) и 280 ГПа (Pb). Проведем аналогичные расчеты для раствора Li₁₇Pb₈₃. С этой целью следует найти зависимости энергии и давления раствора от объема и температуры до степеней сжатия Y = V/V₀ = 0.4–0.5, где V – мольный объем в сжатом состоянии, а V₀ – нормальный объем в исходном состоянии. При обычных условиях сплав со средним составом Li₁₇Pb₈₃ включает в себя две фазы: LiPb и почти чистый свинец. Для упрощения расчетов в качестве исходного состояния при 298 К принята аморфная фаза состава Li₁₇Pb₈₃. Эта модель имеет плотность 10.096 г/см³ (V₀ = 17.152 см³/моль, аддитивный объем – 20.29) и энергию –186.79 кДж/моль (аддитивная энергия –181.5). Фактическая плотность сплава Li₁₇Pb₈₃ при 298 К равна 10.22 г/см³ [40]. При расчетах электронных вкладов в энергию E_el и давление p_el применена методика из [29–31]. Электронные вклады в энергию рассчитаны при средней электронной плотности 3.49 эл/ат. Они приведены в табл. 5. Электронные вклады в давление рассчитаны по формуле p_elV = (2/3)E_el.

Далее построены модели некристаллического сплава Li₁₇Pb₈₃ при различных температурах и объемах с использованием потенциала ЕАМ-2. Заряды ионов считались не зависящими от температуры. Значения энергии и давления этих моделей с добавлением электронных поправок приведены в табл. 6 и 7, с помощью которых можно рассчитать ударную адиабату (рис. 13) методом, описанным в [29–33]. До давления 100 ГПа она почти совпадает с адиабатой чистого свинца, а при более высоких давлениях уходит вверх. В верхней точке адиабаты сплава температура достигает 21 400 К. На адиабате свинца при той же степени сжатия (Y = 0.55) температура составляет ∼14 200 К.

ОБСУЖДЕНИЕ

Выбором потенциалов ЕАМ для пар 1–1 и 2–2 и потенциала “8–4” для пар 1–2, а также включением экранированного кулоновского взаимодействия компонентов удалось получить хорошее согласие ряда термодинамических свойств моделей системы Li–Pb (плотности, энергии, коэффициента теплового расширения, модуля всестороннего сжатия, скорости звука) с опытными данными при температурах до 1000 К. Получено также согласие с дифракционными данными о структуре этих расплавов при 600–1073 К. Дополнительного улучшения согласия с экспериментом по упругим свойствам расплавов можно достичь, немного уменьшая степень жесткости парного потенциала для пар 1–2, но не допуская сближения частиц в этих парах ближе 1.9–2.0 Å. Применение потенциала Леннард-Джонса для пар 1–2 способствует росту модуля всестороннего сжатия моделей, а мягкий сплайновый потенциал пар 1–2 приводит к слишком тесному сближению частиц и к расхождению с дифракционными данными о структуре.

Структурные характеристики расплава Li₄Pb рассчитаны в [16] на основе средне-сферического и гиперцепного приближений теории жидкостей, разработанных для систем с парным потенциалом. В расчетах [16] использованы дифракционные данные [13, 14]. Первый максимум ППКФ пар 1–2 в [16] при 1075 К с высотой около 2.6 расположен при r = 2.6 Å, а в случае 1073 К (рис. 7) он имеет высоту 2.75 и находится при r = 2.70 Å. Для пар 1–1 получено согласие и по положению, и по высоте пика. Однако для пар 2–2 расхождения гораздо больше: r = 4.8 Å в [16] и 5.00 Å (данная работа), а высоты этих пиков составляют соответственно 1.67 и 1.30. Эти расхождения указывают на недостаточную адекватность приближения парного взаимодействия при описании систем с сильным взаимодействием компонентов.

Значительные расхождения получены и с результатами МД-расчетов [27]. Так, при 1000 К первый пик ППКФ пар 1–2 расплава Li₄Pb в [27] имеет высоту 2.00 и расположен при r = 3.00 Å, в данной работе: высота – 2.85 и r = 2.75 Å. Второй пик этой ППКФ в [27] расположен при r = 5.46 Å (высота – 1.45) и r = 5.57 Å (высота – 1.17) соответственно. Эта разница обусловлена различием потенциалов.

Ударная адиабата сплава Li₁₇Pb₈₃ очень близка к адиабате чистого свинца при давлениях до 100 ГПа, а при бóльших давлениях идет выше последней. Следовательно, нельзя рассчитать адиабату сплава, используя только адиабаты чистых компонентов, без расчетов межчастичных потенциалов. Для проверки адекватности потенциала ЕАМ-2 полезно построить экспериментально адиабату Гюгонио реального сплава Li₁₇Pb₈₃ или любого другого.

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

ЗАКЛЮЧЕНИЕ

Построены две серии моделей жидких растворов Li–Pb: 1) в трансферабельном варианте с потенциалами ЕАМ и 2) с дополнительным включением экранированного кулоновского потенциала. Первый вариант хорошо описывает свойства растворов при содержании лития меньше 60 ат. %. Второй позволяет получить правильное положение минимумов на изотермах теплоты смешения и плотности, а также достичь согласия с опытом по совокупности термодинамических свойств и по структуре расплавов. Рассчитаны свойства раствора Li₁₇Pb₈₃ при температурах до 25 000 К и давлениях до 470 ГПа, и построена адиабата Гюгонио. Вариант с включением кулоновского взаимодействия расширяет набор используемых потенциалов и позволяет моделировать растворы с сильным взаимодействием компонентов. При использовании потенциалов ЕАМ без учета переноса заряда с менее электроотрицательного компонента на более электроотрицательный вряд ли можно добиться хорошего согласия расчетных свойств моделей растворов типа Li–Pb с экспериментальными данными.