Оценка формы цепных углеводородных молекул методом Монте-Карло

Журкин Дмитрий Викторович; Рабинович Александр Львович

УЧЕНЫЕ ЗАПИСКИ ПЕТРОЗАВОДСКОГО ГОСУДАРСТВЕННОГО УНИВЕРСИТЕТА

Сентябрь, № 6 Физико-математические науки 2014

УДК 541.64:539.199

ДМИТРИЙ ВИКТОРОВИЧ ЖУРКИН

аспирант кафедры физики твердого тела физико-технического факультета, Петрозаводский государственный университет (Петрозаводск, Российская Федерация) zhurkindm @mail. ru

АЛЕКСАНДР ЛЬВОВИЧ РАБИНОВИЧ

доктор физико-математических наук, старший научный сотрудник лаборатории экологической биохимии, Институт биологии Карельского научного центра РАН (Петрозаводск, Российская Федерация) rabinov@krc. karelia. ru

ОЦЕНКА ФОРМЫ ЦЕПНЫХ УГЛЕВОДОРОДНЫХ МОЛЕКУЛ МЕТОДОМ МОНТЕ-КАРЛО*

Изучение физических и химических свойств углеводородных цепных молекул очень важно для развития многих технологических областей и для углубления понимания структурных особенностей и функционирования различных биологических систем. При описании разных явлений (вязкого потока молекул, адсорбции на поверхностях, абсорбции и диффузии в пористых системах и т. д.) важно учитывать форму молекул. В настоящей работе методом Монте-Карло при температурах T = 293, 303 и 313 K исследовано конформационное поведение 65 цепных углеводородных молекул следующего строения: CH3- (CH2) a- (CH=CH-CH2)d- (CH2) b-CH3 (где a, b, d - целые); количество атомов углерода в остове цепей N = 16, 18, 20, 22, количество cis-двойных связей d = 0, 1, ..., 6. Прототипы молекул - остатки жирных кислот молекул фосфолипидов. Моделирование всех молекул проводили единообразно, в невозмущенном состоянии, использовали «классическую гибкую модель» цепи. При генерировании значений торсионных углов использовали метод существенной выборки в диапазоне 0-360°, с учетом взаимозависимости каждых трех углов вдоль по цепи. Проведена аппроксимация конформаций каждой молекулы прямоугольными параллелепипедами и вычислены средние характеристики формы («спэны»). Изучены зависимости рассчитанных характеристик от количества атомов углерода в цепи, количества двойных связей и их положения в молекулах (зависимости между структурой и свойствами). В их основе лежат характерные особенности внутренних вращений в цепях рассмотренного химического строения.

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

ВВЕДЕНИЕ

Углеводородные цепные молекулы (молекулы, состоящие из групп CHn (n = 1, 2, 3), валентно связанных в цепь) являются объектами исследования многих областей науки: физики конденсированного состояния, физики высокомолекулярных соединений (полимеров), молекулярной физики, биофизики и т. д. Кроме того, такие молекулы широко применяются в производстве, например, поверхностно-активных веществ [6], синтетических волокон [22], [26], [35], эластомеров [23], [24], [32], лакокрасочных материалов, изделий [4]. Углеводородные цепи остатков жирных кислот входят в состав молекул фосфолипидов [2],

[12], [25], которые в природе образуют бислойные и другие структурные образования (мембраны, мицеллы и т. п.). Таким образом, изучение физических и химических свойств углеводородных цепных молекул имеет большое значение как для развития многих технологических областей, так и для углубления понимания свойств различных биологических систем.

Важную роль в теоретических интерпретациях поведения и свойств цепных молекул в разно© Журкин Д. В., Рабинович А. Л., 2014

образных молекулярных системах играет учет их формы [1], [19], [28], [29], [31] (последняя была изучена впервые в работе [21], было показано, что форма усредненной по конформациям полимерной цепи не является сферической). Учет формы важен, например, при описании вязкого потока [33] молекул, при адсорбции на поверхностях, при абсорбции и диффузии в пористых системах [30], [34]. Информация о форме отдельных углеводородных цепей липидов позволяет оценивать степень их совместимости, комплементарности при сосуществовании в бислое при разных условиях [8], [10], [11]. Для исследования формы молекул удобно применять методы компьютерного моделирования: метод молекулярной динамики (МД), метод Монте-Карло (МК), поскольку они позволяют извлечь любую детальную информацию об исследуемых системах.

Жирнокислотные цепи в липидах биомембран могут различаться по химическому строению

[2], [12], [25] (мембраны, как правило, гетеро-генны, в них сосуществуют десятки липидных молекул разной структуры). В настоящей работе рассмотрены лишь неразветвленные углеводо-

110

Д. В. Журкин, А. Л. Рабинович

родные цепные молекулы. Это могут быть цепи без двойных связей (насыщенные цепи), а также таковые с двойными связями (ненасыщенные цепи). Методом МК моделировали конформаци-онное поведение молекул общего строения CH3 -(CH2)a - (CH = CH - CH2)d - (CH2)b - CH3 (где a, b, d - целые) при температурах T = 293, 303, 313

K. Прототипами указанных молекул являлись остатки жирных кислот молекул фосфолипидов, в которых проведена замена концевой группы COO на группу CH3. Исследованы цепи с количеством атомов углерода в остове N = 16, 18, 20, 22, количеством двойных связей d = 0, 1, ..., 6 (все -в конфигурации cis). Моделируемые молекулы отличались также положением двойных связей вдоль по цепи. Всего исследовано 65 различных молекул, изучены характеристики их формы.

МОДЕЛЬ

Для сравнения формы углеводородных цепных молекул разного химического строения сле-

2п 2п 2п

11...| Н(^1,^2,...,^n-1)• exp[-

^ Н \ = 0 0 0______________________

^ П ^ ~ 2п2п 2п

И ... I exp[- Ushort (^ 0 0 0

довало избрать некоторое однотипное состояние молекул и провести расчет свойств единым методом в идентичных условиях. В качестве такового избрано невозмущенное состояние (©-условия) [9]. Средние характеристики молекул в ©-условиях приблизительно соответствуют таковым в жидком или аморфном состояниях вещества [9], [13], [36]. В соответствии с определением невозмущенного состояния [9] моделировали одиночную цепь, в которой учитывали взаимодействия только тех ее атомов, которые находились в пределах некоторых фрагментов вдоль по цепи (то есть ближние взаимодействия). Равновесные свойства такой одиночной цепи вычисляли в рамках так называемой классической гибкой модели [18], согласно которой цепная молекула - это система взаимодействующих материальных точек (разработано несколько типов «гибких моделей» [14], [16], [18]). Для расчета среднего значения <Н> какой-либо характеристики Н цепной молекулы в NVT-ансамбле использовали выражение

Ushort ^.- ^N -1) / квТ] d^d^2... d^N-i

---------------------------------------, (1)

^2,..., &N-1 ) / kBT ] dP1 dty 2... d^N -1

где N - количество атомов углерода в остове цепи; ф1, ф2, ..., ^N-1 - углы внутреннего вращения (торсионные углы) вокруг C-C связей основной цепи; Н (ф1, ф2, ., ф№1) - заданная характеристика и Ushort (ф1, ф2, ..., ^N-1) - энергия ближних взаимодействий цепной молекулы в конформации с торсионными углами ф1, ф2, ..., ф№1; валентные углы и валентные связи предполагаются зафиксированными вблизи их равновесных значений; kB - постоянная Больцмана, Т - температура. В использованной модели цепи углы внутреннего вращения могли принимать любые значения из диапазона 0-360°.

В данной работе учитывали влияние взаимозависимости изменения каждых трех последовательных торсионных углов молекулы на ее конформационные свойства. Для этого разбивали молекулы на фрагменты, в рамках которых вычисляли ближние взаимодействия. Это осуществляли таким образом, что каждый у-й фрагмент (при отсчете вдоль по цепи начиная от одного из ее концов) содержал три соседних угла внутреннего вращения ф, фу+1, фу+2. Смежные фрагменты у-й и (у+1)-й при этом имели два общих торсионных угла - фу+1, фу+2. Для конструирования любой из рассмотренных молекул по этим правилам

требуется несколько фрагментов из 16 различных типов, которые указаны ниже:

CH3 - CH2 - CH2 - CH2CH2, CH3CH2 - CH2 -CH2 - CH2CH2, CH2CH2 - CH2 - CH2 - CH2CH2, CH2CH2 - CH2 - CH2 - CH2CH, CH2CH2 - CH2 -CH2 - CHCHCH2, CH2CH2 - CH2 - CHCH -CH2CH2, CH3CH2 - CH2 - CH2 - CH2CH, CH2CH2 - CH2 - CHCH - CH2CH, CH2CH2 -CHCH - CH2 - CHCHCH2, CH3 - CH2 - CHCH -CH2CH, CH3CH2 - CHCH - CH2 - CHCHCH2, CH2CHCH - CH2 - CHCH - CH2CH, CH2CHCH -CH22 - CH2 - CH22CH3, CHCH2 - 2CH2 - CH2 2 - CH3,

chch2 - chch - ch2 - ch2ch3, ch2chch -

CH2 - CH2 - CH3.

Символом «-» здесь отмечены валентные связи во фрагменте, взаимозависимость вращений вокруг которых учитывали при моделировании. Вычисление заданной величины <Н> из (1) проводили методом МК (метод описан, например, в [17]), то есть определяли ее оценку Н . Для этого использовали ранее разработанный алгоритм МК [7] генерирования случайных точек ($ ,ф\,..., ^N-1) в пространстве торсионных углов с плотностью вероятностей

(2)

P(^1,^2,...,^N-1)

eXP[- Ushort fob ^.^ ^N-1) / kBT]

1Ы exp[-

Ushort (^1,^2,...,^N-1)/квТ ] <М<7ф2... d^N-1

Оценка формы цепных углеводородных молекул методом Монте-Карло

111

В данном алгоритме конфигурационное пространство трех торсионных углов каждого из фрагментов, из которых состоит рассматриваемая цепь, разбивают на 1 000 000 состояний («параллелепипедов») так, что вероятности их

реализации одинаковы, а затем генерируют конформации этой цепи с помощью полученных разбиений. Оценку H величины <H> при использовании этого алгоритма МК следует вычислять согласно выражению [7]

. (3)

Здесь v - номер текущей конформации; а -

объем выборки конформаций; {ф}У - координаты случайной точки в пространстве торсионных углов, выбираемой с плотностью вероятностей p (ф1, ф2, ..., ^N-1) из выражения (2); величины

(Lj) , (L2) , (L3) - длины сторон,

в угловых единицах, того параллелепипеда, который был выбран случайно в процедуре МК в конформации номер v во фрагменте номер у типа т. Для вычисления энергии ближних взаимодействий исследуемых молекул использовали силовое поле CHARMM27 [15] с поправками [20]. Для большинства исследованных цепных молекул при вычислении оценок средних характеристик (по формуле (3)) с заданной точностью были использованы выборки в 129 600 000 000 конформаций. В случае пентаеновых и гексаеновых цепей выборки были увеличены в 10 раз, поскольку поверхности потенциальной энергии почти всех их фрагментов характеризуются наиболее пологими минимумами [8].

Для оценки формы молекул в литературе используется несколько подходов (см., например, [19]). В настоящей работе использовали метод аппроксимации формы параллелепипедом. Для этого в каждой конформации данной цепи находили систему координат O'x'y'z' с началом O' в центре масс этой конформации и осями O'x', O'y', O'z', совпадающими с главными осями инерции данной конформации (рис. 1). Алгоритм поиска системы координат O'x'y'z' состоял в следующем. В исходной системе координат Oxyz (см. рис. 1) вычисляли координаты всех атомов текущей конформации цепи, затем координаты ее центра масс O и осуществляли параллельный перенос исходной системы Oxyz, совмещая начало с центром масс. В этой системе вычисляли компоненты тензора инерции [5] рассматриваемой конформации, проводили его диагонализацию и в итоге определяли направления главных осей инерции данной конформации. Надлежащим поворотом осей координат в точке O совмещали их с осями инерции, переходя в итоге в искомую систему O 'x'y 'z '. Обозначения O 'x', O 'y ', O 'z ' ее осей проводили одинаково в каждой конформации: ось O 'x' всегда отвечала наибольшему главному моменту инерции (то есть O x - ось минимальной протяженности конформации), а ось O z - наименьше-

му главному моменту инерции (то есть O z - ось максимальной протяженности конформации). Затем в текущей конформации цепи вычисляли разности между максимальными и минимальными проекциями атомов углерода и водорода на каждую ось O x , O y и O z . Обозначали полученные разности через g1, g2, g3 соответственно. Их можно рассматривать как минимальные размеры ребер параллелепипеда (ориентированного параллельно главным осям инерции), в который конформация молекулы может поместиться целиком (см. рис. 1). В итоге компьютерного моделирования вычисляли средние величины <g1>, <g2>, <g3> и средние произведения двух размеров <gi‘g2> (в определенном смысле - это площадь поперечного сечения молекулы).

Рис. 1. Расположение «молекулярной» системы координат O 'x'y 'z ', которая определяется главными осями инерции O x , O y , O z данной конформации цепной молекулы, представленной на рисунке только остовом - ломаной линией, соединяющей атомы углерода (•). Переход в систему O 'x'y 'z ' осуществлялся в каждой конформации из некоторой исходной системы координат (Oxyz). Пунктирные линии -трехосный эллипсоид инерции данной конформации цепи, O ' - ее центр масс. Главные оси инерции обозначали так, чтобы ось инерции O 'z ' отвечала направлению наибольшего протяжения конформации цепи (то есть наименьшему собственному значению ее тензора инерции), а ось O x -наименьшему протяжению. Величины g1, g2, g3 - минимальные размеры ребер параллелепипеда, построенного по следующему правилу: его грани выбраны перпендикулярно главным осям инерции данной конформации цепной молекулы (ребро g1 параллельно оси O 'x', g2 - оси O 'y ', g3 - оси O 'z ') таким образом, чтобы центры всех атомов углерода и водорода молекулы в данной конформации поместились внутри параллелепипеда

112

Д. В. Журкин, А. Л. Рабинович

Для удобства представления результатов расчетов целесообразно указывать местоположение нескольких двойных связей в цепи с помощью лишь одного параметра. В качестве такового было использовано (аналогично работе [27]) среднее арифметическое номеров всех атомов углерода, участвующих в образовании двойных связей в цепи (атомы углерода нумеровали вдоль по цепи, начиная от одного из концевых). Этот параметр, обозначенный в работе [27] через X, характеризует положение «центра» двойных связей в цепи. Наряду с этим в ряде случаев вместо параметра X оказалось удобным использовать параметр А [27] - номер атома углерода, ближайшего к заданному концу цепи и участвующего в образовании первой (от данного конца цепи) двойной связи. Соответствие местоположений всех двойных связей в молекулах параметру X или параметру А указано в таблице.

РЕЗУЛЬТАТЫ И ОБСУЖДЕНИЕ

Химическое строение любой из рассмотренных молекул идентифицируется тремя параметрами - N, d и X (или А), и поэтому сравнение между собой любых характеристик для разных молекул требует изучения конкуренции этих трех параметров. Учет влияния всех трех на заданную характеристику эквивалентен задаче рассмотрения 4-мерной гиперповерхности; это весьма затруднительно. Исследование 3-мерных зависимостей (поверхностей) данного свойства от отдельных пар параметров не всегда эффективно, поскольку поверхности могут пересекаться

в пространстве. Поэтому наиболее целесообразно проводить анализ зависимостей свойств от одного из параметров при фиксировании двух других. На рис. 2, 3 приведены полученные в настоящей работе по итогам компьютерных экспериментов при T = 293 K зависимости оценок средних размеров <gj>, <g2>, <g3> и средней «площади поперечного сечения» <gx'g2> от параметра X, а на рис. 4, 5 - от количества d двойных связей в цепи. Доверительные интервалы на всех представленных графиках (они меньше размеров символов) отвечают 95 %-ной надежности согласно t-распределению Стьюдента. Анализ полученных кривых позволил выявить ряд закономерностей:

I. При смещении центра X группы из d двойных связей от концов к середине цепи при N = Const, d = Const для всех исследованных молекул (N = 16, 18, 20, 22) с количеством d = 1-4, а также для молекул N = 22, d = 5 средний размер <g3> уменьшается (см. рис. 2). В молекулах N = 20, d =5 величина <g3> меняется незначительно.

II. Зависимости величин <g1'g2> и <g2> от X для молекул с количеством d = 1-3 симбатны (см. рис. 2, 3): при смещении центра X от концов к середине цепи обе величины растут. С увеличением количества d зависимости от X ослабевают: для тетраеновых (d = 4) молекул с N = 20, 22 величина <g1'g2> при смещении X от концов к середине еще растет, но для N = 18 уже почти не меняется. Величина размера <g2> тетраенов почти не зависит от X для всех изученных N (N = 16, 18, 20, 22). В пентаеновых (d = 5) молекулах зависимость от X очень слаба и для величин <g1'g2>,

Соответствие положений метиленпрерывающихся двойных связей (d - их количество) в углеводородных цепях значениям параметров А и X

Номера атомов C, участвующих в образовании двойных связей Д X Номера атомов C, участвующих в образовании двойных связей Д X Номера атомов C, участвующих в образовании двойных связей Д X

3, 4 (d=1) 3 3.5 3, 4, 6, 7 (d=2) 3 5 3, 4, 6, 7, 9, 10 (d=3) 3 6.5

4, 5 (d=1) 4 4.5 4, 5, 7, 8 (d=2) 4 6 4, 5, 7, 8, 10, 11 (d=3) 4 7.5

5, 6 (d=1) 5 5.5 5, 6, 8, 9 (d=2) 5 7 5, 6, 8, 9, 11, 12 (d=3) 5 8.5

6, 7 (d=1) 6 6.5 6, 7, 9, 10 (d=2) 6 8 6, 7, 9, 10, 12, 13 (d=3) 6 9.5

7, 8 (d=1) 7 7.5 7, 8, 10, 11 (d=2) 7 9 7, 8, 10, 11, 13, 14 (d=3) 7 10.5

8, 9 (d=1) 8 8.5 8, 9, 11, 12 (d=2) 8 10

9, 10 (d=1) 9 9.5 9, 10, 12, 13 (d=2) 9 11 3, 4, 6, 7, 9, 10, 12, 13 (d=4) 3 8

10, 11 (d=1) 10 10.5 10, 11, 13, 14 (d=2) 10 12 4, 5, 7, 8, 10, 11, 13, 14 (d=4) 4 9

11, 12 (d=1) 11 11.5 11, 12, 14, 15 (d=2) 11 13 5, 6, 8, 9, 11, 12, 14, 15 (d=4) 5 10

12, 13 (d=1) 12 12.5 12, 13, 15, 16 (d=2) 12 14

13, 14 (d=1) 13 13.5 13, 14, 16, 17 (d=2) 13 15 3, 4, 6, 7, 9, 10, 12, 13, 15, 16 (d=5) 3 9.5

14, 15 (d=1) 14 14.5 14, 15, 17, 18 (d=2) 14 16 4, 5, 7, 8, 10, 11, 13, 14, 16, 17 (d=5) 4 10.5

15, 16 (d=1) 15 15.5 15, 16, 18, 19 (d=2) 15 17 5, 6, 8, 9, 11, 12, 14, 15, 17, 18 (d=5) 5 11.5

16, 17 (d=1) 16 16.5 16, 17, 19, 20 (d=2) 16 18 6, 7, 9, 10, 12, 13, 15, 16, 18, 19 (d=5) 6 12.5

17, 18 (d=1) 17 17.5

18, 19 (d=1) 18 18.5 3, 4, 6, 7, 9, 10, 12, 13, 15, 16, 18, 19 (d=6) 3 11

19, 20 (d=1) 19 19.5 4, 5, 7, 8, 10, 11, 13, 14, 16, 17, 19, 20 (d=6) 4 12

Оценка формы цепных углеводородных молекул методом Монте-Карло

113

Рис. 2. Зависимости средних размеров <gj>, <g2>, <g3> ребер прямоугольного параллелепипеда для невозмущенных углеводородных молекул с количеством N атомов углерода в цепи 16 (^^), 18 (^^), 20 и 22 (^^), по группам

молекул с одинаковым количеством d двойных связей cis (насыщенные, d = 0; моноены, d = 1; диены, d = 2; триены, d = 3; тетраены, d = 4; пентаены, d = 5; гексаены, d = 6), от местоположения в цепи «центра» метиленпрерывающихся двойных связей - параметра X. Соответствие значений параметра X положениям всех двойных связей в молекулах указано в таблице. Расчет методом МК при температуре T = 293 K

Рис. 3. Зависимости величин средней «площади поперечного сечения» <gpg2> невозмущенных углеводородных молекул от параметра X. Условия и обозначения, как на рис. 2

и для <g2>. Средний размер <gj> для молекул со всеми рассмотренными N и d почти не зависит от X.

III. С ростом d при N = Const, А = Const (для всех рассмотренных N и А) величина <g3> уменьшается, <gj> увеличивается (хотя и слабо), а величины <g2>, <gi'g2> изменяются немонотонно (см. рис. 4, 5), и результат зависит как от N, так и от А. Немонотонность в большей степени проявляет зависимость от d величины <gpg2>, а для цепей с А = 3 она “S-образна” (имеет участки убывания и возрастания). В каждой совокупнос-

ти молекул с одинаковым N, которые содержат разное количество d, от нуля до максимально возможного (d = 0, 1, 2, ..., dmax), метиленпрерывающихся двойных связей с одинаковым расположением А = 3, существует молекула с количеством двойных связей d1 (не обязательно равным нулю), площадь <g1'g2> которой минимальна, и молекула с количеством двойных связей d2 (не равным dmax), площадь <g1'g2> которой максимальна; конкретные значения d1 и d2 зависят от N (см. рис. 5).

IV. С ростом N при d = Const, А = Const средние значения величин <g1>, <g2>, <g3>, <g1'g2>

114

Д. В. Журкин, А. Л. Рабинович

Рис. 4. Зависимости средних размеров <gj>, <g2>, <g3> (спэ-нов, рис. 1) от количества d двойных связей cis для невозмущенных углеводородных молекул с параметром А = 3 (-9-), 4 (-*-), 5 (-*-), 6 (-в-), 7 (-#-), 8 (-*-), 9 11 (~■“) по группам молекул с одинаковым количеством N атомов углерода в остове цепи. Параметр А - это номер атома углерода, участвующего в образовании первой двойной связи от конца цепи. Соответствие значений параметра А положениям всех двойных связей в молекулах указано в табл. 1. Расчет методом МК при температуре T = 293 K

Рис. 5. Зависимости величин средней «площади поперечного сечения» <g1'g2> от количества d двойных связей cis для невозмущенных углеводородных молекул. Условия и обозначения, как на рис. 4

растут для всех рассмотренных вариантов цепей с параметрами d = 0-5 и А = 3-11 (данные не представлены), что и ожидалось для коротких цепей в невозмущенном состоянии [9].

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

Наиболее вероятной в каноническом ансамбле конформацией цепи является конформация, отвечающая глобальному минимуму ее потенци-

альной энергии. Такая конформация насыщенной цепи - это «транс-зигзаг» [3], [9], когда цепь предельно вытянута, а все атомы углерода лежат в одной плоскости. Из соображений симметрии следует, что ось O'z' конформации «транс-зигзаг» лежит в плоскости углеродного остова цепи и направлена от атома углерода одного ее конца к атому углерода другого; ось O'y' также принадлежит плоскости остова и перпендикулярна оси O'z', а ось O'x' перпендикулярна плоскости остова; размер g3 соответствующего параллелепипеда будет максимальным, а g1<<g3, g2<<g3. Конформации, отвечающие минимумам потенциальной энергии ненасыщенных цепей с метиленпреры-вающимися двойными связями cis, являются более свернутыми, поскольку атомы углерода участков с двойными связями одной плоскости не принадлежат [3], [7], [9]. Чем больше таких двойных связей при N = Const содержит молекула, тем короче ее насыщенные участки (низко -энергетические конформации которых являются вытянутыми) с одной и/или другой стороны от группы двойных связей, поэтому меньше становится общий размер g3, меньше разница между g3 и размерами g1, g2. Эти закономерности сохранятся и для средних величин <g1>, <g2>, <g3>, что согласуется с данными рис. 4. Зависимость расчетных характеристик от местоположения метиленпрерывающихся двойных связей в цепях также объяснима. Действительно, пусть d = 1 и двойная связь в молекуле расположена у одного из концов цепи. Тогда большинство атомов углерода приходится на насыщенный участок, а в наиболее низкоэнергетической конформации они лежат в одной плоскости, влияние «свернутости» малого концевого участка на конформацию цепи минимально. Если же d = 1 и двойная связь в молекуле расположена вблизи середины цепи, то тенденция к «свернутости» конформации, наоборот, проявляется максимально: молекула имеет два наиболее протяженных насыщенных участка (по обе стороны от двойной связи), их конформации «транс-зигзаг» принадлежат двум разным плоскостям, под углом друг к другу. Очевидно, размер g3 второй моноеновой молекулы будет в итоге существенно меньше (а размеры g2, g1 - больше), чем соответствующие размеры первой. Такие же соотношения следует ожидать и для средних значений <g3>, <g2>, <g1>, что подтверждается результатами расчетов (см. рис. 2). Можно также ожидать, что обсуждаемое влияние двойной связи будет возрастать по мере перехода от предельного варианта 1 к варианту 2, то есть при смещении двойной связи от любого конца цепи к ее середине (см. рис. 2). Эти же соображения объясняют сходные зависимости расчетных характеристик от X, наблюдающиеся в полиено-вых молекулах (d > 1), с учетом того, что с увеличением d в цепи влияние местоположения двойных связей ослабевает, поскольку участок цепи,

Оценка формы цепных углеводородных молекул методом Монте-Карло

iНе можете найти то, что вам нужно? Попробуйте сервис подбора литературы.

115

на котором можно изменить положения двойных связей, сокращается.

Среди всех изученных в настоящей работе молекул в ©-условиях существуют совокупности, молекулы в которых, будучи различными по строению (то есть имея различные N, d, Л), обладают тем не менее соизмеримыми, сравнительно близкими «продольными» размерами <g3>. Если задавать различные интервалы значений <g3> (по оси ординат на рис. 2), то можно выделять совокупности молекул, размеры <g3> которых попадают в заданные интервалы. Зададим, например, отклонения ±5 % от величины <g3> молекулы 16:0. Можно убедиться, что в этот интервал значений попадают величины <g3> молекул 18:1©9cis, 18:2©6cis, 18:3©3cis, 20:3©6cis, 20:4©6cis, 20:5©3cis, 22:5©3cis, 22:6©3cis и ряда других (см. рис. 2); эти молекулы образуют одну из упомянутых совокупностей. Параметры N, d, Л молекул с метиленпрерывающимися двойными связями cis указаны здесь согласно «©-номенклатуре», через сокращенную формулу N: dюЛcis. Важно отметить, что все молекулы в приведенном примере, как оказалось, входят в список таких цепей, которые встречаются в биологических мембранах в качестве компонентов молекул фосфолипидов наиболее часто [2], [12], [25]. Очевидно, в бислоях в жидкокристаллическом состоянии продольные размеры этих цепей (проекции на нормаль к поверхности бислоя) должны быть согласованы друг с другом. Результаты наших расчетов свидетельствуют об их «предрасположенности» к аналогичной комплементарности и в невозмущенном состоянии.

Подобный результат ранее уже был отмечен в литературе при МК-моделировании и расчетах величин <g3> молекул [27]. При этом использованная в [27] модель цепи и алгоритм МК существенно отличались от таковых настоящей работы

[7]. В частности, при генерировании конформаций цепей в [27] учитывалась взаимозависимость каждых 2 последовательных торсионных углов основной цепи, а не 3, как в настоящей модели [7]; поверхность потенциальной энергии молекулярных фрагментов, зависящих от пар углов, разбивалась в [27] на 1 600 равновероятных состояний, тогда как энергии фрагментов, зависящие от троек углов, были разделены в [7] на 1 000 000 равновероятных состояний; различны также и параметры силовых полей, использованные в [27] и [7]; наконец, расчеты <g1>, <g2>, <g3> в [27] осуществлялись с учетом только атомов углерода, а в настоящей работе по модели [7] - с учетом атомов углерода и водорода. В итоге и численные значения соответственных величин <g1>, <g2>, <g3> одних и тех же молекул, рассчитанные в [27] и в настоящей работе, различаются существенно (на ~20-35 %). Однако перечисленные выше закономерности качественно, в главном сохраняются, и в этом состоит их основная ценность. Устойчи -вость обсуждаемых зависимостей и тенденций к деталям используемой модели и к параметрам силового поля свидетельствует о том, что они отражают фундаментальные соотношения между структурой и свойствами данных молекул, в основе которых лежат характерные особенности внутренних вращений в цепях рассмотренного химического строения.

Совокупности цепей, отвечающих разным диапазонам размеров <g3>, могут, по-видимому, примерно соответствовать основным жирнокислотным составам липидов в мембранах разных биологических объектов, разных органов и тканей.

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

* Работа выполнена при поддержке программ президента РФ «Ведущие научные школы» (гранты НШ-1642.2012.4, НШ-1410.2014.4).

СПИСОК ЛИТЕРАТУРЫ

1. Волькенштейн М. В. Конфигурационная статистика полимерных цепей. М.; Л.: АН СССР, 1959. 468 с.

2. Геннис Р. Биомембраны: Молекулярная структура и функции. М.: Мир, 1997. 624 с.

3. Дашевский В. Г. Конформации органических молекул. М.: Химия, 1974. 428 с.

4. Лакокрасочные материалы и покрытия. Теория и практика / Пер. с англ. Л. Н. Машляковского и др.; Под ред. Р. Лам-бурна. СПб.: Химия, 1991. 512 с.

5. Ландау Л. Д., Лифшиц Е. М. Теоретическая физика. Т. 1: Механика: Учеб. пособие для вузов. 5-е изд., стереотип. М.: ФИЗМАТЛИТ, 2004. 224 с.

6. Поверхностно-активные вещества и полимеры в водных растворах / К. Холмберг, Б. Йёнссон, Б. Кронберг, Б. Линдман; Пер. с англ. Г. П. Ямпольской; Под ред. Б. Д. Сумма. М.: Бином. Лаборатория знаний, 2007. 528 с.

7. Рабинович А. Л., Жу ркин Д. В. Существенная выборка при моделировании непрерывного спектра конформаций макромолекул методом Монте-Карло // Труды КарНЦ РАН. Сер. «Математическое моделирование и информационные технологии». 2013. Вып. 4. С. 96-111.

8. Рабинович А. Л. Цепные молекулы как компоненты мембранных систем: компьютерное моделирование // Методы компьютерного моделирования для исследования полимеров и биополимеров / Отв. ред.: В. А. Иванов, А. Л. Рабинович, А. Р. Хохлов. М.: Книжный дом ЛИБРОКОМ, 2009. С. 410-454.

9. Флори П. Статистическая механика цепных молекул. М.: Мир, 1971. 440 с.

10. Applegate K. R., Glomset J. A. Effect of acyl chain unsaturation on the packing of model diacylglycerols in simulated monolayers // J. Lipid Res. 1991. Vol. 32. P. 1645-1655.

11. Brockman H. L., Applegate K. R., Momsen M. M., King W. C., Glomset J. A. Packing and Electrostatic Behavior of s«-2-Docosahexaenoyl and -Arachidonoyl Phosphoglycerides // Biophys. J. 2003. Vol. 85. P. 2384-2396.

116

Д. В. Журкин, А. Л. Рабинович

12. Cook H. W., McMaster C. R. Fatty acid desaturation and chain elongation in eukaryotes // New comprehensive biochemistry / G. Bernardi (gen. ed.). Vol. 36: Biochemistry of Lipids, Lipoproteins and Membranes / Eds.: D. E. Vance, J. E. Vance. 4th ed. S. l.: Elsevier, 2002. P. 181-204.

13. Dettenmaier M. Conformation of n-alkane molecules in the melt and in cyclohexane solution studied by small-angle neutron scattering // J. Chem. Phys. 1978. Vol. 68, № 5. P. 2319-2322.

14. Echenique P., Cavasotto C. N., Garcia-Risueno P. The canonical equilibrium of constrained molecular models // Eur. Phys. J. Special Topics. 2011. Vol. 200. P. 5-54.

15. Feller S. E., Gawrisch G., MacKerell Jr. A. D. Polyunsaturated Fatty Acids in Lipid Bilayers: Intrinsic and Environmental Contributions to Their Unique Physical Properties // J. Am. Chem. Soc. 2002. Vol. 124, № 2. P. 318-326.

16. F l o r y P. J. Foundations of Rotational Isomeric State Theory and General Methods for Generating Configurational Averages // Macromolecules. 1974. Vol. 7, № 3. P. 381-392.

17. Gould H., Tobochnik J., Christian W. An Introduction to Computer Simulation Methods: Applications to Physical Systems. San Francisco etc.: Pearson, 2007. 796 p.

18. G o N., Scheraga H. A. On the Use of Classical Statistical Mechanics in the Treatment of Polymer Chain Conformation // Macromolecules. 1976. Vol. 9, № 4. P. 535-542.

19. Hadizadeh S., Linhananta A., Plotkin S. S. Improved Measures for the Shape of a Disordered Polymer To Test a Mean-Field Theory of Collapse // Macromolecules. 2011. Vol. 44. P. 6182-6197.

20. Hogberg C. J., Nikitin A. M., Lyubartsev A. P. Modification of the CHARMM Force Field for DMPC Lipid Bilayer // J. Comput. Chem. 2008. Vol. 29. P. 2359-2369.

21. Kuhn W. Uber die Gestalt fadenformiger Molekule in Losungen // Kolloid-Z. 1934. Vol. 68. P. 2-15.

22. Malpass D. B. Introduction to Industrial Polyethylene. Properties, Catalysts, Processes. New Jersey: Wiley-Scrivener, 2010. 156 p.

23. Mark J. E., Erman B., Eirich F. R. Science and Technology of Rubber. 3d ed. Amsterdam etc.: Elsevier, 2005. 762 p.

24. Mark J. E., Erman B. Rubberlike elasticity. A molecular primer. 2d ed. Cambridge etc.: Cambridge University Press, 2007. 272 p.

25. Nelson D. L., Cox M. M. Lehninger Principles of Biochemistry. 5th ed. N. Y.: Freeman W. H. and Co., 2008. Ch. 10. P. 343-370.

26. Peacock A. J. Handbook of polyethylene. Structures, properties, and applications. N. Y.; Basel: Marcel Dekker, Inc., 2000. 539 p.

27. Rabinovich A. L., Ripatti P. O. Monte Carlo simulations of hydrocarbon oligomeric chains. Shape and dimension characteristics // Proc. SPIE. 2001. Vol. 4348. P. 225-236.

28. Rubin R. J., Mazur J. Ordered spans of unrestricted and self-avoiding random-walk models of polymer chains. I. Space-fixed axes // J. Chem. Phys. 1975. Vol. 63, № 12. P. 5362-5374.

29. Rubin R. J., Mazur J. Spans of polymer chains measured with respect to chain-fixed axes // Macromolecules. 1977. Vol. 10, № 1. P. 139-149.

30. S v e 11 e S., B j 0 r g e n M. Selective Catalysts for Petrochemical Industry. Shape Selectivity in Microporous Materials // Selective Nanocatalysts and Nanoscience: Concepts for Heterogeneous and Homogeneous Catalysis / Eds.: A. Zecchina,

S. Bordiga, E. Groppo. Weinheim: Wiley-VCH Verlag GmbH & Co. KGaA., 2011. P 237-269.

31. S o l c К. Shape of flexible polymer molecules // Polymer News. 1977. Vol. 4, № 2. P 67-74.

32. T r e l o a r L. R. G. The physics of rubber elasticity. 3d ed. Oxford; N. Y.: Oxford University Press, 2005. 323 p.

33. Tsolou G., Stratikis N., Baig C., Stephanou P S., Mavrantzas V. G. Melt Structure and Dynamics of Unentangled Polyethylene Rings: Rouse Theory, Atomistic Molecular Dynamics Simulation, and Comparison with the Linear Analogues // Macromolecules. 2010. Vol. 43. P 10692-10713.

34. Van Santen R. A., Neurock M. Molecular Heterogeneous Catalysis: A Conceptual and Computational Approach. Weinheim: Wiley-VCH Verlag GmbH & Co. KGaA, 2006. 488 p.

35. White J. L., Choi D. D. Polyolefins: Processing, Structure Development and Properties. Munich: Carl Hanser Verlag, 2004. 271 p.

36. Yoon D. Y., F l o r y P J. Small angle neutron scattering by n-alkane chains // J. Chem. Phys. 1978. Vol. 69, № 6. P 2536-2538.

Zhurkin D. V., Petrozavodsk State University (Petrozavodsk, Russian Federation) Rabinovich A. L., Institute of Biology of Karelian Research Centre of RAS (Petrozavodsk, Russian Federation)

MONTE CARLO EVALUATION OF HYDROCARBON CHAIN MOLECULES’ SHAPE

The study of physical and chemical properties of hydrocarbon chains is important for development of many technological fields and understanding of the structure and function features of biological systems. The shape of the molecules should be taken into account in the description of various phenomena, such as viscous flow of molecules, adsorption on surfaces, absorption and diffusion in porous systems. In this paper, conformational properties of 65 hydrocarbon chains CH3-(CH2)a-(CH=CH-CH2)d-(CH2)b-CH3 (where a, b, d are integer) were studied at temperatures T = 293, 303 and 313 K, using Monte Carlo simulations. The number of carbons was 16, 18, 20, and 22, the number of cis double bonds was 0, 1, ..., 6. Prototypes of the molecules are residues of fatty acids of phospholipids. The simulations of the molecules were carried out in an unperturbed state; the “classical flexible model” of chain was used. The torsion angles were generated in 0-360 deg range using important sampling technique based on three-wise interdependence of torsions along the chain. All generated conformations were approximated by rectangular parallelepipeds, and average shape characteristics (“spans”) were calculated. Dependencies of the calculated characteristics on the carbon number, the number of double bonds, and their location in the molecules (i. e., interdependencies between structure and properties) were studied. They are caused by certain features of internal rotations in chains of the considered structure.

Key words: chain molecules, unsaturated hydrocarbons, conformational analysis, Monte Carlo simulation, unperturbed state, molecular shape

REFERENCES

1. Volkenshtein M. V. Configuration Statistics of Polymer chains. N. Y.: Wiley-Interscience, 1963. 562 p.

2. G e n n i s R. Biomembranes. Molecular structure and functions. N. Y., etc: Springer, 1989. 533 p.

Оценка формы цепных углеводородных молекул методом Монте-Карло

117

3. Dashevskiy V G. Konformatsii organicheskikh molekul [Conformations of Organic Molecules]. Moscow, Khimiya Publ., 1974. 428 p.

4. Paint and surface coatings: theory and practice / R. Lambourne, eds. N. Y. etc.: Wiley, 1987. 696 p.

5. Landau L. D., Lifshitz E. M. Mechanics. Course of theoretical physics. Vol. 1. Oxford etc.: Butterworth-Heinemann, 1976. 197 p.

6. Surfactants and Polymers in Aqueous Solution / K. Holmberg, B. Jonsson, B. Kronberg, B. Lindman, eds. Chichester: Wiley,

2003. 562 p.

7. Rabinovich A. L., Zhurkin D. V. . Significant sampling for Monte Carlo simulation of a continuous spectrum of conformations of macromolecules [Sushchestvennaya vyborka pri modelirovanii nepreryvnogo spektra konformatsiy makromolekul metodom Monte-Karlo]. Trudy KarNTs RAN. Seriya "Matematicheskoe modelirovanie i informatsionnye tekhnologii” [Mathematical modeling and information technologies]. 2013. Vol. 4. P. 96-111.

8. Rabinovich A. L. Chain molecules as constituents of membrane systems: computer simulations [Tsepnye molekuly kak komponenty membrannykh sistem: komp’yuternoe modelirovanie]. Metody komp’yuternogo modelirovaniya dlya issledovaniya polimerov i biopolimerov [Computer Simulation Methods of Polymer and Biopolimer Researches] / Otv. red.: V. A. Ivanov, A. L. Rabinovich, A. R. Khokhlov. Moscow, Knizhnyy dom LIBROKOM Publ., 2009. P. 410-454.

9. F l o r y P. J. Statistical Mechanics of Chain Molecules. N. Y. etc.: Interscience Publ., 1969. 432 p.

10. Applegate K. R., Glomset J. A. Effect of acyl chain unsaturation on the packing of model diacylglycerols in simulated monolayers // J. Lipid Res. 1991. Vol. 32. P. 1645-1655.

11. Brockman H. L., Applegate K. R., Mo m s en M. M., King W. C., Glomset J. A. Packing and Electrostatic Behavior of sn-2-Docosahexaenoyl and -Arachidonoyl Phosphoglycerides // Biophys. J. 2003. Vol. 85. P. 2384-2396.

12. Cook H. W., M c M a s t e r C. R. Fatty acid desaturation and chain elongation in eukaryotes // New comprehensive biochemistry / G. Bernardi (gen. ed.). Vol. 36: Biochemistry of Lipids, Lipoproteins and Membranes / Eds.: D. E. Vance, J. E. Vance. 4th ed. S. l.: Elsevier, 2002. P. 181-204.