Содержание к диссертации
Введение
Глава I Изучение поведения разбавленного раствора полимерных макромолекул в различных внепіних условиях: обзор литературы 9
1. Теория и компьютерное моделирование разбавленного полимерного раствора 11
1.1. Теория фазовых переходов в одиночной полимерной цепи 11
1.1.1. Переход клубок — глобула 11
1.1.2. Несферические глобулярные структуры 13
1.2. Компьютерное моделирование различных внутримолекулярных структур и переходов между ними 15
2. Исследования полимерных цепей вблизи поверхности 19
2.1. Теория, компьююрное моделирование 19
2.1.1. Адсорбция гибкой цепи 22
2.1.2. Адсорбция жесткоцепных полимеров 26
2.2. Экспериментальные результаты 28
3. Описание состояния полимерной цепи в компьютерном эксперименте 30
3.1. Основные величины, характеризующие состояние полимерной цепи 30
3.2. Определение точек перехода между структурами 34
4. Заключение 35
Глава II Одиночная жесткоцепная макромолекула в пространстве 37
1. Модели 37
1.1. Модель цепи с флуктуирующей длиной векторов связи 37
1.2. Раеппфснный ансамбль, алгоритм Ландау - Ванга по внешнему полю 39
1.3. Алгоритм Ландау - Вапга по энергии для моделирования цепи в объеме 42
2. Диаграмма состояний одиночной жесткоцегаюй гомополимерной молекулы 44
2.1. Диаграмма состояний 44
2.2. Переходы клубок — жидкая глобула — твердая глобула . 47
2.3. Переход клубок — цилиндрическая глобула 51
2.4. Переход сферическая — цилиндрическая глобула 53
3. Переход ктубок - глобула, жидкая — твердая глобула в термодина мическом пределе 57
3.1. Сравнение результатов, полученных с помощью алгоритма Метрополией и алгоритма Ландау Ванга 57
3.2. Определение температур переходов для бесконечно длинной цепи методом конечномерного масштабирования 58
4. Общие замечания 64
Глава III Одиночная жесткоцепная макромолекула вблизи адсор бирующей поверхности 67
1. Описание модели и алгоритма 67
1.1. Модель и техника моделирования макромолекулы вблизи адсорбирующей поверхности 67
1.2. Сравнение результатов, полученных с помощью алгоритма Ландау - Ванга и обычного Монте-Карло 72
2. Диаграммы состояний гибких цеиеіі вблизи адсорбирующей поверхности 74
3. Диаграммы состояний полужестких цепей вблизи адсорбирующей поверхности 85
4. Форма адсорбированной гибкой и жесткой макромолекул 94
5. Заключение 103
Выводы 104
Литература 106
- Теория фазовых переходов в одиночной полимерной цепи
- Основные величины, характеризующие состояние полимерной цепи
- Раеппфснный ансамбль, алгоритм Ландау - Ванга по внешнему полю
- Модель и техника моделирования макромолекулы вблизи адсорбирующей поверхности
Введение к работе
Введение
Диссертационная работа посвящена компьютерному моделированию одиночной жесткоцепной макромолекулы в объеме и вблизи адсорбирующей поверхности. Основной задачей работы является изучение равновесных состояний макромолекулы и переходов между ними для цепей различной длины и жесткости в различных внешних условиях.
Особенностью полимерных молекул является способность принимать огромное число разнообразных пространственных структур - конформаций. Вид образуемых макромолекулой конформаций зависит как от внешних условий (качество растворителя, температуру, наличие поверхности и т.д.), так и от свойств самого полимера (длина, заряд, жесткость, состав и т.д.). В работе изучаются свойства нейтральной гомополимерной молекулы с различной внутрицепной жесткостью.
При изменении внешних условий конформация цепи меняется. Теоретическая фазовая диаграмма для жесткоцепных макромолекул в координатах "жесткость - температура" содержит область, в которой тороидальная конформация является стабильной для цепей конечной длины. Большой интерес вызывает вопрос о характере возможного внутримолекулярного жидкокристаллического (ЖК) упорядочения в плотной глобуле, а также образование различных несферических глобулярных конформаций, таких как конформация "ракетки", цилиндрической и тороидальной глобул. Нерешенной проблемой остается вопрос о существовании тех или и-іьіх структур в термодинамическом пределе, когда длина макромолекулы N —> со, в частности, задача о существовании и стабильности тороидальных и цилиндрических глобул, фазы жидкой глобулы для бесконечно длинных макромолекул. Поэтому важным шагом в исследовании переходов является переход к асимптотике N —> со на основе данных, полученных для цепей конечной длины.
Целью диссертационной работы является исследование конформационных переходов и построение полной диаграммы состояний для цепей различной длины и жесткости в объеме и вблизи адсорбирующей поверхности с помощью методов компьютерного моделирования.
Актуальность обусловлена важностью задачи установить, как влияют локальные внутримолекулярные характеристики, в частности, жесткость полимерной цепи, и наличие пространственных ограничений, в частности, плоской поверхности, на макроскопические свойства полимерного раствора. Несмотря на значительные успехи, многие детали процессов переходов между структу-
Введение
рами остаются неясными. В эксперименте сложно изучать свойства одиночной макромолекулы в объеме. Многие экспериментальные работы проводятся с молекулами, находящимися на поверхности (например, сканирующая зондовая микроскопия), поэтому вопрос о том как меняются конформационные свойства полимера и линии переходов между структурами при уменьшении размерности пространства, то есть при адсорбции на плоскую поверхность, остается важным и актуальным для теории и компьютерного моделирования. Далее, большинство биополимеров являются жесткоцепными макромолекулами, которые взаимодействуют с мембранами; учитывая, что молекула ДНК способна образовывать тороидальную и цилиндрическую глобулы, проблема существования этих глобул на поверхности имеет большое значение и требует дополнительных исследований. Интересен вопрос о конкуренции между адсорбцией и коллапсом цепи, то есть о том какой из этих переходов происходит первым при ухудшении качества растворителя. Актуальным является систематическое исследование влияния жесткости на оба эти перехода.
Практическая ценность работы заключается в согласии полученных с помощью компьютерного моделирования результатов с имеющимися экспериментальными данными, что подтверждает применимость выбранной модели для описания исследуемых полимерных систем и открывает возможность для последующей направленной постановки новых экспериментов. С другой стороны методическая ценность работы состоит в разработке нового алгоритма моделирования методом Монте-Карло с использованием расширенного ансамбля.
Научная новизна работы состоит в следующем:
Впервые построена полная диаграмма-состояний одиночной-жесткоцеп-ной макромолекулы в объеме, включающая в себя области стабильности клубка, жидкой сферической и твердой (кристаллической) сферической глобул, а также область цилиндрической и тороидальной глобул.
Разработан метод Монте-Карло для расширенного ансамбля с использованием четырехмерного пространства и алгоритма Ландау-Ванга для ускорения уравновешивания плотных глобулярных структур.
Впервые построена полная диаграмма состояний для жесткоцепной макромолекулы вблизи адсорбирующей поверхности в координатах температура - сила притяжения к поверхности. Диаграмма содержит области клубка, жидкой и твердой сферических глобул, адсорбированного клубка и вытянутой двумерной (жидко)кристаллической глобулы.
Введение
Впервые показано, что адсорбирующая поверхность способствует образованию цилиндрической анизотропной глобулы для случая малой внутри-цепной жесткости, в то время как та же макромолекула в объеме образует сферическую изотропную глобулу.
На защиту выносятся следующие положения:
Полная диаграмма состояний одиночной жесткоцепной макромолекулы в объеме содержит области существования клубка, жидкой глобулы, твердой глобулы, области цилиндрической и тороидальной глобул. Переход клубок — глобула для жестких цепей происходит по типу фазового перехода 1 рода. Переходы жидкая — твердая глобула и сферическая твердая - (жидко)кристаллическая цилиндрическая глобула также являются переходами 1 рода. С увеличением длины цепи линия перехода жидкая — твердая глобула сдвигается в область больших температур.
Метод расширенного ансамбля с использованием алгоритма Ландау-Ванга в четырехмерном пространстве эффективен для моделирования переходов между плотными структурами (переход жидкая — твердая глобула), в то время как для моделирования неплотных структур этот метод дает неправильные оценки точек переходов.
Температура перехода жидкая — твердая глобула в случае бесконечно длинной цепи сдвигается в сторону меньших температур с увеличением жесткости. Фаза жидкой гл^оулы для полимеров с небольшой жесткостью исчезает в случае бесконечно длинной цепи, также как и для гибких цепей.
Диаграммы состояний для гибких полимерных макромолекул вблизи адсорбирующей поверхности для цепей длиной 64 и 128 мономерных звеньев содержат области стабильности клубка, жидкой и твердой глобул, адсорбированного клубка, а также адсорбированной кристаллической глобулы.
Диаграммы состояний для полужестких полимерных макромолекул вблизи адсорбирующей поверхности для цепей длиной 64 и 128 мономерных звеньев содержат области стабильности клубка, жидкой и твердой глобул, адсорбированного клубка, а также адсорбированной (жидко)кристаллической глобулы. Линий переходов сдвигаются в область более высоких температур с увеличением длины цепи.
Введение
6. Результаты моделирования с использованием алгоритма Ландау-Ванга хорошо совпадают с результатами, полученными при моделировании с помощью алгоритма Метрополиса.
Работа имеет следующую структуру:
первая глава содержит обзор литературы, касающейся теории фазовых переходов, методов моделирования полимерных систем, результатов моделирования различными методами; описание основных величин, которые характеризуют изучаемую систему; рассматривается современное сосюя-ние исследований свойств одиночной цепи теоретическими, экспериментальными методами и с помощью компьютерного моделирования.
вторая глава посвящена исследованию с помощью компьютерного моделирования состояний цепи в трехмерном пространстве с использованием расширенного ансамбля для изучения плотных глобулярных структур; переходов между конформациями цепей различной длины и жесткости, изучение перехода клубок — глобула в термодинамическом пределе.
в третьей главе представлены результаты моделирования цепи вблизи адсорбирующей поверхности. Проводится сравнение методов моделирования. Исследуются переходы для двух длин цепей (N = 64,128), делается попытка изучить адсорбционный переход в термодинамическом пределе. Рассматривается зависимость линий переходов от длины и жесткости цепи при изменении силы притяжения к поверхности.
Апробация работы проводилась на следующих конференциях:
Школа "Forces, Growth and Form in Soft Condensed Matter: At the Interface between Physics and Biology" (Гейло, Норвегия, 23 марта - 2 апреля, 2003).
Конференция студентов и аспирантов по химии и физике полимеров и тонких органических пленок (г. Тверь, 28 - 30 мая 2003).
Зимняя школа "Computational Soft Matter: From Synthetic Polymers to Proteins" (Бонн, Германия, 29 февраля - 6 марта, 2004).
Третья Всероссийская Каргинская конференция "Полимеры — 2004" (Москва, 27 января - 1 февраля, 2004).
Конференция студентов и аспирантов по химии и физике полимеров и тонких органических пленок (г. Солнечногорск, 16-17 сентября 2004).
Введение
Школа "Computer Simulations in Condensed Matter: from Materials to Chemical Biology" (Эриче, Италия, 20 июля - 1 августа 2005).
Международный симпозиум "Statistical Mechanics of Polymers: New Developments" (г. Москва, 6-11 июня 2006)
Всероссийская Школа по математическим методам для исследования полимеров и биополимеров (г. Петрозаводск, база отдыха "Деревня Алек-сандровка", 13 - 17 июня 2006).
Четвертая Всероссийская Каргинская конференция "Наука о полимерах 21-му веку" (Москва, 29 января - 2 февраля, 2007).
Результаты опубликованы в следующих работах:
J. A. Martemyanova, М. R. Stukan, V. A. Ivanov, М. Muller, W. Paul, К. Binder, Dense orientationally ordered states of a single semiflexible macromolecule: An expanded ensemble Monte Carlo simulation//J. Chem. Phys. - 2005.-Vol. 122, p. 174907.
Ю.А. Мартемьянова, M.P. Стукан, В.А. Иванов, Изучение конформаций одиночной жесткоцепной макромолекулы методом компьютерного моделирования//Вестник МГУ, Серия 3. Физика. Астрономия - 2005.- Том. 3, р. 58-60.
V.A. Ivanov, J.A. Martemyanova, Monte Carlo Computer Simulation of a Single Semi-Flexible Macromolecule at a Plane Surface//Macromol. Symp. - 2007.- Vol.252, с 12-23.'
Теория фазовых переходов в одиночной полимерной цепи
На данный момент существует хорошо развитая теория перехода клубок — глобула как для гибких, так и для жестких цепей. Хорошо разработана терминология для описания переходов, происходящих в макромолекуле конечной длины и в термодинамическом пределе, когда N — со [30-34]. Конформационный переход называется фазовым, если его ширина AT стремится к нулю при JV — со. Фазовый переход является фазовым переходом I рода, если в области перехода существует два-минимума свободной энергии в-пространстве-макропеременных и, соответственно, два устойчивых состояния, каждое из которых является термодинамически стабильным по одну сторону от точки перехода и метаста-бильным по другую ее сторону. Фазовый переход является фазовым переходом II рода, если в области перехода существует лишь один минимум свободной энергии в пространстве макропеременных, т.е. не существует метастабильных состояний по другую сторону от точки перехода.
В работе [30] было показано, что переход клубок — глобула является фазовым переходом, так как его ширина AT —» 0 при N —» со. В случае гибкой цепи переход клубок — глобула происходит как плавный фазовый переход II рода, и температурный интервал перехода растянут на всю -область. Для жестких цепей переход клубок — глобула происходит как резкий фазовый переход I рода с конечным скачком плотности в относительно узком интервале температур, точно определенном и заметно удаленном от #-точки. Однако по некоторым своим чертам этот переход близок к переходу II рода, в частности, из-за малости теплоты перехода.
В серии работ А.Ю. Гросберга и Д.В. Кузнецова [31-34] было проведено детальное рассмотрение перехода клубок — глобула с учетом жесткости цепи. Все результаты получены для линейного незаряженного гомополимера. В этих работах были теоретически рассчитаны радиус энерции, гидродинамический радиус и их флуктуации. Были детально изучены бинодальный и спинодальный распад. Важным отличием глобулы от клубка является то, что в малом объеме глобулы содержатся некоррелированные участки цепи, что напоминает раствор или смесь отдельных цепей.
В клубке (по Флори) осмотическое давление компенсируется энтропийной упругостью полимера. Более того, это давление равно нулю во всем объеме глобулы, плотное ядро глобулы имеет однородную плотность По, то есть концентрация звеньев в глобуле должна установится такой, чтобы на соответствующем ей характерном расстоянии между звеньями притягивающие и отталкивающие объемные взаимодействия уравновешивали друг друга - в отличие от клубка, где объемное отталкивание звеньев уравновешивается цепными связями.
Существует два приближения для описания флуктуации радиуса инерции. Для жестких цепей-(С1/,2/а3- 0.05, где С - третий вириальный коэффициент; a - расстояние между бусинками в стандартной Гауссовой модели) переход клубок — глобула происходит по типу фазового перехода первого рода, при этом наблюдается присутствие двух минимумов свободной энергии и бимодальное распределение по состояниям. Для гибких цепей (С1//2/а3 0.05) переход относительно плавный и выглядит как переход второго рода.
Свободная энергия в области низких температур может иметь несколько минимумов, это связано, в частности, с существованием различных фаз в системе разорванных мономерных звеньев [30]. Уравнение (1) может иметь несколько решений по, которые соответствуют глобулам с различными плотностями в ядрах. Если за счет изменения температуры разность энергий двух состояний изменит знак, то произойдет фазовый переход I рода глобула — глобула, связанный с перестройкой структуры ядра.
Флуктуации в случае жестких цепей обусловлены спонтанными переходами системы из одного состояния в другое, с другой стороны для гибких цепей флуктуации связаны с увеличением ширины пика функции распределения. Несмотря на существенную разницу для двух случаев, оба приближения предсказывают существование максимума флуктуации в районе перехода, которые наблюдаются экспериментально.
Кроме того, в переходной области существенную роль играют тройные взаимодействия [35]. Длинные и короткие цепи ведут себя по разному в б-точке. Одни и те же тройные взаимодействия приводят к набуханию длинной макромолекулы и к сжатию клубка короткой цепи. Однако следует отметить, что эти выводы были получены для приближения жесткой цепи. Стоит отметить тот факт, что клубок, образованный жесткой цепью, не является истинно гауссовым клубком, что связано с присутствием в клубке на малых масштабах почти прямолинейных участков цепи. Таким образом, гауссово распределение звеньев цепи достигается только на масштабах значительно больших размеров кунов-ского сегмента макромолекулы [36].
Впервые вопрос об устойчивости несферических глобулярных структур, в том числе имеющих внутриглобулярный ориентационный порядок звеньев, был теоретически рассмотрен в работе [37]. В работах [38,39] теоретически изучался процесс коллапса жесткоцепной макромолекулы в плохом растворителе (моделировался процесс коллапса цепи ДНК) и проводилось сравнение полученных данных с экспериментом. В ходе экспериментальных наблюдений при помощи электронного микроскопа было обнаружено, что длинная одиночная молекула ДНК может образовывать тороидальную структуру в качестве устойчивой глобулярной конформации. Это происходит вследствие того, что в процессе коллапса жесткая макромолекула ДНК принимает конформацию, которая удовлетворяет двум условиям: каждое мономерное звено имеет максимально возможное число соседей, и, одновременно, в цепи нет резких изломов. В исследовании была предпринята попытка определить условия, при которых тороидальная конформация является стабильной.
Для длинных молекул ДНК область существования тороидальной конформации достаточно узкая, и ее ширина уменьшается с ростом длины цепи как AT l/(L//)2, где L - контурная длина цепи, а I - персистентная длина цепи. Для длинных макромолекул радиус тора может быть меньше, чем персистент,-ная длина цепи (R I). Напротив, для коротких цепей возможна ситуация (R L). Было показано, что переход от конформации набухшего клубка к конформации плотной глобулы происходит по типу фазового перехода первого рода. Профиль свободной энергии в районе перехода имеет два минимума, поэтому в районе перехода сосуществуют конформации клубка и глобулы.
Для сравнения использовались экспериментальные данные полученные при наблюдении под электронным микроскопом бактериофага Т4-ДНК и Л-ДНК. Область стабильности тороидальной конформации была также получена при построении диаграммы состояний жесткоцепной макромолекулы путем численного решения уравнений в приближении самосогласованного поля [40,41].
Основные величины, характеризующие состояние полимерной цепи
Положение границ перехода и вид конформаций определялись с помощью величин, которые рассчитывались в ходе моделирования. Прежде всего определялись размеры цепи как целого, а именно, средний квадрат расстояния между концами цепи (Щ) = ((гдг — Т\)2 ) и среднеквадратичный радиус инерции цепи N х N Rl==N (ГІ"" Лсм)2 ЙСМ = N Гг - (?) г—\ г=1 Вычислялся тензор инерции N XaP = luYl &iSi & = fi- &см\ a,0 = x, y, z, i=i (8) здесь fi означает положение г-го мономерного звена вдоль по цепи, ( Щ ) = (Ххх ) + (ХуУ ) + {Xzz ). При моделировании цепи вблизи поверхности удобно было рассматривать компоненты среднеквадратичного радиуса инерции цепи ( Щх ) = (Ххх ),(R2gy) = (Хуу), {Щж ) = {Xzz), которые являются тремя главными моментами тензора инерции. По z - компоненте, например, хорошо определялся адсорбционный переход. При использовании метода Монте-Карло под средним подразумевается усреднение по всем значениям величины, полученной в течение моделирования для каждой конформации, и усреднение по всем различным конформациям.
Соответствие между значением параметров и формой цепи следующее: идеальная сфера К\ = K i = 1; бесконечно тонкий стержень К\ = Q,K% = 1; бесконечно тонкий диск К\ = К2 = 1/2 (рисунок 10). Как можно видеть эти параметры очень полезны при определении формы конформации. Для того чтобы определить ориентационный порядок векторов связей внутри цепи, вычислялся ориентационных тензор: Q = j ilUs($e!)-sJ)- (Ю) 1=1 V / Собственные значения этого тензора есть параметры ориентационного порядка 1, 2, 3- И, наконец, вычислялся локальный параметр порядка, характеризующий среднее ориентационное упорядочение вокруг мономерных звеньев: = 3 2%) (11) Здесь $ij - угол между векторами связи і и j, находящимися на близком расстоянии друг от друга в пространстве, а внутренняя сумма по j берется по N векторам связи, являющимися соседями в пространстве для вектора і (был задан радиус для определния критерия попадания в список соседей).
Кроме того строились гистограммы и температурные зависимости для следующих параметров: Щ, В , 771, Щ, Щ, Цос, Кі, if2, полная энергия цепи (Е/иц), угловая энергия (Еа), энергия объемного взаимодействия [Ец), энергия притяжения к поверхности {Eadsorb) и ее вклады: пристеночный {Ewau) и дально-действующий (Ее[г) (см. определения ниже). Очень информативными являются температурные зависимости флуктуации квадрата радиуса инерции, энергетических составляющих и полной энергии и других величин.
Для анализа внутренней структуры цепи рассчитывался профиль плотности. Профиль плотности р(г) проще всего получить путем усреднения локальной плотности в сферических слоях dr, находящихся на расстоянии г от центра масс. Профиль также строится после усреднения по всем независимым конфор-мациям.
Вводились характеристики, позволяющие определить направления векторов связи в сферическом слое dr на расстоянии г от центра масс для описания внутренней структуры глобулы, а именно, первый и второй полиномы Лежанд-ра между векторами связи и вектором, проведенным из центра масс в точку приложения данного вектора связи, cos# и (3cos2# — 1)/2, соответственно.
Краевой угол а равен углу, образованному радиусом, проведенным в точку, где граничат три фазы, и радиусом, проведенным перпендикулярно поверхности. Касательная к окружности в точке, где граничат 3 фазы, перпендикулярна радиусу, проведенному в точку касания, а линия поверхности перпендикулярна Я, поэтому краевой угол а равен углу, образованному радиусом, проведенным в точку, где граничат три фазы, и радиусом, проведенным перпендикулярно поверхности. Принимая во внимание, что Щ = Щх + Щу + R2gz — R2gxy + R2gz и предполагая, что Н = 2y/Rjz, г — у/Щх = уЩу = у "ipS рассчитываем значение краевого угла для случая частичного смачивания поверхности (рис. 11а) следующим образом.
Для перехода, который происходит по типу фазового перехода первого рода, характерна бимодальная структура на гистограмме параметра порядка. Точка перехода между структурами определяется, исходя из равенства площадей под максимумами на гистограмме параметра порядка. Если переход идет по типу фазового перехода второго рода, то происходит плавный сдвиг максимума на гистограмме параметра порядка. Можно определить точку перехода и как среднее положение максимума гистограмм параметра порядка, так как температурный интервал, в котором происходило это смещение максимума в процессе перехода, был не очень широким. В этом случае точка перехода определялась из максимумов на температурных зависимостях флуктуации энергии или параметре порядка.
Соответствие этих максимумов той или иной структуре определялось путем визуализации и анализа отдельных максимумов. Для различных переходов выбирались определенные параметры порядка, например, как известно, переход клубок — глобула, а также переход жидкая — твердая глобула, хорошо характеризуются такими параметрами порядка как радиус инерции или число контактов между мономерными звеньями. Переход сферическая — цилиндрическая глобула хорошо определяется с помощью ориентационного параметра порядка Г1\(Ъ), локального ориентационного параметра порядка r]ioc(b), структурных параметров, формы глобулы К\ и Кч В случае перехода клубок — цилиндрическая глобула наблюдался гистерезис на зависимости энергии контактов от температуры. Поэтому определение температуры перехода свелось к определению средней температуры при переходе из клубка в цилиндрическую глобулу и обратно. Эти температуры различны, поэтому при таком определении температуры перехода ошибкой является раз брос значений температуры при переходе из клубка в цилиндрическую глобулу и обратно.
При моделировании с помощью алгоритма Ландау - Ванга точки перехода определялись с помощью производных параметров порядка по температуре или с помощью температурных зависимостей флуктуации параметров порядка. Температура, соответствующая максимуму или минимуму производной параметра порядка, который характеризует определенный переход, является температурой перехода. Погрешность определялась как полуширина на 0.7 высоте максимума или минимума. Также проводился пересчет и строились гистограммы параметра порядка при температуре перехода. Если переход первого рода, например, жидкая — твердая глобула, то наблюдалась бимодальная гистограмма, как было сказано выше. В случае перехода 2 рода происходил сдвиг максимума, и температура перехода определялась средним положением максимумов. Значение температур перехода, определенные выше описанными способами, совпадали в пределах погрешности.
Как следует из приведенного обзора литературы, на данный момент в физике полимеров существует достаточно хорошо развитая теоретическая база, позволяющая описывать поведение одиночных макромолекул. Также к настоящему моменту разработано большое количество методов, позволяющих производить численное моделирование и определять значения различных параметров полимерных систем. Однако, многие вопросы, связанные с поведением одиночных цепей в объеме и вблизи адсорбирующей поверхности, до сих-пор остаются открытыми.
К настоящему моменту было проведено большое количество исследований, которые показали, что для достаточно длинных жесткоцепных макромолекул тороидальная и цилиндрическая глобулы стабильны, однако до конца не ясно, при каких условиях стабильна та или иная кон формация, происходит ли переход жидкая — твердая глобула для жесткоцепных полимеров, и как меняются температуры переходов с увеличением жесткости цепи. Поведение жесткоцепных макромолекул в области низких температур также вызывает большой интерес.
Раеппфснный ансамбль, алгоритм Ландау - Ванга по внешнему полю
Для уменьшения времени уравновешивания плотных глобулярных конформации полимерной цепи применяется расширенный ансамбль в четырехмерном пространстве. Вводится четвертая координата, то есть рассматриваются две трехмерные подрешетки, в одной из которой четвертая координата для мономерных звеньев равна 0, в другой 1. Мономерные звенья с разными значениями четвертой координаты не взаимодействуют друг с другом. Гамильтониан такой системы представлен уравнением = + Х « (19) где Но гамильтониан обычной трехмерной системы, ж4(г) - четвертая координата г-ой частицы, принимающая значения 0 или 1. Внешнее поле h контролирует распределение мономерных звеньев между двумя подрешетками. Для h — О мономерные звенья распределены равномерно между двумя подрешетками. Для h—»00 все мономерные звенья имеют четвертую координату равную ХА{І) = 0, то есть находятся в обычном трехмерном пространстве. При уменьшении внешнего поля некоторые звенья оказываются во второй подрешетке, в то время как в первой подрешетке (трехмерие) уменьшается эффективная плотность структуры. Взаимодействие между мономерными звеньями в разных подрешетках отсутствует. Таким образом различные части полимерной глобулы лучше перемешиваются, и укоряется процесс получения независимых конформаций.
В расширенном ансамбле элементарным шагом Монте-Карло являются следующие шаги: а) локальное перемещение звена, б) рептация цепи, в) переключение случайно выбранного мономерного звена между подрешетками (т.е. изменение значения 4-ой координаты), г) случайное изменение значения внешнего поля h. Значение поля h меняется только на значения, которые являются соседними к данному, чтобы не нарушать равновесия системы. Изначально значения поля выбираются таким образом, чтобы происходило перекрывание гистограмм числа мономерных звеньев Пы с четвертой координатой х± = 1 (рис. 13), полученных при h = const. Все шаги Монте-Карло принимаются или отвергаются в соответствии с критерием Метрополиса [124].
Весовые коэффициенты g(h) рассчитываются с помощью алгоритма Ландау - Ванга (см. Приложение 1). Алгоритм устроен таким образом, что получающиеся в результате значения g(h) позволяют системе посещать с равной частотой конформаций, отвечающие различным значениям поля h. Вероятность изменения значения "внешнего поля" равна 9(h) , 9(h) (23) p(h\ — h-i) — vain Функция g(h) "регулирует" количество мономерных звеньев, находящихся в "четвертом измерении". Из всех получающихся конформаций выбираются только трехмерные и по ним проводи .ся усреднение наблюдаемых величин.
В качестве стартовых использовались конформаций, полученные с использованием алгоритма блуждания цепи без самопересечения. Эти конформаций приводились к равновесию при фиксированных значениях параметров жесткости и температуры, после чего проводилось усреднение интересующих нас величин по всем конформациям. Выбор момента времени для начала усреднения (и всего времени моделирования) проводился с таким расчетом, чтобы во второй половине временного интервала не наблюдалось систематического отклонения (в пределах статистической ошибки) среднего значения интересующих параметров от значений, полученных в первой половине. Для получения уравновешенного состояния системы гистограммы наблюдаемых величин не должны изменять своей формы на протяжении последних 106 MCS.
Еще раз отметим, что во второй главе при построении диаграммы состояний для цепи длиной 256 мономерных звеньев в качестве параметра жесткости был использован параметр b = еа/квТ. При фиксированном значении b не происходит эффективного увеличения жесткости цепи при понижении температуры. То есть в отсутствие потенциала объемного взаимодействия распределение цепей по конформациям будет одинаковым при любых значениях температуры.
Для сравнения полученных нами результатов с уже существующими необходимо отметить, что параметр жесткости b связан с используемым в литературе [4] параметром р — Ік/d (1к - длина сегмента Куна, a d - диаметр цепи) как pw 2.16+ 5.0 [48].
Для моделирования цепи в объеме и определения температур перехода клубок — жидкая глобула — твердая глобула в термодинамическом пределе для жест-коцепных макромолекул мы использовали алгоритм Ландау - Ванга по полной энергии. В данном случае полная энергия является суммой энергии контактов и энергии жесткости. Потенциалы остаются такими же как и в случае, описанном выше (уравнения (К)), (1 ))- Однако вместо фиксирования параметра жесткости b = єа/квТ, мы фиксировали параметр энергии изгиба цепи єа. Это озна-чает, что длина сегмента Куна меняется, а энергия жесткости явно не зависит от температуры, то есть при понижении температуры происходит эффективное ожестчение цепи. Для изменения конформации цепи мы использовали локальные шаги и шаги Монте - Карло с конфигурационным смещением выборки (cortfigurational bias 1 (СВ)): выбиралось случайным образом вдоль по цепи мономерное звено, уда , лялась случайно выбранная часть цепи, начиная от выбранного мономерного , звена, и достраивалась недостающая часть цепи. То есть мы вычисляли вес конформации с удаленной частью цепи и)0ц и вес новой конформации с постро енной частью цепи w w В нашем случае один шаг Монте - Карло — это 10 f локальных шагов и 1 шаг с конфигурационным смещением выборки (СВ). Моделирование проводилось с, использованием алгоритма Ландау - Ванга (см. Приложение 1), который основан на получении плоской гистограммы пол ной энергии системы. Мы считали функцию плотности состояний д(Е) отдельно в нескольких интервалах во всем диапазоне возможных значений полной энергии. В течение моделирования в отдельных интервалах мы накапливали двумерную гистограмму полной энергии и радиуса инерции Н(Е,Щ), а также других величин, которые описывают нашу систему, например, Щ, Econt, Estljf и т.д. После того как в каждом энергетическом интервале мы получили функцию плотности состояний и накопшш необходимую статистику в двумерных гистограммах, мы сшивали функции плотности состояний для всех интервалов и получали эту функцию на полном интервале энергий. Далее сшивали файлы гистограмм и пересчитывали температурные зависимости и гистограммы величин для NVT ансамбля (см. Приложение 2). Гистограмма считалась равномерной, если отклонения от среднего значения гистограммы было не больше 20%.
Рассмотрим вопрос о принятии шага для локального и СВ шагов в алгоритме Ландау - Ванга. Выражение для принятия локального шага в алгоритме, который мы использовали, отличается от выражения в методе обычного Монте - Карло. В данном случае при накоплении равномерной гистограммы по полной энергии использовать логарифм функции плотности состояний \п(д(Е)) удобнее нежели весовую функцию д{Е), так как значения д(Е) достаточно большие, а взяв логарифм от этой функции, в памяти компьютера нужно хранить не такие большие числа. Кроме функции плотности состояний больше ничего не влияет на принятие конформации после совершение локального шага в этом алгоритме, т.е. вероятность принятия шага является отношением функции плотности состояний для энергии старой конформации (E0id) к функции плотности состояний для энергии новой конформации (Enew) и записывается следующим образом: pacc = min {1, д{ЕоЫ)/g(Enew)} . (24) А при получении новой конформации с помощью СВ шага, надо вычислять добавочные весовые функции при удалении и построении части цепи. Они вычисляются следующим образом: выбирается случайным образом мономерное звено і вдоль по цепи и удаляется конец цепи, начиная с конца цепи и заканчивая выбранным мономерным звеном. Считаем число к\ свободных (не занятых мономерными звеньями) узлов решетки, вокруг мономерного звена, в которых это мономерное звено могло бы находиться.
Модель и техника моделирования макромолекулы вблизи адсорбирующей поверхности
Как и для моделирования цепи в пространстве, поведение одиночной жестко-цепной макромолекулы вблизи адсорбирующей поверхности мы изучали с использованием решеточной модели с флуктуирующей длиной связи между мономерными звеньями [14-16,29,48,50,51,122,127]. Один конец цепи закреплён на поверхности (z = 0). Для цепи длиной 32 мономерных звена координаты точки закрепления цепи на поверхности (x,y,z)=(30,30,0) размер ячейки моделирования 60 х 60 х 80, для цепи N = 64 точка закрепления цепи (x,y,z)=(40,40,0) размер ячейки моделирования 80 х 80 х 150, для цепи длиной N = 128 точка закрепления (x,y,z)=(75,75,0) в ячейке размером 150 х 150 х 300). Размер ячейки вдоль осей х к у был выбран одинаковым Lx = Ly = 60,80 или 150 (использовались периодические граничные условия), размер ячейки вдоль оси z был взят достаточно большим, чтобы избежать влияния второй плоскости z = L 2N, так как мы изучаем адсорбцию цепи на поверхность, а не состояние цепи в слое между двумя плоскостями. Потенциал взаимодействия мономерных звеньев между собой и потенциал жесткости такие же как и в предыдущей главе. ТТ t\ j -є, r = 2, VE, л/6, .... lWr) = 0) г (29) W/(0O = ee(l-cos0t), (ЗО) где вг это угол между последовательными векторами связи 1г и Іг+і (рис. 12), єа - это параметр энергии изгиба. Потенциал взаимодействия с поверхностью выражается следующим образом: _4ш Uad3(z) = J z Z (31) " "Ш) — - ) J- где z расстояние от мономерного звена до поверхности (выраженное в единицах длины ребра элементарного куба решетки), то есть адсорбирующая поверхность - это плоскость z = 0. Мономерные звенья, имеющие z координату равную 0 или 1, имеют наименьшее значение энергии взаимодействия с поверхностью. Мы выбрали притягивающую часть потенциала адсорбции такую же как в статье [63]. Этот вид потенциала является результатом интегрирования стандартного потенциала Леннард-Джонса между мономерными звеньями и отдельными атомами под поверхностью в трехмерном полупространстве.
Параметр взаимодействия мо:юмерных звеньев между собой мы положили равным единице є = 1, таким образом мы зафиксировали единицы, в которых вычисляется энергия в нашей системе. Температура рассматривается нами, как параметр, по которому мы определяем переход клубок — глобула, то есть этот параметр отвечает за качество растворителя. Также вводим энергетический параметр притяжения мономерного звена к поверхности sw, как параметр, по которому определяем адсорбционный переход. Жесткость определяется параметром єа. В этом случае энергия жесткости не зависит явно от температуры, а сегмент Куна меняется при понижении или повышении температуры (во 2 главе, мы оставляли постоянным параметр b = єа/квТ, то есть энергия жесткости зависела от температуры, а сегмент Куна оставался постоянным при изменении температуры [42,48,50,51]). Таким образом при моделировании нашей системы у нас было 3 входных параметра: температура Т, параметр энергии изгиба еа и параметр притяжения мономерных звеньев к поверхности ew.
При выборе конкретных значений входных параметров мы опирались на ранее опубликованные работы и полученные данные для следующих физических величин: 9 точка [29, 52,125,128] или точка перехода клубок — глобула и зависимость этой точки от жесткости [48,50,51,128]; точка перехода жидкая — твердая глобула [29,52,128] и ее зависимость от жесткости [128], параметры внутриглобулярного ориентационного перехода [50,51,128]; параметры перехода из изотропного в нематически упорядоченное состояние в разбавленных растворах [42]. Значения входных параметров мы выбирали так, чтобы примерно воспроизводились все вышеупомянутые результаты.
Моделирование проводилось с использованием алгоритма Ландау - Ванга (см. Приложение 1). В отличие от второй главы, где полная энергия являлась суммой энергии контактов и энергии жесткости, в данном случае полная энергия состоит из трех вкладов - еще добавляется энергия взаимодействия с поверхностью. В течение моделирования в отдельных интервалах мы накапливали двумерную гистограмму полной энергии и радиуса инерции Н(Е,Щ), а также других величин, которые описывают нашу систему, например, Rl, Econt, Estiff, Ewaii и т.д. Потом мы пересчитывали данные для NVT ансамбля (см. Приложение 2). Для изменения конформации цепи мы использовали локальные шаги и шаги Монте - Карло с конфигурационным смещением выборки (configurational bias (СВ)) как и при исследовании свободной цепи алгоритмом Ландау - Ванга, шаги принимались как в разделе 1 .. . В нашем случае один шаг Монте - Карло это 10 локальных шагов и один СВ шаг. В файл гистограмм значения величин выписывались после каждых 10 тыс. шагов.
На рис. 42 и рис. 4Н приведены графики для зависимостей энергии контактов и квадрата радиуса инерции от температуры для двух независимых случаев при моделировании алгоритмом Ландау-Ванга гибкой (еа — 0) цепи длиной 128 мономерных звеньев, притяжение к поверхности отсутствует. Видно, что зависимости энергии контактов хорошо накладываются друг на друга, в то время как зависимости квадрата радиуса инерции немного расходятся при достаточно больших іемпературах. То есть при моделировании состояний с маленькой энергией (глобулярных состояний при низкой температуре) совпадение хорошее. При повышении температуры цепь находится в состоянии клубка с большой энтропией и видно небольшое расхождение в значениях квадрата радиуса инерции (порядка 10-20% на различных вкладах в квадрат радиуса инерции), однако на определение температуры перехода и на характер поведения цепи при повышении температуры это не влияет.
Для цепей с большей жесткостью и цепей большей длины моделирование проводилось отдельно в нескольких интервалах по полной энергии, как пра вило ширина интервала была равна 200 единицам полной энергии. Соседние интервалы пересекались, ширина пересечения была около 20 единиц. Это делалось для того, чтобы убедиться, что функция плотности состояний гладкая и у нее нет изломов. На рисунках 44 и 15 показаны логарифмы весовых функций для гибкой цепи и цепи с небольшой жесткостью ea = 4, соответственно. Эти функции гладкие, с их помощью достигается равномерное посещение состояний с различной полной энергией. Так как на графиках представлен логарифм весовой функции, любой график можно сдвигать он оси у вверх или вниз, учитывая формулы пересчета в приложений 2.
Проводилось сравнение результатов, полученных с помощью алгоритма Ландау - Ванга, с результатами, полученными моделированием методом Монте -Карло с использованием критерия Метрополиса. На рис. 46 и рис. 47 приведены зависимости энергии контактов при моделировании гибкой цепи N = 128, еа = 0, еш — 2 и энергии жесткости при моделировании цепи iV = 128, еа = 4, Для гибких цепей длиной N = 64,128 моно мерных звеньев были построены диаграммы состояний на плоскости значений параметров квТ и ew(a) и г и (б), где ew - энергия взаимодействия мономерного звена и поверхности в первых двух слоях у поверхности, кв - постоянная Больцмана, Т - температура (рис. 51). Сплошные черные линии соответствуют линиям перехода для цепи длиной 64 мономерных звена, штриховые красные линии - цепи длиной 128 мономерных звеньев, кружками обозначены точки переходов между состояниями клубка, жидкой и твердой глобул макромолекулы, квадратами обозначен адсорбционный переход. Наблюдаются области существования трехмерного (область I на диаграмме) и адсорбированного(IV) клубков, жидкой(Н), твердой(Ш) сферических глобул, область V соответствует области адсорбированной кристаллической глобулы. Для цепи с одним концом, закрепленным на поверхности, при ew 0 как и для свободной цепи при увеличении ее длины температура переходов жидкая - твердая глобула, клубок - жидкая глобула повышается. Как и в 3D случае переход жидкая — твердая глобула происходит по типу фазового перехода 1 рода, а переход клубок — жидкая глобула по типу фазового перехода 2 рода. Так как цепь гибкая (ea = 0), то форма глобул сферическая в объеме и диско-подобная на поверхности.
