Моделирование и анализ виртуального растяжения эластичного и пластичного кристаллов дигалогенфенолов

Квантово-механическое моделирование с учетом периодических граничных условий, основанное на теории функционала плотности (DFT) [1, 2], представляет собой эффективное средство для предсказания и характеристики свойств кристаллических структур как в атомистическом масштабе, так и на уровне свойств электронной плотности. Детальное химическое и физическое понимание, обеспечиваемое первопринципным теоретическим моделированием, достигается с помощью квантово-химических расчетов структуры и свойств кристаллов, которые уже можно относить к материалам. Вычислительный подход является стандартным инструментом для экспериментальной проверки, интерпретации и прогнозирования во многих областях исследований [3, 4]. В настоящее время существует потребность в новых методах прогнозирования механических свойств синтезируемых кристаллических материалов для фармацевтической промышленности. К этому располагает как накопленный практический опыт, так и подходы, основанные на варьировании заместителей в биологически активных соединениях для получения новых более эффективных и менее токсичных лекарств. Однако это может повлиять на механические свойства кристаллических материалов, что важно учитывать при разработке технологических этапов таблетирования и измельчения [5].

Имеющиеся компьютерные исследования в основном полагаются на изучение влияния свойств добавок путем моделирования методом Монте-Карло [6–8]. Часть исследований посвящена исследованию физических свойств и фазовых переходов методами классической молекулярной динамики [9], при анализе работ, посвященных деформациям в кристаллах, можно заметить, что большинство работ посвящено исследованию деформаций под нагрузкой [10–12], в то время как работ по растягивающим деформациям крайне мало [13]. Мотивацией нашего исследования является заполнение пробела в области знания о том, что происходит с кристаллами на уровне укладки отдельных молекул и химических связей, и как при этом изменяются их свойства. Полученные результаты позволят в будущем создавать материалы с заранее заданными свойствами, что имеет большое значение для их применения в медицине и технике.

Существенной особенностью кристаллических тел является их структурная анизотропия, которая приводит к зависимости ряда свойств от направления в кристалле. Кроме того, известно, что многие свойства кристаллов определяются их структурной организацией, которая, в свою очередь, может включать дефекты. Поэтому большой интерес представляют знания о деформации кристаллов и отклик структуры на механические воздействия, а именно, изучение нековалентных взаимодействий и эффектов, вызываемых растягивающими деформациями, что позволит прогнозировать влияние структурных дефектов в кристаллах на свойства состоящих из них материалов. Определенный интерес вызывает поведение галогенных связей [14, 15] и их сравнение с водородными связями в контексте структурной химии. Это важно для развития инженерии и проектировании материалов [16], поскольку совместное присутствие, конкурентное и синергетическое поведение водородных и галогенных связей оказывает сильное влияние на свойство кристалла, подвергающийся деформации.

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

Методы расчетов

Моделирование всех рассматриваемых кристаллических структур было проведено методом Кона-Шэма (DFT) с периодическими граничными условиями при использовании пакета программ CRYSTAL17 (v.1.0.2) [17, 18]. Начальные данные (пространственная группа симметрии, параметры элементарной ячейки и координаты атомов) взяты из Кембриджской кристаллографической базы данных [19]. Для расчётов был использован скорректированный метод Хартри – Фока HF-3c [20]. Этот метод включает в себя базисный набор MINIX [21] и три полуклассиче-ские поправки: дисперсионная поправка Гримме D3 [22], атомное попарное приближение коррекции противовеса gCP [23, 24] и ближняя поправка базиса SRB. Для каждой структуры предварительно проводилась полная оптимизация параметров ячейки и координат атомов (опция FULLOPT) с сохранением пространственной группы, то есть поиск минимума на поверхности потенциальной энергии с полной релаксацией кристаллической структуры. При этом использовались следующие параметры сходимости: по энергии в SCF-цикле 10–10 aт. ед., по градиенту 10–5 aт. ед. и по смещению 3·10–5 aт. ед.

Виртуальный тест растяжения кристаллических структур дигалогенфенолов проводился согласно алгоритму, предложенному в работе [13]. Для фиксации определённой деформации кристалла по осям a, b и c, по-отдельности для каждой из них, применялись опции CRYDEF и FIXDEF программы CRYSTAL17, с сохранением пространственной группы и всех операторов симметрии. Эти опции позволяют проводить оптимизацию координат атомов и параметров кристаллической ячейки, когда один из этих параметров сохраняется постоянным во время оптимизации. Чтобы смоделировать постоянную деформацию растяжения, мы взяли предварительно полностью оптимизированную, равновесную в заданных условиях расчёта, кристаллическую структуру, увеличили один из параметров ячейки на 3 % и зафиксировали его при оптимизации новой, растянутой структуры. После того, как остальные параметры ячейки и координаты атомов были вновь оптимизированы, величина этого параметра увеличивалась ещё на 3 % относительно начального значения (в равновесной структуре), фиксировалась, и далее проводилась новая оптимизация. Такую пошаговую оптимизацию с растяжением кристалла по одной из осей мы повторяли десять раз, до 30 % удлинения кристалла по выбранной оси.

В аналогичных условиях был проведён повторный тест растяжения с уменьшенным до 1 % шагом увеличения каждого параметра ячейки. Это было сделано для выяснения, зависят ли тренды изменения структуры и её рассматриваемых свойств от величины шага, то есть от «плавности» или «равномерности» последовательного увеличения растягивающей деформации. Следует отметить, что в силу симметрии для кристаллических структур 3,4-дихлорфенола и 4-бром-3-хлорфенола, кристаллизующихся в тетрагональной пространственной группе, растяжение проводилось одновременно по осям а и b . Фактически, кристалл растягивался изотропно в любом направлении, параллельном плоскости ab . Далее, после выполнения расчётов, были построены зависимости изменения объёма и полной энергии кристаллографической ячейки относительно величин исходной оптимизированной структуры от процента увеличения каждого параметра ячейки и проанализированы изменения конфигурации связей и структурных особенностей полученных кристаллических форм.

Результаты и обсуждение

Нами были проанализированы результаты моделирования структуры и растягивающих деформаций пластичного 3,4-дихлорфенола (I 4 1 /a), эластичного 4-бром-3-хлорфенола (I 4 1 /a) (рис. 1). Эти кристаллы ранее были исследованы экспериментально в работе [25] на предмет сопоставления их механических свойств с особенностями межмолекулярных взаимодействий.

Рис. 1. Расположение молекул, входящих в состав кристаллов: а – 3,4-дихлорфенола, б – 4-бром-3-хлорфенола

В кристалле сформированы стопки из синтонов, сформированных четверками молекул, связанных в слое водородными связями OI...H4, между кислородом гидроксильной группы одной молекулы и водородом соседней молекулы. Следует сказать, что эти водородные связи повторяются вдоль оси c , и четыре связи образуют виток спирали в стопке высотой в 3,691 Å. При одновременном растяжении по кристаллографическим осям a и b одна из молекул приподнимается и встраивается между молекулами соседнего слоя. При этой деформации синтоны не разрушаются, но молекулы в стопках немного поворачиваются относительно друг друга.

Интересна смена типа нековалентного взаимодействия галоген-галоген: с типа I - ван-дер-ваальсового взаимодействия, при котором а-дырка одного хлора ориентирована с небольшим смещением на σ-дырку другого хлора, на тип II – галогенную связь, при которой σ-дырка ориентируется на неподеленную пару электронов (рис. 2). При 15 % растяжении имеется взаимодействие Сl1…Cl1 с длиной 3,553 Å, соответствующее типу I, которое при растяжении уменьшает свою длину, хотя связь Cl1…Cl1 все ещё остается I типа. При 18% растяжении происходит «переключение» взаимодействий, при котором вместо одного контакта Cl1…Cl1 возникают две галогенные связи Cl1…Cl2, причем они становятся более длинными, чем предшествующее взаимодействие типа I.

Рис. 2. «Переключение» в кристалле 3,4-дихлорфенола нековалентной связи Cl…Cl с типа I (ван-дер-ваальсового взаимодействия) на тип II (галогенную связь)

При растяжении по оси c происходит заметное укорочение контакта Cl1…Cl1 с 3,257 Å в равновесной форме до 3,236 Å при 18 % растяжении с дальнейшим ростом до 3,250 Å при максимальном растяжении. При этом расстояние между соседними слоями, которое можно оценить по расстоянию между атомами Cl1 и H2, равному 3,045 Å в равновесной структуре, уменьшается при 30 % растяжении – указанное расстояние составляет 2,972 Å. Это свидетельствует о сближении молекул по оси с за счет их разворота в одной из стопок. Расстояние между атомами Cl2 и H2, наоборот, растет: с 3,139 до 3,216 Å соответственно. Таким образом, при растяжении по оси с изменение межслоевых расстояний проявляется заметно, но при этом не появляется никаких полостей или трещин.

При рассмотрении графиков на рис. 3а–б видно, что изменение энергии ячейки ΔE перестает резко расти и выходит на небольшое плато, начиная с 15 % растяжения по осям а и b . Это наблюдение можно объяснить, во-первых, «переключением» типа нековалентной связи Cl1…Cl1, описанным выше, во-вторых, появлением зазоров в двух ортогональных направлениях. При рассмотрении графика растяжения кристалла по оси с (рис. 3в-г), выхода на плато не наблюдается.

Из сравнения графиков на рис. 3а-г, демонстрирующих зависимости объёма и изменения полной энергии кристаллографической ячейки 3,4-дихлорфенола от процента удлинения каждого параметра ячейки, следует, что выбор шага не повлиял на поведение кристаллической структуры при ее растяжении. Изменения энергии и объема ячейки происходят симбатно в обоих режимах растяжения. При растяжении с шагом в 1 и 3 % относительно исходной структуры не появляется значимых различий как в достигнутом объеме кристаллографической ячейки, так и в увеличении энергии ячейки. При сравнении нековалентных взаимодействий и межатомных расстояний, получаемых на одинаковом проценте растяжения, можно заметить, что самая большая разница наблюдается для взаимодействия Cl1…Cl1 на 12 % растяжения, и составляет всего 0,022 Å. В этой области мы видим не только изменение энергии, но и небольшие различия в относительных ориентациях молекул. При дальнейших шагах растяжения такое расхождение нивелируется. Следует также отметить, что в режиме растяжения с малым шагом в 1 % изменение энергии ячейки при одинаковых растяжениях оказывается немного меньше, чем при растяжении с большим шагом в 3 %, что видно на рис. 4.

а)

в)

г)

Рис. 3. Зависимости объёма и изменения полной энергии кристаллической ячейки относительно исходной оптимизированной    структуры    от   процента   удлинения    каждого    параметра   ячейки для 3,4-дихлорфенола: а – при растяжении вдоль осей a и b с шагом в 3 %; б – при растяжении вдоль осей a и b с шагом в 1 %; в – при растяжении вдоль оси c с шагом в 3 %; г – при растяжении вдоль оси c с шагом в 1 %

Рис. 4. Разница в изменении энергии кристаллической ячейки при виртуальном растяжении с разным шагом 1 и 3 % кристалла 3,4-дихлорфенола

Эластичный кристалл 4-бром-3-хлорфенола имеет сходную структуру кристаллических мотивов, состоящих из синтонов, сформированных четырьмя ближайшими молекулами, связанными водородными связями (рис. 5). Эти четверки водородно-связанных молекул имеют тенденцию формировать стопки вдоль оси с , которые сохраняются при растягивающих деформациях. При одновременном растяжении вдоль осей а и b с шагом 3 %, при 12 % удлинении появляются кластеры, состоящие из отделяющихся друг от друга стопок. Этот процесс сопровождается началом формирования полостей, возникающих при удлинении межатомных расстояний. Следует заметить, что галогенные связи Br1…Br1 при этом утрачиваются, позволяя сформироваться пустотам между стопками.

Рис. 5. «Переключение» в кристалле 4-бром-3-хлорфенола нековалентной связи Cl…Cl с типа I (ван-дер-ваальсового взаимодействия) на тип II (галогенную связь Cl…Br)

Из-за появляющихся полостей в кристаллической упаковке возникает возможность существенного синхронного поворота всех четырех молекул в стопке. В результате, на 15 % растяжении происходит смена типа структурообразующих нековалентных взаимодействий, реализуемых посредством галогенов: взаимодействие типа I Cl…Cl сменяется парными галогенными связями Сl…Br и Br…Cl, формирующимися между соседними фенольными циклами.

При растяжении кристалла 4-бром-3-хлорфенола по оси с формирования полостей не наблюдается. Поведение нековалентных взаимодействий между молекулами сильно отличается: с растяжением межстопочное взаимодействие Cl1…Cl1 становится короче, при этом галогенная связь Br1…Br1 не теряется, но становится длиннее.

При рассмотрении графика (рис. 6а), видно, что изменение энергии ячейки ΔE перестает резко расти и осциллирует, начиная с 12 % растяжения по осям а и b . Это явление сопровождает тенденцию к возникновению полостей. При рассмотрении графика растяжения кристалла по оси с (рис. 6б), выхода на плато функции изменения энергии ячейки ΔE не наблюдается. Таким образом, основное отличие в поведении объема V и изменения энергии ячейки ΔE для пластичного 3,4-дихлорфенола и эластичного 4-бром-3-хлорфенола состоят в том, что в случае пластичного кристалла его структура, начиная с определенного критического уровня растяжения, сравнительно долго остается мало изменяющейся, а эластичный кристалл продолжает оказывать сопротивление деформации, что вызывает осцилляции энергии, но при большом растяжении структура в нем все же перестраивается окончательно.

Следует отметить, что оба кристалла ведут себя практически одинаково при среднем и высоком растяжении, отличаясь только скоростью роста энергии ячейки и некоторыми особенностями нековалентных взаимодействий. Если сравнивать графики изменения энергии ячейки, то видно, что область с выходом на «плато» у пластичного кристалла меньше (см. рис. 3а), и в дальнейшем происходит рост ΔE, в то время, как у эластичного кристалла, во-первых, наблюдаются сравнительно длительные осцилляции ΔE, во-вторых, на начальной стадии выше скорость роста ΔE (см. рис. 6а). Это можно заметить, исходя из линейной аппроксимации начального участка (до 12 %) зависимости роста энергии ячейки от растяжения: коэффициент пропорциональности линейной зависимости ΔE от процента растяжения для 3,4-дихлорфенола ниже, чем для 4-бром-3-хлорфенола. Мы видим более сильное относительное увеличение энергии ячейки именно

а)

б)

Рис. 6. Зависимости объёма и изменения полной энергии кристаллографической ячейки относительно исходной оптимизированной структуры от процента удлинения каждого параметра ячейки 4-бром-3-хлорфенола: а – при растяжении вдоль осей a и b с шагом в 3 %; б – при растяжении вдоль оси с с шагом в 3 % для эластичного кристалла, что предполагает для него более вероятную релаксацию структуры до исходного состояния при снятии растягивающего напряжения. Более пологий тренд изменения энергии ячейки пластичного кристалла при растяжении, наоборот, означает вероятность сохранения структуры кристалла при снятии напряжения в некотором деформированном состоянии или её перехода в другую равновесную структуру, отличающуюся от начальной. Рассматривая же идентичное поведение этих кристаллов на протяжении большей части диапазона растяжения, мы можем утверждать о достижении эластичным кристаллом предела эластической деформации примерно на уровне 10 % растяжения. Предполагается, что при дальнейшем его растяжении вероятно его поведение уже как пластичной структуры.

Выводы

В результате применения теста на виртуальное растяжение, представляющего собой моделирование равновесных, но деформированных структур галогенфенолов, были даны объяснения пластичного и эластичного поведения этих кристаллических структур. В растянутой кристаллической структуре пластичного 3,4-дихлорфенола формируются сравнительно небольшие зазоры между стопками. Стопки, в свою очередь, с каждым новым шагом растяжения получают возможность разворачиваться друг относительно друга, сохраняя при этом синтоны с четырьмя водородными связями и «переключая» менее прочные межстопочные взаимодействия галогенов с типа I - ван-дер-ваальсова взаимодействия на более прочные - галогенные связи. С появлением зазоров между стопками изменение энергии кристаллической ячейки замедляется с выходом на короткое «плато», затем продолжает расти. Объем ячейки показывает монотонный рост во всем диапазоне растягивающих деформаций.

Изменение энергии и объема кристаллической ячейки при ее растяжении для эластичного кристалла 4-бром-3-хлорфенола происходит сходным образом, но лишь до определённого предела. Затем, с появлением значительных по размеру полостей между стопками в двух ортогональных направлениях, изменение энергии ячейки начинает осциллировать, а потом продолжает расти. Область замедления роста энергии в случае эластичного кристалла 4-бром-3-хлорфенола значительно больше, чем для пластичного 3,4-дихлорфенола, и она сопровождается осцилляциями. Это наблюдение можно интерпретировать как стремление эластичной структуры сохранять возможность возврата к предыдущему состоянию без разрушения.

Показано, что изменения энергии и объема ячейки кристаллов галогенфенолов при моделировании растягивающих деформаций мало зависит от режима растяжения - шага, изменяющегося в пределах от 1 до 3 %.

Таким образом, виртуальный тест на растяжение показал себя подходящим инструментом для объяснения механических свойств молекулярных кристаллов, в том числе, при сравнении эластичных и пластичных форм изоструктурных кристаллов, каковыми являются 3,4-дихлорфенол и 4-бром-3-хлорфенол.

Работа выполнена по гранту РФФИ №20-03-00240.