МОЛЕКУЛЯРНАЯ ДИНАМИКА

Практическая молекулярная динамика. Часть 1

Серия статей, в которых я постараюсь научить людей работать с методом классическоймолекулярной динамики и программным пакетом GROMACS в частности.

Они основаны намоих приватных лекциях студентам, моем собственном опыте и обсуждениях с коллегами.

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

Часть 1 содержит введение и некоторые основы метода молекулярной динамики.

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

Определения метода молекулярной динамики

Метод молекулярной динамики — метод, в котором временная эволюция системывзаимодействующих атомов или частиц отслеживается интегрированием их уравненийдвижения (Wikipedia). Напрактике это означает примерно следующее

• Для описания движения используются положения классической ньютоновской механики
• Силы взаимодействия между частицами задаются в форме каких либо классическихпотенциальных функций.

Давайте рассмотрим эти положения подробнее.

Надеюсь, вы знаете что такое ньютоновская механика, это что то, строящееся навтором законе Ньютона:

$$mddot{x} = F$$

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

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

$$frac{d}{dt} frac{partial L}{partial dot{x_i}} — frac{partial L}{partial x_i} = Qn_i$$

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

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

$$
U = frac{k}{2}(r-r_0)2 $$

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

Таким образом, энергия взаимодействия одной частицы со всей системой описываетсяследующим уравнением:

$$U_{total} = U_{12} + U_{13} + U_{14} + U_{vdw} + U_{coul}$$

где (U_{12}) и (U_{13}) — гармоники, описывающие взаимодействие междусоседними частицами (оно же растяжение связей) и через одну частицу (оно же деформацияплоских углов) соответственно, (U_{14}) — деформация торсионных углов,(U_{vdw}) и (U_{coul}) — энергии ван-дер-ваальсового и кулоновскоговзаимодействия соответственно.

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

$$U_i = f(r_i, {r_j}N_0)$$

Ну а полная энергия системы будет суммой потенциальных и кинетических энергийвсех частиц в системе, что-то вроде:

$$E = sum_{i}{U_i} + sum_{i}{frac{m_i {dot x}_i2}{2}}$$

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

Интегрирование уравнений движения

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

Методы численного интегрирования

Для численного (в отличие от аналитического с первообразными функциями, знакомогопо школьной программе) интегрирования используется геометрическое определениеинтегрирования — определенный интеграл от заданной функции есть площадь фигуры подней.

Для того, чтобы рассчитать площадь фигуры, можно разбить ее на маленькие части,например, прямоугольные трапеции с заданной шириной (Delta x). Площадь такойтрапеции можно рассчитать по формуле

$$S_i = Delta x frac{f(x_i) + f(x_i + Delta x)}{2}$$

тогда площадь всей фигуры, если ее разбить на (N) равных по ширине трапеций,будет равна

$$int f(x) dx = sum_{i=1}{N} {S_i} = sum_{i=0}{N-1} { Delta x frac{f(x + i Delta x) + f(x + (i + 1) Delta x)}{2} }$$

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

Алгоритм Верле

Допустим, у нас задано время (t) и малый шаг интегрирования (Delta t).Тогда в моменты времени (t + Delta t) и (t — Delta t) координаты можнозадать как нечто, что раскладывается в ряд Тейлора в окрестности (t)(пренебрегая членами порядка 3 и выше, поскольку мы предполагаем, что шагинтегрирования бесконечно мал):

$$x(t + Delta t) = x(t) + dot x (t) Delta t + frac{ddot x(t) {Delta t}2}{2}$$

$$x(t — Delta t) = x(t) — dot x (t) Delta t + frac{ddot x(t) {Delta t}2}{2}$$

Ничто нам не мешает сложить два уравнения друг с другом и получить следующее:

$$x(t + Delta t) + x(t — Delta t) = 2x(t) + ddot x(t) {Delta t}2$$

или

$$x(t + Delta t) = 2x(t) — x(t — Delta t) + ddot x(t) {Delta t}2$$

Если перевести на человеческий язык. то, зная координаты частицы в момент времени(t) и немного (на (Delta t)) раньше, а также зная ускорение, мы можемрассчитать координаты частицы в будущем на (Delta t).

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

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

Допустим, мы проводим моделирование при заданной температуре(T); воспользовавшись аппаратом статистической физики, мы можем описатьраспределение частиц по скоростям (например,Максвелловское распределение).

Таким образом, можно случайным образом задать частице такуюскорость, чтобы температура всей системы была равна заданной. Далее, рассчитатьускорение частицы в нулевой момент времени и, в дальнейшем, собственно координатычастиц в момент времени (t — Delta t).

inb4: Следует понимать, что указанный алгоритм имеет свои недостатки.

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

Силовые поля

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

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

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

Термостат, баростат и все-все-все

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

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

• канонический ансамбль, в котором фиксируется число частиц, объем и температура(он же (NVT) ансамбль);
• изотермо-изобарический ансамбль, в котором фиксируется число частиц, температураи давление (он же (NpT) ансамбль).

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

Тогда задача сводится к — так или иначе — заданию определенного давления итемпературы. Для поддержания определенной температуры в системе используется алгоритм,который описывает термостат (обмен температуры системы с окружающей средой).

Длядавления, соответственно, баростат — «обмен» давления с окружающей средой.

Молекулярно-динамические траектории

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

Молекулярная динамика с жесткими внутримолекулярными связями

Page 3

Чаще всего в молекулярной динамике со связями используется процедура, основанная на алгоритме RATTLE, естественно встраиваемом в скоростную форму алгоритма Верлета интегрирования уравнений движения.

Запишем схему для интегрирования уравнений движения:

с условиями связи:

Скоростная форма алгоритма Верлета для уравнений движения записывается следующим образом:

С учетом связей этот алгоритм назывется алгоритмом RATTLE. Опишем алгоритм RATTLE. Обозначим силы связей do

через gi —

Шаг 1.

ar,

Шаг 2.

Шаг 3. Выберем Л так, что ст(г”1″1) = |г112 — d!у = 0.

Шаг 4.

Шаг 5. Выберем неопределенный множитель Лагранжа р так, чтобы ' Pijn+1) = 0.

затем решаем нелинейное уравнение для А:

и выполняем обратную подстановку

Шаг 5 может быть выполнен определением

Далее решаем уравнение для второго множителя Лагранжа р

do,

которое представляет линейное соотношение, так как —(г, р) =

dt

= —g(r) • р/т. Выполняем обратную подстановку

Рис. 12.2. Иллюстрация определений моментов сил одной молекулы относительно другой

Алгоритмом SHAKE вычисляют силы связей Agj так, что положения удовлетворяют условиям связи. Для многих связей необходимы итерации.

Page 4

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

Если внутримолекулярные силы выражены как суммы атом- атомных членов, то преобразование их в силы в системе координат центра масс и моментов вращения относительно центра масс позволяют просто сформулировать уравнения вращения. Рассмотрим две молекулы А и В, изображенные на рис. 12.2.

Пусть векторы и R/j — положения молекул А и В в системе координат центра масс. Определим внутримолекулярный вектор Rb = R/j — R и предположим, что потенциал взаимодействия может быть представлен в виде

где i и j — номера атомов соответствующих молекул. Теперь вычислим силу, действующую на молекулу А со стороны В:

и вращающий момент, действующий на А со стороны В,

и вращающий момент, действующий на В со стороны А,

где гiA = Ti — R/i — положение атома i относительно центра молекулы А и аналогично для г3в-

Заметим, что Ид# ф — Nдд. как это предполагается в учебниках, но, как следует из уравнений, прямо приводит к закону сохранения момента:

откуда следует, что силы удовлетворяют соотношению fij = = —fji и действуют вдоль атом-атомного вектора Тц. Это выражение сохранения локального углового момента. Оно прямо следует из инвариантности потенциальной энергии Еав во вращающейся системе координат.

За оси xi,yi,zi, принимают главные оси инерции молекулы. Уравнения движения тела (молекулы), отнесенные к осям xi,yi,zi, называются динамическими уравнениями Эйлера:

или в координатном представлении:

где Nx,Ny,Nz — главные моменты сил относительно осей х 1, уь ZЭти уравнения интегрируются совместно с кинетическими уравнениями Эйлера. Проекции угловой скорости си на оси xi,yi,zi, связанные с телом, выражаются через углы Эйлера формулами, называемыми кинематическими уравнениями Эйлера:

Рис. 12.3. Углы Эйлера. Две системы координат с одним центром: x,y,z— неподвижная система координат; X,yi,zi — система координат, связанная с телом где точками обозначены производные по времени. Скорость любой точки равна моменту вектора угловой скорости относительно этой точки: v1 = [си1,г].

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

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

Метод классической молекулярной динамики

Метод молекулярной динамики (метод МД) — метод, в котором временная эволюция системы взаимодействующих атомов или частиц отслеживается интегрированием их уравнений движения

• 1 Основные положения
• 2 Ограничения применимости метода
• 3 Временные и пространственные параметры исследуемых систем
• 4 История развития метода
• 5 Применение
• 6 Литература
• 7 Примечания
• 8 Ссылки

Основные положения

• Для описания движения атомов или частиц применяется классическая механика. Закон движения частиц находят при помощи аналитической механики.

• Силы межатомного взаимодействия можно представить в форме классических потенциальных сил (как градиент потенциальной энергии системы).

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

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

Ограничения применимости метода

Метод молекулярной динамики применим, если длина волны Де Бройля атома (или частицы) много меньше, чем межатомное расстояние.
Также классическая молекулярная динамика не применима для моделирования систем, состоящих из легких атомов, таких как гелий или водород.

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

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

Временные и пространственные параметры исследуемых систем

Метод классической (полноатомной) молекулярной динамики позволяет с использованием современных ЭВМ рассматривать системы, состоящие из нескольких миллионов атомов на временах порядка нескольких пикосекунд.

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

Для решения таких задач все чаще требуются большие вычислительные мощности, которыми обладают суперкомпьютеры.

История развития метода

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

Однако, как хорошо известно, даже уже для трех частиц появляются препятствия, затрудняющие аналитическое решение. Примером может служить простая химическая реакция H + H2 = H2 + H. Для такой реакции Hirschfelder, Eyring, Topley в 1936 году провели попытку расчета нескольких шагов вдоль одной из траекторий.

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

Следует отметить несколько ключевых работ, определивших развитие метода молекулярной динамики. Первая работа, посвященная моделированию методом молекулярной динамики, вышла в 1957 году. Её авторами были Alder и Waingwright.

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

Вычисления проводились на компьютерах UNIVAC и на IBM 704.

Статья Dynamics of radiation damage , J.B. Gibson, A. N. Goland, M.Milgram, G.H. Vineyard выполненная в Брукхейвенской национальной лаборатории и появившаяся в 1960 году была, возможно, первым примером моделирования с непрерывным потенциалом.

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

Aneesur Rahman из Аргоннской национальной лаборатории в своей статье 1964 года Correlation in the motion of atoms in liquid argon изучил свойства жидкого аргона, используя потенциал Леннард-Джонса. Система состояла из 864 атомов. Результаты были получены на компьютере CDC 3600 (англ.).

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

Loup Verlet вычислил в 1967 году фазовую диаграмму аргона, используя потенциал Леннард-Джонса и смоделировал корреляционные функции, чтобы проверить теорию жидкого состояния. В своей работе он разработал процедуру сохранения вычислительных ресурсов, ныне известную как Verlet neighbor list, а также предложил новый метод численного интегрирования уравнений движения.

Применение

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

также кристаллография, ЯМР).

Взаимодействие между объектами может быть описано силовым полем (классическая молекулярная динамика), квантовохимической моделью или смешанной теорией, содержащей элементы двух предыдущих (QM/MM (quantum mechanics/molecular mechanics, QMMM (англ.)).

Наиболее популярными пакетами программного обеспечения для моделирования динамики биологических молекул являются: AMBER, CHARMM (и коммерческая версия CHARMm), GROMACS, GROMOS,Lammps, HOOMD-blue и NAMD.

Литература

• Рапапорт Д. К. Искусство молекулярной динамики. — Ижевск: ИКИ, 2012. — 632 с. — ISBN 978-5-4344-0083-1.
• J. A. McCammon, S. C. Harvey (1987) Dynamics of Proteins and Nucleic Acids. Cambridge University Press. ISBN 0-521-35652-0 (paperback); ISBN 0-521-30750-3 (hardback).
• Daan Frenkel, Berend Smit (2001) Understanding Molecular Simulation. Academic Press. ISBN 0-12-267351-4.
• Oren M. Becker, Alexander D. Mackerell Jr, Benoît Roux, Masakatsu Watanabe (2001) Computational Biochemistry and Biophysics. Marcel Dekker. ISBN 0-8247-0455-X.
• Tamar Schlick (2002) Molecular Modeling and Simulation. Springer. ISBN 0-387-95404-X.

Примечания

1. ↑ 1. J.M. Haile, Molecular dynamics simulation, Wiley, 1992.
2. ↑ M. P. Allen, D. J. D. C. Rapaport The Art of Molecular Dynamics Simulation, 1996.
3. ↑ Tildesley Computer simulation of liquids. Oxford University Press, 1989.
4. ↑ G.C Schatz, A Kopperman // J. Chem. Phys., v.

   62 , p.2502, (1975)

5. ↑ B.J. Alder, T.E. Waingwright// J. Chem. Phys. v. 27, p.1208, (1957)
6. ↑ J.B. Gibson, A. N. Goland, M.Milgram, G.H. Vineyard // Phys Rev, v.120, p.1229, (1960)
7. ↑ A Rahman // Phys. Rev. v.136A, p.405, (1964)
8. ↑ L. Verlet // Phys Rev, v.159, p.

   98, (1967)

Ссылки

• The Blue Gene Project (IBM)
• BioSimGrid, a database for biomolecular simulations
• http://www.moldyn.org/library/md/
• https://moldyngrid.org/ (виртуальная лаборатория MolDynGrid)
• Рудяк В.Я. и другие. Моделирование течений в наноканалах методом молекулярной динамики. НАНОСИСТЕМЫ: Ф, Х, М, — 2011, 2 (4), С. 100–112, УДК 532.137.2; 532.5.013.12

Метод классической молекулярной динамики Информацию О

Метод классической молекулярной динамики

Метод классической молекулярной динамики
Метод классической молекулярной динамики Вы просматриваете субъект
Метод классической молекулярной динамики что, Метод классической молекулярной динамики кто, Метод классической молекулярной динамики описание

There are excerpts from wikipedia on this article and video

Наш сайт имеет систему в функции поисковой системы. Выше: «что вы искали?»вы можете запросить все в системе с коробкой. Добро пожаловать в нашу простую, стильную и быструю поисковую систему, которую мы подготовили, чтобы предоставить вам самую точную и актуальную информацию.

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

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

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

Глава 1

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

Сейчас МД применяется в квантовой химии и физике твердого тела для [95]:

·        изучения дефектов в кристаллах, варьирующихся от точечных (вакансии, дефекты внедрения) до линейных (дислокации) и плоских (межфазные, междоменные границы и т. д.). Для расчета подобных структур нужно применять метод суперячеек, требующих использования большого количества атомов в системе [96];

·        реконструкции поверхности кристалла, связанной с перестройкой координат множества атомов на поверхности. Ее можно описать с помощью МД, причем можно проследить изменение поверхности в зависимости от температуры [97];

·        изучения кластеров, величина которых варьируется от нескольких атомов до нескольких тысяч, С помощью МД в основном проводится процедура численного отжига при оптимизации геометрии [98, 99];

·        изучения биологических молекул (белки, ДНК, РНК и др.), обладающих низкой симметрией и содержащих огромное количество атомов [100].

1.3.3.1       Классическая молекулярная динамика

В классической молекулярной динамике используются законы механики Ньютона. Сила, действующая на атом I, записывается, соответственно, в виде

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

Рассмотрим системы N частиц, двигающихся под влиянием потенциальной функции U [101]. Движения частиц можно описать с помощью их координат   и импульсов . Обозначим и   как и   соответственно. Потенциальная функция U зависит только от координат частиц: .

Гамильтониан такой системы имеет следующий вид:

Таким образом, уравнения движения такой системы можно записать в следующем виде:

Из этих уравнений легко получить второй закон Ньютона:

Уравнения движения также можно получить, используя формализм Лагранжа:

Последнее выражение применяется для получения ab initio МД-уравнений.

Интегрирования уравнений движения, алгоритм Верле

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

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

Все они основаны на методе конечных разностей, где время изменяется дискретно с некоторым шагом Δt.

Среди наиболее известных методов интегрирования уравнений движения можно выделить алгоритм Верле [102], алгоритм «прыжков лягушки» (leap frog) [103], метод предиктора-корректора [104] (рис. 1-14). В данной книге будет описан первый из этих алгоритмов.

Алгоритм Верле в литературе по численному анализу часто называется неявной симметричной разностной схемой.

Основная идея алгоритма Верле состоит в записи разложения положения частицы и :

Складывая (1.200) и (1.201), можно получить

где    ускорение частицы ( ), которое, однако, можно вывести также из выражения

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

Кроме вышеописанной версии алгоритма Верле, разработан скоростной алгоритм Верле, где положения, скорости и ускорения на шаге t + Δt вычисляются следующим образом:

Рис. 1-14. Схематическое описание различных форм алгоритма Верле: а   основной алгоритм Верле (1.202); b   алгоритм «прыжков лягушки»; c   скоростная форма алгоритма Верле (1.204) [105].

Метод неэмпирической молекулярной динамики, впервые предложенный в работе Кара Паринелло (Car Parinello) [106], имеет широкое применение в современной физике и химии.

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

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

При этом поведение системы строится исходя из лагранжиана системы, где вводятся фиктивное время и фиктивная масса электронов μ.

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

где    функционал полной энергии, получаемый в любом квантово-химическом методе, набор {αν}   любые возможные внешние параметры   температура, давление, объем и пр.

, μ   фиктивная масса или параметр инерции, связанный с орбитальными степенями свободы, μν   произвольный параметр соответствующей размерности. Матрица   является набором множителей Лагранжа, которые обеспечивают ортонормированность .

Из этих уравнений можно получить уравнения движения Эйлера Лагранжа [107]:

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

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

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

Основные положения

• Для описания движения атомов или частиц применяется классическая механика. Закон движения частиц находят при помощи аналитической механики.

• Силы межатомного взаимодействия можно представить в форме классических потенциальных сил (как градиент потенциальной энергии системы).

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

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

Ограничения применимости метода

Метод молекулярной динамики применим, если длина волны Де Бройля атома (или частицы) много меньше чем межатомное расстояние.
Также классическая молекулярная динамика не применима для моделирования систем состоящих из легких атомов, таких как гелий или водород.

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

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

Временные и пространнственные параметры исследуемых систем

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

История развития метода

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

Однако, как хорошо известно, даже уже для трех частиц появляются трудности, затрудняющие аналитическое решение. Примером может служить простая химическая реакция H + H2 = H2 + H. Для такой реакции Hirschfelder, Eyring, Topley в 1936 году провели попытку расчета нескольких шагов вдоль одной из траекторий.

Это было за 30 лет до того, как возможности такого расчета стали возможны на компьютере. Позднее классический подход был подкреплен полуклассическими и квантовохимическими расчетами в тех областях, где влияние квантовых эффектов становилось значимым [4].

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

Следует отметить несколько ключевых работ, определивших развитие метода молекулярной динамики. Первая работа, посвященная моделированию методом молекулярной динамики, вышла в 1957 году. Ее авторами были Alder и Waingwright [5].

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

Вычисления проводились на компьютерах UNIVAC и на IBM 704.

Статья Dynamics of radiation damage , J.B. Gibson, A. N. Goland, M.Milgram, G.H. Vinhyard [6] выполненная в Брукхейвенской национальной лаборатории и появившаяся в 1960 году была возможно первым примером моделирования с непрерывным потенциалом.

Эта одна из самых лучших работ по данной тематике!
Aneesur Rahman из Аргонской национальной лаборатории в своей статье 1964 года Correlation in the motion of atoms in luquid argon [7] изучил свойства жидкого аргона, используя потенциал Леннарда-Джонса. Система состояла из 864 атомов.

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

Loup Verlet вычислил в 1967 [8] году фазовую диаграмму аргона, используя потенциал Леннарда-Джонса и смоделировал корреляционные функции, чтобы проверить теорию жидкого состояния. В своей работе он разработал процедуру сохранения вычислительных ресурсов, ныне известную, как Verlet neighboor list, а также предложил новый метод численного интегрирования уравнений движения.

Применение

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

Он играет важную роль в определении структуры белка и уточнении его свойств (см. также кристаллография,ЯМР).

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

Наиболее популярными пакетами программного обеспечения для моделирования динамики биологических молекул являются: AMBER, CHARMM (и коммерческая версия CHARMm), GROMACS, GROMOS, и NAMD.

Литература

• J. A. McCammon, S. C. Harvey (1987) Dynamics of Proteins and Nucleic Acids. Cambridge University Press. ISBN 0—52—135652—0 (paperback); ISBN 0—52—130750 (hardback).
• D. C. Rapaport (1996) The Art of Molecular Dynamics Simulation. ISBN 0521445612.

• Daan Frenkel, Berend Smit (2001) Understanding Molecular Simulation. Academic Press. ISBN 0122673514.
• Oren M. Becker, Alexander D. Mackerell Jr, Benoît Roux, Masakatsu Watanabe (2001) Computational Biochemistry and Biophysics. Marcel Dekker. ISBN 082470455X.

• Tamar Schlick (2002) Molecular Modeling and Simulation. Springer. ISBN 038795404X.

Источники

1. ↑ 1. J.M. Haile, Molecular dynamics simulation, Wiley, 1992.
2. ↑ M. P. Allen, D. J. D. C. Rapaport The Art of Molecular Dynamics Simulation, 1996.
3. ↑ Tildesley Computer simulation of liquids. Oxford University Press, 1989.
4. ↑ G.C Schatz, A Kopperman // J. Chem. Phys., v.

   62 , p.2502, (1975)

5. ↑ B.J. Alder, T.E. Waingwright// J. Chem. Phys. v. 27, p.1208, (1957)
6. ↑ J.B. Gibson, A. N. Goland, M.Milgram, G.H. Vinhyard // Phys Rev, v.120, p.1229, (1960)
7. ↑ A Rahman // Phys. Rev. v.136A, p.405, (1964)
8. ↑ L. Verlet // Phys Rev, v.159, p.

   98, (1967)

Внешние ссылки