Физики научились считать распределение плотности в шесть раз быстрее

Двое немецких физиков придумали новый способ, с помощью которого пространственное распределение плотности системы можно рассчитывать в шесть раз быстрее, чем обычно. Для этого ученые предложили рассчитывать сначала среднюю силу, действующую на частицу, а потом восстанавливать из нее распределение плотности. Статья опубликована в Physical Review Letters, препринт работы выложен на сайте arXiv.org.

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

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

Физики Даниэл де ля Герас (Daniel de las Heras) и Маттиас Шмидт (Matthias Schmidt) предлагают альтернативный способ для численного расчета плотности распределения частиц. Этот способ полагается на тот факт, что в равновесии средняя сила, действующая на частицу, и распределение плотности в ее окрестности связаны простым дифференциальным уравнением (F(r)−kT∇ρ(r) = 0) — следовательно, первым делом можно численно рассчитать силу, а потом найти плотность с помощью интегрирования по пространству. При этом в расчетах будет возникать константа интегрирования, которую можно зафиксировать, выбирая нормировку — вычисляя полное число частиц с помощью найденного распределения и требуя, чтобы оно совпало с заданным числом. Определить силу можно тем же способом, что и плотность, то есть усредняя ее по всем частицам системы.

Затем ученые проверили, как предложенный способ будет работать на простейшей системе из N = 25 молекул, взаимодействие которых описывается потенциалом Леннард-Джонса. Частицы помещались в кубическую коробку с периодическими граничными условиями и со стороной, равной десяти характерным длинам потенциала. Кроме того, на частицы накладывалась слабая внешняя сила, которая удерживала их в равновесии. Наконец, для расчета средней силы или плотности ученые использовали три разных метода — метод Монте-Карло, броуновскую динамику и молекулярную динамику, — а потом сравнивали их работу. В результате оказалось, что при одинаковом числе шагов и одинаковом шаге сетки, на которую разбивается пространство параметров, новый метод дает в десять раз меньшую ошибку, чем расчет плотности «в лоб». Это эквивалентно шестикратному увеличению скорости вычислений, обеспечивающих заданную погрешность. Расчеты с бо́льшим числом частиц (N = 1000) только подтвердили полученные результаты.

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

Тем не менее, у предложенного учеными метода есть и недостатки — из-за того, что константу интегрирования в нем нужно восстанавливать заново из нормировки, он дает несимметричные распределения плотности для симметричных систем. Например, при наложении на систему параболического внешнего потенциала плотность при положительных значениях координаты оказывалась немного занижена по сравнению с точным значением. Конечно, при увеличении числа шагов эта разница становится практически незаметной, однако ее нужно иметь в виду при аккуратных расчетах термодинамических параметров, поскольку она может привести к появлению нежелательных артефактов. Кроме того, новый способ работает только для равновесных распределений, поскольку в неравновесном случае основное уравнение модифицируется (F(r)−kT∇ρ(r) = J(r)), и распределение тока нужно искать отдельно.

В настоящее время численное моделирование применяется во многих областях науки — гораздо дешевле проверить поведение потенциально интересной системы на компьютере, чем изготавливать ее на практике. Поэтому важно разрабатывать новые и оптимизировать старые методы вычислений, чтобы уменьшить время, затрачиваемое на расчеты. Например, мы писали о том, как российские математики разработали подобные методы для решения уравнений Навье-Стокса или Кортевега — де Фриза, описывающих движение вязкой жидкости. Однако чаще всего для оптимизации схемы вычислений используют машинное обучение — например, в течение последнего года его применяли для расчета потенциала атомов бора, определения химических свойств произвольных молекул, нахождения топологического инварианта одномерной системы и даже для вычисления функциональных интегралов.

Дмитрий Трунин