Алгоритм Ванга и Ландау

Алгоритм Ванга и Ландау , предложенный Фугао Вангом и Дэвидом П. Ландау , ^[1] — это метод Монте-Карло, предназначенный для оценки плотности состояний системы. Метод выполняет немарковское случайное блуждание для построения плотности состояний, быстро посещая весь доступный энергетический спектр. Алгоритм Ванга и Ландау — важный метод получения плотности состояний, необходимой для выполнения мультиканонического моделирования .

Алгоритм Ванга – Ландау можно применить к любой системе, которая характеризуется функцией стоимости (или энергии). Например, он был применен для решения числовых интегралов ^[2] и сворачивание белков. ^[3]^[4] Выборка Ванга – Ландау связана с алгоритмом метадинамики . ^[5]

Обзор

Алгоритм Ванга и Ландау используется для получения оценки плотности состояний системы, характеризуемой функцией стоимости. Он использует немарковский случайный процесс , который асимптотически сходится к мультиканоническому ансамблю . ^[1] (То есть к алгоритму Метрополиса – Гастингса с распределением выборки, обратным плотности состояний). Основным следствием является то, что такое распределение выборки приводит к моделированию, в котором энергетические барьеры невидимы. Это означает, что алгоритм посещает все доступные состояния (благоприятные и менее благоприятные) гораздо быстрее, чем алгоритм Метрополиса. ^[6]

Алгоритм

Рассмотрим систему, определенную в фазовом пространстве $\Omega$ и функция стоимости E (например, энергия), ограниченная спектром $E\in \Gamma =[E_{\min },E_{\max }]$ , которому соответствует плотность состояний $\rho (E)$ , который подлежит оценке. Оценщик ${\hat {\rho }}(E)\equiv \exp(S(E))$ . Поскольку алгоритм Ванга и Ландау работает в дискретных спектрах, ^[1] спектр $\Gamma$ делится на N дискретных значений с разницей между ними в $\Delta$ , такой, что

N={\frac {E_{\max }-E_{\min }}{\Delta }}

.

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

обнуление всех записей микроканонической энтропии, $S(E_{i})=0\ \ i=1,2,...,N$
инициализация $f=1$ и
случайная инициализация системы путем ввода случайной конфигурации ${\boldsymbol {r}}\in \Omega$ .

Затем алгоритм выполняет мультиканоническое ансамблевое моделирование: ^[1] случайное блуждание Метрополиса – Гастингса в фазовом пространстве системы с распределением вероятностей, заданным выражением $P({\boldsymbol {r}})=1/{\hat {\rho }}(E({\boldsymbol {r}}))=\exp(-S(E({\boldsymbol {r}})))$ и вероятность предложения нового состояния, заданная распределением вероятностей $g({\boldsymbol {r}}\rightarrow {\boldsymbol {r}}')$ . Гистограмма $H(E)$ посещенных энергий сохраняется. Как и в алгоритме Метрополиса–Гастингса, выполняется этап принятия предложения и заключается в (см. обзор алгоритма Метрополиса–Гастингса ):

предлагая государство ${\boldsymbol {r}}'\in \Omega$ по произвольному распределению предложений $g({\boldsymbol {r}}\rightarrow {\boldsymbol {r}}')$
принять/отклонить предложенное состояние в соответствии с

A({\boldsymbol {r}}\rightarrow {\boldsymbol {r}}')=\min \left(1,e^{S-S'}{\frac {g({\boldsymbol {r}}'\rightarrow {\boldsymbol {r}})}{g({\boldsymbol {r}}\rightarrow {\boldsymbol {r}}')}}\right)

где

S=S(E({\boldsymbol {r}}))

и

S'=S(E({\boldsymbol {r}}'))

.

После каждого шага принятия предложения система переходит к некоторому значению $E_{i}$ , $H(E_{i})$ увеличивается на единицу и выполняется следующее обновление:

S(E_{i})\leftarrow S(E_{i})+f

.

Это решающий шаг алгоритма, и именно он делает алгоритм Ванга и Ландау немарковским: теперь случайный процесс зависит от истории процесса. Значит в следующий раз будет предложение государству именно с такой энергией $E_{i}$ , это предложение теперь, скорее всего, будет отклонено; в этом смысле алгоритм заставляет систему одинаково посещать весь спектр. ^[1] Следствием этого является то, что гистограмма $H(E)$ становится все более плоским. Однако эта неравномерность зависит от того, насколько хорошо вычисленная энтропия приближена к точной энтропии, которая, естественно, зависит от значения f. ^[7] Чтобы все лучше и лучше аппроксимировать точную энтропию (и, следовательно, плоскостность гистограммы), f уменьшается после M шагов принятия предложения:

f\leftarrow f/2

.

Позже было показано, что обновление f путем постоянного деления на два может привести к ошибкам насыщения. ^[7] Небольшая модификация метода Ванга и Ландау, позволяющая избежать этой проблемы, заключается в использовании коэффициента f, пропорционального $1/t$ , где $t$ пропорционален количеству шагов моделирования. ^[7]

Тестовая система

Мы хотим получить DOS для потенциала гармонического осциллятора .

E(x)=x^{2},\,

Аналитический DOS определяется выражением:

\rho (E)=\int \delta (E(x)-E)\,dx=\int \delta (x^{2}-E)\,dx,

выполнив последний интеграл, получим

\rho (E)\propto {\begin{cases}E^{-1/2},{\text{for }}x\in \mathbb {R} ^{1}\\{\text{const}},{\text{for }}x\in \mathbb {R} ^{2}\\E^{1/2},{\text{for }}x\in \mathbb {R} ^{3}\\\end{cases}}

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

Следовательно, мы можем использовать простой потенциал гармонического осциллятора для проверки точности алгоритма Ванга – Ландау, поскольку мы уже знаем аналитическую форму плотности состояний. Поэтому сравним оцененную плотность состояний ${\hat {\rho }}(E)$ полученный алгоритмом Ванга–Ландау с $\rho (E)$ .

Пример кода

Ниже приведен пример кода алгоритма Ванга-Ландау на Python , где мы предполагаем, что используется симметричное распределение предложений g:

g({\boldsymbol {x}}'\rightarrow {\boldsymbol {x}})=g({\boldsymbol {x}}\rightarrow {\boldsymbol {x}}')

Код рассматривает «систему», которая является базовой изучаемой системой.

currentEnergy = system.randomConfiguration()  # A random initial configuration

while f > epsilon:
    system.proposeConfiguration()  # A proposed configuration is proposed
    proposedEnergy = system.proposedEnergy()  # The energy of the proposed configuration computed

    if random() < exp(entropy[currentEnergy] - entropy[proposedEnergy]):
        # If accepted, update the energy and the system:
        currentEnergy = proposedEnergy
        system.acceptProposedConfiguration()
    else:
        # If rejected
        system.rejectProposedConfiguration()

    H[currentEnergy] += 1
    entropy[currentEnergy] += f

    if isFlat(H):  # isFlat tests whether the histogram is flat (e.g. 95% flatness)
        H[:] = 0
        f *= 0.5  # Refine the f parameter

Молекулярная динамика Ванга и Ландау: статистическая температурная молекулярная динамика (STMD)

Молекулярная динамика (МД) обычно предпочтительнее Монте-Карло (МК), поэтому желательно иметь алгоритм МД, включающий основную идею WL для выборки с плоской энергией. Этот алгоритм — статистическая температурная молекулярная динамика (STMD), разработанный ^[8] Джегил Ким и др. из Бостонского университета.

Важный первый шаг был сделан с помощью алгоритма статистической температуры Монте-Карло (STMC). WLMC требует значительного увеличения количества энергетических ячеек в зависимости от размера системы, что вызвано непосредственной работой с плотностью состояний. STMC сосредоточен на интенсивной величине, статистической температуре, $T(E)=1/(dS(E)/dE)$ , где E — потенциальная энергия. В сочетании с отношением $\Omega (E)=e^{S(E)}$ , где мы установили $k_{B}=1$ , правило WL для обновления плотности состояний дает правило для обновления дискретной статистической температуры,

{\tilde {T}}_{j\pm 1}^{\prime }=\alpha _{j\pm 1}{\tilde {T}}_{j\pm 1},

где где $\alpha _{j\pm 1}=1/(1\mp \delta f{\tilde {T}}_{j\pm 1}),\delta f=(lnf/2\Delta E),\Delta E$ - размер энергетического бункера, и ${\tilde {T}}$ обозначает текущую оценку. Мы определяем f следующим образом: ^[1] коэффициент >1, который умножает оценку DOS для i-го бункера энергии, когда система посещает энергию в этом бункере.

Подробности приведены в Ref. ^[8] С первоначальным предположением о $T(E)$ и диапазон ограничен диапазоном между $T_{L}$ и $T_{U}$ , моделирование происходит так же, как в WLMC, со значительными численными различиями. Интерполяция ${\tilde {T}}(E)$ дает континуальное выражение оцененного $S(E)$ после интегрирования обратного, что позволяет использовать более крупные энергетические элементы, чем в WL. Различные значения $S(E)$ доступны в одном и том же энергетическом интервале при оценке вероятности принятия. Когда колебания гистограммы составляют менее 20% от среднего значения, $f$ уменьшается в соответствии с $f\rightarrow {\sqrt {f}}$ .

STMC сравнивали с WL для модели Изинга и жидкости Леннарда-Джонса. При увеличении размера энергетического бункера STMC получает те же результаты в значительном диапазоне, в то время как производительность WL быстро ухудшается. STMD может использовать меньшие начальные значения $f_{d}=f-1$ для более быстрого сходимости. В целом, STMC требуется меньше шагов для получения результатов того же качества.

Теперь рассмотрим основной результат — STMD. Он основан на наблюдении, что при стандартном МД-моделировании при температуре $T_{0}$ с силами, производными от потенциальной энергии $E([x])$ , где $[x]$ обозначает все позиции, вес выборки для конфигурации равен $e^{-E([x])/T_{0}}$ . Более того, если силы получены из функции $W(E)$ , вес выборки $e^{-W(E([x]))/T_{0}}$ .

Для выборки плоской энергии пусть эффективный потенциал будет $T_{0}S(E)$ - энтропийная молекулярная динамика. Тогда вес $e^{-S(E)}$ . Поскольку плотность состояний $e^{+S(E)}$ , их продукт дает плоскую выборку энергии.

Силы рассчитываются как

F=(-d/dx)T_{0}S(E)=T_{0}(dS/dE)(-d/dx)E([x])=(T_{0}/T(E))F^{0}

где $F^{0}$ обозначает обычную силу, производную от потенциальной энергии. Масштабирование обычных сил с помощью коэффициента $(T_{0}/T(E))$ производит плоскую выборку энергии.

STMD начинается с обычного алгоритма MD при постоянной $T_{0}$ и V. Силы масштабируются, как указано, а статистическая температура обновляется на каждом временном шаге с использованием той же процедуры, что и в STMC. Поскольку моделирование сходится к выборке с фиксированной энергией, текущая оценка ${\tilde {T}}(E)$ сходится к истине $T(E)$ . Технические подробности, включая шаги по ускорению конвергенции, описаны в ^[8] и. ^[9]

В СТМД $T_{0}$ называется кинетической температурой, поскольку она, как обычно, контролирует скорости, но не входит в конфигурационную выборку, что необычно. Таким образом, STMD может исследовать низкие энергии с помощью быстрых частиц. Любое каноническое среднее значение можно рассчитать с помощью повторного взвешивания, но статистическую температуру, $T(E)$ , доступен сразу же без дополнительного анализа. Это чрезвычайно ценно для изучения фазовых переходов. В конечных наносистемах $T(E)$ имеет особенность, соответствующую каждому «субфазовому переходу». Для достаточно сильного перехода равновеликая конструкция на S-петле в $1/T(E)$ дает температуру перехода.

STMD был доработан группой BU, ^[9] и применялся ими и другими к нескольким системам. Д. Стелтер признал, что, несмотря на наш упор на работу с интенсивными объемами, $ln(f)$ является обширным. Однако $\delta f=(ln(f)/2\Delta E)$ является интенсивным, и процедура $f\rightarrow {\sqrt {f}}$ на основе гистограммы плоскостность заменяется обрезкой $\delta f$ пополам через каждое фиксированное количество временных шагов. Это простое изменение делает STMD полностью интенсивным и существенно повышает производительность больших систем. ^[9] Кроме того, конечная стоимость интенсивного $\delta f$ — константа, определяющая величину ошибки в сходившейся $T(E)$ и не зависит от размера системы. STMD реализован в LAMMPS как исправление stmd.

STMD особенно полезен для фазовых переходов. Информацию о равновесии невозможно получить с помощью канонического моделирования, поскольку для возникновения перехода необходимо переохлаждение или перегрев. Однако прогон STMD позволяет получить плоскую выборку энергии с естественным развитием нагрева и охлаждения, не попадая в состояние низкой или высокой энергии. Совсем недавно его стали применять для перехода жидкость/гель. ^[9] в наночастицах, завернутых в липиды.

Обмен репликами STMD ^[10] также был представлен группой BU.

Ссылки

^ Jump up to: ^а ^б ^с ^д ^и ^ж Ван, Фугао и Ландау, ДП (март 2001 г.). «Эффективный многодиапазонный алгоритм случайного блуждания для расчета плотности состояний». Физ. Преподобный Летт . 86 (10): 2050–2053. arXiv : cond-mat/0011174 . Бибкод : 2001PhRvL..86.2050W . doi : 10.1103/PhysRevLett.86.2050 . ПМИД 11289852 . S2CID 2941153 .
^ Р.Э. Белардинелли, С. Манзи и В.Д. Перейра (декабрь 2008 г.). «Анализ сходимости алгоритмов 1/t и Ванга – Ландау при вычислении многомерных интегралов». Физ. Преподобный Е. 78 (6): 067701. arXiv : 0806.0268 . Бибкод : 2008PhRvE..78f7701B . дои : 10.1103/PhysRevE.78.067701 . ПМИД 19256982 . S2CID 8645288 .
^ П. Охеда, М. Гарсия, А. Лондоно и Н.Ю. Чен (февраль 2009 г.). «Моделирование белков в клетках методом Монте-Карло: влияние заключения на стабильность промежуточных состояний» . Биофиз. Дж . 96 (3): 1076–1082. arXiv : 0711.0916 . Бибкод : 2009BpJ....96.1076O . дои : 10.1529/biophysj.107.125369 . ПМК 2716574 . ПМИД 18849410 .
^ П. Охеда и М. Гарсия (июль 2010 г.). «Управляемое электрическим полем разрушение нативной конформации белка бета-листа и образование альфа-спиральной структуры» . Биофиз. Дж . 99 (2): 595–599. Бибкод : 2010BpJ....99..595O . дои : 10.1016/j.bpj.2010.04.040 . ПМК 2905109 . ПМИД 20643079 .
^ Кристоф Юнгханс, Дэнни Перес и Томас Фогель. «Молекулярная динамика в мультиканоническом ансамбле: эквивалентность выборки Ванга – Ландау, статистическая температурная молекулярная динамика и метадинамика». Журнал химической теории и вычислений 10.5 (2014): 1843–1847. дои : 10.1021/ct500077d
^ Берг, Б.; Нейхаус, Т. (1992). «Мультиканонический ансамбль: новый подход к моделированию фазовых переходов первого рода». Письма о физических отзывах . 68 (1): 9–12. arXiv : hep-lat/9202004 . Бибкод : 1992PhRvL..68....9B . дои : 10.1103/PhysRevLett.68.9 . ПМИД 10045099 . S2CID 19478641 .
^ Jump up to: ^а ^б ^с Белардинелли, Р.Э. и Перейра, В.Д. (2007). «Алгоритм Ванга – Ландау: теоретический анализ насыщения ошибки». Журнал химической физики . 127 (18): 184105. arXiv : cond-mat/0702414 . Бибкод : 2007JChPh.127r4105B . дои : 10.1063/1.2803061 . ПМИД 18020628 . S2CID 25162388 .
^ Jump up to: ^а ^б ^с Ким, Джегил; Штрауб, Джон и Киз, Том (август 2006 г.). «Статистически-температурный метод Монте-Карло и алгоритмы молекулярной динамики». Физ. Преподобный Летт . 97 (5): 50601–50604. doi : 10.1103/PhysRevLett.97.050601 .
^ Jump up to: ^а ^б ^с ^д Стелтер, Дэвид и Киз, Том (2019). «Моделирование фазового равновесия жидкость/гель в липидных везикулах». Мягкая материя . 15 : 8102–8112. дои : 10.1039/c9sm00854c .
^ Ким, Джегил; Штрауб, Джон и Киз, Том (апрель 2012 г.). «Алгоритм статистико-температурной молекулярной динамики репликационного обмена» . Журнал физической химии Б. 116 : 8646–8653. дои : 10.1021/jp300366j . ПМК 11240102 .

[WangLandau-1] Jump up to: ^а ^б ^с ^д ^и ^ж Ван, Фугао и Ландау, ДП (март 2001 г.). «Эффективный многодиапазонный алгоритм случайного блуждания для расчета плотности состояний». Физ. Преподобный Летт . 86 (10): 2050–2053. arXiv : cond-mat/0011174 . Бибкод : 2001PhRvL..86.2050W . doi : 10.1103/PhysRevLett.86.2050 . ПМИД 11289852 . S2CID 2941153 .

[Belardinelli_Integrals-2] Р.Э. Белардинелли, С. Манзи и В.Д. Перейра (декабрь 2008 г.). «Анализ сходимости алгоритмов 1/t и Ванга – Ландау при вычислении многомерных интегралов». Физ. Преподобный Е. 78 (6): 067701. arXiv : 0806.0268 . Бибкод : 2008PhRvE..78f7701B . дои : 10.1103/PhysRevE.78.067701 . ПМИД 19256982 . S2CID 8645288 .

[Ojeda1-3] П. Охеда, М. Гарсия, А. Лондоно и Н.Ю. Чен (февраль 2009 г.). «Моделирование белков в клетках методом Монте-Карло: влияние заключения на стабильность промежуточных состояний» . Биофиз. Дж . 96 (3): 1076–1082. arXiv : 0711.0916 . Бибкод : 2009BpJ....96.1076O . дои : 10.1529/biophysj.107.125369 . ПМК 2716574 . ПМИД 18849410 .

[Ojeda2-4] П. Охеда и М. Гарсия (июль 2010 г.). «Управляемое электрическим полем разрушение нативной конформации белка бета-листа и образование альфа-спиральной структуры» . Биофиз. Дж . 99 (2): 595–599. Бибкод : 2010BpJ....99..595O . дои : 10.1016/j.bpj.2010.04.040 . ПМК 2905109 . ПМИД 20643079 .

[5] Кристоф Юнгханс, Дэнни Перес и Томас Фогель. «Молекулярная динамика в мультиканоническом ансамбле: эквивалентность выборки Ванга – Ландау, статистическая температурная молекулярная динамика и метадинамика». Журнал химической теории и вычислений 10.5 (2014): 1843–1847. дои : 10.1021/ct500077d

[Berg-6] Берг, Б.; Нейхаус, Т. (1992). «Мультиканонический ансамбль: новый подход к моделированию фазовых переходов первого рода». Письма о физических отзывах . 68 (1): 9–12. arXiv : hep-lat/9202004 . Бибкод : 1992PhRvL..68....9B . дои : 10.1103/PhysRevLett.68.9 . ПМИД 10045099 . S2CID 19478641 .

[Belardinelli_Saturation-7] Jump up to: ^а ^б ^с Белардинелли, Р.Э. и Перейра, В.Д. (2007). «Алгоритм Ванга – Ландау: теоретический анализ насыщения ошибки». Журнал химической физики . 127 (18): 184105. arXiv : cond-mat/0702414 . Бибкод : 2007JChPh.127r4105B . дои : 10.1063/1.2803061 . ПМИД 18020628 . S2CID 25162388 .

[STMD-8] Jump up to: ^а ^б ^с Ким, Джегил; Штрауб, Джон и Киз, Том (август 2006 г.). «Статистически-температурный метод Монте-Карло и алгоритмы молекулярной динамики». Физ. Преподобный Летт . 97 (5): 50601–50604. doi : 10.1103/PhysRevLett.97.050601 .

[fluidgel-9] Jump up to: ^а ^б ^с ^д Стелтер, Дэвид и Киз, Том (2019). «Моделирование фазового равновесия жидкость/гель в липидных везикулах». Мягкая материя . 15 : 8102–8112. дои : 10.1039/c9sm00854c .

[RESTMD-10] Ким, Джегил; Штрауб, Джон и Киз, Том (апрель 2012 г.). «Алгоритм статистико-температурной молекулярной динамики репликационного обмена» . Журнал физической химии Б. 116 : 8646–8653. дои : 10.1021/jp300366j . ПМК 11240102 .

[1]

[2]

[3]

[4]

[5]

[6]

[7]

[8]

[9]

[10]