УДК 544.77.022.54
Силовые постоянные двумерного коллоидного кристалла с квадратной статической решеткой
Ю.В.Александров, П.Е.Дышловенко Ульяновский государственный технический университет
Определены силовые постоянные двумерного коллоидного кристалла с квадратной решеткой в широком диапазоне плотностей. Электростатическое взаимодействие в кристалле описано в рамках модели уравнения Пуассона-Больцмана. Учет свойств симметрии системы позволяет значительно уменьшить объем вычислений и заранее определить структуру матриц силовых постоянных. Показано, что эффективное взаимодействие в рассматриваемом кристалле исчерпывается взаимодействием ближайших соседей только первого и второго порядков.
Ключевые слова: коллоидные кристаллы, силовые постоянные, уравнение Пуассона-Больцмана.
Коллоидные кристаллы - это коллоидные суспензии, в которых, несмотря на наличие жидкой среды, имеет место регулярное расположение коллоидных частиц в пространстве. В этом случае частицы твердой фазы расположены в узлах кристаллической решетки того или иного типа [1]. В технологическом отношении коллоидные кристаллы вызывают большой интерес как основа для создания оптических фильтров, оптических переключателей, а также материалов с фотонной запрещенной зоной [2, 3].
Важным классом коллоидных систем, обладающих способностью образовывать коллоидные кристаллы, являются электрически стабилизированные системы, в которых частицы твердой фазы (макроионы) электрически заряжены [1]. Одно из направлений исследования электрически стабилизированных коллоидных кристаллов - изучение их упругих и структурных свойств, в частности силовых постоянных. Следует отметить, что коллоидные кристаллы описанного типа, в отличие от обычных кристаллов, есть системы с начальным напряжением. Теория упругости таких систем применительно к коллоидам в настоящее время не разработана.
Один из возможных подходов к описанию упругих и структурных свойств электрически стабилизированных коллоидных кристаллов заключается в использовании модели на основе нелинейного уравнения Пуассона-Больцмана [1, 4], которое учитывает силы электростатического взаимодействия и гидростатическое давление ионов в растворе. Другие возможные типы взаимодействий, такие как ван-дер-ваальсовы, сольватационные, энтропийные и прочие [4], проявляют себя на сравнительно коротких расстояниях, не превышающих нескольких сотен ангстрем, и в данной модели не учитываются.
Таким образом, предлагаемая модель позволяет учесть вклад электростатического и гидростатического взаимодействий во всем диапазоне межчастичных расстояний, при этом с уменьшением плотности кристалла (увеличением расстояния между частицами) этот вклад становится доминирующим. Уравнение Пуассона-Больцмана широко используется для описания обычных коллоидных систем [5]. В кристаллических системах уравнение Пуассона-Больцмана использовалось в [6] при моделировании плавления. В то же время для исследования свойств упругости электрически стабилизированных коллоидных кристаллов такой подход, насколько известно, ранее не применялся.
© Ю.В.Александров, П.Е.Дышловенко, 2011
В настоящей работе вычисляются силовые постоянные двумерного электрически стабилизированного коллоидного кристалла с квадратной решеткой с постоянным потенциалом частиц в рамках модели уравнения Пуассона-Больцмана. Описывается процедура определения силовых постоянных коллоидного кристалла на основе численного решения уравнения Пуассона-Больцмана, при этом используются свойства симметрии кристалла. Использование соображений симметрии, которые в случае квадратной решетки весьма наглядны, позволяет существенно упростить процедуру и значительно сократить объем вычислений, а также предсказать структуру матриц силовых постоянных. Информация о силовых постоянных позволяет оценить вклад ближайших соседей различных порядков в общее силовое взаимодействие в кристалле. Для этого предложено использовать специально введенные параметры.
Описание задачи и основные уравнения. Рассматриваемая система схематично показана на рис.1. Она включает в себя бесконечно длинные твердые цилиндрические частицы радиусом Я, расположенные в узлах квадратной решетки с постоянной й . Векторы а1 и а2 - векторы примитивных трансляций решетки, |а1| = |а2| = d. Частицы имеют постоянный потенциал поверхности ф/;.
Пространство вне частиц заполнено симметричным бинарным одновалентным электролитом. В ходе компьютерных экспериментов центральная частица испытывает смещения в положительном направлении оси х.
Для приведения задачи к безразмерному виду вводятся: для единицы длины
-1 ( 2 / 7 VI2 (длины Дебая) к = ^2пде/ ев 0кТ) ,
для электрического потенциала кТ^е .
Единицы для других величин могут быть при необходимости получены из этих двух, в частности силовые постоянные измеряются в единицах 2щ1 /ее 0. Здесь п - равновесная концентрация одной (любой) из двух компонент в объеме электролита; де - элементарный заряд; е - относительная диэлектрическая проницаемость электролита; е0 - электрическая постоянная; к - постоянная Больцмана; Т - абсолютная температура. В этих обозначениях уравнение Пу-ассона-Больцмана для электрического потенциала ф рассматриваемой системы имеет вид [7, 8]
Рис.1. Поперечное сечение системы бесконечно длинных цилиндрических частиц, образующих двумерную квадратную решетку
V ф = бЬ ф.
(1)
Стандартные электростатические граничные условия на поверхности частицы принимают вид
ф = ф г
(2)
где фр - потенциал на поверхности частицы.
Область определения численной задачи включает в себя ячейку Вигнера-Зейтца произвольной частицы, а также ячейки Вигнера-Зейтца ее ближайших соседей до пятого порядка включительно. В силу симметрии задачи в вычислениях использовалась толь-
ко половина указанной области. На внешней границе области определения задано однородное граничное условие Неймана:
Еп = 0, (3)
где п - внешняя единичная нормаль; Еп - нормальная компонента электрического поля Е = -Уф.
Решение уравнения (1) совместно с граничными условиями (2), (3) полностью описывает распределение электрического потенциала ф = ф(х, у) в области вне частиц. Система уравнений решалась численно методом, описанным в [7, 8]. Метод включает в себя конечно-элементное решение дифференциального уравнения в сочетании с адаптивным перестроением сеток конечных элементов. Используя решение ф(х, у), с помощью тензора напряжений, ассоциированного с уравнением (1), вычислялись силы, действующие на частицы системы, при различных смещениях центральной частицы. Затем путем дифференцирования сил по смещениям находились силовые постоянные.
Определение силовых постоянных. Рассмотрим кристалл с моноатомной решеткой, все частицы которого находятся в своих положениях равновесия. Потенциальная энергия кристалла за счет внутренних сил взаимодействия при малых смещениях частиц из положения равновесия может быть представлена следующим разложением [9]:
V=1 I
2 а, Р, N М
( 2 ^ а2 V
Ха, N ^ Хр, ^ М )
Ха, N Хр, N+М + К , (4)
где X - а-компонента смещения Z частицы из ее положения равновесия, задаваемого вектором N N и М - векторы решетки Бравэ; а = х, у; Р = х, у. Все производные в (4) вычисляются в нуле; для упрощения обозначений в дальнейшем это специально не отмечается. Коэффициенты СМ квадратичной формы (4), называемые
силовыми постоянными, зависят только от относительного положения частиц, определяемого вектором М. (В дальнейшем вектор М используется для обозначения как узла решетки Бравэ, так и частицы, центр которой в положении равновесия находится в этом узле.) Это позволяет без ограничения общности положить N = 0 :
СМр=-—-. (5)
ар г) 7 г) 7
а ха, 0 а ХР, М
В случае двух измерений для данного вектора М имеются четыре силовые постоянные, записываемые в виде матрицы (2 х 2) :
СМ =
^м ^МЛ
С хх Сху ММ
V Сух Суу У
Силовые постоянные (5) обладают рядом общих свойств [10]:
ГМ _ ,^-М _ ,^-М — П (£Л
Сар = Сра , Сар = Сар , IСар = 0 . (6)
М
Из первых двух свойств следует свойство симметричности
гы _ гм
Сар = Сра
(7)
матриц См для любого вектора М.
В (5) производная энергии V по смещению 2р м является с точностью до знака
компонентой Fр м силы, действующей на частицу М :
F
д V
р, м
д 2,
(8)
р,м
Компоненты сил можно вычислить непосредственно, без использования численного дифференцирования, посредством интегрирования тензора напряжений:
м ={
'Р
Е®Е-I 1Е2 + соб ф-111
•п ё!, Р = х, у.
(9)
Здесь Е = —Уф; I - единичная матрица; Ь - длина границы Г ячейки Вигнера-Зейтца, выбранной в качестве контура интегрирования; п - вектор внешней единичной нормали к элементу контура интегрирования; ер - соответствующий единичный базисный вектор декартовой системы координат.
Из (5) и (8) следует, что для получения силовых постоянных силы, определяемые из (9), необходимо однократно продифференцировать по смещению 2а 0 :
гм _ _д м
Са,р = "
д 2,
(10)
а, 0
В соответствии с (10) силовая постоянная Сар интерпретируется как взятая с обратным знаком Р-компонента силы, действующей на частицу м, рассчитанная на единицу смещения, при бесконечно малом смещении частицы 0 в а-направлении.
Для системы на рис. 1 вектор м принимает значение (0, 0) для частицы, произвольно выбранной в качестве центральной, (1, 0) (0,1) (—1, 0) (0, —1) - для соседей 1-го порядка и (1,1) (—1,1) (—1, —1) (1, —1) - для соседей 2-го порядка. При этом используется базис (я1, а2) примитивных векторов решетки. Используя выражение (10) и вытекающую из него интерпретацию силовых постоянных, а также свойства симметрии квадратной решетки, можно сделать заключение о структуре матриц силовых постоянных. Для центральной частицы и ее ближайших соседей до 2-го порядка включительно результаты приведены ниже.
Матрица силовых постоянных центральной частицы определяется одной независимой постоянной А и имеет вид
А(0'0) =
ГА 0 >! 0 А
Матрицы силовых постоянных для соседей 1-го порядка определяются двумя независимыми постоянными А и В и имеют вид
Г
(А 0 ^ ( В 0 ^
А(1'0) = А(_1'0) = А1 0 а(од) = А(0'_1) = 0
о в
0А
Матрицы силовых постоянных для соседей 2-го порядка тоже определяются двумя независимыми постоянными А и В2 и имеют вид
А(1Л) = А(-1-1) =
Г а2
V В2
в2 Л А2 )
А("11) = А(1-1) =
Г А,
V
2 - В,
2
- В2 Л А2 )
Из приведенных выражений следует, что для определения силовых постоянных квадратной решетки до 2-го порядка включительно требуются следующие пять функциональных зависимостей компонент сил, действующих на различные частицы, от смещения X центральной частицы в направлении оси х в окрестности нуля: FХ(00)(X), Ех,(10)(X) , ^.д0д)(X) , ^сд1д)(Х),
Еу (1 1)(Х) . Производные этих зависимостей по X, вычисленные в нуле, дадут искомые силовые постоянные.
Функциональные зависимости сил от смещения определены для набора равноотстоящих значений X, включающих нулевое значение, десять шагов в положительном х-направлении и десять таких же шагов в отрицательном х-направлении. Наибольшее смещение частицы из ее положения равновесия составляло 10% от удаления ближайших соседей. Вычисления выполнялись только для положительных смещений. Данные для отрицательных смещений получены из данных для положительных на основе свойств симметрии. Например, х-компонента силы, действующей на частицу (1, 0) при отрицательных смещениях, определялась как взятая ^ _о,2 с обратным знаком х-компонента силы, дей- -0,3 ствующей на частицу (-1, 0) при положитель- _п 4 ных смещениях. Пример зависимостей сил от смещения, полученных в результате компьютерных экспериментов, для квадратной решетки с периодом ё = 4,0 показан на рис.2.
0,5-.
0,1 0 0,1 0,2 Смещение X центральной частицы
Рис.2. Зависимости компонент сил, действующих на Численно дифференцирование согласно различные частицы, от смещения X центральной час-пугем аппроксимации тицы в направлении оси х в окрестности положения
(10) выполнялось
функций м (X) полиномами и взятием р^швестя X =0 (сила ^(щСЮ не шказ^а коэффициента при линейном члене. поск°льку блгока к ^Х, (ц) да)
Результаты вычисления силовых постоянных. На рис.3 показаны результаты вычислений силовых постоянных для коллоидных кристаллов с частицами радиусом Я = 1,0, поверхностным потенциалом частиц ф = 2,0 и параметром решетки ё, изменяющимся в диапазоне от 2Я до 8,0. Поведение силовых постоянных при больших значениях параметра решетки при ё > 5,0 близко к экспоненциальному, однако на меньших расстояниях наблюдаются существенные отклонения от этого закона. Особенности поведения кривых на рис.3 при ё < 2,5, т.е. при малых удалениях частиц, объясняются спецификой модели постоянного потенциала: при тесном сближении одноименно заряженных частиц заряд на их поверхности может перемещаться, приводя к уменьшению величины заряда в области наибольшего сближения и даже к локальному изменению знака заряда.
Полученные результаты позволяют провести оценку относительного вклада соседей различных порядков в общее силовое взаимодействие в кристалле. Оценка основана на третьем из свойств (6) силовых постоянных. Это свойство может быть переписано в виде
С0Р--У С (11)
м
для любой пары индексов (а, р).
Таким образом, силовые постоянные центральной частицы с точностью до знака равны сумме соответствующих силовых постоянных всех остальных частиц кристалла. Равенство (11) выполняется строго только тогда, когда суммирование справа осуществляется по всем ненулевым векторам М решетки. Если же ограничиться суммированием только по ближайшим соседям до некоторого определенного порядка, то левая и правая части в (11) будут отличаться друг от друга. Удобно ввести параметры ^ и £2 :
У С М У С М
/—I XX XX
П _ М_ ? — _ М 2
XX
.^0 ' 2 ^0
у у у. у
где М - вектор решетки Бравэ для ближайших соседей только 1-го порядка; ММ2 - для ближайших соседей 1-го и 2-го порядка вместе. Близость этих параметров к единице является мерой влияния ближайших соседей указанных порядков.
Зависимости параметров Sj и S2 от параметра решетки d показаны на рис. 4. Для параметра Sj наблюдается некоторое отклонение от единицы, особенно при высоких плотностях. В то же время значения параметра S2 очень близки к единице во всем исследованном диапазоне плотностей. Из приведенных данных следует, что силовое взаимодействие в двумерном коллоидном кристалле с квадратной решеткой практически полностью исчерпывается взаимодействием ближайших соседей только 1-го и 2-го порядка; вклады соседей более высоких порядков пренебрежимо малы для любых плотностей вплоть до полного контакта макроионов.
Работа выполнена при финансовой поддержке Российского фонда фундаментальных исследований (грант № 09-01-97012).
Литература
1. Дерягин Б.В., Чураев Н.В., Муллер В.М. Поверхностные силы. - М.: Наука, 1987.- 398 с.
2. Van Blaaderen A., RuelR., Wiltzius P. Template-directed colloidal crystallization // Nature. - 1997. -Vol. 385. - P. 321-324.
3. Joannopoulos J.D., Villeneuve P.R., Fan S.H. Photonic crystals: putting a new twist on light // Nature. -1997. - Vol. 386. - P. 143-149.
4. Israelachvili J.N. Intermolecular and Surface Forces. - 2-nd ed. - London: Academic Press, 1991.
5. BelloniL. Colloidal interactions // J. Phys.: Condens. Matter. - 2000. - Vol. 12. - P. R549-R587.
6. Dobnikar J., Chen Y., R. Rzehak, H.H. von Grunberg. Many-body interactions and the melting of colloidal crystals // J. Chem. Phys. - 2003. - Vol. 119, № 9. - P. 4971-4985.
7. Dyshlovenko P.E. Adaptive mesh enrichment for the Poisson-Boltzmann equation // J. Comp. Phys. -2001. - Vol. 172. - P. 198-208.
8. Dyshlovenko P.E. Adaptive numerical method for Poisson-Boltzmann equation and its application // Comp. Phys. Commun. - 2002. - Vol. 147. - P. 335-338.
9. Фейнман Р.Ф. Статистическая механика: Пер с англ. / Под ред. Д. Н. Зубарева. - М.: Мир, 1978. -407 с.
10. Ашкрофт Н., Мермин Н. Физика твердого тела. Т. 1: Пер. с англ. / Под ред. М.И. Каганова. -М.: Мир, 1979. - 399 с.
Статья поступила 10 августа 2010 г.
Александров Юрий Владимирович - аспирант Ульяновского государственного технического университета. Область научных интересов: компьютерное моделирование коллоидных систем.
Дышловенко Павел Евгеньевич - кандидат физико-математических наук, доцент Ульяновского государственного технического университета. Область научных интересов: физика конденсированных сред, компьютерное моделирование, высокопроизводительные вычислительные системы. E-mail: pavel@ulstu.ru
М
c-j
1,0
s
И
3 0,9
и
1 0,
G 0.7Н
: S2
:
III! 2 3 ¡■■i 4 ■ ¡■i iiii б IIII iiii 8
Постоянная кристаллической решетки с1
Рис.4. Параметры ^ и Б2 как функции параметра решетки ё