авторефераты диссертаций БЕСПЛАТНАЯ БИБЛИОТЕКА РОССИИ

КОНФЕРЕНЦИИ, КНИГИ, ПОСОБИЯ, НАУЧНЫЕ ИЗДАНИЯ

<< ГЛАВНАЯ
АГРОИНЖЕНЕРИЯ
АСТРОНОМИЯ
БЕЗОПАСНОСТЬ
БИОЛОГИЯ
ЗЕМЛЯ
ИНФОРМАТИКА
ИСКУССТВОВЕДЕНИЕ
ИСТОРИЯ
КУЛЬТУРОЛОГИЯ
МАШИНОСТРОЕНИЕ
МЕДИЦИНА
МЕТАЛЛУРГИЯ
МЕХАНИКА
ПЕДАГОГИКА
ПОЛИТИКА
ПРИБОРОСТРОЕНИЕ
ПРОДОВОЛЬСТВИЕ
ПСИХОЛОГИЯ
РАДИОТЕХНИКА
СЕЛЬСКОЕ ХОЗЯЙСТВО
СОЦИОЛОГИЯ
СТРОИТЕЛЬСТВО
ТЕХНИЧЕСКИЕ НАУКИ
ТРАНСПОРТ
ФАРМАЦЕВТИКА
ФИЗИКА
ФИЗИОЛОГИЯ
ФИЛОЛОГИЯ
ФИЛОСОФИЯ
ХИМИЯ
ЭКОНОМИКА
ЭЛЕКТРОТЕХНИКА
ЭНЕРГЕТИКА
ЮРИСПРУДЕНЦИЯ
ЯЗЫКОЗНАНИЕ
РАЗНОЕ
КОНТАКТЫ


Pages:     | 1 || 3 | 4 |   ...   | 9 |

«Федеральное государственное образовательное учреждение высшего профессионального образования «МОСКОВСКИЙ ГОСУДАРСТВЕННЫЙ УНИВЕРСИТЕТ ИМЕНИ М.В. ...»

-- [ Страница 2 ] --

Для простых моделей низкомолекулярных жидких кристаллов (растворы абсолютно жестких стержней) капиллярная нематизация в тонких слоях была исследована с помощью компьютерного эксперимента в работах [212, 213]. Само это явление капиллярной нематизации является аналогом хорошо известного явления капиллярной конденсации для перехода пар – жидкость в тонких слоях или капиллярах [214-216]. Для асимптотического случая полубесконечной системы, т.е. в случае очень большой ширины плоского слоя, результаты для индуцированного поверхностью ЖК упорядочения в растворах жесткоцепных макромолекул должны согласоваться с теорией «индуцированного поверхностью упорядочения» [217] для фазовых переходов 1-го рода.

Для простых решеточных моделей антиферромагнетиков в присутствии поверхностей возникновение приповерхностного упорядочения уже было исследовано с помощью моделирования методом МК [218], но для сложных жидкостей, в частности, для полуразбавленных и концентрированных растворов жесткоцепных полимеров, такие исследования ранее не проводились. С помощью компьютерного моделирования поверхности и межфазные границы в растворах гибкоцепных полимеров исследовались в работах [219-222], а в нематических ЖК в коллоидных системах жестких стержней – в работах [223-230]. Растворы жестких стержней в капсулах изучались в работе [231], а фазовые переходы в двумерных жидких кристаллах – в работе [232].

1.5. Некоторые другие актуальные направления исследования жесткоцепных полимеров 1.5.1. Упругое поведение отдельной макромолекулы В связи с появлением нового научного направления – наномеханики полимеров – актуальными становятся начатые много лет назад теоретические исследования и компьютерное моделирование упругого отклика отдельной макромолекулы на растяжение внешним полем. К этой казалось бы хорошо понятой много лет назад задаче приходится возвращаться по двум причинам: с одной стороны, нет единства в понимании применимости того или иного статистического ансамбля при анализе экспериментальных данных, а с другой стороны, в малых системах (отдельная макромолекула) есть сильная зависимость упругих свойств от жесткости цепи, которая, к тому же, по-разному проявляется при различной постановке эксперимента (в разных статистических ансамблях). Все релевантные ссылки можно найти в двух недавно вышедших работах [394, 395]. На упругое поведение отдельной цепи большое влияние может оказывать присутствие адсорбирующей поверхности, а также пространственная структура самой жесткоцепной макромолекулы, которая может принимать конформации тора или цилиндра, что существенно влияет на процесс десорбции при приложении силы [379].

1.5.2. Сетки из жесткоцепных макромолекул Интерес к проблеме неаффинных (нелинейных) деформаций и повышения жесткости при растяжении полимерных волокон не ослабевает в течение многих лет [396 400], что связано с важностью промышленных применений полимерных волокон, а также с изучением морфологических [408], механических [405-407] и реологических [403, 404] свойств цитоскелета клеток, который может состоять из микротрубочек или из микрофиламентов, представляющих собой, в свою очередь, жесткие супрамолекулярные структуры на основе белков тубулина и актина. Теоретически и с помощью компьютерного моделирования изучались огрубленные модели сеток с жесткоцепными субцепями между узлами [396-400], в которых наблюдается ориентационное упорядочение [396] и разные деформационные моды [398]. Были построены различные модели для описания повышения жесткости волокон при растяжении [397, 399].

Результаты моделирования подтверждают предположения, что при малых растяжениях сшитой сетки из жесткоцепных субцепей упругий отклик определяется изгибом отдельных сегментов субцепей, а при больших растяжениях основной вклад в упругие свойства дает растяжение сегментов субцепей. Локальные перестройки, которые обусловливают переход между этим двумя режимами, вызывают неаффинную деформацию сеток [401, 402]. Для жесткоцепных полимерных сеток характерно сложное многомодовое релаксационное поведение [403, 404]. Структура жесткоцепных сеток сильно зависит от концентрации сшивающего агента, и комбинирование жесткоцепных филаментов, в частности, F-актина с различными сшивающими агентами является перспективными направлением разработки новых сетчатых биополимерных материалов с иерархической структурой [408].

1.5.3. Многомасштабное моделирование сопряженных полимеров Сопряженные полимеры представляют собой важный класс полимерных материалов для органической электроники [409]. Как правило, такие полимеры имеют довольно большую внутрицепную жесткость за счет наличия двойных связей и бензольных колец в основной цепи. Для моделирования свойств таких материалов важно, с одной стороны, правильно учесть специфические взаимодействия, а именно, взаимодействия -электронов, которые отвечают за электронные свойства, а с другой стороны, важно получить правильную супрамолекулярную структуру, которая образуется в таких системах, в том числе, за счет -стэкинга, и по которой может перемещаться заряд. Для правильного учета электронных свойств нужны квантово-химические расчеты, а для получения супрамолекулярной структуры нужны крупнозернистые модели, позволяющие исследовать большие пространственные и временные масштабы. Поэтому тут никак не обойтись без многомасштабного моделирования, что представляет собой нетривиальную задачу, которая весьма актуальна в настоящее время [410, 411].

1.6. Выводы по 1-ой главе В области компьютерного моделирования систем жесткоцепных полимеров на момент начала работ в рамках настоящей диссертации было много нерешенных проблем.

Вошедшие в диссертацию результаты были опубликованы в период с 1988 по 2013 год. В течении этого времени работы автора данной диссертации соответствовали мировому научному уровню, всегда являлись новыми и актуальными, а иногда даже пионерскими. В 1-ой главе представлен обзор литературы, включающий в себя как достаточно старые работы, так и более новые, которые были выполнены в других научных группах примерно в одно время с работами автора данной диссертации. В эту главу не вошел обзор непосредственно самих методов компьютерного моделирования, который приводится в главе 2, причем там к нему сразу добавлено описание методов и алгоритмов, существенно развитых или впервые разработанных в рамках настоящей диссертации. Некоторые из наиболее актуальных направлений исследований систем жесткоцепных полимеров, где все еще остается много нерешенных проблем, обсуждаются в разделе 1.5 этой главы, причем у автора настоящей диссертации есть опубликованные работы и по этим направлениям тоже, но эти результаты не были включены в диссертацию, чтобы не расширять без надобности объем и тему диссертации.

ГЛАВА 2. ОПИСАНИЕ МОДЕЛЕЙ И АЛГОРИТМОВ, ИСПОЛЬЗОВАННЫХ И РАЗРАБОТАННЫХ В ДИССЕРТАЦИИ Суть же метода, мной примененного тут, Объяснить я подробней готов, Если есть у вас пара свободных минут И хотя бы крупица мозгов.

Льюис Кэрролл «Охота на Снарка»

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

Огрубленные модели и атомистические модели имеют разные сферы применения.

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

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

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

Главное, конечно, при упрощениях вместе с «водой» не выплеснуть и «ребенка», то есть не потерять самого существенного свойства исследуемой системы. В настоящее время из решеточных моделей наиболее часто используется модель цепи с флуктуирующей длиной связи [15,16,17], которая является квазиконтинуальной и позволяет достаточно хорошо описывать многие процессы в полимерных системах.

2.1. Решеточная модель цепи с флуктуирующей длиной связи «На решетке думается легче»

(высказывание Т.М.Бирштейн) Решеточная модель цепи с флуктуирующей длиной связей [15, 16, 17] относится к классу огрубленных моделей. В таких моделях одно эффективное мономерное звено заменяет некоторую группу атомов реальной полимерной цепи. Модель цепи с флуктуирующей длиной связей для случая трехмерного пространства строится на простой кубической решетке (элементарный шаг решетки выбирается за единицу длины). Базовая ячейка моделирования имеет размер Lx Ly Lz единиц длины. Мономерное звено представляет собой элементарный куб решетки, т.е. занимает 8 узлов решетки (рис.3).

Мономерные звенья связаны в полимерные цепи с помощью специально заданного набора векторов связи, который получается из базового набора (2,0,0), (2,1,0), (2,1,1), (2,2,1), (3,0,0), (3,1,0) с помощью всех возможных перестановок и изменений знаков координат.

Полный набор возможных векторов связей состоит из 108 векторов. Длина векторов связей не является фиксированной, а может принимать значения 2, 5, 6, 3, единиц длины, что и дало название алгоритму. Элементарный пробный шаг изменения конформации цепи состоит в локальном смещении случайно выбранного мономерного звена на единичный шаг решетки в случайно выбранном направлении (одном из шести возможных ± x, ± y, ± z ). Такой шаг принимается, если выполнены три условия: (1) не нарушается условие исключенного объема, т.е. четыре узла решетки, соседние с этим мономером в направлении его смещения, являются не занятыми другими мономерными звеньями;

(2) цепь при таком смещении не рвется, т.е. новые вектора связи от предыдущего по цепи и к последующему по цепи мономерным звеньям также принадлежат к разрешенному набору векторов связи;

(3) выполняется критерий Метрополиса [19] (в случае, когда в системе имеется потенциал взаимодействия).

Рис. 3. Модель цепи с флуктуирующей длиной связей.

Основная идея алгоритма состоит во взаимно согласованном выборе набора векторов связи и направлений смещения мономерных звеньев таким образом, чтобы не надо было уже дополнительно к этим трем условиям проверять, не пересекают ли друг друга вектора связи звеньев вдоль по цепи при смещении мономерных звеньев. И именно по этой причине в базовом наборе нет вектора (2,2,0), хотя его длина и меньше максимально допустимой, а также не разрешены элементарные смещения по диагоналям элементарного куба и по диагоналям его граней. Всего в модели допустимы 108 векторов связи, между которыми возможны 67 различных углов, т.е. такая решеточная модель является по сути квазиконтинуальной.

Рассмотрим теперь более конкретно случай раствора жесткоцепных макромолекул длиной N мономерных звеньев. Для моделирования движения цепей с целью приведения системы к равновесию могут использоваться как шаги локального смещения мономерных звеньев [15, 16, 17], так и рептационные движения по алгоритму «скользящей змеи» [18, 20]. Если моделирование производится в большом каноническом ансамбле, используются элементарные шаги встраивания/удаления цепей (см. ниже). Все шаги МК принимаются или отвергаются в соответствии с критерием Метрополиса [19,20].

Качество растворителя, который не учитывается в данной модели явным образом, описывается введением эффективного потенциала притяжения между мономерными звеньями (объемные взаимодействия по типу ван-дер-ваальсовых):

U thermal ( r ) = = J для r = 2, 5, (1) kBT 0 для других r где T – температура, kB – постоянная Больцмана, r – расстояние на решетке между рассматриваемыми мономерными звеньями, – энергия парного взаимодействия мономерных звеньев, =1/kBT, J=.

Жесткость полимерной цепи задается обычно с помощью потенциала, зависящего от угла между соседними вдоль по цепи связями (рис.3). Этот потенциал может быть выбран, например, в виде U stiff ( ) cos = b cos (2) = k BT k BT где – энергетический параметр жесткости, b=/kBT. Параметр жесткости b связан с обычно используемым в литературе [21] параметром жесткости p = lK/d (отношением сегмента Куна к диаметру) линейной зависимостью [22]. В зависимости от модели можно фиксировать либо b, либо. При постоянном параметре b не происходит эффективного увеличения жесткости цепи при понижении температуры. При постоянном параметре энергия жесткости явно не зависит от температуры.

Если рассматривается полимерная система вблизи плоской адсорбирующей поверхности, расположенной в плоскости z=0, потенциал адсорбции может быть выбран, например, в виде w, z = 0, z = k BT U ads ( z ) (3) = w 4 k BT k BT, z z где w – энергетический параметр притяжения к поверхности, а z – аппликата мономерного звена.

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

2.2. Методы Монте-Карло для моделирования фазового равновесия и вычисления свободной энергии Алгоритм Метрополиса [19] позволяет вычислять средние значения наблюдаемых физических величин, но не позволяет рассчитывать статистический интеграл системы и соответствующий нужному ансамблю термодинамический потенциал и, следовательно, он плохо подходит для исследования фазового равновесия и фазовых переходов. Вместе с тем, метод МК предоставляет практически неограниченные возможности и для этих целей, если использовать другие МК-алгоритмы [238-264,283]. Построение фазовой диаграммы крайне важно для понимания поведения любой молекулярной системы.

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

раздел 2.2.3 и более подробно – главы 3 и 4 в книге [20,303]), мультиканонического моделирования, методы расчета функции плотности состояний (см. главы 3 и 4 в книге [20,303]).

Общее описание идеологии, основ метода расширенных ансамблей содержится в [303], два из четырех авторов которой являются основоположниками этого метода. Там же приведен пример ансамбля, расширенного по энергии (энтропическое моделирование с использованием алгоритма Ванга-Ландау). Эти методы обсуждаются также и в [20], которая содержит достаточно полный обзор литературы по этим методам.

2.2.1. Моделирование различных статистических ансамблей Классическими статистическими ансамблями являются микроканонический, канонический, изотермически-изобарический и большой канонический. Однако, на практике используется много других ансамблей, специфических для конкретных систем.

Если при этом накопленная выборка была бы бесконечно большой, то усреднение значений наблюдаемых физических величин по любому ансамблю привело бы к строго одинаковым результатам [192]. Но реально используемые выборки всегда конечны, поэтому оказывается, что для обеспечения наиболее быстрой сходимости конкретной системы (или какого-либо свойства системы) к равновесию тот или иной ансамбль более предпочтителен. Обзор применимости различных ансамблей обсуждается во многих книгах и обзорах ([35,192,238];

стр.95-108 в [262]).

Например, отметим, что процедура Метрополиса может быть использована для расчета разности свободных энергий двух состояний в большом каноническом ансамбле, что было впервые показано при моделировании простых жидкостей в работе [284], в которой, кстати, был впервые разработан и собственно алгоритм МК для большого канонического ансамбля. Перечислим некоторые другие методы, используемые при расчете термодинамических потенциалов или их разностей: термодинамическое интегрирование, зонтичная выборка, интегрирование по Гиббсу-Дюгему вдоль линий фазовых переходов (см., например, [35, 239, 258, 262]). Следует упомянуть также метод ансамбля Гиббса [285,286], который оказывается весьма эффективным при моделировании сосуществования фаз в плотных многокомпонентных жидкостях.

2.2.2. Стратегии выборки макросостояний в фазовом пространстве При моделировании фазовых переходов с помощью методов МК ключевым моментом является построение траектории (пути) в фазовом пространстве от одного микросостояния к другому таким образом, чтобы можно было установить правильные вероятностные соотношения между макросостояниями системы, отвечающими разным фазам. Приведенный в этом подразделе обзор методов построения таких траекторий основан на работах [239,262,264,287]. Модельная молекулярная система характеризуется некоторым набором «контрольных параметров», включающих в себя как макроскопические термодинамические переменные (например, температуру), так и микроскопические параметры модели (например, параметры потенциала взаимодействия между частицами). Следуя работам [239,264], будем обозначать в данном разделе этот набор контрольных параметров буквой c или, если в этот набор включена еще метка фазового состояния системы, большой буквой C. От набора C контрольных параметров зависит макросостояние системы. Микросостояние системы в конфигурационном пространстве определяется набором (обобщенных) координат {q}. Равновесная функция плотности вероятности определяется подходящей (для каждого конкретного ансамбля – своей) безразмерной (т.е. отнесенной к величине k BT ) конфигурационной энергией E ({q}, C ) согласно следующему соотношению:

1 P0 ({q} | C ) = exp E ({q}, C ), Z (C ) Z ( C ) dqi exp E ({q}, C ). (4) i Задача построения представительной выборки всех макросостояний может быть решена разными способами, три из которых перечислены в следующих подразделах.

Последовательная выборка. Самый простой способ сбора информации о наборе макросостояний {Cj} (j=1,…,) состоит в использовании больцмановской выборки для каждого из макросостояний по очереди:

PS ({q}) = P0 ({q} | C j ), j = 1,...,. (5) Затем, чтобы построить траекторию (путь) между этими макросостояниями, необходимо объединить информацию, полученную в независимых расчетах, с помощью, например, методов термодинамического интегрирования [288].

Параллельная выборка. Вместо накопления последовательной выборки по макросостояниям можно проводить процедуру параллельно. При таком подходе проводится моделирование набора реплик (копий) физической системы, причем j-тая реплика имеет набор контрольных параметров Cj. Выборка в составном конфигурационном пространстве определяется формулой ({q} ) ( ) (1) () ( j) PS...{q} = P0 {q} | C j j =. (6) Состояния в таком составном ансамбле могут изменяться путем обмена наборами координат частиц системы между соседними макросостояниями, например, Cj и Cj+1, так что ( j) ( j +1) ( j +1) ( j) {q} = {q}, {q} = {q}. (7) Близкие алгоритмы такого типа были независимо предложены для разных задач и под разными названиями [287], включая алгоритм «парных цепочек Метрополиса»

(Metropolis coupled chain) [289], «обменные» алгоритмы МК (exchange Monte Carlo) [290], метод «параллельного регулирования» (parallel tempering) [291] (см. также [20]).

Расширенная выборка (расширенные ансамбли). Существует возможность исследовать статистические свойства полного набора макросостояний C1,…, C в рамках одного расширенного моделирования. Общая формула для выборки, которая приводит к достижению этой цели, записывается в виде суперпозиции больцмановских выборок:

PS ({q}) = W0 P0 ({q} | C j ), (8) j = где W0 – нормировочная константа. С начала 90-х годов прошлого века был разработан широкий класс алгоритмов МК для «расширенных ансамблей» (extended ensemble Monte Carlo;

другой вариант объединяющего названия для методов этого класса – «обобщенные ансамбли», generalized ensembles [287]. Методы, использующие выборку в виде суперпозиции больмановских выборок для набора макросостояний, включают в себя собственно «расширенные ансамбли» (expanded ensembles, исторически первый метод этого класса) [292], «условное регулирование» (simulated tempering) [293], «масштабирование температуры» (temperature scaling) [294], «адаптивную зонтичную выборку» (adaptive umbrella sampling) [295], «мультиканоническое» моделирование (multicanonical ensemble) [296], и ряд других. Оба термина – extended и expanded – переводятся на русский язык одинаково («расширенный»), но в современной литературе первый из них принято использовать для обозначения всего класса подобных алгоритмов, а второй закрепился за исторически первым вариантом [292].

2.2.3. Метод расширенных ансамблей Компьютерная имитация макромолекулярной (полимерной) системы методом МК в рамках канонического ансамбля с использованием, например, процедуры Метрополиса [19], имеет ряд недостатков. Один из них – невозможность прямого вычисления свободной энергии системы, хотя некоторые способы ее расчета и были предложены [299,300,301] (в том числе с использованием упомянутой выше “зонтичной” выборки [282]). Идеи этих способов сходны: предлагается вычислить разность свободных энергий для изучаемой и некоторой эталонной системы (например, идеального газа или свободносочлененной полимерной цепи), абсолютное значение свободной энергии которой вычисляется строго. Однако, прямая реализация подобных способов возможна лишь тогда, когда изучаемая система не намного отличается от эталонной. В противном же случае все вычисления надлежит проводить в несколько этапов: конструировать ряд промежуточных систем и рассчитывать разности свободных энергий для всей цепочки пар соседних систем, постепенно приближаясь к изучаемой.

Идея метода расширенных ансамблей состоит в рассмотрении не одной системы, а сразу целой совокупности подсистем (“расширенного ансамбля”). Каждая из этих подсистем является промежуточной между изучаемой системой и эталонной, и все они различаются между собой некоторым “параметром расширения” (контрольным параметром). Роль этого параметра расширения могут играть как обычные термодинамические величины, например, температура или давление, так и вновь введенный для промежуточных систем особый «параметр порядка». При этом в общем случае промежуточные подсистемы включают состояния, обычно запрещенные в основном ансамбле, над которым проводится процедура расширения. Например, основной ансамбль может быть дополнен (расширен) состояниями с большими или меньшими размерами атомов, так что при переключении между подсистемами происходит постепенное увеличение/уменьшение размеров атомов. Другая возможность: подсистемы могут включать состояния с частичным перекрыванием «абсолютно жестких» сфер или состояния с другим значением потенциала межчастичного взаимодействия (более «мягкий» потенциал, вплоть до полного его исчезновения, то есть в пределе вплоть до идеальной системы без взаимодействия). Эти состояния характеризуются некоторым значением отдельного нового «параметра порядка» (например, «силой» потенциала взаимодействия, степенью перекрытия атомов). Заданием последовательных значений “параметра расширения” задается искомый переход от эталонной системы к изучаемой.

При этом нет необходимости проводить отдельные расчеты для каждой из систем с промежуточными параметрами. Блуждание в расширенном ансамбле проводится с помощью процедуры Метрополиса с двумя видами МК-шагов: (а) смещениями частиц в заданной подсистеме и (б) изменениями параметра расширения на соседнее значение при фиксированной конфигурации (или переключение между подсистемами может состоять непосредственно в переключении между конформациями, как это делается в методе «параллельного регулирования», который подробнее описан ниже). В процессе МК эксперимента учитывается количество посещений каждой подсистемы и вычисляется вероятность этого состояния, что в итоге дает возможность вычислить разность свободных энергий для любой пары подсистем. И, поскольку свободная энергия исходной (эталонной) подсистемы известна, можно вычислить искомое абсолютное значение свободной энергии изучаемой системы. Типом параметра “расширения” задается тип расширенного ансамбля.

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

M wy = Z ye y =1 (9) Zy где – статистическая сумма подсистемы со значением «параметра расширения», w равным y, а y – вес этой подсистемы. Веса состояний (подсистем) выбираются так, чтобы все эти состояния (подсистемы) посещались с равной вероятностью. Как уже было сказано выше, каждое состояние (подсистема) уравновешивается с помощью обычных элементарных шагов смещения с вероятностью принятия пробного смещения, которая отвечает соответствующему «головному» ансамблю, но дополнительно вводится новый элементарный шаг, состоящий в переходе между этими состояниями (подсистемами) и принимаемый с вероятностью:

exp [ wx U x ] pacc ( y x) = min 1, exp wy U y (10) энергия состояния в подсистеме i, = 1/ kBT, и вероятность принятия шага Ui – Здесь переключения между двумя состояниями в разных подсистемах записана для случая, когда головной ансамбль – канонический, а расширение проводится не по температуре.

При этом функция плотности вероятности подсистемы и разность свободной энергии записывается в виде py = exp( wy ) exp( U y ) (11) p y Z y Fy Fx = log = log + wy wx Z x px (12) 2.2.4. Метод моделирования в расширенном ансамбле в четырехмерном пространстве Для уменьшения времени уравновешивания плотных глобулярных конформаций полимерной цепи в работе [23] был предложен алгоритм для расширенного ансамбля в четырехмерном пространстве, который был в дальнейшем разработан в настоящей диссертации для моделирования внутримолекулярных структур в одиночной жесткоцепной макромолекуле (см. работу [24]). Общая идеология расширенных ансамблей обсуждается в обзоре [20]. Основная идея данного конкретного алгоритма состоит в расширении пространства конформаций на четырехмерное пространство и разработке специальной процедуры в некотором смысле «равномерного» блуждания цепи в конформационном пространстве между чисто трехмерными и четырехмерными конформациями с разной «степенью выхода» из трехмерного в четырехмерное пространство. Четырехмерные конформации могут рассматриваться как трехмерные со снятым условием исключенного объема, т.е. с допущением перекрывания мономерных звеньев. С этой целью вводится четвертая координата для мономерных звеньев, то есть рассматриваются две трехмерные подрешетки, в одной из которой четвертая координата для мономерных звеньев равна 0, в другой 1. Это, конечно, не есть четырехмерное пространство в истинном значении этого слова, а фактически две «параллельные»

решетки. Мономерные звенья с разными значениями четвертой координаты не взаимодействуют друг с другом (см. рис.4).

(x,y,z,1) (x,y,z,0) Рис.4. Схематическое изображение двух трехмерных подрешеток в квази четырехмерном пространстве.

Рассмотрим далее этот алгоритм на примере одиночной жесткоцепной макромолекулы. Гамильтониан такой системы представлен уравнением N H = H 0 + hx4 ( i ) (13) i = где H0 – гамильтониан обычной трехмерной системы, x4(i) - четвертая координата i ой частицы, принимающая значения 0 или 1. «Внешнее поле» h контролирует распределение мономерных звеньев между двумя подрешетками. Для h=0 мономерные звенья распределены равномерно между двумя подрешетками. Для все h мономерные звенья имеют четвертую координату равную x4(i)=0, то есть находятся в обычном трехмерном пространстве. При уменьшении значений «внешнего поля»

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

« », ( ) h.

h.

, z 0, z w. x4(i)=0, k BT.x4(i)=0, ).

(3). «, » записывается в виде (см.

4 w Гамильтониан одиночной полужесткой полимерной цепи k BT « », ( ),z2.

( ).

уравнения (1)-(2) выше):

z,.

.

, z–,. N 2, U stiff (i ) H.

= Jn 0.

N,2U H. +,i (14) stiff k BT i =1 Jn BT k (5) k BT k BT,, i где n – полное число контактов между мономерными звеньями H0полимерной цепи. В i N 2U, n– stiff Jn N 2U H0.

общем случае параметр взаимодействия мономерных звеньев J kявляется функцией BT k BT stiff i i Jn [23] (5) J h, J(h). n – параметра h, J(h). Тогда гамильтонианkрасширенного ансамбля (13) выражается формулой k BT BT, i (4).

n– J h, J(h).

H [24]. N. (4) (6) J h n b J(h).

h, cos hn4 d J i H k(4) BT N, J h n b cos (15) 4 d i hn i k BT n4d – x4(i)=1 i H – число мономерных звеньев с четвертой координатой x4(i)=1 (это, N ( где n4d J, h n b cos i,hn4 d, n4d – (6) k BT (,степень выхода конформаций из,, ), i собственно, и есть параметр, характеризующий n4d – x4(i)= ».

(, h h k BT.,, h h k BT. Статсумма цепи в таком.

трехмерного пространства в четырехмерное),,.. ), N 1 N 1 (7).

hансамбле.выражается формулой exp J h n b cos i Z hn4 d Z exp J h n b cos h k BT, h g ( h) c i h g ( h) c i1 i N hn4 d h, – » g(h) (7) « » h, Z g(h) – exp J h n « b, 1.,, cos (16) i h g ( h), c i « ».

. « g(h) – ». « » h,.

где g(h) – весовая функция «внешнего поля» h, а суммы берутся соответственно по :) :), « всем конформациям цепи и по всем значениям «внешнего поля».

». ) «.

) « » [18, 20], ) В расширенном ансамбле элементарным шагом МК являются следующие шаги: а) (,.

. 4 :) N H0 hx4 i ) h.

« » [18, 20], ) (4) (.. 4- ), ) локальное перемещение случайно выбранного звена, б) билокальные перемещения по, i h. h, x (i) алгоритму (.. 4- ), ) «скользящей 4 змеи», [18, 20], в) переключение случайно. выбранного i-, 0 1. h. h, x =1).

h мономерного звена между подрешетками (т.е. изменение значения 4-ой координаты этого,.

( ), ( ) () n4d (. h=0, звена), г) случайное изменение значения внешнего=1). (h. Значение поля h x4 поля меняется n4d ( ), ( ) () x4=1).

только на значения, которые являются соседними к данному. Изначально дискретный ( ), ( ) () набор значений поля выбирается таким образом, чтобы происходило перекрывание гистограмм величины n4d (рис.5а). Элементарные шаги (а), (б) и (в) изменения конформации принимаются в соответствии с критерием Метрополиса [19], а шаги изменения значения внешнего поля (г) принимаются в соответствии с алгоритмом Ванга Ландау (см. следующий подраздел 2.3, уравнение (18в)).

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

Весовые коэффициенты g(h) в формуле (16) рассчитываются с помощью алгоритма Ванга-Ландау (см. следующий подраздел 2.3). Алгоритм устроен таким образом, что получающиеся в результате значения g(h) позволяют системе посещать с равной частотой (рис.5b, 5c) все значения из дискретного набора для поля h, т. е. функция g(h) «регулирует» количество мономерных звеньев, находящихся в «четвертом измерении». Из всех получающихся конформаций выбираются только чисто трехмерные, и по ним проводится усреднение наблюдаемых величин.

a) 0, 0, 0, 0, 0, Hist Hist 0, 0, 0, 0, 0, 192 0 64 -16 -14 -12 -10 -8 -6 -4 -2 N µ b) c) Рис.5. (a) Зависимость гистограммы числа звеньев, имеющих четвертую координату, равную 1, от величины внешнего поля. (b) Гистограмма величины внешнего поля. (c) Гистограмма числа звеньев, имеющих четвертую координату, равную 1.

С помощью этого нового метода, разработанного в рамках настоящей диссертации, были исследованы плотные конформации жидкой и твердой (кристаллической) глобулы и определена точка перехода между этими структурами (см. разделы 3.1.7 – 3.1.9).

2.2.5. Метод «параллельного регулирования»

Метод «параллельного регулирования» (parallel tempering technique) идейно относится к категории обобщенных (расширенных) ансамблей. В системах с большим числом минимумов на профиле свободной энергии (как говорят, в системах с сильно изрезанным ландшафтом свободной энергии) возникает проблема: каким образом преодолеть барьеры между состояниями и помочь системе найти глобальный минимум? В методе «параллельного регулирования» рассматривается набор невзаимодействующих реплик системы. Одна реплика – это отдельное состояние (конфигурация) при некотором значении какого-либо контрольного параметра, например, температуры, или параметра потенциала межчастичного взаимодействия и т.п. Статистическая сумма полной системы записывается в виде r = Z ( N,V, Ti ) i = где Z ( N,V, Ti ) – статистическая сумма i-той реплики. В дополнение к обычным шагам МК разрешается элементарный шаг, состоящий во взаимном обмене конформациями между двумя репликами i и j с вероятностью принятия шага { }, pacc (i j ) = min 1, exp ( ij U ij ) (17) где учитывается не только разность энергий исходного и пробного состояний, но и разность значений контрольного параметра у двух этих реплик.

Метод параллельного регулирования в данной диссертации не использовался.

2.2.6. Метод пересчета гистограмм Нельзя не упомянуть о важном техническом приеме – методе пересчета гистограмм [297]. Исторически это была одна из первых методик расчета функции плотности состояний (а через нее и статистического интеграла системы) на основе использования того факта, что гистограммы энергии и параметра порядка в окрестности критической точки являются очень широкими и включают в себя достаточно много информации о большом числе состояний системы. Эти гистограммы могут быть использованы для пересчета данных на другие значения контрольных макропараметров системы, например, температуры;

о взаимосвязи метода пересчета гистограмм и мультиканонического моделирования см. [262], стр.111-157. Метод пересчета гистограмм использовался в данной диссертации в разделах 3.1 и 5.1.

2.2.7. Метод конечномерного масштабирования Еще одна из основополагающих методик, без которой не обойтись при исследовании фазовых переходов с помощью машинного моделирования, – это метод конечномерного масштабирования [14,238,256,298]. Этот метод позволяет экстраполировать результаты компьютерного моделирования, полученные для систем достаточно малых размеров, на системы бесконечно больших размеров, для которых, собственно, только и можно говорить о фазовых переходах. Данный метод использовался в диссертации в разделах 3.1 и 4.2.

2.3. Алгоритм Ванга-Ландау К методам расширенных ансамблей близок метод мультиканонического (или энтропического) моделирования [296,304]. Этот метод представляет собой фактически моделирование в микроканоническом ансамбле, “расширенном” по энергии. Чуть более 10 лет назад был предложен весьма эффективный его вариант (вариант с автоматической настройкой параметров) – алгоритм Ванга-Ландау [25,26], который используется в диссертации как в оригинальном виде (в разделах 3.1, 3.2, 5.2), так и в модификации для расширенного ансамбля (в подразделе 3.1.7).

Алгоритм Ванга-Ландау (ВЛ) позволяет рассчитать плотность состояний g(A) для некоторой наблюдаемой физической величины A, т.е. число (микро)состояний системы с данным значением величины A. Изначально этот алгоритм [25,26] был предложен для расчета плотности состояний g(E) для полной энергии E, но в общем случае в качестве наблюдаемой величины A может фигурировать не только полная энергия, но и, например, один из вкладов в полную энергию, или характерный линейный размер системы, или величина характерного для данной системы параметра порядка и т.д.

Более того, в качестве величины A может выбираться даже не значение некоторой внутренней физической величины, которая может быть рассчитана для данной системы, но и значение какого-либо внешнего входного параметра моделируемой системы, если само моделирование проводится не при фиксированном значении, а сразу в некотором допустимом интервале значений этого параметра. Обычно здесь имеется в виду «параметр расширения» в некотором расширенном ансамбле. В этом случае терминология другая:

алгоритм ВЛ позволяет рассчитать весовую функцию g(A) для внешнего входного параметра A, которая есть функция, обратно пропорциональная функции плотности вероятности выбрать данное значение входного параметра A для моделирования системы.

В данном разделе обсуждаются обе такие реализации алгоритма ВЛ, и соответственно под величиной A мы будем подразумевать либо полную энергию системы E, либо значение параметра расширения («внешнее поле») h расширенного ансамбля в четырехмерном пространстве. Для удобства будем здесь называть эти две реализации соответственно алгоритмами ВЛ по внутреннему и по внешнему параметру.

Суть алгоритма ВЛ по внутреннему параметру состоит в использовании небольцмановской выборки с целью получить равномерную гистограмму для энергии (или другой внутренней наблюдаемой величины A). В вероятность принятия пробного шага входит очередная (актуальная на данный момент) итерация функции плотности состояний. Функция плотности состояний изменяется в ходе специальной итерационной процедуры в ходе моделирования до тех пор, пока ее изменения не станут меньше заранее заданной точности, т.е. пока итерационная процедура не сойдется к правильной функции плотности состояний. Зная функцию плотности состояний и накопив гистограмму значений наблюдаемых физических величин, можно рассчитать статсумму и пересчитать равномерные гистограммы на функции распределения наблюдаемых величин в нужном статистическом ансамбле.

Суть алгоритма ВЛ по внешнему параметру немного другая. При моделировании в расширенном ансамбле по параметру расширения надо обычно делать равномерное блуждание, т.е. обеспечить примерно одинаковую частоту использования всех значений внешнего параметра расширения в ходе моделирования, чтобы дойти от одной предельной системы (идеальной в некотором смысле) до нужной нам реальной системы.

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

Алгоритм Ванга-Ландау реализуется следующим образом:

1) в самом начале необходимо задать стартовую конформацию системы, значения гистограммы величины A приравнять нулю, а значения функции g приравнять единице, H(A)=0, g(A)=1 для всех значений A, задать начальное f0 и конечное ffinal значения параметра f, который является множителем для функции g, и максимальное значение max параметра, характеризующего равномерность гистограммы величины A;

таким параметром может служить, например, относительная разность максимального и минимального значений гистограммы, = ( H max - H min ) H max ;

значения этих параметров, как правило, выбираются такими: f0 = 2 или f0 = e (основание натурального логарифма), ffinal = 1.0000001 (некоторые пояснения относительно выбора даются ниже), 0.1 max 0.2;

2) производится попытка изменения конформации системы (попытка перехода из состояния 1 в новое состояние 2);

для этого используются любые элементарные пробные шаги смещения, которые допускаются данным конкретным алгоритмом МК, с помощью которого производится моделирование системы (локальные шаги смещения, алгоритм «скользящей змеи», алгоритм с конформационным смещением выборки, алгоритмы в 2) ( расширенных ансамблях и т.д.);

но вот только вероятность перехода в новую 1 2);

конформацию определяется теперь не критерием Метрополиса [19], а другими, f, формулами;

для ВЛ по внутреннему параметру вероятность перехода в новую (, конформацию есть», «,..);

g ( A1 ), (18а) p ( A1 A2 ) = min 1, g ( A2 ) [19], ;

где A1 и A2 есть значения величины A соответственно в исходном (старом) и g A p A1 а для ВЛ min 1, A2, (8 ) пробном (новом) состоянии, по внешнему параметру вероятность принятия шага g A переключения значения внешнего параметра A из A1 в A2 (шаг (г) в алгоритме в A1 A2 A предыдущем подразделе 2.2) есть ( ) ( ), g ( A1 ) {H ( A22.2) )} A A1 A2 ( (), (18б) p ( A1 A2 ) = min 1, 1 g ( A2 ) {H ( A1 )} H A gA p A1 A2 min 1,, (8 ) H A g A где – плотность вероятности данного микросостояния с гамильтонианом H, зависящим от – значения параметра A (само микросостояние не изменяется при шаге H, A( переключения значения внешнего параметра A);

в частности, для конкретного случая расширенного ансамбля A);

четырехмерном пространстве для гибкой цепи (параметр в, жесткости b=0) эта формула принимает вид (см. уравнение (15)):

(. (6) (7)) g h1 exp J h2 n h2 n4 d p h1 h2 min 1,.. (8 ) (18в) g h2 exp J h1 n h1n4 d 3) 3) если g(A )=g(A )f, H(Aреализовалась и система перешла в новое состояние, то эта вероятность )=H(A )+1,, 1 1 1 g(A2)=g(A2)f, H(A2)=H(A2)+1, а если попытка оказалась,не принятой и 2система осталась в g(A2)=g(A )f, H(A2)=H(A2)+1, g старой конформации, то g(A1)=g(A1)f, H(A1)=H(A1)+1, то есть соответствующий элемент f, 1;

массива g умножается на параметр f, а соответствующий элемент массива гистограммы 4) H(A) : ( увеличивается на 1;

, max, 4) производится проверка гистограммы H(A) на равномерность: рассчитывается, 2 3;

величина параметра, и если, то есть гистограмма не удовлетворяет критерию max равномерности, повторяются пункты 2 и 3;

5) если же критерий равномерности гистограммы H(A) выполняется, то есть max, то производится изменение параметра f по формуле f = f (если fffinal), обнуляется гистограмма величины A, H(A)=0 для всех значений A, и алгоритм продолжает работать опять с пункта 2;

следует подчеркнуть, что значения функции g на этом шаге никак не меняются, то есть функция g продолжает накапливать изменения на протяжении всего итерационного процесса (пункты 2-5), до тех пор пока множитель f не станет равным ffinal;

после этого считается, что итерационная процедура сошлась и в результате ее работы получена с заранее заданной точностью функция g (функция плотности состояний величины A или весовая функция внешнего входного параметра A), то есть можно продолжать расчет согласно пунктам 2-4, только теперь не надо производить никаких изменений функции g в пункте 2.

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

(a) (b) Рис.6. (a) Логарифм функции плотности состояний для гибких цепей, привитых к плоской поверхности одним концом, для разных значений длины цепи и параметра притяжения к поверхности, указанных в легенде (см. раздел 3.2 ниже). (b) Логарифм функции плотности состояний для полужестких цепей (значение параметра жесткости b=4), привитых к плоской поверхности одним концом, для разных значений длины цепи и параметра притяжения к поверхности, указанных в легенде (см. раздел 3.2 ниже).

Рисунки взяты из работы [234].

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

Z (T ) = exp(E / kBT ) g ( E)dE, A A( E ) exp( E / kBT ) g ( E )dE = Z (T ) T (19) где A – среднее по всем конформациям с данной энергией Е. Для ВЛ по внешнему параметру расширения используется простое арифметическое среднее наблюдаемых физических величин, так как усреднение требуется только в одной из подсистем расширенного ансамбля, а именно в той, где находится нужная нам система (например, только для чисто трехмерных конформаций), а все переходы между состояниями в этой подсистеме осуществляются с помощью обычного критерия Метрополиса [19].

Осталось сделать еще несколько замечаний. Первое из них касается условия детального баланса (микроскопической обратимости) [19]. Это условие состоит в том, что произведение функции плотности вероятности найти какое-либо состояние на вероятность перехода (за один шаг алгоритма) в другое (доступное за один шаг) состояние должно равняться произведению функции плотности вероятности найти это другое состояние на вероятность обратного перехода. Это условие является достаточным для того, чтобы у марковского процесса (каковым является последовательность состояний системы, получаемая в ходе моделирования) существовало бы стационарное (инвариантное) предельное распределение вероятности состояний. Как видно из уравнения (18), условие детального баланса в алгоритме Ванга-Ландау имеет вид 1 p ( A2 A1 ), (20) p ( A1 A2 ) = g ( A1 ) g ( A2 ) где 1/g(A1) есть вероятность находится в состоянии A1, а p ( A1 A2 ) есть вероятность перехода из состояния A1 в состояние A2 (за один шаг алгоритма). Так вот, это условие выполняется в алгоритме Ванга-Ландау только в предельном случае, когда множитель f стремится к 1. По этой причине конечное значение ffinal для множителя f выбирается достаточно близким к 1 (например, с точностью до седьмого знака после запятой, то есть конечное значение множителя f равно 1 в пределах точности обычного действительного числа на компьютере). Таким образом, производить набор статистики для вычисляемых в ходе компьютерного эксперимента физических величин имеет смысл только после окончания итерационной процедуры для расчета функции плотности состояний, то есть когда моделирование проводится уже с неизменяемой функцией g (в том числе, в пункте 2). Для алгоритма ВЛ по внешнему параметру условие детального баланса должно фактически выполняться только при одном значении внешнего параметра расширения (при том, который соответствует реальной физической системе), а это должен обеспечить «основной» алгоритм МК для данного расширенного ансамбля (лежащий в его основе еще до процедуры расширения).

Второе замечание – чисто техническое. При расчете функции g необходимо производить все операции с действительными числами с двойной точностью, а также накапливать не саму весовую функцию g, а ее логарифм, ибо иначе слишком быстро можно достигнуть максимально допустимого (на данном компьютере и для данного компилятора) значения действительного числа. Кстати, именно поэтому начальное значение множителя f для весовой функции выбирается равным f0 = e (на первом шаге итерационной процедуры к логарифму g прибавляется 1 при каждом посещении данного значения величины A).

Третье замечание касается разбиения полного интервала исследуемых значений параметра на «окна», или подынтервалы. Если интервал доступных значений энергии слишком широк, приходится считать функцию плотности состояний g(E) отдельно в нескольких интервалах во всем диапазоне возможных значений полной энергии. После этого производится «сшивка» функции плотности состояний и «сшивка» гистограмм энергии и других наблюдаемых физических величин (радиуса инерции, параметра порядка и т.д.) для всех интервалов. Далее пересчет температурных зависимостей и гистограмм физических величин для канонического ансамбля производится стандартным способом.

Четвертое замечание: могут возникнуть проблемы со сходимостью алгоритма ВЛ, особенно, в сложных полимерных системах.


Пятое замечание: в литературе есть примеры использования двумерного алгоритма ВЛ, когда строится функция плотности состояний по двум параметрам [27], однако такая процедура требует намного больше компьютерного времени.

2.4. Методы расчета давления в решеточных моделях Монте-Карло для расчета уравнения состояния полимерного раствора В этом разделе изложены методы измерения давления в компьютерном моделировании с помощью решеточных алгоритмов Монте-Карло, часть из которых была впервые разработана в рамках данной диссертации. Было также детально исследовано влияния конечных размеров системы на величину рассчитанного давления (результаты подробно изложены в работах [193] и [233].

Проблема точного вычисления давления в компьютерном моделировании полимерных растворов и расплавов является очень важной при изучении фазовых переходов первого рода (например, нематического упорядочения в растворе жесткоцепных макромолекул). Температура перехода и соответствующие плотности фаз (изотропной и нематической) в точке перехода могут быть определены из условия равенства химических потенциалов и давлений в разных фазах. Так, например, в главе настоящей диссертации (в работах [188, 190]) изучался переход изотроп-нематик в растворах жесткоцепных макромолекул с потенциалом притяжения между мономерными звеньями в рамках решеточной модели и большого канонического ансамбля. Важный методический аспект этих работ связан с методом определения перехода и фазовой диаграммы. Из трех условий равновесия [191] между фазами 1 и 2: T1 = T2, µ1 = µ2, P = P (где T – температура, µ - химический потенциал, P – давление), первые два выполняются автоматически в большом каноническом ансамбле. Задача, таким образом, сводится к нахождению значения химического потенциала µ, µmin µ µmax, при котором значения P ( µ ) = P2 ( µ ). Таким образом, возникает давления в обеих фазах равны, т.е. необходимость рассчитать давление в полимерных растворах различной плотности, т.е.

рассчитать уравнение состояния растворов макромолекул.

Рассчитать давление при моделировании методами МД и МК с использованием континуальных моделей (и вычислить уравнение состояния, т.е. зависимость давления от плотности) можно несколькими способами: 1) при помощи использования теоремы вириала, 2) путем моделирования в изотермо-изобарическом ансамбле, 3) методом встраивания пробной частицы (цепи), 4) на основе межмолекулярных корреляционных функций (детальное описание этих методов можно найти в [35,192]). Например, согласно теореме вириала [192], давление в каноническом ансамбле может быть получено как, (21) P = k BT + f (r ) rij ij dV i j где =N/V – плотность, N – число частиц, V – объем, d – размерность пространства, f ( rij ) - сила, действующая между частицами i и j, находящимися на расстоянии rij, сумма берется по всем парам частиц, а усреднение проводится по полному времени моделирования после достижения системой равновесия. Этот метод неприменим в решеточных моделях в силу их дискретной природы.

Детальная разработка методов расчета давления в полимерных растворах для решеточных моделей МК (включая сравнение различных методик) была выполнена в рамках настоящей диссертации в работах [193, 205]. Эти исследования позволили изучить зависимость уравнения состояния раствора жесткоцепных макромолекул от жесткости цепей, а также изучить поверхности раздела изотропной и нематической фаз и эффекты смачивания плоской поверхности нематическим раствором полимера [205].

В случае решеточных моделей наиболее часто используются метод встраивания пробной цепи (МВПЦ, test chain insertion method, TCIM) [194, 195] и метод отталкивающей поверхности (МОП, repulsive wall method, RWM) [196-200]. Используя эти методы для расчета давления при разных плотностях системы, мы можем получить уравнение состояния раствора макромолекул. В конце 80-х годов Р.Дикманом с соавторами был опубликован ряд работ [194-204], посвященных вычислению давления в атермических полимерных растворах (где присутствуют взаимодействия, обусловленные только исключенным объемом). Были представлены теоретические расчеты и результаты компьютерного моделирования методом МК в каноническом ансамбле для решеточной и континуальной моделей в двух- и трехмерном пространстве с помощью этих двух методов.

2.4.1. Метод термодинамического интегрирования Метод вставки пробной цепи является методом термодинамического интегрирования, при котором требуется определить вероятность вставки пробной цепи в P раствор с объемной долей полимера. Осмотическое давление ( ) для системы k BT n-меров рассчитывается затем по формуле ( ) = [1 ln p(, n)] + ln p(, n)d, (22) n n где p(, n) есть вероятность вставить n-мер в раствор n-меров с объемной долей.

Следует отметить, что метод вставки пробной цепи трудоемок и применим лишь для разбавленных и полуразбавленных растворов и не слишком длинных цепей, т.к. при больших длинах цепей вероятность вставки n-мера становится ничтожно малой. Это ограничение можно обойти, если использовать алгоритм МК с конформационным смещением выборки (CBMC, [28-34]).

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

2.4.2. Метод отталкивающей поверхности В то же время метод отталкивающей поверхности (МОП) остается эффективным и при больших плотностях растворов, а также в расплавах. В ячейке моделирования размером L L H вводится отталкивающий потенциал = exp{1 / kBT}, 0 1, действующий на частицы, находящиеся в пристеночном слое z = H. Тогда осмотическое давление может быть получено по формуле d H ( ), (23) ( ) = где H ( ) - доля занятых ячеек в слое z = H при данном потенциале отталкивания от стенки. Таким образом, этот метод предполагает серию запусков программы для различных значений отталкивающего потенциала для последующего интегрирования и расчета давления раствора макромолекул с данным значением плотности.

Отметим, что при использовании любого из вышеописанных методов для того, чтобы получить уравнение состояния раствора макромолекул (т.е. зависимость давления от плотности) необходимо проводить серии запусков программы для ряда значений плотности системы. Это, в свою очередь, приводит к большим затратам по времени. Более того, в работе [193] было показано, что методу отталкивающей стенки присущи значительные эффекты конечного размера (в особенности при больших значениях плотности), которых можно избежать, лишь проводя моделирование в большом каноническом ансамбле (более подробно этот материал изложен в разделе 4.1.4 ниже).

2.4.3. Метод седиментационного равновесия В работе [206] был предложен новый метод получения уравнения состояния растворов макромолекул, основанный на уравнении седиментационного (гидростатического) равновесия (рис.7). В частности, было показано, что уравнение состояния разбавленных и полуразбавленных растворов можно гораздо эффективнее рассчитать, поместив систему в гравитационное поле (в дальнейшем считается, что внешнее поле действует вдоль оси z). Из установившегося профиля плотности можно получить, используя условие гидростатического равновесия, искомую зависимость осмотического давления от плотности (или, соответственно, от объемной доли), проведя лишь один запуск программы. Авторами изучалась система из 1600 цепей длиной N = на простой кубической решетке. Во внешнем поле потенциальная энергия каждого мономерного звена mgz, где z – аппликата мономерного звена, m – его масса, g – ускорение «свободного падения». В результате действия гравитационного поля на раствор макромолекул устанавливаются равновесный профиль численной плотности мономерных звеньев (z ) и соответствующий ему профиль численной плотности центров масс cm (z ).

Для данной системы седиментационная длина равна m = k BT / mg, или cm = kBT / Mg = m / N, где M = Nm есть масса полимерной цепи. Используется также безразмерный параметр =g=a/m, который характеризует силу гравитационного поля (здесь a есть шаг решетки или какая-либо другая единица длины). Описываемый метод состоит в том, что для достаточно медленно меняющегося профиля плотности центров масс цепей cm (z ) (т.е. для достаточно слабого гравитационного поля) можно записать макроскопическое уравнение гидростатического равновесия:

dP ( z ) = Mg cm ( z ), (24) dz где P(z) –давление на высоте z. Интегрируя уравнение (24), получим для осмотического давления ( z)dz. (25) ( z) = cm cm z Формула (25) допускает наглядную интерпретацию: в достигшей равновесия системе осмотическое давление ( z ) компенсирует вес раствора выше уровня z. Таким образом, известны зависимость осмотического давления ( z ) и численной плотности мономерных звеньев m (z ) от координаты z. Исключая ее, получаем зависимость давления от плотности ( ), т.е. искомое уравнение состояния.

z 80 * 80 * 180, = 0. Nmonomers a 0 20 40 60 80 100 120 140 160 180 z П роф ил ь пл от нос т и x Рис.7. Метод седиментационного (гидростатического) равновесия.

Аналогичный метод был успешно применен для изучения систем твердых коллоидных частиц (сферических [207, 208] и плоских [209, 210]). В работе [210] было в частности показано, что во внешнем гравитационном или седиментационном поле в растворе плоских коллоидных частиц наблюдается фазовое равновесие между тремя фазами (изотропной, нематической и колоннообразной). Эти фазы наблюдались как в реальном (с использованием гидраргиллита [Al(OH)3], так и в компьютерном эксперименте. С помощью метода, аналогичного описанному выше, в работе было получено уравнение состояния раствора плоских коллоидных частиц. При этом переходы между фазами отчетливо видны на графике полученной зависимости давления от плотности. Для растворов жестких стержней в отсутствии внешнего гравитационного или седиментационного поля известно, что в них могут сосуществовать различные фазы:


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

Метод гидростатического равновесия был использован для получения уравнения состояния раствора жесткоцепных полимерных молекул в рамках решеточной модели с флуктуирующей длиной связи в главе 4 данной диссертации (и в нашей работе [205]), где проводилось его сравнение с методом термодинамического интегрирования и методом отталкивающей стенки. При изучении растворов гибких макромолекул мы убедились, что рассматриваемый метод приводит к тем же результатам, что и другие методы, и смогли сформулировать критерии его применимости. Сформулированные критерии были использованы для подбора параметров при моделировании растворов жесткоцепных макромолекул. При этом было обнаружено, что уравнение состояния, рассчитанное с помощью метода гидростатического равновесия, согласуется с результатами, полученными с помощью метода термодинамического интегрирования, лишь для значений объемной доли полимера ниже точки перехода изотроп-нематик. Для больших плотностей метод не дает правильных результатов, во всяком случае при тех размерах системы и величинах внешнего поля, которые использовались в главе 4 (и в работе [205]).

Метод гидростатического равновесия может быть эффективным при исследовании поведения систем жесткоцепных макромолекул вблизи плоской поверхности, а также для изучения фазового равновесия и межфазной границы. Проблема вычисления анизотропного давления в ЖК фазе в решеточных алгоритмах МК остается актуальной. В главе 4 диссертации (и в работе [205]) показано, что как минимум для случая слабого перехода первого рода три разных метода, изначально верные только для изотропных систем, дают правильный результат для давления в анизотропной фазе.

2.5. Моделирование в большом каноническом ансамбле с использованием алгоритма с конформационным смещением выборки В главе 4 моделирование производится в большом каноническом ансамбле. Для реализации этого ансамбля используются элементарные шаги встраивания/удаления цепей с использованием считающегося сегодня наиболее эффективным метода с конформационным смещением выборки (Configurational Bias Monte Carlo, далее CBMC) [28-34]. Алгоритм CBMC был изначально предложен для решеточной модели в каноническом ансамбле [28], а потом расширен для континуальных моделей [29-32]. В принципе, алгоритм CBMC основан на алгоритме Розенблютов [36]. Элементарный шаг алгоритма CBMC в каноническом ансамбле состоит в перемещении целой цепи в растворе и реализуется с помощью процедуры удаления одной случайно выбранной цепи и вставки новой цепи в новом месте (наиболее подробно этот алгоритм изложен в оригинальной работе [28] и на стр.331-340 книги [35]). Для большого канонического ансамбля (см.

рис.8), который хорошо подходит для моделирования многих фазовых переходов [20], 290 Глава 13. Схемы Монте-Карло со смещением выборки алгоритм CBMC был предложен в работах [33, 34].

Рис.8. Алгоритм CBMC в большом каноническом ансамбле (схематическое изображение шагов добавления и удаления цепи) [35].

Рис. 13.3. Эскиз схемы Монте-Карло с конфигурационным смещением выборки. Левы рисунок показывает создание новой конфигурации, а рисунок справа показывает обратны Можно использовать CBMC-алгоритм и для перестроения не всей цепи, а только ее ход по старой конформации. Стрелки указывают три пробных положения части. При изучении адсорбции одиночной привитой жесткоцепной макромолекулы [37, 38] для изменения конформации цепи мы использовали локальные шаги и шаги МК с конкретной конформации с помощью этого методаалгоритма CBMC производились больц конфигурационным смещением выборки (CBMC). Шаги не пропорциональна ее мановскому весу. Казалось бы, это смещение выборки полимерных конформаци следующим образом: на цепи мы выбирали случайным образом мономерное звено, было исправлено в схемеприкрепленную к поверхности часть цепи, начиная от выбранного удаляли свободную, не Розенблютов с помощью введения конформационн зависимых весовых звена, вычисляя W. этой удаленной было цепи wold, и Батулисом и Кр мономерного множителей вес Однако, как части показано достраивали мером [300], эта поправка, вычисляя веса новой построенной части цепи wnew. Таким же хотя и является правильной по сути, на практик недостающую часть цепи, работает образом можно относительно коротких цепей (см. пример 13). только только для менять и конформацию отдельной свободной цепи в растворе, Решение этой проблемы заключается в смещении выборки Розенблютов таки тогда надо еще случайным образом выбирать, какую именно часть цепи от выбранного на образом, чтобы в последовательности Монте-Карло восстанавливалось бы пр ней мономерного звена надо перестраивать (рис.9).

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

13.2.2 Решеточные модели Рис.9. Алгоритм CBMC для перестроения некоторой части одиночной макромолекулы.

Весовые функции старой и новой конформаций вычисляются следующим образом.

Вначале выбирается случайным образом мономерное звено i* вдоль по цепи и удаляется часть цепи, начиная с конца цепи и заканчивая выбранным мономерным звеном;

при этом подсчитывается число kiold свободных (не занятых мономерными звеньями) узлов решетки, вокруг (i-1)- го мономерного звена, в которых i-тое мономерное звено могло бы находиться, если бы мы пытались пошагово встроить эту цепь в раствор. Затем вместо удаленной части строится новая конформация для этой части цепи и подсчитывается число свободных узлов kinew, куда можно поместить i-тое мономерное звено. Весовые функции старой и новой конформаций при удалении и построении части цепи вычисляются по формулам i + N new, wold = kiold (26) wnew = k i i=N i =i + Наконец, такой элементарный пробный шаг CBMC-алгоритма принимается с вероятностью (если при этом еще используется алгоритм ВЛ):

g ( Eold ) wnew (27) p ( old new) = min 1, g ( Enew ) wold где Eold и Enew энергии старой и новой конформаций соответственно.

2.6. Выводы по 2-ой главе Если бы строители строили здания так же, как программисты пишут программы, первый залетевший дятел разрушил бы цивилизацию.

(Второй закон Вейнберга, цит. по книге "Законы Мерфи") В главе 2 рассмотрены все методы, которые были использованы в диссертации.

Некоторые из этих методов были впервые адаптированы для систем жесткоцепных полимеров, а некоторые – впервые разработаны. Впервые разработан метод расширенного ансамбля в четырехмерном пространстве (раздел 2.2.4), в рамках которого алгоритм Ванга-Ландау (раздел 2.3) впервые адаптирован для равновероятного выбора значений внешнего параметра. Впервые выполнено систематическое сравнение различных методов расчета давления в ходе компьютерного моделирования полимерных систем методом Монте-Карло (раздел 2.4).

Детали новых и усовершенствованных методов, изложенные в этой главе, опубликованы в работах [6, 12, 16, 22, 26, 27] из списка основных публикаций по теме данной диссертации.

Общее примечание к цитированию литературы в диссертации: ссылки в разделах 2.6, 3.3, 4.5, 5.5, где сформулированы выводы по каждой главе, приводятся на список основных публикаций автора по теме данной диссертации, а во всех остальных местах текста ссылки приводятся на общий список литературы (который включает, в том числе, и статьи автора диссертации).

ГЛАВА 3. КОМПЬЮТЕРНОЕ МОДЕЛИРОВАНИЕ ОДИНОЧНОЙ ЖЕСТКОЦЕПНОЙ МАКРОМОЛЕКУЛЫ В этом разделе обсуждаются исследования одиночной макромолекулы в предельно разбавленном растворе. Рассматриваются диаграммы состояний (структур) цепи в объеме (свободной цепи) и вблизи адсорбирующей поверхности (в условиях пространственных ограничений). Для цепей конечной длины говорят именно о «диаграммах состояний» или «структур», а не о «фазовых диаграммах», так как само понятие фазы применимо только к случаю термодинамического (ТД) предела бесконечно длинной цепи. При этом диаграммы состояний для конечной цепи оказываются богаче, так как не все структуры «выживают» в ТД пределе и становятся фазами.

Физические величины, которые обычно рассчитываются в ходе компьютерного эксперимента, включают в себя средний квадрат расстояния между концами цепи, lb2, средний квадрат радиуса инерции цепи, средний квадрат длины векторов связи тензор инерции отдельной макромолекулы и параметры асимметрии (именно они, а не тензор, усредняются по всем цепям системы), концентрация полимера в системе (в случае моделирования в большом каноническом ансамбле), полная энергия E и отдельные вклады в нее: энергия изгиба E, энергия растяжения связей El, энергия объемного взаимодействия Ethermal, E = E+El+Ethermal. Весьма полезно накапливать зависимости всех этих величин от времени, как для получения информации о степени уравновешенности системы и оценки времен релаксации, так и для расчета гистограмм всех этих физических величин.

Поскольку основное внимание уделялось ориентационному упорядочению в системе, вычислялся ориентационный тензор:

( N 1) 1 2 (3e ) (28) Q = e, ii (N 1) i= где !! – -ая компонента единичного вектора вдоль направления вектора связи i-го мономерного звена (вектор, соединяющий i-е и (i+1)-е мономерные звенья). Три собственных значения этого тензора являются параметрами ориентационного упорядочения, причем максимальное собственное значение соответствует параметру ЖК порядка S, стандартно определяемому как среднее значение второго полинома Лежандра от углов между векторами связи и направлением директора.

Чтобы охарактеризовать форму возникающих структур, удобно использовать собственные значения тензора инерции полимерной цепи S,.

, N X Si S i ;

Si ri RCM ;

, x, y, z, (29) (13) N i r length near the coil-globule transition... i- i-того мономерного звена вдоль по цепи, а R есть rii обозначает положение где CM RCM, -.

радиус-вектор центра масс.

es of finite chain length near the coil-globule transition...

3.1.

depending on the angle c between subsequent bond vec r stiff macro- Одиночная жесткоцепная макромолекула в предельно разбавленном растворе 3.1.

Noguchi and tors, al structures for stiff macro- depending on the angle c between subsequent bond vec A very inter simulations by Noguchi –,.

and U tors, a c 2, evious work. the very inter- системы методика cos c stood, is 3.1.1. Модель и bcos c компьютерного эксперимента [29] A T Ua c id-crystalline [39-43].настоящей работе для bcos c cos c0, yet fully understood, is the проведения компьютерного эксперимента В использовалась T, T molecular liquid-crystalline T is the temperature, b, is the stiffness parameter other 8,models Here ules.[5, 30] For решеточная модель цепи is сthe temperature, b is длинами векторов связей между флуктуирующими the stiffness parameter.

other models Here T, alsoand differ- and(that we consider as a variable independent of temperature цепи были d complicated I differ- (that we consider 2.1 variable independent of temperature мономерными звеньями (см. раздел as a и рис.3). Для моделирования движения, in our локальныеsimulation), c0,мономерных most preferable angle simulation), and and is0 thethe звеньев, билокальные движения most preferable angle [31–36], in our eported. перемещения c is использованы ethe collapsethe collapse between subsequent bond vectors,which is taken as as between subsequent bond vectors, which is taken have studied.

рептации цепи по алгоритму «ползущей змеи» [18, 273], а также шаги алгоритма с c0 =II 8 here. To model the quality the solvent, we intro olecule with variable chain 8 here.0To model the quality ofof thesolvent, we intro ariable chain c0 = 0, ain lengths from N = 20 to duce an interaction that перестроения части, цепи (см. раздел 2.5).

конформационным смещением выборки для acts only between monomeric..

m N = 20 to duceofan units that are not nearest-neighborsbetween chain (non interaction that acts only along the monomeric idence for the existence Потенциал взаимодействия мономерных звеньев состоит из двух вкладов –.

existence of lengths.that are not nearest-neighbors along the chain (non units intermediate chain bonded attractive interaction), and is taken in a discrete потенциала объемного взаимодействия между мономерными звеньями [29] потенциала и hain lengths. bonded “quasi-Lennard-Jones” form, as is taken in work, ere somewhat preliminary attractive interaction), and in our previous a discrete s strong fluctuations “quasi-Lennard-Jones” растворителя our previous work,[29] между изгибной жесткости. Качество описывается притяжением preliminary in their so that the treatment of the non-bonded interactions is form, as in d the identification of parti- звеньями, in spirit to off-lattice work, e. g.[11–13] от друга, и для различных tions in their so that the treatment of the расстоянии r другinteractions is мономерными similar находящимися на non-bonded исследованных в диссертации цепей этот потенциал задавался либо в виде потенциальной [11–13] tion of parti take this problem upsimilar in ULJ r to off-lattice work, e. g.

spirit again, gth N = 80 but ямы, похожей по виду на потенциал Леннард-Джонса:

performing a T ( p p p emparticular, we record LJ r In up again, U his- 3 b2r 2 3r 2 1;

r 2;

5;

6;

8;

performing a the chain various orientational order 0;

other r T (30) different states of ( e record his- show up виде прямоугольной 22 1;

rямы2;

p;

p;

pmono oices of parameters либоb2r 2 is 3r distance between the considered 8;

контактов между where r the потенциальной (энергия 6 5 парных в ational these distributions. meric units, and b = 1/T. This potential vanishes for the some of order to select subensembles of звеньями): r 0;

3 and has zero derivatives at both r = 2 and other r мономерными sof the chain distance = ns which belong to “pure eters showaverage where r ris 3. The standardbetween the considered mono = “random hopping” algorithm is up proper- the distance, to extract implemented by attempting a displacement of a radomly sdistributions.

separately. In this manner, units, and b = 1/T.unit by potential vanishes for the meric chosen monomeric This one lattice site in one (ran ensembles of distance r = 3 and has the possible lattice both r = 2 and (31) cation of the possible states domly chosen) ofzero6 derivatives at directions. One ng can “pure rДля3. The standard “randomмодель вводится algorithmзависящий от моделирования жесткости цепи в потенциал, hains to exist near the coil- Monte Carlo step corresponds to N such attempted mono = hopping” is угла между последовательными связями (рис.3), only if the new state mer moves (which are accepted который был либо линейным по erage proper- describe fine the model and implemented by with excluded a displacementlength restric attempting volume and bond of a radomly complies косинусу этого угла:

this manner, chosen tions and passes the usual Metropolis site in one (ran ile Section III presents the monomeric unit by one lattice acceptance criter ossible states parameters, chosen)However, the possible lattice properties of the ntational order domly ion[37]). of the 6 equilibration of the directions. One near theone path varies the Carlo step corresponds becomes very slow whenmono diagram: coil Monte chain with this algorithm to N such attempted the (32) ess, while the other path var- globular state of the polymer is realized, because the emperature constant. mer moves (which are accepted only ifmoves new state Section acceptance probability of the attempted the becomes and describe Section mble analysis, while compliesvery small. This fact makesand computationally restric with excluded volume it bond length very presents the tions and passes the generate Metropolis acceptance criter s. demanding to usual enough statistically independent r parameters, ion[37]). However, the equilibration of theregion where col-the chain configurations in the parameter properties of bond length l bond can take the values 2, A5, A6, 3, A10 and only interested in stat there are 108 different bond vectors and 87 different angles make use of the repta between bonds. To model the chain stiffness we introduced made to remove one m into the model a potential depending on the angles q be- and to add it to the othe либо квадратичным tween adjacent bond vectors ~per monomer unit!

moves we performed m U a~ q ! 5b•~cos q2cos q 0 !, ~1! the chain. Typically w T (33) our simulation depend where T is the(30)-(33) энергетические параметрыparameter ~b В формулах temperature, b is the stiffness и измеряются type of motion makes o в единицах kBT, 5e a /T where e a is the energetical parameter of the angle cal properties of the co b= /kBT, =/kBT, где kB - постоянная Больцмана, T - температура. Потенциал на длину potential!, and q 0 is the most preferable angle between ad- The simulations w связи не учитывается. Параметр взаимодействия 0°. To model the мономерных звеньев между собой from a self We started мы jacent bond vectors which is taken equal to положили равным единице =1, таким образом мы зафиксировали единицы, в которых вычисляется энергия в нашей системе. Угол соответствует наиболее предпочтительному значению угла между векторами связей соседних вдоль по цепи мономерных звеньев (в нашем случае 0= 180). В данной модели параметр b является входным параметром, который поддерживался постоянным. Это означает, что энергия жесткости зависит от температуры, а сегмент Куна остается постоянным при изменении температуры (энергетический параметр жесткости b прямо пропорционален традиционному параметру жесткости p=lK/d). Далее в работе для конкретных систем указывается, какие именно потенциалы использовались.

Движение цепи осуществляется в соответствии с методом Монте-Карло путем случайного изменения конформации за счет локального движения мономерных звеньев и рептации по так называемому алгоритму «ползущей змеи» [18, 273]. Локальные перемещения выполняются путем смещения случайно выбранного мономерного звена на одну ячейку решетки в одном из 6 направлений вдоль осей решетки, которое также выбирается случайным образом. Рептация цепи моделируется так: случайно выбирается один из двух концов цепи, удаляется одно мономерное звено и присоединяется к другому концу цепи с использованием случайного вектора из набора разрешенных векторов связей. Кроме этого использовался также алгоритм с конформационным смещением выборки для перестроения части цепи (см. выше гл.2).



Pages:     | 1 || 3 | 4 |   ...   | 9 |
 





 
© 2013 www.libed.ru - «Бесплатная библиотека научно-практических конференций»

Материалы этого сайта размещены для ознакомления, все права принадлежат их авторам.
Если Вы не согласны с тем, что Ваш материал размещён на этом сайте, пожалуйста, напишите нам, мы в течении 1-2 рабочих дней удалим его.