Свойства стекла, как и любого другого вещества, напрямую зависят от его структуры, поэтому, зная структуру, мы можем прогнозировать свойства. Экспериментально структуру стекла определяют на основании спектрографического анализа (также существуют другие методы анализа) полученного образца и знания структуры соответствующего кристалла. Многие методы анализа являются очень трудоемкими, многие дорогостоящими, в любом случае создание принципиально нового образца ? это большие затраты.
Компьютерное моделирование значительно упрощает решение данной задачи. В частности метод молекулярной динамики позволяет рассчитывать траектории движения частиц моделируемого вещества. Знание траекторий каждой частицы в модели объекта позволяет рассчитать практически любое свойство системы - как термодинамическое (например, энергию, давление, энтропию), так и кинетическое (коэффициенты диффузии, частоты колебания атомов). Размер модели определяется быстродействием и оперативной памятью компьютера, что в настоящее время не является большой проблемой.
Основой для создания программного комплекса была выбрана программа, позволяющая выполнить основной расчет на основе метода молекулярной динамики, - Moldy [1]. Данная программа предоставляет возможность обрабатывать любую компоновку твердых многоатомных молекул, атомов или ионов, а также любую их смесь. Для вычисления короткодействующих сил используется метод связной сетки и метод сумм Эвальда для обработки дальнодействующих электростатических сил. Моделирование может быть выполнено в различных ансамблях. Поскольку МД ячейка не обязательно может быть кубической, программа подходит как для моделирования твердого вещества, так и для моделирования жидкостей.
Большинство существующих программ имеет ограничения в их применимости, например, поддерживает только одну функцию потенциала, или требуется молекулярная симметрия, или число молекул ограничено. Moldy (насколько это возможно) не имеет таких ограничений. Система определяется перед запуском, и ее размер может быть ограничен только количеством памяти доступной программе: для обработки большей системы требуется больше памяти. В программе заложены различные функции потенциалов межчастичного взаимодействия.
На рис. 1 представлена схема разработанного программного комплекса для моделирования стекол методом молекулярной динамики. На вход системы подаются исходные данные, например, количество атомов каждого типа в расчетной ячейке, плотность системы для NVT ансамбля или давление системы для NPT ансамбля, шаг по времени и количество шагов по времени для алгоритма метода молекулярной динамики, начальная температура системы, конечная температура системы, шаг по температуре, возможно задание начальной структуры и др. Также указываются характеристики, которые требуется рассчитать (если сохранить данные основного расчета, то характеристики, не указанные ранее, можно рассчитать в любой момент).
Рис. 1. Схема программного комплекса для моделирования германатных стекол
Полученные исходные данные обрабатываются для начала проведения основного расчета методом молекулярной динамики. Основной расчет представляет собой цикл расчетов системы с заданными характеристиками в программе Moldy при одной температуре, расчет повторяется столько раз, сколько необходимо для достижения конечной температуры. Все данные основного расчета сохраняются для дальнейшей обработки в модулях системы (см. схему рис. 1).
С помощью данного комплекса произведено моделирование стеклообразного оксида германия GeO2, нескольких составов свинцово-германатных стекол - (x)PbO-(1-x)GeO2 (с содержанием оксида свинца от 20 до 62,5%) и щелочно-германатных стекол - (x)Na2O-(1-x)GeO2 (с содержанием оксида натрия от 15 до 50%). Для всех систем сначала были получены модельные образцы стекла, для чего было произведено моделирование охлаждения системы от 6000 до 300 К. Количество шагов каждого МД прогона составляло 200000, шаг по времени 0,001 пс. В начальной конфигурации не задавалось начального расположения атомов, они были вычислены автоматически.
Функции радиального распределения являются одной из важнейших структурных характеристик стекол. На рис. 2 представлены результаты работы модуля получения функций радиального распределения для щелочно-германатных стекол (x)Na2O-(1-x)GeO2 (с содержанием оксида натрия от 0 до 50%). Первый пик полученных функций радиального распределения рассмотрен в табл. 1. В сравнении с литературными данными, среди которых присутствуют как расчетные (полученные в результате моделирования), так и экспериментальные значения, отличие не превышает 5%, что позволяет говорить об адекватности полученных результатов.
Рис. 2. Функции радиального распределения для щелочно-германатных стекол составов SGXX, где XX -% содержания оксида натрия
Модуль расчета координационных чисел позволяет получить координационные числа любого типа атомов системы, причем существует возможность получить как средние координационные числа, так и процентное содержание атомов определенной координации. В табл. 2 представлены расчетные координационные числа германия для щелочно-германатных стекол и данные из литературных источников. Литературные данные подтверждают адекватность результатов - ошибка минимальна.
Таблица 1
Расстояния между парами атомов для щелочно-германатных стекол составов SGXX,
где XX -% содержания оксида натрия
|
Пара ионов |
Состав стекла |
Первое координационное расстояние, Å |
|
|
Расчетные |
Литературные данные |
||
|
O-O |
SG00 |
2,69 |
2,80 [4], 2,83-2,85 [5], 2,75[3] |
|
|
SG15 |
2,61 |
2,73[3] |
|
|
SG20 |
2,66 |
2,80 (SG18)[2] |
|
|
SG30 |
2,60 |
2,74[3] |
|
O-Na |
SG15 |
2,49 |
2,25[3] |
|
|
SG20 |
2,53 |
2,40 (SG18) [2] |
|
|
SG30 |
2,47 |
2,17[3] |
|
O-Ge |
SG00 |
1,68 |
1,60-1,80[2], 1,60 [4], 1,70-1,74[5], 1,72[3] |
|
|
SG15 |
1,69 |
1,73[3] |
|
|
SG20 |
1,69 |
1,60-1,80 (SG18) [2] |
|
|
SG30 |
1,69 |
1,73[3] |
|
Na-Na |
SG15 |
3,33 |
3,20[3] |
|
|
SG20 |
3,07 |
|
|
|
SG30 |
3,20 |
3,20[3] |
|
Na-Ge |
SG15 |
3,30 |
3,30[3] |
|
|
SG20 |
3,40 |
3,40 (SG18) [2] |
|
|
SG30 |
3,40 |
3,36[3] |
|
Ge-Ge |
SG00 |
3,31 |
3,10-3,20 [4], 3,15-3,45 [5], 3,26[3] |
|
|
SG15 |
3,34 |
3,24[3] |
|
|
SG20 |
3,36 |
3,00-3,20 (SG18) [2] |
|
|
SG30 |
3,38 |
3,32[3] |
Таблица 2
Координационные числа атомов германия для щелочно-германатных стекол составов SGXX, где XX -% содержания оксида натрия
|
Состав стекла |
Координационные числа атомов германия |
Среднее знач. к.ч. |
Литературные данные |
|||
|
3 коорд |
4 коорд |
5 коорд |
6 коорд |
|||
|
SG00 |
1,50% |
95,00% |
3,50% |
0,00% |
4,02 |
3,95[2], 4,03[3] |
|
SG15 |
0,30% |
63,50% |
36,20% |
0,00% |
4,36 |
4,36[3] |
|
SG20 |
0,00% |
56,90% |
42,80% |
0,30% |
4,44 |
4,45-4,63 (SG18) [2] |
|
SG30 |
1,80% |
66,40% |
31,80% |
0,00% |
4,30 |
4,27[3] |
|
SG40 |
1,70% |
67,10% |
31,20% |
0,00% |
4,30 |
- |
|
SG50 |
3,50% |
72,00% |
24,50% |
0,00% |
4,21 |
- |
Наличие 5-координированного германия в германатных стеклах также подтверждено последними экспериментальными исследованиями [6] (ранее считалось, что в стеклах присутствует только 4-х и 6-координированный германий).
Так как в моделировании методом молекулярной динамики имитируются реальные движения атомов, он позволяет рассчитывать динамические характеристики материала. Примером одной из самых важных таких характеристик является колебательный спектр F(ω). Для его получения вычисляется автокорреляционная функция скоростей γ(t) и затем производится интегрирование в заданном интервале частот. Вибрационная плотность состояний (VDOS) представляет собой зависимость квадрата интенсивности колебательного спектра от частоты.
На рис. 3 представлены расчетная вибрационная плотность состояний для модели стекла GeO2 при 300 К (модуль расчета колебательного спектра) и экспериментальный спектр Рамана стекла GeO2 [7]. Наблюдается хорошее согласие по проявляющимся частотам с экспериментальным спектром Рамана. Частоты в районе 600 см-1 наблюдаются как в спектре Рамана и ИК [7], так и расчетном спектре. Из этого следует, что построенная модель стекла близка к оригинальному стеклу GeO2 по колебательным характеристикам.
Рис. 3 Вибрационная плотность состояний модели стекла GeO2 (SG00) при 300 К и экспериментальный спектр Рамана стекла GeO2 [7] (жирная линия)
Наиболее распространенным методом исследования фазовых превращений и определения температуры стеклования является дилатометрия - наблюдение за изменением объема аморфного вещества при изменении температуры. Для определения температуры стеклования было проведено моделирование нагревания образцов свинцово-германатных стекол, полученных ранее. Шаг по температуре составлял 50К, для обработки результатов было задействовано 2 модуля комплекса - модуль получения физических характеристик системы и собственно модуль определения температуры стеклования, который основан на дилатометрии. На рис. 4 представлены расчетная и экспериментальные зависимости температуры стеклования от содержания оксида свинца в свинцово-германатном стекле. Как видно из приведенных данных, даже при шаге по температуре в 50 К зависимости практически совпадают, из чего делаем вывод, что данный программный комплекс подходит для определения температуры стеклования.
Рис. 4. Зависимость температуры стеклования моделей свинцово-германатных стекол от состава и экспериментальная зависимость [8]
Таким образом, на основе метода молекулярной динамики и программы Moldy был разработан программный комплекс для моделирования стекол, позволяющий определять основные структурные и физические свойства. Исходя из анализа данных, полученных в результате расчета и литературных источников, можно говорить об их адекватности и, следовательно, о применимости данного комплекса для моделирования других оксидных стекол.
Данный расчет показал возможность применения метода молекулярной динамики в составе данного программного комплекса для изучения структуры и свойств стекла. При этом моделирование дает больше представлений о системе, например, визуальное представление (модуль визуализации). Например, можно производить моделирование для новых разрабатываемых стекол. Моделирование поможет сократить материальные затраты на некоторые виды дорогостоящих методов исследования стекла, поскольку для новых стекол необходимо проводить множество экспериментов. Моделирование также поможет скорректировать начальные условия для проведения экспериментов, что также сокращает затраты ресурсов.
Работа выполнена в рамках государственного контракта № 11.519.11.4004.
Рецензенты:
-
Куркина Е.С., д.ф-м.н., в.н.с. факультета ВМК МГУ им. М.В. Ломоносова, г. Москва;
-
Носов Г.А., д.т.н., профессор, зав. кафедрой «Процессы и аппараты химической технологии» имени Н.И. Гельперина Московского государственного университета тонких химических технологий им. М.В. Ло- моносова, г. Москва.
Работа поступила в редакцию 12.03.2012.